1 00:00:00,140 --> 00:00:05,900 In the current chapter of the numerical methods course, we are concerned with single step 2 00:00:06,060 --> 00:00:14,040 methods meant to provide approximate numerical solutions of initial value problems for ordinary 3 00:00:14,200 --> 00:00:15,600 differential equations. 4 00:00:15,760 --> 00:00:21,600 In the previous course unit, I presented a few simple numerical methods that could solve 5 00:00:21,760 --> 00:00:26,720 initial value problems, for instance the explicit and implicit Euler method. 6 00:00:26,880 --> 00:00:32,880 And I told you that they were representatives of the class of single step methods. 7 00:00:33,040 --> 00:00:38,160 Well, now it is time to introduce that class in full generality. 8 00:00:38,320 --> 00:00:43,720 Without loss of generality, we only consider autonomous ODEs in this unit. 9 00:00:43,880 --> 00:00:48,220 Now from the previous section, recall the explicit Euler method. 10 00:00:48,380 --> 00:00:54,200 To define it, we need a temporal mesh, time points tk, covering the time interval of interest 11 00:00:54,360 --> 00:00:58,780 that for an autonomous ODE will invariably start at zero. 12 00:00:58,940 --> 00:01:04,520 Then the explicit Euler method produces a sequence of states according to the recursion 13 00:01:04,680 --> 00:01:13,000 that yk plus one is obtained from yk by adding hkf of yk in the autonomous case, where hk 14 00:01:13,160 --> 00:01:18,380 was defined as the current time step difference of tk plus one and tk. 15 00:01:18,540 --> 00:01:25,120 Now it is clear that the next state yk plus one is a function of current time step and 16 00:01:25,280 --> 00:01:26,920 current state yk. 17 00:01:27,080 --> 00:01:34,980 So we can simply write yk plus one psi of hk and yk, where psi is a suitable function. 18 00:01:35,140 --> 00:01:42,740 And now instead of writing psi as a function of these two arguments, often the time step 19 00:01:42,900 --> 00:01:47,660 argument is introduced as a superscript. 20 00:01:47,820 --> 00:01:48,820 This is just notation. 21 00:01:48,980 --> 00:01:54,020 Well, of course, psi k has to have suitable mapping properties. 22 00:01:54,180 --> 00:01:59,900 Its first argument is, well, a real number, time step, and the second argument must be 23 00:02:00,060 --> 00:02:01,060 a state. 24 00:02:01,220 --> 00:02:08,140 And it's supposed to return another state, say, and generically a vector in Rn. 25 00:02:08,300 --> 00:02:11,900 We hope that it's another valid state in the set D. 26 00:02:12,060 --> 00:02:17,180 So, this is an abstract way to write the explicit Euler method. 27 00:02:17,340 --> 00:02:23,660 In fact, it's easy to specify this mapping psi for the explicit Euler method. 28 00:02:23,820 --> 00:02:28,140 And I would like to ask you to pause the video briefly and write down the formula for the 29 00:02:28,300 --> 00:02:29,340 mapping psi. 30 00:02:29,500 --> 00:02:36,460 Not a big deal here, because psi tau y, just y, plus hf of y. 31 00:02:36,620 --> 00:02:40,940 And you can simply read this off the defining recursion. 32 00:02:41,100 --> 00:02:48,380 Now let's relate this function psi to the concept of an evolution operator that we saw 33 00:02:48,540 --> 00:02:49,540 earlier. 34 00:02:49,700 --> 00:02:56,740 In fact, a notation with a superscript time step index already suggests a relationship. 35 00:02:56,900 --> 00:03:01,060 So we again look at the explicit Euler method as an example. 36 00:03:01,220 --> 00:03:06,140 And now write it in this abstract way using the psi we have just introduced. 37 00:03:06,300 --> 00:03:12,460 And now, well, let's write y of t for the solution of the initial value problem for 38 00:03:12,620 --> 00:03:16,420 the autonomous ODE, of course, and with initial value y naught. 39 00:03:16,580 --> 00:03:20,220 Then we can write y of tk plus 1. 40 00:03:20,380 --> 00:03:26,620 So the state on the trajectory at time point tk plus 1, of course, can be computed from 41 00:03:26,780 --> 00:03:32,860 the state on the trajectory at time point tk by just applying the evolution operator 42 00:03:33,020 --> 00:03:34,420 over time hk. 43 00:03:34,580 --> 00:03:38,980 Well, this phi is this important concept. 44 00:03:39,140 --> 00:03:43,180 It's the evolution operator associated with the ODE. 45 00:03:43,340 --> 00:03:49,360 And recall that it maps from an interval and a state to another state. 46 00:03:49,520 --> 00:03:55,500 And of course, let's write i here as to be a sufficiently small interval around zero. 47 00:03:55,660 --> 00:03:57,340 Delta standing for a small number. 48 00:03:57,500 --> 00:03:59,820 And now we can relate the two mappings. 49 00:03:59,980 --> 00:04:02,260 Well, they perform very similar functions. 50 00:04:02,420 --> 00:04:08,260 You see that psi and phi are somehow close in the sense that if the method, the single 51 00:04:08,420 --> 00:04:14,620 step method provides states that are close to the trajectory, then psi should be an approximation 52 00:04:14,780 --> 00:04:15,780 of phi. 53 00:04:15,940 --> 00:04:20,940 And this is why psi is now called the discrete evolution operator. 54 00:04:21,100 --> 00:04:26,860 And this already suggests, and I've already said this, that this is supposed to approximate 55 00:04:27,020 --> 00:04:33,260 the continuous, the genuine evolution operator in some sense that will be specified in a 56 00:04:33,420 --> 00:04:34,920 moment. 57 00:04:35,080 --> 00:04:37,540 So psi should be close to phi. 58 00:04:37,700 --> 00:04:42,620 And now you see that everything that is important about the single step method is captured by 59 00:04:42,780 --> 00:04:43,780 the psi. 60 00:04:43,940 --> 00:04:50,900 And now this, I hope, shed light on my initial statement that I made in the beginning of 61 00:04:51,060 --> 00:04:57,100 the chapter, namely that the construction of single step methods and numerical integration 62 00:04:57,260 --> 00:05:03,060 in general aims to approximate evolution operators induced by ordinary differential equations. 63 00:05:03,220 --> 00:05:08,580 And now, once we have this concept of a discrete evolution operator at our disposal, we can 64 00:05:08,740 --> 00:05:11,940 give a rigorous definition of what constitutes a single step method. 65 00:05:12,100 --> 00:05:14,920 Here for the autonomous case, for the sake of simplicity. 66 00:05:15,080 --> 00:05:17,180 So we need a discrete evolution. 67 00:05:17,340 --> 00:05:25,600 And then we apply this discrete evolution on a temporal mesh and so generate a sequence 68 00:05:25,760 --> 00:05:26,740 of states. 69 00:05:26,900 --> 00:05:29,660 And this is what is called a single step method. 70 00:05:29,820 --> 00:05:36,040 Keep in mind that the psi may not have such a nice formula as for the explicit Euler method. 71 00:05:36,200 --> 00:05:42,580 So psi may be defined only implicitly, like in the case of implicit Euler method, implicit 72 00:05:42,740 --> 00:05:46,180 midpoint rule, maybe, if a defining equation. 73 00:05:46,340 --> 00:05:52,940 Let's look at the implicit Euler method because there we find that psi H of y is such that 74 00:05:53,100 --> 00:05:55,180 it satisfies this equation. 75 00:05:55,340 --> 00:05:58,380 Now this is an implicit definition of a state. 76 00:05:58,540 --> 00:06:01,340 As you know, of course, only for sufficiently small h. 77 00:06:01,500 --> 00:06:05,900 This is a sort of tested assumption that is made again and again in order to render the 78 00:06:06,060 --> 00:06:09,300 discrete evolution operator well defined, of course. 79 00:06:09,460 --> 00:06:13,220 Now I've told you that psi has to be close to phi. 80 00:06:13,380 --> 00:06:18,320 And so we have to find a way how to express this closeness rigorously. 81 00:06:18,480 --> 00:06:23,400 And one way is the notion of consistency, giving rise to consistent single step methods. 82 00:06:23,560 --> 00:06:26,420 So what does it mean that psi is close to phi? 83 00:06:26,580 --> 00:06:29,780 Let's recall a key property of phi. 84 00:06:29,940 --> 00:06:36,180 When you vary the t argument in phi, then you get a trajectory in state space. 85 00:06:36,340 --> 00:06:39,520 And so this trajectory satisfies the ODE. 86 00:06:39,680 --> 00:06:42,280 So we differentiate the trajectory. 87 00:06:42,440 --> 00:06:44,900 Of course, we get f of the state. 88 00:06:45,060 --> 00:06:48,420 And when we do this at t equals zero, we get f of the initial state. 89 00:06:48,580 --> 00:06:50,860 And here this is just y, y as state. 90 00:06:51,020 --> 00:06:53,660 Well, this is a simple property. 91 00:06:53,820 --> 00:06:59,620 So and in order to connect psi and phi, we want the same to hold for the discrete evolution. 92 00:06:59,780 --> 00:07:06,580 But now I choose h as that scalar independent variable to indicate that it's related to 93 00:07:06,740 --> 00:07:07,740 the step size. 94 00:07:07,900 --> 00:07:13,140 So we want this to be f of y, of course, for any state y. 95 00:07:13,300 --> 00:07:19,900 And we impose this as a minimal requirement, this property, and call it consistency. 96 00:07:20,060 --> 00:07:23,700 The formal definition is slightly different, and I'm giving it now. 97 00:07:23,860 --> 00:07:31,300 The formal definition states that the discrete evolution has to be of a particular form. 98 00:07:31,460 --> 00:07:36,260 It's psi, h, y is y plus h times little psi. 99 00:07:36,420 --> 00:07:41,540 And the little psi, when the first argument is zero, has to be equal to f of y. 100 00:07:41,700 --> 00:07:44,180 y is the second argument. 101 00:07:44,340 --> 00:07:46,060 And it's supposed to be continuous. 102 00:07:46,220 --> 00:07:49,540 Now let's relate the two definitions. 103 00:07:49,700 --> 00:07:52,500 This is not hard. 104 00:07:52,660 --> 00:07:58,060 So now let's compute this derivative for h to zero. 105 00:07:58,220 --> 00:08:01,820 Product rule, this is now, this is the capital psi. 106 00:08:01,980 --> 00:08:06,460 I should now really indicate that this is a capital psi. 107 00:08:06,620 --> 00:08:13,180 Because now here this is the little psi, little psi, the h, h, y, and the whole thing for 108 00:08:13,340 --> 00:08:14,840 h equals zero. 109 00:08:15,000 --> 00:08:19,060 So h equals zero means that this partial derivative is suppressed. 110 00:08:19,220 --> 00:08:22,420 What we end up with is the zero y. 111 00:08:22,580 --> 00:08:25,740 And now by definition, this is f of y. 112 00:08:25,900 --> 00:08:31,620 And so, well, the two definitions for consistency, they are more or less equivalent. 113 00:08:31,780 --> 00:08:37,220 Now, what is the little psi for the simple case we have already seen, namely the explicit 114 00:08:37,380 --> 00:08:38,500 Euler method? 115 00:08:38,660 --> 00:08:39,660 What is this? 116 00:08:39,820 --> 00:08:44,280 Please pause the video briefly and give the formula for the little psi. 117 00:08:44,440 --> 00:08:48,900 Of course, when we go back to the formula for the explicit Euler method, this one, we 118 00:08:49,060 --> 00:08:54,700 copy it, and then we immediately see that here the little psi is really just the f of 119 00:08:54,860 --> 00:08:56,980 y, no h dependence. 120 00:08:57,140 --> 00:09:01,820 And of course, this satisfies the requirements of the definition. 121 00:09:01,980 --> 00:09:05,460 So for explicit Euler, it's straightforward to show consistency. 122 00:09:05,620 --> 00:09:12,300 Now let's look at an implicit method and wonder whether we can verify the consistency. 123 00:09:12,460 --> 00:09:16,580 And we don't do it for the simplest case, it would be implicit Euler, but for the next 124 00:09:16,740 --> 00:09:24,100 simpler case, implicit midpoint method, well, we assume f is smooth, while continuity would 125 00:09:24,260 --> 00:09:25,260 be okay. 126 00:09:25,420 --> 00:09:31,380 And then let's look at the discrete evolution of the implicit midpoint method. 127 00:09:31,540 --> 00:09:40,860 So we have psi hy is implicitly defined to solve hf of y plus, and here it's yk plus 128 00:09:41,020 --> 00:09:46,020 one in the recursion, but when we formulate things in terms of discrete evolution operator, 129 00:09:46,180 --> 00:09:49,240 then replace yk plus one by psi hy. 130 00:09:49,400 --> 00:09:57,820 And now the trick is to, well, use this definition of psi hy and plug it into itself. 131 00:09:57,980 --> 00:10:07,300 y plus hf y, sorry, there's a small error, one half, and so now we have y plus one half 132 00:10:07,460 --> 00:10:14,340 hf of one half y plus psi h of y. 133 00:10:14,500 --> 00:10:17,860 Now it's getting a bit cumbersome, but this is what we have. 134 00:10:18,020 --> 00:10:23,700 Well, look at the division on the left page, and well, there's little choice we have. 135 00:10:23,860 --> 00:10:26,780 This should be the psi, the little psi. 136 00:10:26,940 --> 00:10:33,660 Second requirement that psi zero y is f of y is clear, because well, there's an h here, 137 00:10:33,820 --> 00:10:37,660 setting it to zero, then we are left with just the y argument. 138 00:10:37,820 --> 00:10:38,820 This is clear. 139 00:10:38,980 --> 00:10:40,020 Continuity. 140 00:10:40,180 --> 00:10:44,460 Well, here we can appeal to the implicit function theorem. 141 00:10:44,620 --> 00:10:50,180 Because of the continuity of f, we know that this implicitly defined psi, capital psi to 142 00:10:50,340 --> 00:10:53,220 be precise, is continuous. 143 00:10:53,380 --> 00:10:58,720 And so now we have verified all the properties of psi required for consistency. 144 00:10:58,880 --> 00:11:02,180 So this is how it works in the case of an implicit method. 145 00:11:02,340 --> 00:11:05,860 Finally, a remark on single step methods. 146 00:11:06,020 --> 00:11:10,860 Going back to the generic definition, we see that they generate a sequence of states. 147 00:11:11,020 --> 00:11:16,100 Then one is happy with this, in particular when you only want to know the state at some 148 00:11:16,260 --> 00:11:17,460 final time. 149 00:11:17,620 --> 00:11:24,020 But in other contexts, one would really like to have a trajectory, and one can obtain a 150 00:11:24,180 --> 00:11:25,180 trajectory. 151 00:11:25,340 --> 00:11:31,620 We have already seen this in the case of the explicit implicit Euler method, where I mentioned 152 00:11:31,780 --> 00:11:38,460 the concept of an Euler polygon, the trajectory obtained by piecewise linear interpolation 153 00:11:38,620 --> 00:11:40,280 of the approximate states. 154 00:11:40,440 --> 00:11:41,700 And this can be generalized. 155 00:11:41,860 --> 00:11:49,160 So when you're given a sequence of states, we can obtain a trajectory by interpolation. 156 00:11:49,320 --> 00:11:55,620 And then we get a function of time defined for all times between zero and say final time, 157 00:11:55,780 --> 00:11:59,540 capital T, and the values in the state space. 158 00:11:59,700 --> 00:12:05,540 Well, there's one twist that gives us more possibilities for interpolation here. 159 00:12:05,700 --> 00:12:14,300 Remember that the derivative of the trajectory, say at time tk, of course y is the trajectory. 160 00:12:14,460 --> 00:12:18,380 This is yf applied to ytk for an autonomous ODE. 161 00:12:18,540 --> 00:12:25,840 And we have an approximation of ytk, so we can expect that f of yk is an approximation 162 00:12:26,000 --> 00:12:27,040 for the derivative. 163 00:12:27,200 --> 00:12:30,820 So we can also use derivatives for the interpolation. 164 00:12:30,980 --> 00:12:35,260 So because derivatives are cheap, say on derivative values. 165 00:12:35,420 --> 00:12:41,620 For instance, one can use piecewise cubic Hermite interpolation, which is a piecewise 166 00:12:41,780 --> 00:12:49,040 polynomial interpolation based on cubic polynomials and harnessing both function values and derivative 167 00:12:49,200 --> 00:12:50,180 values. 168 00:12:50,340 --> 00:12:56,660 And this automatically gives us a C1 approximation of the trajectories. 169 00:12:56,820 --> 00:13:03,460 And with this remark, I would like to close this unit and then we move on to more concrete, 170 00:13:03,620 --> 00:13:06,220 but also more complex single step methods.