Page - 126 - in Programming for Computations – Python - A Gentle Introduction to Numerical Simulations with Python

126 4 SolvingOrdinaryDifferentialEquations 4.3.2 NumericalSolution Wehave not looked at numericalmethods for handling second-order derivatives, and suchmethods are an option, butwe knowhow to solve first-order differential equationsandevensystemsoffirst-orderequations.Witha little,yetverycommon, trickwecan rewrite (4.42)as afirst-order systemof twodifferential equations.We introduce u D x and v D x0 D u0 as two new unknown functions. The two correspondingequationsarise from thedefinitionv Du0 and theoriginal equation (4.42): u0 Dv; (4.43) v0 D !2u: (4.44) (Notice thatwecanuseu00 D v0 to remove the second-orderderivative fromNew- ton’s2nd law.) Wecannowapply theForwardEulermethod to(4.43)–(4.44),exactlyaswedid inSect. 4.2.2: unC1 un t Dvn; (4.45) vnC1 vn t D !2un; (4.46) resulting in thecomputational scheme unC1 DunC tvn; (4.47) vnC1 Dvn t!2un: (4.48) 4.3.3 ProgrammingtheNumericalMethod;theSpecialCase Asimpleprogramfor (4.47)–(4.48)follows the same ideasas inSect. 4.2.3: from numpy import zeros, linspace, pi, cos, array import matplotlib.pyplot as plt omega = 2 P = 2*pi/omega dt = P/20 T = 3*P N_t = int(round(T/dt)) t = linspace(0, N_t*dt, N_t+1) u = zeros(N_t+1) v = zeros(N_t+1) # Initial condition X_0 = 2 u[0] = X_0 v[0] = 0