Page - 267 - in Programming for Computations – Python - A Gentle Introduction to Numerical Simulations with Python 3.6, Volume Second Edition

8.4 Oscillating1DSystems:ASecondOrderODE 267 A major problem arises when we want to start the scheme. We know thatu0 = U0, butapplying(8.76) forn=0 tocomputeu1 leads to u1 =2u0−u−1 −Δt2ω2u0, (8.77) wherewedonotknowu−1.Theinitial conditionu′(0)=0canhelpus toeliminate u−1—andthisconditionmustanywaybeincorporatedinsomeway.Tothisend,we discretizeu′(0)=0bya centereddifference, u′(0)≈ u 1−u−1 2Δt =0 . It follows thatu−1 =u1, andwecanuse this relation toeliminateu−1 in (8.77): u1 =u0− 1 2 Δt2ω2u0 . (8.78) Withu0 =U0 andu1 computed from (8.78), we can computeu2,u3, and so forth from(8.76).Exercise8.25 asksyou to explorehow the steps aboveare modified in casewe havea nonzero initial conditionu′(0)=V0. Remarkona simpler methodforcomputingu1 We couldapproximate the initial conditionu′(0)bya forwarddifference: u′(0)≈ u 1−u0 Δt =0, leading to u1 = u0. Then we can use (8.76) for the coming time steps. However, this forward difference has an error proportional to Δt, while the centered difference we used has an error proportional to Δt2, which is compatible with the accuracy (error goes likeΔt2) used in the discretization of the differentialequation. The method for the second-order ODE described above goes under the name Störmer’s method or Verlet integration.7 It turns out that this method is math- ematically equivalent with the Euler-Cromer scheme (!). Or more precisely, the generalformula(8.76) isequivalentwith theEuler-Cromerformula,but thescheme for the first time level (8.78) implements the initial conditionu′(0) slightly more accurately than what is naturally done in the Euler-Cromer scheme. The latter will do v1 =v0 −Δtω2u0, u1 =u0+Δtv1 =u0−Δt2ω2u0, 7 http://en.wikipedia.org/wiki/Verlet_integration.