 We are discussing numerical methods for initial value problems involving first order ordinary differential equations. In the last two classes we have developed Euler's method and studied the convergence analysis also. Although Euler's method is easy to implement, this method is not so efficient in the sense that to get a better approximation, one needs a very small step size. This is mainly because Euler's method is a first order method. One way to get better accuracy with less step size h is to go for higher order methods. Recall that in Euler's method, we considered approximation of the derivative by taking the Taylor expansion with degree 1 and truncated all terms of order 2 and higher. In order to achieve higher order accuracy, we will have to include higher order terms in the approximation of y dash. This will naturally bring in higher order total derivatives of y, which is not preferred computationally. Rangekutta methods are focused to achieve higher order accuracy without involving higher order total derivatives of y. This makes the methods more popular till to date. Our interest in this method is to understand the rigorous derivation of the Rangekutta method of order 2. Let us recall that we are interested in devising numerical method for the initial value problems of the form y dash is equal to f of x comma y for all x in some interval x naught to b, where y of x naught is equal to y naught. Now, let us use the Taylor approximation for y of x plus h above the point x, where we now retain the terms up to order 2 in h and truncate the terms of order 3 and above. Well, we have to write the reminder term at the third order, but it is just comfortable to simply write the big O of h 3, because anyway at the end of the derivation, we are going to neglect this terms. Therefore, it is more comfortable for us to simply write these terms in the big O notation to just indicate that whatever expression we have after this term is going to involve something of order h cube. Recall in Euler method, we have taken only the first order term and made this term with a reminder. That is how we got only the first order accuracy in Euler method. Now, to achieve the second order accuracy, we are also including the second order term now into our expression and we are truncating the Taylor expansion only from order 3. That is the main idea of achieving order 2 now in the numerical method. Let us see how to go ahead with this. You can see that by including the higher order term, we have y double dash of x. If you recall, we have substituted y dash equal to f of x y, when we were deriving the Euler method. Now, we have to similarly handle this second order term also for that what we will do is we will differentiate the given ODE with respect to x and we get y double dash is equal to dou f by dou x into dx by dx that is equal to 1 plus dou f by dou y into dy by dx that is y dash. This is simply the chain rule that is applied on f. Remember f is a function of two variables x and y where y is a function of x. Therefore, we get this expression if we use the chain rule. Therefore, we have the expression for y dash which is coming from the equation and we have y double dash which we obtained by differentiating the ODE once with respect to x. Now, we will substitute these two expressions into the Taylor expansion and see what happens. We get y of x plus h is equal to y of x which is already there plus h into instead of y dash we will put f of x y plus h square by 2 which is already there into instead of y double dash we put this expression here and then finally, we are going to neglect all other terms which are of big O of h cube. So, this is what we have achieved so far using Taylor expansion and then substituting y dash and y double dash. Now, let us write this expression in a slightly different way what we will do is we will split this term into two equal parts h by 2 into f of x comma y plus h by 2 into f of x comma y. One part we will retain separately and merge the other part into this third term and we get this expression now as the third term. We have not done anything we have just split this second term into two parts and merge the second part into the third term that is all we did. Well, what is the reason to do this small adjustment? This will be clear if you look at the Taylor expansion of a two variable function up to degree 1. It will happen to be of this form, but before doing this let us fix our notation at par with the numerical method for that let us replace x by x j and write the above expansion now replacing x by x j. Therefore, x plus h will be x j plus h and wherever x comes we will put x j and we are just writing the same expression which is written here now here with x replaced by x j. So, we have this expression for y of x j plus 1 with us let us now see how to handle this term because this term involves partial derivative of f which are not so comfortable from the implementation point of view. Well, we have to provide all these informations into our code in order to compute the solutions and that is not very comfortable in general. Let us see if we can somehow replace this expression with something which can be handled computationally for this let us first recall the Taylor representation of f about the points xi comma tau. Let us take some point xi comma tau and let us recall how to write the Taylor representation of a two variable function f about this point xi comma tau the value of f at some point st which is very close to xi comma tau can be written as f of xi comma tau plus s minus xi into dou f by dou s evaluated at xi comma tau plus t minus tau. Now, this is coming for the second argument of f times dou f by dou t evaluated at xi comma tau you can observe that this is just a straight forward generalization of 1 d Taylor expansion to two dimensional Taylor expansion. We just want the Taylor expansion up to order 1 therefore, from the second order term onwards we will simply write in the big o notation because at the end anyway we will try to observe this second order term into the term which we already have why it is so we already have a h here. Now, we will try to see that s minus xi is of order h therefore, its square will be of order h square into h will be of order h cube already we made our mind to neglect all the terms which are of order h cube. Therefore, whatever may be the terms which are sitting after this will finally, be neglected in our formula therefore, we will not write them explicitly or write as a reminder, but we will simply write them in the big o notation. Now, let us take xi comma tau equal to x j comma y of x j and look at this expression in the bracket and compare it with the expression here. We can see that if you take s is equal to x j plus 1 and t is equal to y of x j plus h into f of x j comma y of x j then this term will match exactly with the term in the bracket. Therefore, we will take these expressions and put them in this and then replace f of s comma t where s is x j plus 1 and t is this expression we will replace this entire expression by f of this. So, that is the idea this is how you are getting rid of this partial derivatives of f which we saw that they are not very comfortable from the computational point of view and they are just going off from our scene by simply replacing this expression by f of x j plus 1 comma y of x j plus h f right. If you do that then you will have y of x j plus 1 equal to y of x j plus h by 2 into f of x j comma y of x j up to this we already have here. Now, what we are doing is only in the third term we are simply replacing the expression inside the bracket by f of this using this Taylor expansion of a bivariate function right and that gives us the third term as h by 2 into f of x j plus 1 comma y of x j plus h into f of x j comma y of x j ok. Plus of course, we have some higher order terms which we anyway made our mind to neglect them therefore, you just simply freeze them into this notation and keep and finally, we will simply neglect this higher order terms and write the formula and just to denote that by doing this we only get an approximation to this solution. We will use the notation y of j plus 1 in order to denote the approximate value of the exact solution y of x j plus 1. Therefore, we have this formula y j plus 1 which is the approximate value of y of x j plus 1 obtained using this formula where all the exact values are now replaced by the approximate values of the solution otherwise this expression is exactly the same as this expression just we have neglected all the higher order terms and this is precisely the Runge-Kutta method of order 2. It is very nice and interesting to see how we have achieved the higher order without using the partial derivatives of f with respect to x and y right. We have cleverly eliminated this partial derivatives simply using the Taylor expansion of a bivariate function up to order 1. So, that is the idea behind the Runge-Kutta method of order 2. What is the truncation error of this method? Well, the truncation error of this method is something of order 3 that is what we can see. In fact, you can easily write this expression it is simply the reminder term in the Taylor formula right. Now, the question is why this method is of order 2 when the truncation error of the method is of order 3? Well, you can easily answer this question if you have understood why we had Euler method as of order 1. If you recall the truncation error in Euler method is of order 2, but the Euler method is of order 1 right. As far as the formulas that we derived so far for the first order ODs are concerned the order of the method will be 1 less than the order of the truncation error. Let us see why it is so in the present case recall that we started our derivation of the method with the Taylor expansion given like this with the truncation part being order h cube right and this is finally used in obtaining an approximation of y dash. So, let us write this expression in this form where you bring all the terms on the left hand side keeping only the y dash term on the right hand side plus of course the terms which you are going to neglect. Now, in order to get y dash on the right hand side what you have to do you have to divide both sides by h right that is what we are doing. So, in order to get y dash on the right hand side from here we are dividing both sides by h and that takes away one order from your truncated part. So, that is why the method is order 1 less than the truncation error order. This happens with the methods that we have derived so far that is it happened with Euler method also and also it happened with the Runge-Kutta method that we have derived so far because we are approximating y dash and in the Taylor formula y dash appears with h. Therefore, in order to get an approximation for y dash finally, you have to divide both sides by h that will naturally reduce one order from the truncation error. Therefore, you should always remember that in such methods if your truncation error is of order something say n then the order of the method will be 1 less that is it will be n minus 1. Here the truncation error is 3 therefore, the order of the method is 2. Similarly, in the Euler method the truncation error was of order 2 therefore, the method as such is of order 1. Now, let us put the Runge-Kutta method of order 2 in a slightly simple form which is more easy for us to remember generally how to write it recall the Runge-Kutta method of order 2 that we have derived in our previous slide is given like this. Now, what we will do is we will write it in a different form which is easy for us to remember how to do that well you can write this expression as y j plus 1 equal to y j plus 1 by 2 into k 1 plus k 2 where k 1 is h into f of x j comma y j right. Now, you can see that that expression is sitting here therefore, k 2 can be written as h into f of x j plus 1 plus y j is already there I am sorry this should be y j here plus this part is simply now named as k 1 therefore, you can simply plug in k 1 here and get this expression well this is the same as this it is only written in a different form generally it is also easy for us to implement as a computer code and also from the examination point of you this may be easy for us to remember right. Since the method that we have derived is of order 2 means it is truncation error is of order 3 remember that the Euler method is of order 1 and therefore, it is truncation error is of order 2. So, we gained one order in the Runge-Kutta method we derived so far and therefore, for a given h we expect better accuracy from the Runge-Kutta method of order 2 when compared to the Euler method right. Let us see this by an example let us consider the initial value problem y dash equal to y on the interval 0 comma 0.04 with the initial condition as y of 0 equal to 1 it is a very simple linear ordinary differential equation and therefore, we of course, know the exact solution of this initial value problem. Now, let us compute the numerical solution of this initial value problem using the Runge-Kutta method of order 2 and let us present the numerical result in this table well we start with x equal to 0 and x equal to 0 and we are given that at x equal to 0 y takes the value 1. Therefore, we do not do any computation at x equal to 0 we simply put the initial condition here and now we have to find k 1 and k 2 using the formula shown in the previous slide and they are given by 0.01 for k 1 and 0.0101 for k 2 you can check this and once you have k 1 and k 2 you can plug in that into the formula of y 1 and you can get it as 1.01005 and once you have y 1 you can get again k 1 and k 2 and then you can obtain y 2 from there remember you have to use the formula that we have displayed in the previous slide I hope I do not need to show you the calculation you can do it easily right once you have y 2 then again put y 2 and x 2 and get k 1 and k 2 once you have k 1 and k 2 you can go to get y 3 right. So, like that you can keep on going till you reach the right hand limit of the domain that you are interested in here we have taken it as 0.04 and we have taken h as 0.01 from this column you can see that we have taken h equal to 0.01 what is the exact solution the exact solution is y of x equal to e power x and the value of the exact solution at the point x equal to 0.04 is given by this number approximately at least up to the six significant digits and what we obtained as the numerical solution from the rangikutta method of order 2 we obtained this value now what is the error error is something nearby 10 to the power of minus 6 and that is pretty good when compared to the Euler method you may recall that in Euler method we computed the numerical solution of the same initial value problem with the same h is equal to 0.01 and we got the error in the numerical solution of y of 0,04 as something nearby 2 into 10 to the power of minus 4. Therefore, you can see that we have a good improvement in the numerical method when it is of order 2 when compared to a numerical method of order 1 at least in this particular example we can clearly see the advantage of going for the higher order method when compared to the lower order methods. Once you understand the derivation of the rangikutta method of order 2 the same idea can be adopted to derive the rangikutta method of any higher order the idea is to retain adequate number of terms in the Taylor expansion and then you have to use the bivariate Taylor approximation accordingly. For instance the rangikutta method of order 4 is given by this expression where y j plus 1 is equal to y j plus 1 by 6 into k 1 plus 2 k 2 plus 2 k 3 plus k 4 where k 1 is given like this k 2 is given like this k 3 is given like this and k 4 is given like this and k 4 is given like this well the idea of the derivation of the rangikutta method of order 4 is very much similar to what we did with the rangikutta method of order 2. However, the derivation is very lengthy, but it is straight forward we omit the derivation of rangikutta method of order 4 and just take this formula alone for our computation. Let us compute the solution for the same initial value problem that we did in our previous example and see how the accuracy is in the present case. Let us consider the initial value problem y dash equal to y on the interval 0 to 0.004 with the initial condition as y of 0 equal to 1 and now we will use the rangikutta method of order 4. You can see that at the initial point you have the exact value because it is given to us and therefore there is no error involved in the numerical solution when compared to the exact solution at the initial point then when you go to compute the solution at the later points we are expected to accumulate some error right. At x equal to 0.01 remember again we are taking h equal to 0.01 numerical solution is given like this and it is nice to see that exact solution is also the same at least up to the 6 significant digits that is shown here therefore the error is 0. Now let us go to the next grid point 0.02 here also you can see that up to 6 significant digits we do not have any error similarly for x equal to 0.03 and even for x equal to 0.04 we got almost the exact solution at least up to 6 significant digits. This is very nice to see and also gives us a clear understanding that the rangikutta method of order 4 is better than the rangikutta method of order 2. In fact, it is in general true that as you go on for higher order methods you tend to get better results but this may not hold for all the problems for all the methods but this is a general idea. With this note let us finish this lecture. Thank you for your attention.