lec15_ODE2 -         LeapFrog and Verlet method for...

Info iconThis preview shows page 1. Sign up to view the full content.

View Full Document Right Arrow Icon
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview:         LeapFrog and Verlet method for ODEs Example: Simple Pendulum ODE of vectors Example: Projectile motion with air drag       HW#6 due today Get HW#6 off website Reading for Differential Equations in Appendix B.     Consider general ODE with initial condition: We are interested in a solution from t=0 to t=T in n time steps: du = f (u(t)) and u(0) = u0 dt T ∆t = n     Notation for solution at time step k: Using foward difference, ODE becomes Finally solution for time k+1 is uk+1 − uk = f (uk ) ∆t uk+1 = uk + ∆tf (uk ) uk ≡ u(tk )     Need value of u at k=0 time step to get this started (initial condition).   In many systems, the rate at which a quantity changes depends linearly on the value of that quantity. du (chemical reactions, nuclear fission, = −bu dt heat transfer,…) Initial condition (t=0):   u(0) = A   Solution is an exponential: u(t) = Ae −bt     Euler is based on forward difference approximation with error of order τ 2 . Big problem is it is NOT energy conserving. Using center difference method, we get the leap-frog method:       Error is of order τ 3 , and it IS energy conserving. However it is not self starting (need k=-1 step at start). One solution is to take an initial Euler step, then go with LeapFrog after the first step is computed. Why “leap-frog”? Let’s apply to a standard mechanics problem… uk+1 = uk−1 + 2τ f (uk )   By Newton’s laws and basic mechanics force (F), acceleration (a) velocity (v), position (r): F = ma   dv = a(r(t)) dt dr = v (t) dt Apply Leap-Frog method to solve for position and velocity at time step k based on previous values: vk+1 = vk−1 + 2τ a(rk ) rk+2 = rk + 2τ vk+1   Note how v and r values leap-frog each other.   Here we use the central difference formulas for BOTH the first derivative (v) AND second derivative (r). 2 dr dr = a(r) 2 dt dt = v (t) rk+1 − rk−1   So scheme becomes: v k= 2τ rk+1 = 2rk − rk−1 + τ 2 ak     Error for v is of order τ 2 and error for r is of order Calculation of v is not required to compute r. τ 4     Physics – velocity becomes angular velocity and position becomes angular position. Equations of motion are: dφ =ω dt dω g = α(φ) = − sin φ dt l   Or numerically using Verlet (just interested in position) φk+1 = 2φk − φk−1 + τ αk 2 Euler Verlet     Often in science the quantities we are interested in are vectors which require a magnitude and direction. Typically we express a vector by it’s individual components: ￿ = (vx , vy , vz ) v   Computationally, our 1D problems now become 3 1-D problems. We just perform the calculations on each of the components individually.     dvx 1 = Fdx (vx ) dt m dvy 1 = Fdy (vy ) − g dt m Restrict to 2-D, so only need x and y components of r and v. Force of gravity is only in y direction, drag force is always opposite to velocity (x and y components constantly change). drx = vx dt dry = vy dt 1 ￿ Fd = − Cd ρA|￿ |￿ vv 2 ...
View Full Document

This note was uploaded on 10/05/2010 for the course PHYS phy503 taught by Professor Gladden during the Spring '09 term at Ole Miss.

Ask a homework question - tutors are online