NONLINEAR MECHANICAL SYSTEMS CANONICAL TRANSFORMATION S AND NUMERICAL INTEGRATION Jacobi Canonical Transformations A Jacobi canonical transformations yields a Hamiltonian that depends on only one of the conjugate variable sets. Assume dependence on new momentum alone. H(p*,q*) = K(p*) ∂K(p*)/∂q* = 0 Thus dp*/dt = e* dq*/dt = ∂K(p*)/∂p* – f* The simple relation between effort and the rate of change of momentum is recovered in the new coordinates. Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 1EXAMPLE: SIMPLE HARMONIC OSCILLATOR Hamiltonian H(p,q) = 21(p2/I + q2/C) Hamilton's equations dq/dt = ∂H/∂p = p/I dp/dt = –∂H/∂q = q/C Change variables from old (q,p) to new (P,Q) Define Zo = IC and the generating function S(q,Q) = Zo(q2/2) cotQ The transformation equations are p = ∂S/∂q = Zoq cotQ P = –∂S/∂Q = Zo(q2/2)/sin2Q Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 2Express the old variables in terms of the new p = 2P cosQ Zo q = 2P sinQ (1/ Zo ) Define ωo = 1/IC and the new Hamiltonian is H(P,Q) = ωo P = K(P) Hamilton’s equations in the new coordinates dQ/dt = ∂K/∂P = ωo dP/dt = –∂K/∂Q = 0 Their solution is Q(t) = ωo t + constant P(t) = constant Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 3In essence this variable change has integrated the equations. As the product of P and Q has the units of action (energy by time) it is sometimes called a (simple harmonic) actional transformation. PHYSICAL INTERPRETATION: P is proportional to the total system energy. Its square root is proportional to oscillation amplitude. Q is the phase angle of the oscillations. Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 4In general, finding Jacobi canonical transformations requires solving a non-trivial partial differential equation. A practical alternative is to separate the Hamiltonian into two parts, one with a known Jacobi canonical transform. H(p,q) = Hj(p,q) + Hn(p,q) Apply the known Jacobi canonical transformation H*(P,Q) = Hj*(P) + H*n(P,Q) Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 5We may represent the second term as a set of canonical forces e*(P,Q) = –∂H*n/∂Q f*(P,Q) = –∂H*n/∂P The transformed equations become dP/dt = e*(P,Q) dQ/dt = ∂H*j/∂P – f*(P,Q) An advantage of this change of variables is that, in effect, it integrates the fundamental oscillatory mode of the solution. Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 6EXAMPLE: SIMPLE PENDULUM For large amplitudes, the simple pendulum is a nonlinear oscillator. H(η,θ) = η2/2 + 1 – cos θ where θ angle with respect to the vertical η corresponding angular momentum Expand the cosine as a power series H(η,θ) = η2/2 + θ2/2 – θ4/4! + θ6/6! – ... The Hamiltonian is quadratic in momentum and displacement with additional terms in displacement of fourth power and higher. Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 7Until the fourth power of the angle in radians becomes significant, the nonlinear pendulum may be treated as linear system with a Hamiltonian that is quadratic in momentum and displacement. For the quadratic terms have a knownJacobi canonical transformation: the simple harmonic actional. Split the Hamiltonian as follows H(η,θ) = η2/2 + θ2/2 + (1 – cosθ – θ2/2) H(η,θ) = K(η,θ) + N(η,θ) Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 8Apply the simple harmonic actional θ = 2P sinQ η = 2P cosQ The Hamiltonian becomes H*(P,Q) = K*(P) + N*(P,Q) In the original variables, the system equations are dη/dt = –∂H/∂θ dθ/dt = ∂H/∂η In the new variables, the system equations become dP/dt = –∂N*/∂Q dQ/dt = 1 + ∂N*/∂P Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 9Transformation does not change the value of either K or N. Use the chain rule on the original N which depends only on θ. ∂N*/∂Q = (∂N/∂θ) (∂θ/∂Q) ∂N*/∂P = (∂N/∂θ) (∂θ/∂P) ∂N/∂θ = sinθ – θ The transformed equations become dP/dt = [ 2P sinQ – sin( 2P sinQ)][ 2P cosQ] dQ/dt = 1 + [sin( 2P sinQ) – 2P sinQ][sinQ (1/ 2P )] ∂θ/∂Q = 2P cosQ ∂θ/∂P = sinQ (1/ 2P ) Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 10Use the transformation equations to express the rates of change as a function of both old and new variables. dP/dt = (θ – sinθ)η dQ/dt = 1 + (sinθ – θ)θ/2P θ = 2P sinQ η = 2P cosQ What have we gained? The system equations are simpler in the old variables dη/dt = –sinθ dθ/dt = η In the new variables, the solution is far more stable numerically. Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 11Simple Euler integration algorithm starting time 0 seconds, final time 50 seconds, time step 0.1 seconds. start from rest at an angle of 0.1 radians (≈6°) In old coordinates, simulation is unstable. Total system energy grows exponentially. Lagrangian formulation 0 5 10 15 20 25 30 35 40 45 50 time in seconds Lagrangian formulation -1.5 -1 -0.5 0 0.5 1 angle in radians 0 0.2 0.4 0.6 total energy 0 5 10 15 20 25 30 35 40 45 50 time in seconds Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 12In new coodinates, the simulation is stable. Total system energy variation: 5.6x10-7. Hamilton-Jacobi formulation -0.1 -0.05 0 0.05 0.1 angle in radians 0 5 10 15 20 25 30 35 40 45 50 time in seconds x 10-3 Hamilton-Jacobi formulation 4.9954 4.9956 4.9958 4.996 4.9962 total energy 0 5 10 15 20 25 30 35 40 45 time in seconds Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 13 50Perform the same integration using a third-order fixed-step Runge-Kutta algorithm. Lagrangian formulation 0 5 10 15 20 25 30 35 40 45 50 time in seconds x 10-3 Lagrangian formulation 0 5 10 15 20 25 30 35 40 45 50 time in seconds Total energy, declines steadily by 2.1x10-5 over 50 seconds. -0.1 -0.05 0 0.05 0.1 angle in radians 4.975 4.98 4.985 4.99 4.995 5 total energy Mod. Sim. Dyn. Sys. Transforms & Numerical Integration page 14Hamilton-Jacobi formulation 0 5 10 15 20 25 30 35 40 45 time in seconds x 10-3 Hamilton-Jacobi formulation -0.1 -0.05 0 0.05 0.1 angle in radians 4.9958 4.9958 4.9958 4.9958 4.9958 total energy 0 5 10 15 20 25 30 35 40 45 time in seconds Total energy also decreases, but by 2.3x10-10 — a hundred thousand time less. Mod. Sim. Dyn. Sys. Transforms & Numerical Integration