11 30 10 cs412 introduction to numerical analysis Lecture 21 Numerical Solution of Differential Equations Instructor Professor Amos Ron 1 Scribes Yunpeng Li Mark Cowlishaw Nathanael Fillmore Introduction The last topic of this semester is numerical solution of differential equations We will not be so much concerned with the theory of differential equations but rather with algorithms for finding numerical approximations Recall that differential equations are classified by order For this discussion we shall consider only first order differential equations In general when trying to solve a first order differential equation we are searching for a function y t given some constraints expressed in terms of the independent variable t the dependent variable y and the derivative of y y 0 When we say that we are looking for a function y t of course what we really mean is that we need to be able to determine the value of y at some point b i e y b To illustrate consider the following problem Example 1 1 Determine y t given the first order differential equation y 0 t 2t y t Solution The trivial solution y t 0 is obvious and it is not too difficult to find other solutions 2 y t et 2 y t 2 et It turns out that there are infinitely many solutions to example 1 1 each of the form y t k et 2 for all k R We can rewrite example 1 1 in terms of a function f t y f t y 2t y y 0 t f t y t As we noted there are infinitely many solutions to this problem In fact the solutions to this problem set up a one to one correspondence for the entire 2 space defined by t and y In other 1 Y t1 y1 t Figure 1 Solutions to y 0 t 2t y t words for any point t1 y1 there is a unique solution that passes through this point as shown in Figure 1 This means that example 1 1 poses an ill formed problem since infinitely many functions solve this differential equation we have no way to pick one Our analysis above shows the way to make this problem well formed since every point t y has a unique solution that passes through it if we provide a single data point with the differential equation our problem will be well formed In general we will discuss problems of the form Problem 1 1 Initial Value Problem Find y t given y 0 t f t y t y a Ya A problem of this form is called an Initial Value Problem It is interesting to note that first order differential equations which are inherently difficult in terms of the theory of differential equations such as f t y 2t2 y 10 3t3 y do not pose any special difficulty to numerical methods since f is a polynomial from a numerical perspective this problem is not much different from a simpler problem such as example 1 1 2 Solution of Initial Value Problems To solve initial value problems we will use numerical integration Note that we can determine the difference between values of y using the definite integral of y 0 Z b y 0 t dt y b y a a 2 So given an initial value problem with y a defined and a differential equation for y 0 t we can find y b using Z b y b y a y 0 t dt a We will use this idea to define an iterative method for estimating y b given y a As usual we will split the interval between a and b into many small equal sized subintervals with partition points a t0 t1 tN b Using numerical integration we will generate values for y and y 0 at each partition point until we reach the final value y tn Yn This process is depicted in Figure 2 h a t0 t1 t2 Y0 Y1 Y2 Y 0 Y 1 Y 2 b tn Yn Y n Figure 2 An Iterative Method for Solving Initial Value Problems A generic step in this method consists of 0 Given Y0 Y1 Yi 1 and Y00 Y10 Yi 1 Find 1 Yi and 2 Yi0 Notice that given Yi Yi0 can be easily computed using Yi0 f ti Yi Thus our method really consists of finding each Yi then using this value to calculate Yi0 One way to solve this problem is to calculate each value Yi 1 based only on the previous point Yi that is Z y ti 1 y ti ti 1 ti y 0 t dt z but an obvious problem is how to calculate the definite integral when we only know the value of the function at the left endpoint y 0 ti Since the intervals are small we may be able to use simple rules for approximating the definite integral to estimate this value If we use the rectangle rule to approximate the definite integral then the resulting iterative process is known as Euler s method 3 3 Euler s Method Recall that the rectangle rule uses the value of the function at the left end point of an interval to approximate the definite integral Z ti 1 y 0 t dt Yi0 h ti Euler s method applies the above idea to generate each Yi for 0 i N 1 Yi 1 Yi Yi0 h for i 0 1 N 1 0 using the given differential equation then generates each Yi 1 0 Yi 1 f ti 1 Yi 1 To illustrate this process consider the following example Example 3 1 Given the following initial value problem y 0 t 3t2 y t y 0 1 find y 0 4 3 Of course it is easy to see that y t e t is the solution For the purpose of comparing methods of solving the IVP we will solve this problem using an interval width of 0 1 and we will only compare the solution over the last interval that is given Y0 y 0 Y1 y 0 1 Y2 y 0 2 3 and Y3 y 0 3 find Y4 y 0 4 We can calculate the initial values using y t e t and 3 y 0 t 3t2 e t Y00 0 Y0 1 Y1 0 999 Y10 0 03 Y2 0 992 Y20 0 119 Y3 0 973 Y30 0 263 Y4 0 938 Solution Using Euler s method we calculate y 0 4 as Z 0 4 y 0 4 y 0 3 y 0 t 0 3 Y3 Y30 0 1 0 973 0 263 0 1 0 947 This yields an error of about 0 009 which is rather large for such a small single step 4 4 Modified Euler s Method How can we improve this estimate Our knowledge of approximating the definite integral tells us that the midpoint rule does a far better job than the rectangle rule so we might try applying the midpoint rule However there is a problem The midpoint rule requires knowing the value of the …