Numerical Approximation of Differential Equations with Midpoint Method

Question

I want to proof that the local truncation error of the Midpoint Method is

$d_{k+1}=O\left(h^{3}\right)$

---

Approach

The local truncation error is defined as:

$d_{k+1}=u(t_{k+1})-\underset{u_{k+1}}{\underbrace{\left[u(t_{k})+\triangle t\cdot\phi\right]}}$

• $u(t_{k+1})\rightarrow$ exact solution of the differential equation $u'(t)=f(u,t)$, which is being approximated
• $\Phi\rightarrow f\left(u(t_{k+1})+\frac{1}{2}\left(\triangle t\right)\cdot f\left(u(t_{k}),t_{k}\right),t_{k}+\frac{1}{2}\left(\triangle t\right)\right)$ incremental Function called the Midpoint Method

The first step is to write the term $u(t_{k+1})$ as its Taylor expansion: $u(t_{k+1})=u(t_{k})+\frac{1}{1!}\cdot h^{1}\cdot u^{(1)}(t_{k})+\frac{1}{2!}\cdot h^{2}\cdot u^{(2)}(t_{k})+\frac{1}{3!}\cdot h^{3}\cdot u^{(3)}(t_{k})+O\left(h^{4}\right)$

Then the second step is to write the incremental function in a fashion that the elements can cancel out with the previously introduced Taylor expansion of $u(t_{k+1})$.

I found a solution to the second step, but i do not understand how it is formed. And i would be happy if someone could explain in detail and write out the rules which lead to it. The solution goes as follows:

$\Phi\rightarrow f\left(u(t_{k+1})+\frac{1}{2}\left(\triangle t\right)\cdot f\left(u(t_{k}),t_{k}\right),t_{k}+\frac{1}{2}\left(\triangle t\right)\right)$

What happened here?

$=f\left(u\left(t_{k}\right),t_{k}\right)+\frac{1}{2}\cdot h\cdot f_{t}\left(u\left(t_{k}\right),t_{k}\right)+\frac{1}{2}\cdot h\cdot f\left(u\left(t_{k}\right),t_{k}\right)\cdot f_{y}\left(u\left(t_{k}\right),t_{k}\right)$

$+\frac{1}{2}\cdot\left(\frac{1}{2}\cdot h\right)^{2}\cdot f_{tt}+\left(\frac{1}{2}\cdot h\right)^{2}f\cdot f_{tu}+\frac{1}{2}\left(\frac{1}{2}\cdot h\right)^{2}\cdot f^{2}f_{uu}+O\left(h^{3}\right)$

I understand that there is some kind of differentiating going on, but i cant figure it out. Can someone provide the steps involved?

---

Edits

I found out that total derivatives are used somehow, but i do not have written out the exact transformation steps.

Answer 1

It is easier to do this for an autonomous DE and then only for the scalar case, the error terms for the other cases have a similar structure. Using $u'(t)=f(u(t))$ and thus $u''(t)=f'(u(t))u'(t)=f'(u(t))f(u(t))$ etc. gives \begin{align} u(t+h)&=u(t)+u'(t)·h+\frac12·u''(t)·h^2+\frac16·u'''(t)·h^3+\frac1{24}·u^{(4)}(t)·h^4+… \\ &=u(t)+f·h+\frac12·f'·f·h^2+\frac16·(f''·f^2+f'^2·f)·h^3 \\ &\qquad\qquad+\frac1{24}·(f'''·f^3+4·f''f'f^2+f'^3·f)·h^4+… \end{align} and for the improved Euler or midpoint method at $u=u(t)$ \begin{align} u_+&=u+h·f\left(u+\frac h2·f\right) \\ &=u+h·\left(f+f'·\frac h2·f+\frac12·f''·\left(\frac h2·f\right)^2+\frac16·f'''·\left(\frac h2·f\right)^3\right)+… \\ &=u+f·h+\frac12·f'·f·h^2+\frac18·f''·f^2·h^3+\frac1{48}·f'''·f^3·h^4+… \end{align} where we see several differences in the coefficient of $h^3$.

---

See also the (Heun-)Ralston method. https://en.wikipedia.org/wiki/Runge%E2%80%93Kutta_methods#Usage

---

For a comprehensive methodology to treat the vector case, see Butcher's rooted trees as presented in the introduction here

---

In the vector case, the terms in the expansion of the midpoint step all have the form $f^{(k)}(u)[f(u),…,f(u)]_{k\text{ repetitions}}$. The terms in the Taylor series get more structure as the derivatives now are tensors, $f'$ is a linear operator, $f''$ a vector valued symmetric 2-form, $f'''$ a vector-valued symmetric 3-form etc. \begin{align} u'&=f(u)\\ u''&=f'(u)[f(u)]\\ u'''&=f''(u)[f(u),f(u)]+f'(u)[f'(u)[f(u)]]\\ u^{(4)}&=f'''(u)[f(u),f(u),f(u)]+3·f''(u)[f'(u)[f(u)],f(u)]+f'(u)[f''(u)[f(u),f(u)]]+f'(u)[f'(u)[f'(u)[f(u)]]] \end{align} From this one can see why rooted trees are a good idea to describe the single terms.

---

As to the last point of the question, this is a Taylor expansion that IMO hides more of the structure of the problem than helps to expose it. The structure is $$ u_+=u+h·f(u+h·a,t+h·b),\text{ where } a=\frac12·f(u,t)\text{ and }b=\frac12 $$ For the error analysis, the first $f$ gets replaced by its Taylor expansion in $(u,t)$. Note that here $a,b$ are the variables, $u,t,f(u,t)$ are constants. Thus $$ u_+=u+h·\Biggl(f(u,t)+h·\Bigl(f_u(u,t)·a+f_t(u,t)·b\Bigr)\\ +\frac{h^2}{2}·\Bigl(f_{uu}(u,t)·a^2+2f_{ut}(u,t)·ab+f_{tt}(u,t)·b^2\Bigr)\\ +\frac{h^3}6·\Bigl( f_{uuu}(u,t)·a^3+3f_{uut}(u,t)·a^2b+3f_{utt}(u,t)·ab^2+f_{ttt}(u,t)·b^3\Bigr) +…\Biggr) $$ Now insert the terms for $a$ and $b$ to get the full form of the expansion.

Answer 2

Too much going on. Keep it simple to not get lost. Just look at the Taylor expansion about the midpoint. Then since one side is larger than the midpoint and one side is smaller, each with a distance of $\frac{\Delta t}{2}$, you get

$$ u(t+\Delta t) = u(t+\frac{\Delta t}{2}) + (\frac{\Delta t}{2})u'(t+ \frac{\Delta t}{2}) + (\frac{\Delta t}{2})^2 u''(t + \frac{\Delta t}{2}) + \mathcal{O}(\Delta t)^3 $$

Here I use shorthand $f(t)=f(u(t),t)$ and plugging in $u' = f$ you get

$$u(t+\Delta t) = u(t+\frac{\Delta t}{2}) + (\frac{\Delta t}{2})f(t+ \frac{\Delta t}{2}) + (\frac{\Delta t}{2})^2 f'(t + \frac{\Delta t}{2}) + \mathcal{O}(\Delta t)^3 $$

Call that equation 1. The same steps for the other side gives

$$u(t) = u(t+\frac{\Delta t}{2}) - (\frac{\Delta t}{2})f(t + \frac{\Delta t}{2}) + (\frac{\Delta t}{2})^2 f'(t + \frac{\Delta t}{2}) + \mathcal{O}(\Delta t)^3 $$

Call that equation 2. This was the only hard part. Now subtract equation 1 from equation 2 to make equation 3:

$$u(t+\Delta t) - u(t) = (\Delta t) f(t + \frac{\Delta t}{2}) + \mathcal{O}(\Delta t)^3 $$

What this tells us is that if you know the solution at the midpoint and use it to calculate the next step then you have an error of $\mathcal{O}(\Delta t)^3$ which is what we want. But how do we find out $u$ at the midpoint, i.e. $u(t+\frac{\Delta t}{2})$? Add equation 1 to equation 2 and divide by 2:

$$\frac{u(t+\Delta t) +u(t)}{2} = u(t+\frac{\Delta t}{2}) + \mathcal{O}(\Delta t)^2$$

Notice the $\mathcal{O}(\Delta t)$ term cancelled. Plugging this into $f$ as our approximation for $u$ at the midpoint, we get

$$ f(t+\frac{\Delta t}{2},u(t+\frac{\Delta t}{2})) = f(t+\frac{\Delta t}{2},\frac{u(t+\Delta t) +u(t)}{2}) + \mathcal{O}(\Delta t)^2 $$

Notice that this is an order 2 approximation to $f$ at the midpoint, but since $f$ is multiplied by $\Delta t$, we get that the error is order 3, that is by substitution into equation 3 we get:

$$u(t+\Delta t) - u(t) = (\Delta t) f(t+\frac{\Delta t}{2},\frac{u(t+\Delta t) +u(t)}{2}) + \mathcal{O}(\Delta t)^3 $$

which is the midpoint method with the error. I hope this illuminates "why" it works: the Taylor expansion at the center makes the order 1 term flip signs and cancel. This trick is used a lot in numerical methods.