Numerical Solution of Differential Equations

1 Numerical Solution of Differential EquationsLiz BradleyDepartment of Computer ScienceUniversity of ColoradoBoulder, Colorado, USA 80309-0430c 1998 Revised versionc 2002, Report on Curricula and Teaching CT003-981 Ordinary Differential EquationsA differential equation expresses a set of constraints amongthe derivatives of an unknown function. Here s asimple example:ddtx(t) =ax(t)(1)What this means is that the derivative of the unknown function is equal toatimes the unknown function compactness, the independent variable (t, here) is often omitted and the notationdxdtis often abbreviatedwith a dot, like so: xor with a prime, like so:x.

2 Using these two notations, equation (1) becomes x=axorx =ax, solve a differential equation, you (generally) can t justintegrate both sides. Rather, you have to findsome function that satisfies the constraints expressed in that equation in the example above, some functionx(t) whose derivativex (t) is equal to a constant multiple of the function itself1. There is no general, foolproofway to do this unless the equation is very simple. Differential Equations textbooks are cookbooks that giveyou lots of suggestions about approaches, but there are lotsof differential Equations (DEs) that simply don thaveanalytic solutions that is, solutions that you can write down.

3 These Equations can only be solvednumerically, using the kinds of methods that are described in these Equations are interesting and useful to scientists and engineers because they model thephysical world: that is, they capture the physics of a system, and their solutions emulate the behavior of thatsystem. The class of differential Equations that have no analytic solutions has a highly specific and interestinganalog in the physical world: they model so-calledchaoticsystems. This means that Numerical methods forsolving ODEs are an essential tool for anyone (like me) who studies equation (ODE) has only one independent variable, and all derivatives in it aretaken with respect to that variable.

4 Most often, this variable is timet, but some books usexas the independentvariable; pay careful attention to what the derivative is taken with respect to so you don t get s an example of where ODEs come from. Consider a massmon a spring with spring constantk:1 The answer isx(t) =beat+c, wherea, bandcare you definexas the position of the end of the spring, as shown in the figure,and setx= 0 at the point wherethe end of the spring would naturally rest if it were unloaded( , ifm= 0), you can write the following forcebalance at the mass: F=mamg kx=ma(2)Accelerationais the second derivative of position:a=x , so the equation becomes:mg kx=mx Themgforce is gravity;kxis Hooke s law for the force exerted by a spring.

5 The signs ofmgandkxareopposite because gravity pulls in the direction of positivexand the spring pulls in the direction of negativex. Gathering terms and dividing through bymgives you the following ODE for the spring-mass system:x (t) +kmx(t) g= 0 Note that this ODE model does not include the effects of friction; if friction plays an important role in thesystem, this model is independent variabletonly appearsimplicitlyin this equation that is, hidden in the time dependenceof the functionx(t). This means that thephysics(the governing Equations ) of the system is not a functionof time.

6 If, on the other hand, someone were holding the top ofthe spring and moving it up and down in asinusoidal pattern, the apparent gravity acting on the masswould change sinusoidally (much as the gravitychanges in an elevator that is starting or stopping), and theequation would look like this:x (t) +kmx(t) g0sint= 0 The variabletappearsexplicitlyin this ODE by itself, and not wrapped in the brackets of a function likex(t). Systems where this occurs are an ODE is the degree of the highest derivative in that =axis a first-orderODE,x tanx = 2 is a third-order ODE, and the spring-mass equation above is second order.

7 Annth-order ODE can be transformed intonfirst-order ODEs (annth-order ODE system ), and vice versa2. Thistransformation requires the introduction of helper variables. I ll illustrate with the third-order (n= 3)examplex 1(t) + 36 logx 1(t) x21(t) + sin 2t= 141. The first step is to rewrite the ODE with the highest-order term by itself on the left-hand side:x 1= 14 +x21 36 logx 1 sin 2t2 This means that the order of an ODE system is equal to the number of first-order ODEs in the corresponding ODE The second step is to definen 1 helper variables, like so:x 1=x2, x 2=x33.

8 The third step is to rewrite the equation from step 1 using the helper variables, with no derivative signsat all on the right-hand side and only one on the left-hand side:x 3= 14 +x21 36 logx2 sin 2t4. Finally, one appends the Equations that define the helper variables to that rewritten equation, obtainingthenth-order system:x 1=x2x 2=x3x 3= 14 +x21 36 logx2 sin 2tThis set of first-order ODEs is equivalent tox 1= 14 +x21 36 logx 1 sin 2t, as you can see by substitutingthe first two Equations into the third. The variables that appear on the left-hand side of an ODE system aretermed thestate variablesof the system.

9 Thestate vector~xof this system is (x1x2x3)Tand the ODE systemis of the form~x =~f(~x, t)I will use this standard form throughout these notes. For autonomous systems,~fdoesn t actually depend ontat all, so the standard form collapses to~x =~f(~x)Here s a famous third-order ODE system, derived by Edward Lorenz, that models the physics of a fluidthat is being heated from below:~x =~f(~x) = x y z = a(y x)rx y xzxy bz (3)This system of ODEs is autonomous (no explicitts on the right-hand side). It is alsononlinear, as you cansee by looking at equation (3) above; if any of right-hand sides of any of these Equations contain any termsthat are nonlinear ( , products, transcendentals, powers.)

10 In the state variables, then the whole systemis nonlinear. Here, the product termsxzandxyare the culprits;a, r, andbare constants. The simpleexamplex =ax, in contrast, islinear;x 1+ 36 logx 1 x21+ sin 2t= 14 is nonlinear because of the log andx21terms. This terminology can be a little confusing; the individual Equations in the system (3) like thex =a(y x) one may look linear by themselves, but you need to rememberthat the whole thing is apackage: asingle(albeit vector-valued) function of a vector-valued argument.