# Numerical Methods Qualification Exam Problems and Solutions (University of Maryland)/Aug07 667

## Problem 4

 Consider the problem of solving a nonlinear system of ODE ${\displaystyle y'=f(t,y)\!\,}$  by an implicit method. The ${\displaystyle n\!\,}$ -th step consists of solving for the unknown ${\displaystyle y\!\,}$  a non-linear algebraic system of the form ${\displaystyle y=\alpha hf(t_{n},y)+g_{n-1}\quad (1)\!\,}$  where ${\displaystyle g_{n-1}\in \mathbb {R} ^{n}\!\,}$  is known, ${\displaystyle \alpha >0\!\,}$  and ${\displaystyle h\!\,}$  is the stepsize. Let ${\displaystyle f\in C^{1}\!\,}$

## Problem 4a

 Write ${\displaystyle (1)\!\,}$  as a fixed point iteration and find conditions in ${\displaystyle h\!\,}$  and ${\displaystyle f(t_{n},y)\!\,}$  that local convergence for this iteration

ononon

## Solution 4a

### Fixed point iteration

Equation ${\displaystyle (1)\!\,}$  is conveniently in fixed point iteration form.

${\displaystyle y_{n}^{j+1}=\underbrace {\alpha hf(t_{n},y_{n}^{j})+g_{n-1}} _{\phi (y_{n}^{j})}\!\,}$

Notice that the right hand side is only a function of ${\displaystyle y_{n}^{j}\!\,}$  since ${\displaystyle \alpha ,h,t_{n},g_{n-1}\!\,}$  are fixed when solving for the fixed point ${\displaystyle y_{n}^{*}\!\,}$  where ${\displaystyle \phi (y_{n}^{*})=y_{n}^{*}\!\,}$

Also note that ${\displaystyle j\!\,}$  is the fixed point iteration index.

### Conditions for local convergence

The fixed point iteration will converge when the norm of the Jacobian of ${\displaystyle \phi \!\,}$  is less than 1 i.e.

${\displaystyle \|D(\phi )\|<1\!\,}$

Since ${\displaystyle \|D(\phi )\|=\|\alpha hD(f)\|\!\,}$ , we equivalently have the condition

${\displaystyle \|D(f)\|<{\frac {1}{\alpha h}}\!\,}$

## Problem 4b

 Write the Newton iteration for ${\displaystyle (1)\!\,}$  and give conditions on ${\displaystyle h\!\,}$  and ${\displaystyle f(t_{n},y)\!\,}$  that guarantee local convergence for this iteration. State precise additional assumptions on ${\displaystyle f\!\,}$  that guarantee quadratic convergence

## Solution 4b

### Newton iteration

The Newton iteration solves ${\displaystyle \psi (y)=0\!\,}$  and the iteration is given by

${\displaystyle y_{i+1}=y_{i}-D(\psi (y_{i}))^{-1}\psi (y_{i})\!\,}$

Let ${\displaystyle \psi (y)=y-\phi (y)\!\,}$

### Conditions for local convergence

If ${\displaystyle D(\psi (y))^{-1}\!\,}$  exists, i.e. ${\displaystyle D(\psi (y))\!\,}$  is invertible or equivalently non-singular, then local convergence is guaranteed.

Note that

${\displaystyle D(\psi (y))=I-D(\phi (y))\!\,}$

If ${\displaystyle D(\psi (y))\!\,}$  is Lipschitz, then we have quadratic convergence and ${\displaystyle \psi (y)\!\,}$  is twice continuously differentiable in a neighborhood of the root

## Problem 5

 This problem is about choosing between a specific single-step and a specific multi-step methods for solving ODE: ${\displaystyle y'=f(t,y)\!\,}$

## Problem 5a

 Write the trapezoid method, define its local truncation error and estimate it.

## Solution 5a

${\displaystyle y_{n+1}=y_{n}+{\tfrac {1}{2}}h{\big (}f(t_{n+1},y_{n+1})+f(t_{n},y_{n}){\big )}\!\,}$

### Define Local Truncation Error

The local truncation error is given as follows:

${\displaystyle {\mbox{error}}=(y(t_{n+1})-y(t_{n}))-{\tfrac {1}{2}}h{\big (}f(t_{n+1},y(t_{n+1}))+f(t_{n},y_{n}){\big )}\!\,}$

### Find Local Truncation Error Using Taylor Expansion

Note that ${\displaystyle y_{i}\approx y(t_{i})\!\,}$ . The uniform step size is ${\displaystyle h\!\,}$ . Hence,

${\displaystyle t_{i+1}=t_{i}+h\!\,}$

Therefore, the given equation may be written as

${\displaystyle y(t_{i}+h)-y(t_{i})\approx {\frac {1}{2}}h(y'(t_{i}+h)+y'(t_{n}))\!\,}$

#### Expand Left Hand Side

Expanding about ${\displaystyle t_{n}\!\,}$ , we get

${\displaystyle {\begin{array}{|c|c|c|c|}{\mbox{Order}}&y(t_{n}+h)&-y(t_{n})&\Sigma \\\hline &&&\\0&y(t_{n})&-y(t_{n})&0\\&&&\\1&y'(t_{n})h&0&y'(t_{n})h\\&&&\\2&{\frac {1}{2!}}y''(t_{n})h^{2}&0&{\frac {1}{2}}y''(t_{n})h^{2}\\&&&\\3&{\frac {1}{3!}}y'''(t_{n})h^{3}&0&{\frac {1}{6}}y'''(t_{n})h^{3}\\&&&\\4&{\mathcal {O}}(h^{4})&0&{\mathcal {O}}(h^{4})\\&&&\\\hline \end{array}}}$

#### Expand Right Hand side

Also expanding about ${\displaystyle t_{n}\!\,}$  gives

${\displaystyle {\begin{array}{|c|c|c|c|}{\mbox{Order }}h&{\frac {1}{2}}hy'(t_{n}+h)&{\frac {1}{2}}hy'(t_{n})&\Sigma \\\hline &&&\\0&0&0&0\\&&&\\1&{\frac {1}{2}}h\cdot y'(t_{n})&{\frac {1}{2}}hy'(t_{n})&y'(t_{n})h\\&&&\\2&{\frac {1}{2}}h\cdot y''(t_{n})h&0&{\frac {1}{2}}y''(t_{n})h^{2}\\&&&\\3&{\frac {1}{2}}h\cdot {\frac {y'''}{2}}(t_{n})h^{2}&0&{\frac {y'''}{4}}(t_{n})h^{3}\\&&&\\4&{\mathcal {O}}(h^{4})&{\mathcal {O}}(h^{4})&{\mathcal {O}}(h^{4})\\\hline \end{array}}}$

### Calculate local truncation error

Since the order 3 terms of ${\displaystyle h\!\,}$  do not agree (${\displaystyle {\frac {y'''(t_{n})}{6}}h^{3}\neq {\frac {y'''(t_{n})}{4}}h^{3}\!\,}$ ), the error is of order ${\displaystyle h^{2}\!\,}$ .

## Problem 5b

 Show that the truncation error for the following multistep method is of the same order as in (a): ${\displaystyle y_{n+1}=2y_{n}-y_{n-1}-hf(t_{n-1},y_{n-1})+hf(t_{n},y_{n})\!\,}$

## Solution 5b

We need to show that ${\displaystyle y(t_{n+1})-2y(t_{n})+y(t_{n-1})+hy'(t_{n-1})-hy'(t_{n})={\mathcal {O}}(h^{3})\!\,}$

Again, note that {\displaystyle {\begin{aligned}t_{n+1}&=t_{n}+h\\t_{n-1}&=t_{n}-h\end{aligned}}\!\,}

${\displaystyle {\begin{array}{|c|c|c|c|c|c|c|}{\mbox{Order of }}h&y(t_{n}+h)&-2y(t_{n})&y(t_{n}-h)&hy'(t_{n}-h)&-hy'(t_{n})&\Sigma \\\hline 0&y(t_{n})&-2y(t_{n})&y(t_{n})&0&0&0\\&&&&&&\\1&y'(t_{n})h&0&-y'(t_{n})h&y'(t_{n})h&-y'(t_{n})&0\\&&&&&&\\2&{\frac {1}{2}}y''(t_{n})h^{2}&0&{\frac {1}{2}}y''(t_{n})h^{2}&-y''(t_{n})h^{2}&0&0\\&&&&&&\\3&{\frac {1}{6}}y'''(t_{n})h^{3}&0&-{\frac {1}{6}}y'''(t_{n})h^{3}&{\frac {1}{2}}y'''(t_{n})h^{3}&0&{\frac {1}{2}}y'''(t_{n})h^{3}\\&&&&&&\\\hline \end{array}}}$

So this method is also consistency order 2.

## Problem 5c

 What could be said about the global convergence rate for these two methods? Justify your conclusions for both methods.

## Solution 5c

The trapezoid is stable because its satisfies the root condition. (The root of the characteristic equation is 1 and has a simple root)

The second method is not stable because the characteristic equation has a double root of 1.

Both the trapezoid method and second method are consistent with order ${\displaystyle h^{2}\!\,}$

Note that convergence occurs if and only the method is both stable and consistent.

Therefore, the trapezoid method converges in general but the second method does not. mkmkmkmlmklml

## Problem 6

 Consider the boundary value problem ${\displaystyle L(u)=-u''+bu=f,\;\;x\in I\equiv [0,1],\;\;u(0)=u(1)=0\qquad (2)\!\,}$  where ${\displaystyle b\geq 0\!\,}$  is constant. Let ${\displaystyle \{0=x_{0}  be a uniform meshsize ${\displaystyle h\!\,}$ . Let ${\displaystyle V_{h}=\{v\in C[0,1]:\left.v\right|_{[x_{i-1},x_{i}]}{\mbox{ is linear for each i, }}v(0)=v(1)=0\}\!\,}$ be the corresponding finite element space, and let ${\displaystyle u_{h}=R_{h}(u)\!\,}$  be the corresponding finite element solution of (2). Note that ${\displaystyle R_{h}\!\,}$  is a projection operator, the Ritz projector, onto the finite dimensional space ${\displaystyle V_{h}\!\,}$  with respect to the element scalar product ${\displaystyle a(\cdot ,\cdot )\!\,}$  induced by problem (2).

## Problem 6a

 Let ${\displaystyle |\cdot |_{1}\!\,}$  be the ${\displaystyle H^{1}\!\,}$ -seminorm, namely ${\displaystyle |v|_{1}^{2}=\int _{I}|v'|^{2}\!\,}$  for all ${\displaystyle v\in H_{0}^{1}(I)\!\,}$ . Find the constant ${\displaystyle \Lambda \!\,}$  in terms of the parameter ${\displaystyle b\!\,}$  such that ${\displaystyle |u_{h}|_{1}\leq \Lambda |u|_{1}\!\,}$  Hint: recall the Poincare inequality ${\displaystyle \|v\|_{0}\leq {\frac {1}{2}}|v|_{1}\!\,}$  for all ${\displaystyle v\in H_{0}^{1}(I)\!\,}$  where ${\displaystyle \|\cdot \|_{0}\!\,}$  denotes the ${\displaystyle L^{2}\!\,}$ -norm

## Solution 6a

### Weak Form

Integrating by parts gives, for all ${\displaystyle v\in V\!\,}$

${\displaystyle \int u'v'+b\int uv=\int fv\!\,}$

Specifically,

${\displaystyle \int u'u_{h}'+b\int uu_{h}=\int fu_{h}\!\,}$

### Discretized Form (Finite Element Formulation)

Similarly, the finite element formulation is find ${\displaystyle u_{h}\in V_{h}\!\,}$  such that for all ${\displaystyle v_{h}\in V_{h}\!\,}$

${\displaystyle \int u_{h}'v_{h}'+b\int u_{h}v_{h}=\int fv_{h}\!\,}$

Specifically,

${\displaystyle \int u_{h}'u_{h}'+b\int u_{h}u_{h}=\int fu_{h}\!\,}$

### Equate Both Sides and Apply Inequalities

{\displaystyle {\begin{aligned}\int _{0}^{1}u_{h}'u_{h}'+b\int _{0}^{1}u_{h}u_{h}&=\int _{0}^{1}u'u_{h}'+b\int _{0}^{1}uu_{h}\\&\leq |u|_{1}|u_{h}|_{1}+b\|u\|_{0}\|u_{h}\|_{0}{\mbox{ Cauchy-Schwartz}}\\&\leq |u|_{1}|u_{h}|_{1}+b{\frac {1}{2}}|u|_{1}{\frac {1}{2}}|u_{h}|_{1}{\mbox{ Poincare}}\\&=|u_{h}|_{1}(1+{\frac {b}{4}})|u|_{1}\end{aligned}}}

Hence we have,

{\displaystyle {\begin{aligned}|u_{h}|_{1}^{2}\leq |u_{h}|_{1}^{2}+b\|u_{h}\|_{0}^{2}&\leq |u_{h}|_{1}(1+{\frac {b}{4}})|u|_{1}\\|u_{h}|_{1}^{2}&\leq |u_{h}|_{1}(1+{\frac {b}{4}})|u|_{1}\\|u_{h}|_{1}&\leq (1+{\frac {b}{4}})|u|_{1}\end{aligned}}}

## Problem 6b

 If ${\displaystyle I_{h}u\in V_{h}\!\,}$  is the Lagrange interpolant of ${\displaystyle u\!\,}$ , then prove ${\displaystyle R_{h}(I_{h}u)=I_{h}u\!\,}$ . Deduce ${\displaystyle |u_{h}-I_{h}u|_{1}\leq \Lambda |u-I_{h}u|_{1}\!\,}$

## Solution 6b

### Prove equality

We have for all ${\displaystyle v\in H_{0}^{1}(0,1)\!\,}$

${\displaystyle a(u,v)=\int u'v'dx+b\int uvdx=\int fv\!\,}$

Specifically, for all ${\displaystyle v_{h}\in V_{h}\!\,}$

${\displaystyle a(u,v_{h})=\int u'v_{h}'dx+b\int uv_{h}dx=\int fv_{h}\quad \quad (1)\!\,}$

The discrete form of the energy scalar product is for all ${\displaystyle u_{h},v_{h}\in V_{h}\!\,}$

${\displaystyle a(u_{h},v_{h})=\int u_{h}'v_{h}'dx+b\int u_{h}v_{h}dx=\int fv_{h}\quad \quad (2)\!\,}$

Subtracting equation (2) from equation (1), we have

${\displaystyle a(u-u_{h},v_{h})=0\!\,}$

Let ${\displaystyle R_{h}(I_{h}u)={\overline {u_{h}}}\in V_{h}\!\,}$ . Note that by hypothesis ${\displaystyle I_{h}u\in V_{h}\!\,}$ . Then,

${\displaystyle a(I_{h}u-{\overline {u_{h}}},I_{h}u-{\overline {u_{h}}})=0\!\,}$

By ellipticity,

${\displaystyle 0=a(I_{h}u-{\overline {u_{h}}},I_{h}u-{\overline {u_{h}}})\geq \alpha \|I_{h}u-{\overline {u_{h}}}\|^{2}\!\,}$

which implies

${\displaystyle I_{h}u-{\overline {u_{h}}}=0\!\,}$

i.e.

${\displaystyle I_{h}u={\overline {u_{h}}}\!\,}$

### Deduce inequality

Hence we have

${\displaystyle R_{h}(u-I_{h}u)=u_{h}-I_{h}u\!\,}$

Arguing as we did in part (a), we have

${\displaystyle |u_{h}-I_{h}u|_{1}\leq \Lambda |u-I_{h}u|_{1}\!\,}$

## Problem 6c

 Use (b) to derive the error estimate ${\displaystyle |u-u_{h}|_{1}\leq (1+\Lambda )|u-I_{h}u|_{1}\!\,}$  and bound the right hand side by suitable power of ${\displaystyle h\!\,}$ . Make explicit the required regularity of ${\displaystyle u\!\,}$

## Solution 6c

### Show inequality

{\displaystyle {\begin{aligned}|u-u_{h}|_{1}&=|u-I_{h}u+I_{h}u-u_{h}|_{1}\\&\leq |u-I_{h}u|_{1}+|I_{h}u-u_{h}|_{1}\\&\leq |u-I_{h}u|_{1}+\Lambda |u-I_{h}u|_{1}\\&=(1+\Lambda )|u-I_{h}u|_{1}\end{aligned}}}

### Bound Right Hand Side

For ${\displaystyle x\in [x_{i-1},x_{i}]\!\,}$ , Newton's polynomial interpolation error gives for some ${\displaystyle \xi _{i}\in [x_{i-1},x_{i}]\!\,}$

{\displaystyle {\begin{aligned}|u-I_{h}u|_{1}^{2}&=\left|{\frac {u^{(2)}(\xi _{i})}{2}}(x-x_{i})(x-x_{i-1})\right|_{1}^{2}dx\\&=\int _{0}^{1}\left(\left[{\frac {u^{(2)}(\xi _{i})}{2}}(x-x_{i})(x-x_{i-1})\right]'\right)^{2}\\&=\left({\frac {u^{(2)}(\xi _{i})}{2}}\right)^{2}\int _{0}^{1}[\underbrace {(x-x_{i-1})} _{\leq h}\underbrace {(x-x_{i})} _{\leq h}]dx\\&\leq (u^{(2)}(\xi _{i})h)^{2}\end{aligned}}\!\,}

Therefore the error on the entire interval is given by

{\displaystyle {\begin{aligned}|u-I_{h}u|_{1}^{2}&\leq \sum _{i=1}^{n}(u^{(2)}(\xi _{i})h)^{2}\\&\leq \max _{0<\xi <1}(u^{(2)}(\xi )h)^{2}\cdot n\end{aligned}}\!\,}

which implies

${\displaystyle |u-I_{h}u|_{1}\leq \max _{0<\xi <1}u^{(2)}(\xi _{i}){\sqrt {n}}h\!\,}$

${\displaystyle u\!\,}$  needs to be twice differentiable.