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

248 8 SolvingOrdinaryDifferentialEquations nonlinear equations, we can approximate or predict un+1 using a Forward Euler step: un+1 =un+Δtf(un,tn). This reasoninggives rise to themethod u∗ =un+Δtf(un,tn), (8.57) un+1 =un+Δt 2 (f(un,tn)+f(u∗,tn+1)). (8.58) Theschemeapplies tobothscalar andvectorODEs. Foranoscillatingsystemwithf = (v,−ω2u) thefileosc_Heun.pyimplements thismethod.Thedemofunctioninthatfilerunsthesimulationfor10periodswith20 time steps per period. The correspondingnumerical and exact solutions are shown in Fig. 8.24. We see that the amplitude grows, but not as much as for the Forward Eulermethod.However, theEuler-Cromermethodperformsbetter! We shouldadd that in problemswhere theForwardEulermethodgives satisfac- toryapproximations,suchasgrowth/decayproblemsor theSIRmodel, thesecond- orderRunge-Kuttamethod (Heun’smethod)usuallyworksconsiderablybetter and produces greater accuracy for the same computational cost. It is therefore a very valuablemethod to be aware of, although it cannotcompetewith the Euler-Cromer scheme for oscillation problems. The derivation of the RK2/Heun scheme is also goodgeneral training in “numerical thinking”. Fig. 8.24 Simulationof 10 periods ofoscillationsby Heun’s method