Numerical Analysis/ODE in vector form Exercises

All of the standard methods for solving ordinary differential equations are intended for first order equations. When you need to solve a higher order differential equation, you first convert it to a system of first order of equations. Then you rewrite as a vector form and solve this ODE using a standard method. On this page we demonstrate how to convert to a system of equations and then apply standard methods in vector form.

Reduction to a first order system
(Based on Reduction of Order and Converting a general higher order equation.)

I want to show how to convert higher order differential equation to a system of the first order differential equation. Any differential equation of order n of the form
 * $$f\left(t, u, u', u'',\ \cdots,\ u^{(n-1)}\right) = u^{(n)}$$

can be written as a system of n first-order differential equations by defining a new family of unknown functions
 * $$y_i = u^{(i-1)}\quad\text{for}\quad i = 1, 2,... n\,.$$

The n-dimensional system of first-order coupled differential equations is then
 * $$\begin{array}{rclcl}

y_1&=&u\\ y_2&=&u'\\ y_3&=&u''\\ &\vdots&\\ y_n&=&u^{(n-1)}.\\ \end{array} $$ Differentiating both sides yields
 * $$\begin{array}{rclclcl}

y_1'&=&u' &=&y_2\\ y_2'&=&u''&=&y_3\\ y_3'&=&u'''&=&y_4\\ &\vdots&\\ y_{n}'&=&u^{(n)}&=&f(t,y_1,\cdots,y_n).\\ \end{array} $$

We can express this more compactly in vector form
 * $$\mathbf{y}'=\mathbf{f}(t,\mathbf{y})$$

where $$\ y_{i+1} = f_i\left(t, \mathbf{y} \right)$$ for $$i < n$$ and $$\ f_n\left(t, \mathbf{y} \right)$$ = $$\ f\left(t, y_1, y_2, \cdots,y_n \right)\,.$$

Exercise
Consider the second order differential equation $$\ u''+u =0$$ with initial conditions $$\ u{(0)}=1 $$ and $$\ u'{(0)}=0 $$. We will use two steps with step size $$\ h = \frac{\pi}{8} $$ and approximate the values of $$\ u{(\frac{\pi}{4})}$$ and $$\ u'{(\frac{\pi}{4})}.$$

Since the exact solution is $$u(t)=\cos(t)$$ we have $$\ u{(\frac{\pi}{4})}=0.707106781$$ and $$u'{(\frac{\pi}{4})}=-0.707106781$$.

Exercise 1: Convert this second order differential equation to a system of first order equations.
Solution: We have second order differential equation $$\ u''{(t)}= -u{(t)}$$  with initial conditions $$\ u{(0)}=1 $$ and $$\ u'{(0)}=0 $$.

Let $$y_1 = u$$ and  $$y_2 = u'$$. Differentiating $$y_1$$ and $$y_2$$ gives
 * $$\ y_1' = u' = y_2 $$,
 * $$\ y_2' = u'' = -u =-y_1.$$

Thus we have a system of first order equation in vector form
 * $$\left[\begin{array}{c}

y_1' \\   y_2' \end{array}\right] = \left[\begin{array}{c} y_2 \\ -y_1 \end{array} \right] = f\left(t, \left[\begin{array}{c}     y_1 \\ y_2    \end{array} \right] \right)$$

Exercise 2: Apply the Euler method twice.
Solution: By the Euler's method, $$y_{n+1} = y_n + h f\left(t_n, y_n \right)$$.

I. First apply to get y1

 * $$y_1 = y_0 + h f\left(t_0, y_0 \right)$$.

Now we convert $$\ y_1$$ as a vector form, where $$\ t_0 = 0 ,h=\frac{\pi}{8} $$ and $$\left[\begin{array}{c} y_{0,1} \\   y_{o,2} \end{array}\right] = \left[\begin{array}{c} 1 \\ 0   \end{array} \right] $$, then we have
 * $$\begin{align}

\left[\begin{array}{c} y_{1,1} \\   y_{1,2} \end{array}\right] & = \left[\begin{array}{c} y_{0,1} \\ y_{0,2} \end{array} \right] + h f\left(t_0, \left[\begin{array}{c}     y_{0,1} \\    y_{0,2}    \end{array}\right]\right) \\   & =  \left[\begin{array}{c} y_{0,1} \\ y_{0,2} \end{array} \right] + h\left[\begin{array}{c} y_{0,2} \\  -y_{0,1} \end{array}\right] \\ & = \left[\begin{array}{c} 1 \\ 0   \end{array} \right] + h \left[\begin{array}{c} 0 \\   -1 \end{array}\right] \\ & = \left[\begin{array}{c} {1+0} \\{ 0-h} \end{array} \right] \\ & = \left[\begin{array}{c} 1 \\ {-h} \end{array} \right] \\ & = \left[\begin{array}{c} 1 \\ {- \frac{\pi}{8}} \end{array} \right] \\ & = \left[\begin{array}{c} 1 \\ {-0.392699082}   \end{array} \right]\,. \end{align}$$

2. Second apply to get y2

 * $$y_2 = y_1 + h f\left(t_1, y_1 \right)$$.

Similarly, we convert $$y_2$$ as a vector form, where $$t_1 = t_0 + h =0+h =h$$, where $$h=\frac{\pi}{8} $$ and $$\left[\begin{array}{c} y_{1,1} \\   y_{1,2} \end{array}\right] = \left[\begin{array}{c} 1 \\ -h \end{array} \right] $$, then we have
 * $$\begin{align}

\left[\begin{array}{c} y_{2,1} \\   y_{2,2} \end{array}\right] & = \left[\begin{array}{c} y_{1,1} \\ y_{1,2} \end{array} \right] + h f\left(t_1, \left[\begin{array}{c}     y_{1,1} \\    y_{1,2}    \end{array}\right]\right) \\   & =  \left[\begin{array}{c} y_{1,1} \\ y_{1,2} \end{array} \right] + h\left[\begin{array}{c} y_{2,1} \\  -y_{1,1} \end{array}\right] \\ & = \left[\begin{array}{c} 1 \\ {-h} \end{array} \right] + h \left[\begin{array}{c} { -h} \\ {  -1} \end{array}\right] \\ & = \left[\begin{array}{c} {1-h^2} \\{ -h-h} \end{array} \right] \\ & = \left[\begin{array}{c} {1-h^2} \\ -2h \end{array} \right] \\ & = \left[\begin{array}{c} {1 - \frac{\pi^2}{8^2}}\\ {- 2 \frac{\pi}{8}} \end{array} \right] \\ & = \left[\begin{array}{c} 0.8457874321 \\ {-0.785398163}   \end{array} \right]\,. \end{align}$$

Exercise 3: Apply the Backward Euler method twice.
Solution: By the Backward Euler's method,$$\ y_{n+1} = y_n + h f\left(t_{n+1}, y_{n+1} \right).$$

I. First apply to get y1

 * $$\ y_1 = y_0 + h f\left(t_1, y_1 \right)$$.

Now we convert $$y_1$$ as a vector form, where $$t_1=t_0+h=0+h=h $$, where $$h=\frac{\pi}{8} $$ and $$\left[\begin{array}{c} y_{0,1} \\   y_{o,2} \end{array}\right] = \left[\begin{array}{c} 1 \\ 0   \end{array} \right] $$, then we have
 * $$\begin{align}

\left[\begin{array}{c} y_{1,1} \\   y_{1,2} \end{array}\right] & = \left[\begin{array}{c} y_{0,1} \\ y_{0,2} \end{array} \right] + h f\left(t_1, \left[\begin{array}{c}     y_{1,1} \\    y_{1,2}    \end{array}\right]\right) \\   & =  \left[\begin{array}{c} y_{0,1} \\ y_{0,2} \end{array} \right] + h\left[\begin{array}{c} y_{1,2} \\  -y_{1,1} \end{array}\right] \\ & = \left[\begin{array}{c} 1 \\ 0   \end{array} \right] + \left[\begin{array}{c} {hy_{1,2}} \\  { -hy_{1,1}} \end{array}\right] \\ & = \left[\begin{array}{c} {1+ hy_{1,2} } \\ { -hy_{1,1} } \end{array} \right]\,. \end{align}$$

Now we have to solve a system of linear equations

\left\{ \begin{array}{l l}   y_{1,1}=1+ hy_{1,2}& \quad \\ y_{1,2}= -hy_{1,1} \quad \\ \end{array} \right. $$ That is,

\left\{ \begin{array}{l l}     y_{1,1}  -hy_{1,2}=1& \quad \\ hy_{1,1}+ y_{1,2}= 0 \quad \\ \end{array} \right. $$ Set up augmented matrix to solve this system,

\left( \begin{array}{cc|c} 1 & -h & 1 \\ h &  1 & 0 \\ \end{array} \right) \Rightarrow \left( \begin{array}{cc|c} 1 & -h & 1 \\ 0 & 1+h^2 & -h \\ \end{array} \right) $$ Thus, the solution is $$\ y_{1,1}=\frac{1}{1+h^2}, y_{1,2}=-\frac{h}{1+h^2}.$$

Plugging in $$\ h=\frac{\pi}{8}$$, then we have
 * $$\left[\begin{array}{c}

y_{1,1} \\   y_{1,2} \end{array}\right] = \left[\begin{array}{c} 0.8457874321 \\ {-0.785398163}   \end{array} \right]\,. $$

2. Second apply to get y2

 * $$\ y_2 = y_1 + h f\left(t_2, y_2 \right)$$.

Similarly, we convert $$\  y_2$$ as a vector form, with $$\ t_2 = t_1 + h =h+h =2h$$, where$$\ h=\frac{\pi}{8} $$ and $$\left[\begin{array}{c} y_{1,1} \\   y_{1,2} \end{array}\right] = \left[\begin{array}{c} {\frac{1}{1+h^2}} \\ {-\frac{h}{1+h^2} } \end{array} \right] $$, then we have
 * $$\begin{align}

\left[\begin{array}{c} y_{2,1} \\   y_{2,2} \end{array}\right] & = \left[\begin{array}{c} y_{1,1} \\ y_{1,2} \end{array} \right] + h f\left(t_2, \left[\begin{array}{c}     y_{2,1} \\    y_{2,2}    \end{array}\right]\right) \\   & =  \left[\begin{array}{c} {\frac{1}{1+h^2}} \\ {-\frac{h}{1+h^2}} \end{array} \right] + h f \left(2h, \left[\begin{array}{c}     y_{2,1} \\   y_{2,2}    \end{array}\right]\right) \\   & =  \left[\begin{array}{c} {\frac{1}{1+h^2}} \\ {-\frac{h}{1+h^2}} \end{array} \right] + h \left[\begin{array}{c} y_{2,2} \\  {- y_{2,1} } \end{array}\right] \\ & = \left[\begin{array}{c} {\frac{1}{1+h^2} + hy_{2,2} } \\ {-\frac{h}{1+h^2} - hy_{2,1} } \end{array} \right]\,. \end{align}$$ Now we have to solve a system of linear equations

\left\{ \begin{array}{l l}   y_{2,1}= \frac{1}{1+h^2} + hy_{2,2}& \quad \\ y_{2,2}= -\frac{h}{1+h^2} - hy_{2,1} \quad \\ \end{array} \right. $$ That is,

\left\{ \begin{array}{l l}     y_{2,1}  -hy_{2,2}=\frac{1}{1+h^2}& \quad \\ hy_{2,1}+ y_{2,2}= -\frac{h}{1+h^2} \quad \\ \end{array} \right. $$ Set up augmented matrix to solve this system,

\left( \begin{array}{cc|c} 1 & -h & \frac{1 }{1+h^2}\\ h &  1 & -\frac{h}{1+h^2} \\ \end{array} \right) \Rightarrow \left( \begin{array}{cc|c} 1 & -h & \frac{1}{1+h^2} \\ 0 & 1+h^2 & -\frac{2h}{1+h^2} \\ \end{array} \right) $$ Thus, the solution is $$\ y_{2,1}=\frac{1-h^2}{(1+h^2)^2}, y_{2,2}=-\frac{2h}{(1+h^2)^2}.$$

Plugging in $$\ h=\frac{\pi}{8}$$, then we have
 * $$\left[\begin{array}{c}

y_{2,1} \\   y_{2,2} \end{array}\right] = \left[\begin{array}{c} 0.63487705 \\ {-0.589546795}   \end{array} \right]\,. $$

Exercise 4: Apply the Midpoint method twice.
Solution: By the Midpoint method, $$y_{n+1} = y_n + h f\left(t_n +\frac{1}{2} h, y_n + \frac{1}{2} h f\left(t_n, y_n\right)\right).$$

I. First apply to get y1

 * $$y_1 = y_0 + h f\left(t_0 +\frac{1}{2} h, y_0 + \frac{1}{2} h f\left(t_0, y_0\right)\right).$$

Now we convert $$\ y_1$$ as a vector form, where $$\ t_0 = 0 ,h=\frac{\pi}{8} $$ and $$\left[\begin{array}{c} y_{0,1} \\   y_{0,2} \end{array}\right] = \left[\begin{array}{c} 1 \\ 0   \end{array} \right] $$, then we have
 * $$\begin{align}

\left[\begin{array}{c} y_{1,1} \\   y_{1,2} \end{array}\right] &= \left[\begin{array}{c} y_{0,1} \\ y_{0,2} \end{array} \right] + h f\left(t_0 +\frac{1}{2} h, \left[\begin{array}{c}     y_{0,1} \\ y_{0,2}   \end{array} \right]    + \frac{1}{2} h f\left(t_0, \left[\begin{array}{c} y_{0,1} \\ y_{0,2} \end{array} \right]  \right)\right) \\ &= \left[\begin{array}{c} y_{0,1} \\ y_{0,2} \end{array} \right] + h f\left(t_0 +\frac{1}{2} h, \left[\begin{array}{c}     y_{0,1} \\ y_{0,2}   \end{array} \right]    + \frac{1}{2} h \left[\begin{array}{c}      y_{0,2} \\ {-y_{0,1}}   \end{array} \right]   \right) \\ &= \left[\begin{array}{c} 1 \\ 0  \end{array} \right] + h f\left(0 +\frac{1}{2} h, \left[\begin{array}{c}     1 \\ 0 \end{array} \right]    + \frac{1}{2} h \left[\begin{array}{c}       0 \\ {-1}   \end{array} \right]   \right) \\ &= \left[\begin{array}{c} 1 \\ 0  \end{array} \right] + h \left[\begin{array}{c} {-\frac{1}{2} h} \\ {-1} \end{array} \right] \\ &= \left[\begin{array}{c} { 1 -\frac{1}{2} h^2}\\{-h} \end{array} \right]\,. \end{align}$$ Plugging in $$\ h=\frac{\pi}{8}$$, then we have
 * $$\left[\begin{array}{c}

y_{1,1} \\   y_{1,2} \end{array}\right] = \left[\begin{array}{c} 0.922893716 \\ {-0.392699082}  \end{array} \right]\,. $$

2. Second apply to get y2

 * $$y_2 = y_1 + h f\left(t_1 +\frac{1}{2} h, y_1 + \frac{1}{2} h f\left(t_1, y_1\right)\right).$$

Now we convert $$\ y_1$$ as a vector form, where $$\ t_1=t_0+h = 0+h=h ,h=\frac{\pi}{8} $$ and $$\left[\begin{array}{c} y_{1,1} \\   y_{1,2} \end{array}\right] = \left[\begin{array}{c} {1-\frac{1}{2}h^2} \\ {-h} \end{array} \right] $$, then we have
 * $$\begin{align}

\left[\begin{array}{c} y_{2,1} \\   y_{2,2} \end{array}\right] &= \left[\begin{array}{c} y_{1,1} \\ y_{1,2} \end{array} \right] + h f\left(t_1 +\frac{1}{2} h, \left[\begin{array}{c}     y_{1,1} \\ y_{1,2}   \end{array} \right]    + \frac{1}{2} h f\left(t_1, \left[\begin{array}{c} y_{1,1} \\ y_{1,2} \end{array} \right]  \right)\right) \\ &= \left[\begin{array}{c} y_{1,1} \\ y_{1,2} \end{array} \right] + h f\left(h +\frac{1}{2} h, \left[\begin{array}{c}     y_{1,1} \\ y_{1,2}   \end{array} \right]    + \frac{1}{2} h \left[\begin{array}{c}      y_{1,2} \\ {-y_{1,1}}   \end{array} \right]   \right) \\ &= \left[\begin{array}{c} {1-\frac{1}{2}h^2} \\ {-h} \end{array} \right] + h f\left(\frac{3}{2} h, \left[\begin{array}{c}    {1-\frac{1}{2}h^2} \\ {-h} \end{array} \right]    + \frac{1}{2} h \left[\begin{array}{c}       {-h}  \\ {-1+\frac{1}{2}h^2}   \end{array} \right]   \right) \\ &= \left[\begin{array}{c} {1-\frac{1}{2}h^2} \\ {-h} \end{array} \right] + h f \left( \frac{3}{2} h, \left[\begin{array}{c}    {1- \frac{1}{2} h^2 - \frac{1}{2} h^2} \\ {-h- \frac{1}{2}h + \frac{1}{4}h^3} \end{array} \right]  \right) \\ &= \left[\begin{array}{c} {1-\frac{1}{2}h^2} \\ {-h} \end{array} \right] + h f\left( \frac{3}{2} h, \left[\begin{array}{c}    {1- h^2} \\ {- \frac{2}{3}h+ \frac{1}{4} h^3} \end{array} \right]  \right) \\ &= \left[\begin{array}{c} {1- \frac{1}{2} h^2} \\ {-h} \end{array} \right] + h \left[\begin{array}{c} {- \frac{2}{3} h+ \frac{1}{4} h^3} \\ {-1+ h^2} \end{array} \right] \\ &= \left[\begin{array}{c} {1-\frac{1}{2}h^2} \\ {-h} \end{array} \right] + \left[\begin{array}{c} {-\frac{2}{3}h^2+\frac{1}{4}h^4} \\ {-h+h^3} \end{array} \right] \\ &= \left[\begin{array}{c} {1- \frac{1}{2}h^2- \frac{2}{3}h^2+ \frac{1}{4}h^4} \\ {-h-h+h^3} \end{array} \right] \\ &= \left[\begin{array}{c} {1- 2h^2+ \frac{1}{4}h^4} \\ {-2h+h^3} \end{array} \right]\,. \end{align}$$ Plugging in $$\ h=\frac{\pi}{8}$$, then we have
 * $$\left[\begin{array}{c}

y_{2,1} \\   y_{2,2} \end{array}\right] = \left[\begin{array}{c} 0.694547552011673 \\ {-0.7248390292562377}  \end{array} \right]\,. $$

Exercise 5: Using the values from the Midpoint method at t = h in exercise3, apply the Two-step Adams-Bashforth method once.
Solution: By the Midpoint method, we found
 * $$\left[\begin{array}{c}

y_{1,1} \\   y_{1,2} \end{array}\right] =\left[\begin{array}{c} {1- \frac{1}{2}h^2} \\ {-h} \end{array} \right]\,. $$

Now apply the Two-step Adams-Bashforth method when $$\ n=0$$, then by formula
 * $$\ y_2 = y_1 + \frac{3}{2}h f\left(t_1, y_1\right) - \frac{1}{2}h f\left(t_0 , y_0\right).$$

Vector form is following by
 * $$\begin{align}

\left[\begin{array}{c} y_{2,1} \\   y_{2,2} \end{array}\right] & = \left[\begin{array}{c} y_{1,1} \\ y_{1,2} \end{array} \right] + \frac{3}{2}h f \left(t_1, \left[\begin{array}{c}     y_{1,1} \\    y_{1,2}    \end{array}\right]\right) - \frac{1}{2}h f \left( t_0, \left[\begin{array}{c}      y_{0,1} \\    y_{0,2}    \end{array}\right]\right) \\ & = \left[\begin{array}{c} {1- \frac{1}{2}h^2 }\\ {-h} \end{array} \right] + \frac{3}{2}h f \left(t_1, \left[\begin{array}{c}     {1- \frac{1}{2}h^2 }\\ {-h}   \end{array}\right]\right) - \frac{1}{2}h f \left( t_0, \left[\begin{array}{c}      1 \\    0    \end{array}\right]\right) \\   & =  \left[\begin{array}{c} {1- \frac{1}{2}h^2} \\ {-h} \end{array} \right] + \frac{3}{2} h \left[\begin{array}{c} {-h}  \\     {-1+ \frac{1}{2}h^2} \end{array}\right] - \frac{1}{2}h   \left[\begin{array}{c} 0 \\   {-1}  \end{array}\right] \\ & = \left[\begin{array}{c} {1-2h^2 } \\ {-2h+ \frac{3}{4}h^3 } \end{array} \right]\,. \end{align}$$

Plugging in $$h=\frac{\pi}{8},$$
 * $$ \left[\begin{array}{c}

y_{2,1} \\  y_{2,2} \end{array}\right] = \left[\begin{array}{c} 0.691574862 \\ {-0.739978812}  \end{array} \right]\,.$$

Reference
http://en.wikipedia.org/wiki/Ordinary_differential_equation

http://www.math.ohiou.edu/courses/math3600/lecture29.pdf

http://www.ohio.edu/people/mohlenka/20131/4600-5600/hw7.pdf