Notes Solving Nonlinear Systems � Most thoroughly explored in the context of optimization � For systems arising in implicit time integration of stiff problems: • Must be more efficient than taking k substeps of an explicit method ruling out e.g. fixed point iteration • But if we have difficulties converging (a solution might not even exist!) we can always reduce time step and try again � Thus Newton � s method is usually chosen cs542g-term1-2006 1 cs542g-term1-2006 2 Newton’s method Variations � Start with initial guess y 0 at solution � Just taking a single step of Newton (e.g. current y) corresponds to “freezing” the coefficients of F(y)=0 in time: sometimes called “semi-implicit” • Just a linear solve, but same stability � Loop until converged: according to linear analysis • Linearize around current guess: • However, usually nonlinear effects cause F(y k + � y) � F(y k ) + dF/dy � y worse problems than for fully implicit methods • Solve linear equations: � In between: keep Jacobian dF/dy constant dF/dy � y = -F(y k ) • Line search along direction � y with initial step but iterate as in Newton size of 1 � And endless variations on inexact Newton cs542g-term1-2006 3 cs542g-term1-2006 4 Second order systems Hamiltonian Systems � One of the most important time differential � For a Hamiltonian function H(p,q), the equations is F=ma, 2nd order in time system: dt = � H dp � Reduction to first order often throws out � q useful structure of the problem dt = � � H dq � In particular, F(x,v) often has special � p properties that may be useful to exploit • E.g. nonlinear in x, but linear in v: mixed � Think q=positions, implicit/explicit methods are natural p=momentum (mass times velocity), � We � ll look at Hamiltonian systems in and H=total energy (kinetic plus potential) particular for a conservative mechanical system cs542g-term1-2006 5 cs542g-term1-2006 6
Conservation The Flow � Take time derivative of Hamiltonian: � For any initial condition (p,q) and any later dt = � H dt + � H time t, can solve to get p(t), q(t) dH dp dq ( ) = p ( t ), q ( t ) ( ) � t p , q � Call the map � p � q dt = � � H � H � q + � H � H the “flow” of the system � p = 0 � p � q � Hamiltonian dynamics possess flows with special properties � Note H is generally like a norm of p and q, so we � re on the edge of stability: solutions neither decay nor grow • Eigenvalues are pure imaginary! cs542g-term1-2006 7 cs542g-term1-2006 8 Area in 2D Area-preserving linear maps � What is the area of a parallelogram with � Let A be a linear map in 2D: x � =Ax vector edges (u 1 ,u 2 ) and (v 1 ,v 2 )? • A is a 2x2 matrix � Then A is area-preserving if the area of area( u , v ) = u � v any parallelogram is equal to the area of = u 1 v 2 � u 2 v 1 the transformed parallelogram: � � � � � u T J � � v = u T Jv � v 1 ) 0 1 ( � = u 1 � u 2 � 1 0 � � � � v 2 u T A T JAv = u T Jv = u T Jv A T JA = J cs542g-term1-2006 9 cs542g-term1-2006 10 Symplectic Matrices Symplectic Maps � y = � ( y ) � We can generalize this to any even dimension � Consider a nonlinear map � � � Let � Assume it � s adequately smooth J = 0 I � � � I 0 � � � At a given point, infinitesimal areas are transformed by Jacobian matrix: � Then matrix A is symplectic if A T JA=J A ( y ) = � � � Note that u T Jv is just the sum of the projected � y areas � Map is symplectic if its Jacobian is everywhere a symplectic matrix • Area (or summed projected area) is preserved cs542g-term1-2006 11 cs542g-term1-2006 12
Hamiltonian Flows Hamiltonian Flows are Symplectic � Let � s look at Jacobian of a Hamiltonian flow � Theorem: for any fixed time t, the flow of a � t ( y 0 ) = y ( t ) Hamiltonian system is symplectic � � � y = p • Note at time 0, the flow map is the identity � � t � � t ( y 0 ) = dy � � dt ( t ) (which is definitely symplectic) q • Differentiate A T JA: dt = J � 1 � H ( ) dy = J � 1 � H � t ( y 0 ) � y � ) = � A JA + A T J � A � t = … = 0 = � J T ( � � t �� t ( y 0 ) = J � 1 �� H � t ( y 0 ) ( ) �� t ( y 0 ) � t A T JA � t � t � A � t = J � 1 �� H A � Important point: �� H, the Hessian, is symmetric. cs542g-term1-2006 13 cs542g-term1-2006 14 Trajectories Symplectic Methods � Volume preservation: � A symplectic method is a numerical • if you start off a set of trajectories occupying some method whose map is symplectic region, that region may get distorted but it will • Note if map from any t n to t n+1 is symplectic, maintain its volume � General ODE � s usually have sources/sinks then composition of maps is symplectic, so full • Trajectories expand away or converge towards a method is symplectic point or a manifold � Example: symplectic Euler • Obviously not area preserving • Goes by many names, e.g. velocity Verlet � General ODE methods don � t respect symplecticity: area not preserved � Also implicit midpoint • In long term, the trajectories have the wrong • Not quite trapezoidal rule, but the two are behaviour essentially equivalent… cs542g-term1-2006 15 cs542g-term1-2006 16 Modified Equations Symplectic Euler � Backwards error analysis for differential � Look at simple example: H=T(p)+V(q) equations � Symplectic Euler is essentially ( ) � Say we are solving p n + 1 = p n � � tH q p n , q n dy ( ) dt = f y ( ) q n + 1 = q n + � tH p p n + 1 , q n � Goal: show that numerical solution {y n } is actually the solution to a modified � Do some Taylor series expansions to find equation: first term in modified equations… dy ( ) , ( ) + hf 2 y ( ) + h 2 f 3 y ( ) + … y n = y t n dt = f y cs542g-term1-2006 17 cs542g-term1-2006 18
Modified Equations are Modified equations in general Hamiltonian! � The first expansion is: � Under some assumptions, any symplectic method is solving a nearby Hamiltonian system dt = � � � q H � � t � � dp exactly 2 H p � H q � + O ( � t 2 ) � � � � Aside: this modified Hamiltonian depends on step size h dt = � � p H � � t � � dq 2 H p � H q � + O ( � t 2 ) � • If you use a variable time step, modified Hamiltonian � � is changing every time step • Numerical flow is still symplectic, but no strong � So numerical solution is, to high order, guarantees on what it represents solving a Hamiltonian system • As a result - may very well see much worse long-time (but with a perturbed H) behaviour with a variable step size than with a fixed • So to high order has the same structure step size! cs542g-term1-2006 19 cs542g-term1-2006 20
Recommend
More recommend