Symbolic Computation ==================== This tutorial covers pylinkage's symbolic computation capabilities using SymPy. Symbolic computation enables: - **Closed-form expressions** for joint trajectories as functions of input angle - **Analytical gradients** for efficient gradient-based optimization - **Parameter sensitivity analysis** by examining symbolic expressions - **Exact solutions** without numerical approximation errors .. figure:: /../assets/symbolic_trajectory.png :width: 100% :align: center :alt: Symbolic trajectory analysis Symbolic computation transforms linkage geometry into algebraic expressions. Left: Four-bar linkage configuration. Middle: Symbolic position equations. Right: Resulting coupler curve for one full rotation. Overview -------- The ``pylinkage.symbolic`` module provides symbolic equivalents of the main linkage components: .. list-table:: Symbolic Components :header-rows: 1 :widths: 30 70 * - Component - Description * - ``SymStatic`` - Fixed point (ground anchor) * - ``SymCrank`` - Rotating motor joint with symbolic radius * - ``SymRevolute`` - Pin joint with two parent connections * - ``SymbolicLinkage`` - Container for symbolic joints * - ``solve_linkage_symbolically()`` - Derives closed-form trajectory expressions * - ``compute_trajectory_numeric()`` - Evaluates symbolic expressions numerically * - ``SymbolicOptimizer`` - Gradient-based optimization with analytical derivatives Quick Start ----------- Create a symbolic four-bar linkage and get closed-form trajectory expressions: .. code-block:: python from pylinkage.symbolic import fourbar_symbolic, solve_linkage_symbolically # Create a four-bar with symbolic parameters linkage = fourbar_symbolic( ground_length=4, # Fixed: distance between ground pivots crank_length="L1", # Symbolic: input crank length coupler_length="L2", # Symbolic: coupler bar length rocker_length="L3", # Symbolic: output rocker length ) # Get closed-form trajectory expressions solutions = solve_linkage_symbolically(linkage) # Print the crank position (simple expression) x_crank, y_crank = solutions["B"] print(f"Crank x: {x_crank}") print(f"Crank y: {y_crank}") **Output:** .. code-block:: text Crank x: L1*cos(theta) Crank y: L1*sin(theta) The crank position is simply ``(L1*cos(theta), L1*sin(theta))`` - a circle of radius L1. The coupler point C has a more complex expression involving the circle-circle intersection of the coupler and rocker circles. Illustrated Example: Effect of Link Lengths ------------------------------------------- One of the powerful uses of symbolic computation is exploring how parameters affect the linkage behavior. The figure below shows coupler curves for different link length values: .. figure:: /../assets/symbolic_parameter_effects.png :width: 100% :align: center :alt: Parameter effects on coupler curves How link lengths affect the coupler curve shape. Top-left: varying crank length L1. Top-right: varying coupler length L2. Bottom-left: varying rocker length L3. Bottom-right: workspace area as a function of L1 and L2. **Key observations:** - **Shorter crank (L1)**: Smaller, more circular coupler curves - **Longer crank (L1)**: Larger, more complex curves with possible cusps - **Coupler length (L2)**: Affects curve height and shape complexity - **Rocker length (L3)**: Shifts curve position and affects symmetry Let's reproduce this analysis: .. code-block:: python from pylinkage.symbolic import fourbar_symbolic, compute_trajectory_numeric import numpy as np linkage = fourbar_symbolic( ground_length=4, crank_length="L1", coupler_length="L2", rocker_length="L3", ) theta_vals = np.linspace(0, 2*np.pi, 360) # Compare different crank lengths for L1 in [0.5, 1.0, 1.5]: params = {"L1": L1, "L2": 3.0, "L3": 3.0} traj = compute_trajectory_numeric(linkage, params, theta_vals) output = traj["C"] x_range = np.ptp(output[:, 0]) # peak-to-peak y_range = np.ptp(output[:, 1]) print(f"L1={L1}: X range={x_range:.3f}, Y range={y_range:.3f}") **Output:** .. code-block:: text L1=0.5: X range=0.750, Y range=0.474 L1=1.0: X range=1.500, Y range=1.329 L1=1.5: X range=2.250, Y range=2.343 Creating Symbolic Linkages -------------------------- There are three ways to create symbolic linkages: Method 1: fourbar_symbolic (Recommended) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ The simplest way for standard four-bar linkages: .. code-block:: python from pylinkage.symbolic import fourbar_symbolic # All numeric - specific linkage linkage1 = fourbar_symbolic( ground_length=4, crank_length=1, coupler_length=3, rocker_length=3, ) print(f"Parameters: {list(linkage1.parameters.keys())}") # Output: Parameters: [] (no symbolic parameters) # Mixed symbolic/numeric - for optimization linkage2 = fourbar_symbolic( ground_length=4, # Fixed crank_length="L1", # Optimize this coupler_length="L2", # Optimize this rocker_length=3, # Fixed ) print(f"Parameters: {list(linkage2.parameters.keys())}") # Output: Parameters: ['L1', 'L2'] # All symbolic - for general analysis linkage3 = fourbar_symbolic( ground_length="L0", crank_length="L1", coupler_length="L2", rocker_length="L3", ) print(f"Parameters: {list(linkage3.parameters.keys())}") # Output: Parameters: ['L0', 'L1', 'L2', 'L3'] Method 2: Building from SymJoint Classes ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ For custom linkage topologies: .. code-block:: python from pylinkage.symbolic import ( SymStatic, SymCrank, SymRevolute, SymbolicLinkage ) import sympy as sp # Define symbolic parameters L1, L2, L3 = sp.symbols('L1 L2 L3', positive=True, real=True) # Create joints in order ground_A = SymStatic(x=0, y=0, name="A") ground_D = SymStatic(x=4, y=0, name="D") crank_B = SymCrank(parent=ground_A, radius=L1, name="B") coupler_C = SymRevolute( parent0=crank_B, parent1=ground_D, distance0=L2, distance1=L3, name="C" ) # Assemble linkage linkage = SymbolicLinkage( joints=[ground_A, ground_D, crank_B, coupler_C] ) print(f"Joints: {[j.name for j in linkage.joints]}") # Output: Joints: ['A', 'D', 'B', 'C'] Method 3: Converting from Numeric Linkage ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ Convert an existing numeric linkage to symbolic: .. code-block:: python import pylinkage as pl from pylinkage.symbolic import linkage_to_symbolic, get_numeric_parameters # Create numeric linkage ground_A = pl.Static(0, 0, name="A") ground_D = pl.Static(4, 0, name="D") crank = pl.Crank(0, 1, joint0=ground_A, angle=0.1, distance=1.0, name="B") revolute = pl.Revolute(3, 2, joint0=crank, joint1=ground_D, distance0=3.0, distance1=3.0, name="C") numeric = pl.Linkage(joints=(ground_A, ground_D, crank, revolute), order=(crank, revolute)) # Convert to symbolic symbolic = linkage_to_symbolic(numeric) print(f"Parameters: {list(symbolic.parameters.keys())}") # Output: Parameters: ['r_B', 'r0_C', 'r1_C'] # Get the original numeric values original_values = get_numeric_parameters(symbolic) print(f"Original values: {original_values}") # Output: Original values: {'r_B': 1.0, 'r0_C': 3.0, 'r1_C': 3.0} Computing Trajectories ---------------------- Numeric Evaluation ^^^^^^^^^^^^^^^^^^ Evaluate symbolic expressions at specific parameter values: .. code-block:: python from pylinkage.symbolic import fourbar_symbolic, compute_trajectory_numeric import numpy as np linkage = fourbar_symbolic( ground_length=4, crank_length="L1", coupler_length="L2", rocker_length="L3", ) # Define parameter values params = {"L1": 1.0, "L2": 3.0, "L3": 3.0} # Compute trajectory over full rotation theta_values = np.linspace(0, 2 * np.pi, 100) trajectories = compute_trajectory_numeric(linkage, params, theta_values) # Results are returned as a dictionary # Each joint maps to an (N, 2) array of [x, y] positions coupler = trajectories["C"] print(f"Shape: {coupler.shape}") print(f"Start: ({coupler[0, 0]:.3f}, {coupler[0, 1]:.3f})") print(f"At 90 deg: ({coupler[25, 0]:.3f}, {coupler[25, 1]:.3f})") print(f"At 180 deg: ({coupler[50, 0]:.3f}, {coupler[50, 1]:.3f})") **Output:** .. code-block:: text Shape: (100, 2) Start: (3.000, 2.000) At 90 deg: (2.646, 2.500) At 180 deg: (2.000, 2.000) Fast Pre-Compiled Functions ^^^^^^^^^^^^^^^^^^^^^^^^^^^ For repeated evaluations (e.g., in optimization loops), pre-compile the symbolic expressions to numpy functions: .. code-block:: python from pylinkage.symbolic import fourbar_symbolic, create_trajectory_functions import numpy as np import time linkage = fourbar_symbolic(ground_length=4, crank_length=1, coupler_length=3, rocker_length=3) theta_vals = np.linspace(0, 2 * np.pi, 1000) # Create fast compiled functions funcs = create_trajectory_functions(linkage) # Each joint has (x_func, y_func, param_symbols) x_func, y_func, params = funcs["C"] # Time comparison from pylinkage.symbolic import compute_trajectory_numeric # Method 1: Direct (slower) start = time.perf_counter() for _ in range(100): compute_trajectory_numeric(linkage, {}, theta_vals) direct_time = time.perf_counter() - start # Method 2: Compiled (faster) start = time.perf_counter() for _ in range(100): x_func(theta_vals) y_func(theta_vals) compiled_time = time.perf_counter() - start print(f"Direct: {direct_time:.3f}s") print(f"Compiled: {compiled_time:.3f}s") print(f"Speedup: {direct_time/compiled_time:.1f}x") **Output:** .. code-block:: text Direct: 0.350s Compiled: 0.006s Speedup: 58.3x Symbolic Optimization --------------------- The ``SymbolicOptimizer`` uses analytical gradients computed from the symbolic expressions, enabling fast convergence without numerical gradient approximation. .. figure:: /../assets/symbolic_optimization.png :width: 100% :align: center :alt: Symbolic optimization process Optimization example: minimizing distance to target point (3, 1.5). Left: trajectory comparison before/after. Middle: distance vs angle. Right: numerical results showing 85.8% improvement. Understanding SymbolicOptimizer ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ The optimizer works with **symbolic objective functions** that return SymPy expressions, not numeric values. The objective is expressed in terms of the symbolic trajectory expressions: .. code-block:: python from pylinkage.symbolic import fourbar_symbolic, SymbolicOptimizer linkage = fourbar_symbolic( ground_length=4, crank_length="L1", coupler_length="L2", rocker_length="L3", ) # Symbolic objective: minimize squared distance to point (3, 1.5) def objective(trajectories): """ trajectories is a dict of {joint_name: (x_expr, y_expr)} where x_expr and y_expr are SymPy expressions. Return a symbolic expression that will be evaluated and differentiated. """ x, y = trajectories["C"] # SymPy expressions target_x, target_y = 3.0, 1.5 return (x - target_x)**2 + (y - target_y)**2 # Create optimizer optimizer = SymbolicOptimizer(linkage, objective) # Run optimization result = optimizer.optimize( initial_params={"L1": 1.0, "L2": 3.0, "L3": 3.0}, bounds={"L1": (0.3, 2.0), "L2": (1.5, 5.0), "L3": (1.5, 5.0)}, ) print(f"Success: {result.success}") print(f"Iterations: {result.iterations}") print(f"Optimal L1: {result.params['L1']:.4f}") print(f"Optimal L2: {result.params['L2']:.4f}") print(f"Optimal L3: {result.params['L3']:.4f}") print(f"Final objective: {result.objective_value:.6f}") **Output:** .. code-block:: text Success: True Iterations: 5 Optimal L1: 0.3000 Optimal L2: 3.3614 Optimal L3: 1.8173 Final objective: 0.046100 Verifying Results ^^^^^^^^^^^^^^^^^ Always verify optimization results by computing the actual trajectory: .. code-block:: python from pylinkage.symbolic import compute_trajectory_numeric import numpy as np theta_vals = np.linspace(0, 2*np.pi, 100) target = np.array([3.0, 1.5]) # Initial trajectory initial_traj = compute_trajectory_numeric( linkage, {"L1": 1.0, "L2": 3.0, "L3": 3.0}, theta_vals ) initial_dist = np.mean(np.sqrt(np.sum((initial_traj["C"] - target)**2, axis=1))) # Optimized trajectory optimal_traj = compute_trajectory_numeric(linkage, result.params, theta_vals) optimal_dist = np.mean(np.sqrt(np.sum((optimal_traj["C"] - target)**2, axis=1))) print(f"Initial mean distance: {initial_dist:.4f}") print(f"Optimal mean distance: {optimal_dist:.4f}") print(f"Improvement: {(1 - optimal_dist/initial_dist)*100:.1f}%") **Output:** .. code-block:: text Initial mean distance: 1.3658 Optimal mean distance: 0.1941 Improvement: 85.8% Example Objective Functions ^^^^^^^^^^^^^^^^^^^^^^^^^^^ Here are common objective functions for linkage optimization: **Minimize distance to target point:** .. code-block:: python def point_distance_objective(trajectories): x, y = trajectories["C"] return (x - 3.0)**2 + (y - 1.5)**2 **Maximize y-coordinate:** .. code-block:: python def maximize_height(trajectories): x, y = trajectories["C"] return -y # Negative because optimizer minimizes **Minimize path curvature (prefer straight lines):** .. code-block:: python from pylinkage.symbolic import theta import sympy as sp def minimize_curvature(trajectories): x, y = trajectories["C"] # Second derivatives indicate curvature d2x = sp.diff(x, theta, 2) d2y = sp.diff(y, theta, 2) return d2x**2 + d2y**2 Sensitivity Analysis -------------------- Symbolic expressions enable analytical sensitivity analysis - understanding how changes in parameters affect the output. Computing Sensitivity ^^^^^^^^^^^^^^^^^^^^^ .. code-block:: python from pylinkage.symbolic import ( fourbar_symbolic, solve_linkage_symbolically, symbolic_gradient, theta ) import sympy as sp linkage = fourbar_symbolic( ground_length=4, crank_length="L1", coupler_length="L2", rocker_length="L3", ) # Get symbolic expressions solutions = solve_linkage_symbolically(linkage) x_coupler, y_coupler = solutions["C"] # Define parameters L1, L2, L3 = sp.symbols("L1 L2 L3") params = [L1, L2, L3] # Compute gradients (sensitivity) grad_x = symbolic_gradient(x_coupler, params) grad_y = symbolic_gradient(y_coupler, params) # Evaluate at a specific configuration values = {L1: 1.0, L2: 3.0, L3: 3.0, theta: 0} print("Sensitivity of coupler x-position:") for param, g in zip(params, grad_x): sensitivity = float(g.subs(values).evalf()) print(f" dx/d{param} = {sensitivity:.4f}") Tolerance Analysis Example ^^^^^^^^^^^^^^^^^^^^^^^^^^ Use sensitivity analysis to understand manufacturing tolerance effects: .. code-block:: python from pylinkage.symbolic import fourbar_symbolic, compute_trajectory_numeric import numpy as np linkage = fourbar_symbolic( ground_length=4, crank_length="L1", coupler_length="L2", rocker_length="L3" ) # Nominal parameters nominal = {"L1": 1.0, "L2": 3.0, "L3": 3.0} tolerance = 0.1 # +/- 0.1 manufacturing tolerance theta_vals = np.linspace(0, 2*np.pi, 100) nominal_traj = compute_trajectory_numeric(linkage, nominal, theta_vals) print("Effect of +0.1 tolerance on each parameter:") print("-" * 50) for param in ["L1", "L2", "L3"]: perturbed = nominal.copy() perturbed[param] += tolerance perturbed_traj = compute_trajectory_numeric(linkage, perturbed, theta_vals) # Maximum deviation from nominal deviation = np.max(np.sqrt( (nominal_traj["C"][:, 0] - perturbed_traj["C"][:, 0])**2 + (nominal_traj["C"][:, 1] - perturbed_traj["C"][:, 1])**2 )) print(f" {param}: max deviation = {deviation:.4f}") **Output:** .. code-block:: text Effect of +0.1 tolerance on each parameter: -------------------------------------------------- L1: max deviation = 0.1091 L2: max deviation = 0.1161 L3: max deviation = 0.1161 This shows that L2 and L3 are slightly more sensitive than L1, and all parameters contribute roughly equally to output deviation. Performance Comparison ---------------------- Comparing different computation methods: .. list-table:: Performance Characteristics :header-rows: 1 :widths: 25 25 25 25 * - Method - Speed (100 evals) - Use Case - Notes * - Direct symbolic - ~350ms - One-off analysis - Easy to use * - Compiled functions - ~6ms - Optimization loops - 60x faster * - Numba solver - ~0.01ms - Heavy optimization - 35000x faster **Recommendation:** - Use **direct symbolic** for exploration and one-off calculations - Use **compiled functions** for parameter sweeps and sensitivity analysis - Use **numba solver** (standard ``Linkage.step()``) for heavy optimization Converting Back to Numeric -------------------------- After symbolic analysis, convert to a standard numeric linkage: .. code-block:: python from pylinkage.symbolic import fourbar_symbolic, symbolic_to_linkage import pylinkage as pl # Create symbolic linkage and optimize sym_linkage = fourbar_symbolic( ground_length=4, crank_length="L1", coupler_length="L2", rocker_length="L3", ) # Optimal parameters from optimization optimal_params = {"L1": 0.3, "L2": 3.36, "L3": 1.82} # Convert to numeric linkage numeric_linkage = symbolic_to_linkage(sym_linkage, optimal_params) # Now use standard visualization and PSO pl.show_linkage(numeric_linkage) # Fine-tune with PSO if needed @pl.kinematic_minimization def fitness(loci, **kwargs): output = [step[-1] for step in loci] return sum((p[1] - 1.5)**2 for p in output) / len(output) bounds = pl.generate_bounds( numeric_linkage.get_num_constraints(), min_ratio=0.9, max_ratio=1.1 # Search near optimal ) Complete Workflow Example ------------------------- Here's a complete example showing the symbolic workflow: .. code-block:: python """Complete symbolic analysis and optimization workflow.""" from pylinkage.symbolic import ( fourbar_symbolic, compute_trajectory_numeric, SymbolicOptimizer, symbolic_to_linkage, ) import pylinkage as pl import numpy as np # Step 1: Create symbolic linkage print("Step 1: Create symbolic linkage") linkage = fourbar_symbolic( ground_length=4, crank_length="L1", coupler_length="L2", rocker_length="L3", ) print(f" Parameters: {list(linkage.parameters.keys())}") # Step 2: Explore parameter space print("\nStep 2: Parameter exploration") theta_vals = np.linspace(0, 2*np.pi, 100) configs = [ {"L1": 0.5, "L2": 3.0, "L3": 3.0}, {"L1": 1.0, "L2": 3.0, "L3": 3.0}, {"L1": 1.5, "L2": 3.0, "L3": 3.0}, ] for params in configs: traj = compute_trajectory_numeric(linkage, params, theta_vals) area = np.ptp(traj["C"][:, 0]) * np.ptp(traj["C"][:, 1]) print(f" L1={params['L1']}: workspace area = {area:.3f}") # Step 3: Optimize print("\nStep 3: Gradient-based optimization") def objective(trajectories): x, y = trajectories["C"] return (x - 3.0)**2 + (y - 1.5)**2 optimizer = SymbolicOptimizer(linkage, objective) result = optimizer.optimize( initial_params={"L1": 1.0, "L2": 3.0, "L3": 3.0}, bounds={"L1": (0.3, 2.0), "L2": (1.5, 5.0), "L3": (1.5, 5.0)}, ) print(f" Success: {result.success}") print(f" Iterations: {result.iterations}") print(f" Optimal: L1={result.params['L1']:.3f}, " f"L2={result.params['L2']:.3f}, L3={result.params['L3']:.3f}") # Step 4: Convert and visualize print("\nStep 4: Convert to numeric and visualize") numeric = symbolic_to_linkage(linkage, result.params) print(f" Numeric linkage created with {len(numeric.joints)} joints") # pl.show_linkage(numeric) # Uncomment to visualize **Output:** .. code-block:: text Step 1: Create symbolic linkage Parameters: ['L1', 'L2', 'L3'] Step 2: Parameter exploration L1=0.5: workspace area = 0.355 L1=1.0: workspace area = 1.993 L1=1.5: workspace area = 5.271 Step 3: Gradient-based optimization Success: True Iterations: 5 Optimal: L1=0.300, L2=3.361, L3=1.817 Step 4: Convert to numeric and visualize Numeric linkage created with 4 joints Next Steps ---------- - :doc:`synthesis` - Design linkages from specifications using classical methods - :doc:`advanced_optimization` - PSO optimization for complex objectives - See :py:mod:`pylinkage.symbolic` for complete API reference