Differential Equations Day 18 -- Runge-Kutta Numerical Method for Systems and Second Order Differential Equations

Overview
In the last lesson we saw how to use Euler's Method to numerically estimate solutions of systemifications of second order differential equations. [br][br]In this lesson we'll see the [b]Runge Kutta (4) Method[/b], or RK4 for short, another numerical tool for estimating solutions, but which performs considerably better than Euler's Method. [br][br]The RK4 method is actually the [i]fourth[/i] member of a class of numerical methods. The first member of the class is actually Euler's Method! [br][br]We won't cover the theory of RK4 in these notes. Instead, in the rest of this lesson we'll study what RK4 is, and how RK4 works in a GeoGebra applet. That said, know that the mathematical practice of RK4 is motivated by the same intuitive principle as Euler's Method: [br][br][quote]The slope field of a differential equation can be used as a guide to move the estimated solution forward. [/quote]Euler's Method implements this intuitive principle in a "linear" fashion. Euler's Method observes the slope field [i]once[/i], and uses this one piece of information to head off in that direction [i]linearly[/i] for a predetermined step size in the independent variable (that predetermined step size was referred to as [code]h[/code] in previous lessons). [br][br]On the other hand, RK4 systematically observes the slope field [i]four times[/i], and uses these four pieces of information to head off in a direction that represents motion of [i]a fourth order polynomial approximation[/i] of the solution function. RK4 also advances for a predetermined step size in the independent variable.[br][br]The details of how exactly to do this in a systematic fashion were worked out at the dawn of the 20th century by German mathematicians Carl Runge and Wilhelm Kutta. It's amazing to think of them doing this work prior to the advent of the computer.
The RK4 Method
Here are the steps for implementing RK4 on a system of first order differential equations[br][br][math]x_1'=f_1(x_1,x_2,t)[/math][br][math]x_2'=f_2(x_1,x_2,t)[/math][br][br]with initial conditions [math]x_1\left(t_0\right)=x_{1_0}[/math] and [math]x_2\left(t_0\right)=x_{2_0}[/math].[br][br]Note that the only difference between RK4 and Euler's method is in step 5. The RK4 formulas are about 4 times as complex as the Euler's Method formulas. [list=1][*]Put the system of first order equations in [b]slope field form [/b](also sometimes called [b]normal form[/b])and identify [math]f_1[/math] and [math]f_2[/math]. What this means in practice is that the equations are as written above with [math]x_i'[/math] on the left, and everything else involving [math]x_1[/math], [math]x_2[/math] and [math]t[/math] on the right. [/*][*]Select a [b]step count[/b], [math]n[/math]. This is how many steps of an estimated solution you need.[/*][*]Select a [b]step size[/b] [math]h[/math] (sometimes also called [math]dt[/math] or [math]\Delta t[/math])[/*][*]The [b]first estimate point[/b] is [math](x_1_{_0},x_{2_0})[/math], the initial conditions on [math]x_1[/math] and [math]x_2[/math]. The [b]first time point [/b]is [math]t_0[/math], the initial time (it's often the case that [math]t_0=0[/math]).[br][/*][*]Apply the [b]RK4 formulas[/b] to calculate [b]the next estimate point [/b]and[b] the next time point [/b](see below for the RK4 formulas)[/*][*][b]Continue [/b]for [math]n[/math] steps.[/*][/list]The RK4 formulas for step 5 are (calculate them in the order they appear, and use results from previous formulas in later formulas):
[math]t_{n+1}=t_n+h[/math][br][br][math]k_{1,1,n+1}=f_1\left(x_{1_n},x_{2_n},t_n\right)\cdot h[/math][br][math]k_{1,2,n+1}=f_2\left(x_{1_n},x_{2_n},t_n\right)\cdot h[/math][br][br][math]k_{2,1,n+1}=f_1\left(x_{1_n}+\frac{k_{1,1,n+1}}{2},x_{2_n}+\frac{k_{1,2,n+1}}{2},t_n+\frac{h}{2}\right)\cdot h[/math][br][math]k_{2,2,n+1}=f_2\left(x_{1_n}+\frac{k_{1,1,n+1}}{2},x_{2_n}+\frac{k_{1,2,n+1}}{2},t_n+\frac{h}{2}\right)\cdot h[/math][br][br][math]k_{3,1,n+1}=f_1\left(x_{1_n}+\frac{k_{2,1,n+1}}{2},x_{2_n}+\frac{k_{2,2,n+1}}{2},t_n+\frac{h}{2}\right)\cdot h[/math][br][math]k_{3,2,n+1}=f_2\left(x_{1_n}+\frac{k_{2,1,n_1}}{2},x_{2_n}+\frac{k_{2,2,n+1}}{2},t_n+\frac{h}{2}\right)\cdot h[/math][br][br][math]k_{4,1,n+1}=f_1\left(x_{1_n}+k_{3,1,n+1},x_{2_n}+k_{3,2,n+1},t_n+h\right)\cdot h[/math][br][math]k_{4,2,n+1}=f_2\left(x_{1_n}+k_{3,1,n+1},x_{2_n}+k_{3,2,n+1},t_n+h\right)\cdot h[/math][br][br][math]x_{1_{n+1}}=x_{1_n}+\frac{1}{6}\left(k_{1,1,n}+2k_{2,1,n}+2k_{3,1,n}+k_{4,1,n}\right)[/math][br][math]x_{2_{n+1}}=x_{2_n}+\frac{1}{6}\left(k_{1,2,n}+2k_{2,2,n}+2k_{3,2,n}+k_{4,2,n}\right)[/math][br][br]Because of the complexities of these RK4 formulas, it's best to think of them as describing a computer program: Job 1 is to compute the next time point, [math]t_{n+1}[/math]. Job 2 is to compute the next pair of intermediate point [math]k_{1,1,n+1}[/math] and [math]k_{1,2,n+1}[/math]. Keep performing the jobs until you get to the final job which is to compute the next estimate point [math]\left(x_{1_{n+1}},x_{2_{n+1}}\right)[/math]. Note that the formulas for the next estimate point can be thought of as a weighted average of the intermediate [math]k_{i,j,n+1}[/math] points from previous jobs.[br][br]I'm not going to lie: the RK4 formulas are not very fun to work with by hand. The good news is that a computer has no trouble with them, and they can be coded into any spreadsheet software, including GeoGebra. In fact, it takes just about as much time to code RK4 into a spreadsheet as it does to implement one full pass through the RK4 formulas. [br][br]This lesson won't go over how precisely to program RK4 into GeoGebra's spreadsheet window, but it's extremely similar to how we saw [url=https://www.geogebra.org/m/cxgtwkqa#material/y9heug3r]how to program Euler's Method for first order equations in GeoGebra earlier[/url]. The only difference is that there are now quite a few more equations than there were then. If you want to look over the details of what I've done to code RK4 in GeoGebra, click the three dots in the top right of this screen, select "Details" and then follow the instructions to download the raw GeoGebra (.ggb) file.
Implementation in GeoGebra
I've implemented RK4 in GeoGebra' spreadsheet window below. [br][br]The [b]bolded column[/b] is the column of estimates of the solution of [code]x_1[/code]. If you are using this tool to estimate solutions of second order differential equations in [math]y(t)[/math], then this bold column represents the estimates of the solution function you seek. [br][br]To interact with the calculator, adjust [code]f_2[/code] to introduce a different second order differential equation. You can also adjust the initial conditions [code]x_{10}[/code] and [code]x_{20}[/code], and the step size, [code]h[/code]. If you want more steps, unfortunately there's no way to directly change [code]n[/code] due to a limitation in GeoGebra (or at least my ability to write code in GeoGebra). [br][br]Instead, to increase the number of steps, select the entire final row of the spreadsheet, and then use the fill handle to drag it down further. Each new row represents one more step.
At the outset, the calculator is set up to use RK4 to compute estimates of the solution of the system [br][br][math]x_1'=x_2[/math][br][math]x_2'=4x_1+3x_2+2\sin(t)[/math][br][br]which we've seen before, and is the systemification of the second order differential equation[br][br][math]y''-3y'-4y=2\sin\left(t\right);y\left(0\right)=2;y'\left(0\right)=3[/math][br][br]The calculator is set up to calculate in steps of [code]h=0.1[/code] out to[code] t=1[/code]. Note the incredible accuracy. The algebraic solution of the second order differential equation (which we saw in the previous lesson) is[br][br][math]y(t)=\frac{4}{5}e^{\left(-1t\right)}+\frac{87}{85}e^{\left(4t\right)}-\frac{5}{17}\sin\left(t\right)+\frac{3}{17}\cos\left(t\right)[/math][br][br]Evaluating this algebraic solution at [math]t=1[/math] yields approximately 56.025 (rounded to three decimal places). In other words: RK4 is only off by about 3 hundredths after 10 steps! In the previous lesson, Euler's Method was off by almost 3 tenths after just two steps. [b][i]Incredible[/i][/b]![br][br]On the next page is a copy of this calculator for quick and each access.

Information: Differential Equations Day 18 -- Runge-Kutta Numerical Method for Systems and Second Order Differential Equations