Numerical Solution of Ordinary Differential Equations

1 Numerical Solution of OrdinaryDifferential EquationsEndre S uliMathematical Institute, University of OxfordMarch 2, 20221 Contents1 Picard s theorem12 One-step Euler s method and its relatives: the -method .. Error analysis of the -method .. General explicit one-step method .. Runge Kutta methods .. absolute stability of Runge Kutta methods .. 193 linear multi-step Construction of linear multi-step methods .. Zero-stability .. Consistency .. Convergence .. Necessary conditions for convergence .. Sufficient conditions for convergence .. Maximum order of a zero-stable linear multi-step method .. absolute stability of linear multistep methods .. General methods for locating the interval of absolute stability.

2 The Schur criterion .. The Routh Hurwitz criterion .. Predictor-corrector methods .. absolute stability of predictor-corrector methods .. The accuracy of predictor-corrector methods .. 524 Stiff Stability of Numerical methods for stiff systems .. Backward differentiation methods for stiff systems .. Gear s method ..585 Nonlinear Stability596 Boundary value Shooting methods .. The method of bisection .. The Newton Raphson method .. Matrix methods .. linear boundary value problem .. Nonlinear boundary value problem .. Collocation method .. 70 Preface. The purpose of these lecture notes is to provide an introduction to compu-tational methods for the approximate Solution of Ordinary differential Equations (ODEs).

3 Only minimal prerequisites in differential and integral calculus, differential equation the-ory, complex analysis and linear algebra are assumed. The notes focus on the constructionof Numerical algorithms for ODEs and the mathematical analysis of their behaviour, cov-ering the material taught in the in Mathematical Modelling and Scientific Compu-tation in the eight-lecture courseNumerical Solution of Ordinary Differential notes begin with a study of well-posedness of initial value problems for a first-order differential Equations and systems of such Equations . The basic ideas of discretisationand error analysis are then introduced in the case of one-step methods. This is followedby an extension of the concepts of stability and accuracy to linear multi-step methods,including predictor corrector methods, and a brief excursion into Numerical methods forstiff systems of ODEs.

4 The final sections are devoted to an overview of classical algorithmsfor the Numerical Solution of two-point boundary value Approximation of initial value problems for Ordinary differential Equations :one-step methods including the explicit and implicit Eulermethods, the trapezium rulemethod, and Runge Kutta methods. linear multi-step methods: consistency, zero-stability and convergence; absolute stability. Predictor-corrector , stability regions, Gear s methods and their implementation. Nonlinear value problems: shooting methods, matrix methodsand List:[1] Keller, Numerical Methods for Two-point Boundary value Problems. SIAM,Philadelphia, 1976.[2] Lambert,Computational Methods in Ordinary Differential Equations . Wiley,Chichester, Reading:[1]E.

5 Hairer, Norsett, and G. Wanner,Solving Ordinary Differential Equa-tions I: Nonstiff Problems. Springer-Verlag, Berlin, 1987.[2]P. Henrici,Discrete Variable Methods in Ordinary Differential Equations . Wiley,New York, 1962.[3] Morton, Numerical Solution of Ordinary Differential Equations . OxfordUniversity Computing Laboratory, 1987.[4] Stuart and Humphries,Dynamical Systems and Numerical University Press, Cambridge, Picard s theoremOrdinary differential Equations frequently occur as mathematical models in many branchesof science, engineering and economy. Unfortunately it is seldom that these Equations havesolutions that can be expressed in closed form, so it is common to seek approximatesolutions by means of Numerical methods; nowadays this can usually be achieved very in-expensively to high accuracy and with a reliable bound on theerror between the analyticalsolution and its Numerical approximation.

6 In this section we shall be concerned with theconstruction and the analysis of Numerical methods for first-order differential Equations ofthe formy =f(x,y)(1)for the real-valued functionyof the real variablex, wherey dy/dx. In order to selecta particular integral from the infinite family of Solution curves that constitute the generalsolution to (1), the differential equation will be consideredin tandem with aninitialcondition: given two real numbersx0andy0, we seek a Solution to (1) forx > x0suchthaty(x0) =y0.(2)The differential equation (1) together with the initial condition (2) is called aninitialvalue general, even iff( , ) is a continuous function, there is no guarantee that theinitial value problem (1 2) possesses a unique solution1.

7 Fortunately, under a further mildcondition on the functionf, the existence and uniqueness of a Solution to (1 2) can beensured: the result is encapsulated in the next 1(Picard s Theorem2)Suppose thatf( , )is a continuous function of itsarguments in a regionUof the(x,y)plane which contains the rectangleR={(x,y) :x0 x XM,|y y0| YM},whereXM> x0andYM>0are constants. Suppose also, that there exists a positiveconstantLsuch that|f(x,y) f(x,z)| L|y z|(3)holds whenever(x,y)and(x,z)lie in the rectangleR. Finally, lettingM= max{|f(x,y)|: (x,y) R},suppose thatM(XM x0) YM. Then there exists a unique continuously differentiablefunctionx7 y(x), defined on the closed interval[x0,XM], which satisfies (1) and (2).The condition (3) is called aLipschitz condition3, andLis called shall not dwell on the proof of Picard s Theorem; for details, the interested readeris referred to any good textbook on the theory of Ordinary differential Equations , or the1 Consider, for example, the initial value problemy =y2/3,y(0) = 0; this has two solutions :y(x) 0andy(x) =x3 Picard (1856 1941)3 Rudolf Lipschitz (1832 1903)1lecture notes of P.

8 J. Collins,Differential and Integral Equations , Part I,Mathematical In-stitute Oxford, 1988 (reprinted 1990). The essence of the proof is to consider the sequenceof functions{yn} n=0, defined recursively through what is known as thePicard Iteration:y0(x) y0,(4)yn(x) =y0+Zxx0f( ,yn 1( )) d , n= 1,2,.. ,and show, using the conditions of the theorem, that{yn} n=0converges uniformly on theinterval [x0,XM] to a functionydefined on [x0,XM] such thaty(x) =y0+Zxx0f( ,y( )) d .This then implies thatyis continuously differentiable on [x0,XM] and it satisfies thedifferential equation (1) and the initial condition (2). The uniqueness of the solutionfollows from the Lipschitz s Theorem has a natural extension to an initial valueproblem for a system ofmdifferential Equations of the formy =f(x,y),y(x0) =y0,(5)wherey0 Rmandf: [x0,XM] Rm Rm.

9 On introducing the Euclidean normk konRmbykvk= mXi=1|vi|2!1/2, v Rm,we can state the following 2(Picard s Theorem)Suppose thatf( , )is a continuous function of itsarguments in a regionUof the(x,y)spaceR1+mwhich contains the parallelepipedR={(x,y) :x0 x XM,ky y0k YM},whereXM> x0andYM>0are constants. Suppose also that there exists a positiveconstantLsuch thatkf(x,y) f(x,z)k Lky zk(6)holds whenever(x,y)and(x,z)lie inR. Finally, lettingM= max{kf(x,y)k: (x,y) R},suppose thatM(XM x0) YM. Then there exists a unique continuously differentiablefunctionx7 y(x), defined on the closed interval[x0,XM], which satisfies (5).A sufficient condition for (6) is thatfis continuous onR, differentiable at each point(x,y) in int(R), the interior ofR, and there existsL >0 such that f y(x,y) Lfor all (x,y) int(R),(7)2where f/ ydenotes them mJacobi matrix ofy Rm7 f(x,y) Rm, andk kis amatrix norm subordinate to the Euclidean vector norm onRm.

10 Indeed, when (7) holds,the Mean value Theorem implies that (6) is also valid. The converse of this statement isnot true; for the functionf(y) = (|y1|,..,|ym|)T, withx0= 0 andy0=0, satisfies (6)but violates (7) becausey7 f(y) is not differentiable at the pointy= the counter-example in the footnote on page 1 indicates, the expression|y z|in (3) andky zkin (6) cannot be replaced by expressions of the form|y z| andky zk , respectively, where 0< <1, for otherwise the uniqueness of the Solution tothe corresponding initial value problem cannot be conclude this section by introducing the notion 1A solutiony=v(x)to (5) is said to bestableon the interval[x0,XM]if for every >0there exists >0such that for allzsatisfyingkv(x0) zk< thesolutiony=w(x)to the differential equationy =f(x,y)satisfying the initial conditionw(x0) =zis defined for allx [x0,XM]and satisfieskv(x) w(x)