1 00:00:04,470 --> 00:00:06,570 PROFESSOR: The most frequently used ODE solver 2 00:00:06,570 --> 00:00:09,865 in MATLAB and Simulink is ODE45. 3 00:00:09,865 --> 00:00:14,700 It is based on method published by British mathematicians JR 4 00:00:14,700 --> 00:00:19,230 Dormand and PJ Prince in 1980. 5 00:00:19,230 --> 00:00:22,760 The basic method is order five. 6 00:00:22,760 --> 00:00:29,350 The error correction uses a companion order four method. 7 00:00:29,350 --> 00:00:33,610 The slope of tn is, first same as last left over 8 00:00:33,610 --> 00:00:37,020 from the previous successful step. 9 00:00:37,020 --> 00:00:44,070 Then there are five more slopes from function values at 1/5 h, 10 00:00:44,070 --> 00:00:52,240 3/10h, 4/5h, 8/9h and then at tn plus 1. 11 00:00:52,240 --> 00:00:55,710 These six slopes, linear combinations of them, 12 00:00:55,710 --> 00:00:59,540 are used to produce yn plus 1. 13 00:00:59,540 --> 00:01:05,900 The function is evaluated at tn plus 1 and yn plus 1 14 00:01:05,900 --> 00:01:07,950 to get a seventh slope. 15 00:01:07,950 --> 00:01:10,940 And then linear combinations of these 16 00:01:10,940 --> 00:01:14,760 are used to produce the error estimate. 17 00:01:14,760 --> 00:01:17,630 Again, if the error estimate is less 18 00:01:17,630 --> 00:01:22,610 than the specified accuracy requirements 19 00:01:22,610 --> 00:01:24,800 the step is successful. 20 00:01:24,800 --> 00:01:31,750 And then that error estimate is used to get the next step size. 21 00:01:31,750 --> 00:01:36,460 If the error is too big, the step is unsuccessful 22 00:01:36,460 --> 00:01:42,390 and that error estimate is used to get the step 23 00:01:42,390 --> 00:01:44,800 size to do the step over again. 24 00:01:44,800 --> 00:01:47,890 If we want to see the actual coefficients that are used, 25 00:01:47,890 --> 00:01:51,550 you can go into the code for ODE45. 26 00:01:51,550 --> 00:01:53,850 There's a table with the coefficients. 27 00:01:53,850 --> 00:01:59,190 Or you go to the Wikipedia page for the Dormand-Prince Method 28 00:01:59,190 --> 00:02:01,160 and there is the same coefficients. 29 00:02:04,590 --> 00:02:08,150 As an aside, here is an interesting fact about higher 30 00:02:08,150 --> 00:02:10,780 order Runge-Kutta methods. 31 00:02:10,780 --> 00:02:15,010 Classical Runge-Kutta required four function evaluations 32 00:02:15,010 --> 00:02:18,800 per step to get order four. 33 00:02:18,800 --> 00:02:24,490 Dormand-Prince requires six function evaluations per step 34 00:02:24,490 --> 00:02:27,060 to get order five. 35 00:02:27,060 --> 00:02:32,760 You can't get order five with just five function evaluations. 36 00:02:32,760 --> 00:02:36,560 And then, if we were to try and achieve higher order, 37 00:02:36,560 --> 00:02:41,690 it would take even more function evaluations per step. 38 00:02:41,690 --> 00:02:47,136 Let's use ODE45 to compute e to the t. 39 00:02:47,136 --> 00:02:50,910 y prime is equal to y. 40 00:02:50,910 --> 00:02:54,660 We can ask for output by supplying 41 00:02:54,660 --> 00:02:58,700 an argument called tspan. 42 00:02:58,700 --> 00:03:02,710 Zero and steps of 0.1 to 1. 43 00:03:02,710 --> 00:03:05,750 If we supply that as the input argument 44 00:03:05,750 --> 00:03:07,700 to solve this differential equation 45 00:03:07,700 --> 00:03:11,700 and get the output at those points, 46 00:03:11,700 --> 00:03:14,220 we get that back as the output. 47 00:03:14,220 --> 00:03:17,860 And now here's the approximations to the solution 48 00:03:17,860 --> 00:03:21,370 to that differential equation at those points. 49 00:03:21,370 --> 00:03:30,500 If we plot it, here's the solution at those points. 50 00:03:30,500 --> 00:03:33,700 And to see how accurate it is, we 51 00:03:33,700 --> 00:03:42,140 see that we're actually getting this result to nine digits. 52 00:03:42,140 --> 00:03:45,435 ODE45 is very accurate. 53 00:03:49,740 --> 00:03:52,520 Let's look at step size choice on our problem 54 00:03:52,520 --> 00:03:57,300 with near singularity, is a quarter. 55 00:03:57,300 --> 00:04:01,270 y0 is close to 16. 56 00:04:01,270 --> 00:04:08,390 The differential equation is y prime is 2(a-t) y squared. 57 00:04:08,390 --> 00:04:15,230 We let ODE45 choose its own step size by indicating we just 58 00:04:15,230 --> 00:04:19,970 want to integrate from 0 to 1. 59 00:04:19,970 --> 00:04:26,000 We capture the output in t and y and plot it. 60 00:04:29,870 --> 00:04:34,360 Now, here, there's a lot of points here, 61 00:04:34,360 --> 00:04:41,310 but this is misleading because ODE45, by default, 62 00:04:41,310 --> 00:04:44,150 is using the refine option. 63 00:04:44,150 --> 00:04:47,530 It's only actually evaluating the function 64 00:04:47,530 --> 00:04:50,440 at every fourth one of these points 65 00:04:50,440 --> 00:04:56,880 and then using the interpolant to fill in in between. 66 00:04:56,880 --> 00:05:03,380 So we actually need a little different plot here. 67 00:05:06,380 --> 00:05:09,240 This plot shows a little better what's going on. 68 00:05:09,240 --> 00:05:13,820 The big dots are the points that ODE45 69 00:05:13,820 --> 00:05:18,180 chose to evaluate the differential equation. 70 00:05:18,180 --> 00:05:23,050 And the little dots are filled in with the interpolant. 71 00:05:23,050 --> 00:05:26,570 So the big dots are every fourth point. 72 00:05:26,570 --> 00:05:33,740 And the refine option says that the big dots are too far apart 73 00:05:33,740 --> 00:05:36,100 and we need to fill it in with the interpolant. 74 00:05:36,100 --> 00:05:41,060 And so this is the continuous interpolant in action. 75 00:05:41,060 --> 00:05:46,600 The big dots are more closely concentrated 76 00:05:46,600 --> 00:05:49,010 as we have to go around the curve. 77 00:05:49,010 --> 00:05:55,430 And then, as we get farther away from the singularity 78 00:05:55,430 --> 00:05:57,970 the step size increases. 79 00:05:57,970 --> 00:06:04,280 So this shows the high accuracy of ODE45 80 00:06:04,280 --> 00:06:07,745 and the automatic step size choice in action. 81 00:06:11,750 --> 00:06:13,870 Here's an exercise. 82 00:06:13,870 --> 00:06:21,700 Compare ODE23 and ODE45 by using each of them to compute pi. 83 00:06:21,700 --> 00:06:28,530 The integral 4 over 1 plus t squared from 0 to 1 is pi. 84 00:06:28,530 --> 00:06:32,700 You can express that as a differential equation, 85 00:06:32,700 --> 00:06:36,040 use each of the routines to integrate that differential 86 00:06:36,040 --> 00:06:40,900 equation and see how close they get to computing pi.