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

Problem 4

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

Problem 4a

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

ononon

Solution 4a

Fixed point iteration

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

$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 $y_{n}^{j}\!\,$ since $\alpha ,h,t_{n},g_{n-1}\!\,$ are fixed when solving for the fixed point $y_{n}^{*}\!\,$ where $\phi (y_{n}^{*})=y_{n}^{*}\!\,$

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

Conditions for local convergence

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

$\|D(\phi )\|<1\!\,$

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

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

Problem 4b

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

Solution 4b

Newton iteration

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

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

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

Conditions for local convergence

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

Note that

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

Conditions for quadratic convergence

If $D(\psi (y))\!\,$ is Lipschitz, then we have quadratic convergence and $\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: $y'=f(t,y)\!\,$

Problem 5a

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

Solution 5a

Trapezoid method (Implicit, Adams-Moulton)

$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:

${\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 $y_{i}\approx y(t_{i})\!\,$. The uniform step size is $h\!\,$. Hence,

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

Therefore, the given equation may be written as

$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 $t_{n}\!\,$, we get

${\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 $t_{n}\!\,$ gives

${\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 $h\!\,$ do not agree (${\frac {y'''(t_{n})}{6}}h^{3}\neq {\frac {y'''(t_{n})}{4}}h^{3}\!\,$), the error is of order $h^{2}\!\,$.

Problem 5b

Show that the truncation error for the following multistep method is of the same order as in (a): $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 $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 ${\begin{aligned}t_{n+1}&=t_{n}+h\\t_{n-1}&=t_{n}-h\end{aligned}}\!\,$

${\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 $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 $L(u)=-u''+bu=f,\;\;x\in I\equiv [0,1],\;\;u(0)=u(1)=0\qquad (2)\!\,$ where $b\geq 0\!\,$ is constant. Let $\{0=x_{0}<x_{1}<\dots <x_{n}=1\}\!\,$ be a uniform meshsize $h\!\,$. Let $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 $u_{h}=R_{h}(u)\!\,$ be the corresponding finite element solution of (2). Note that $R_{h}\!\,$ is a projection operator, the Ritz projector, onto the finite dimensional space $V_{h}\!\,$ with respect to the element scalar product $a(\cdot ,\cdot )\!\,$ induced by problem (2).

Problem 6a

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

Solution 6a

Weak Form

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

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

Specifically,

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

Discretized Form (Finite Element Formulation)

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

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

Specifically,

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

Equate Both Sides and Apply Inequalities

${\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,

${\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 $I_{h}u\in V_{h}\!\,$ is the Lagrange interpolant of $u\!\,$, then prove $R_{h}(I_{h}u)=I_{h}u\!\,$. Deduce $|u_{h}-I_{h}u|_{1}\leq \Lambda |u-I_{h}u|_{1}\!\,$

Solution 6b

Prove equality

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

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

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

$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 $u_{h},v_{h}\in V_{h}\!\,$

$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

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

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

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

By ellipticity,

$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

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

i.e.

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

Deduce inequality

Hence we have

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

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

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

Problem 6c

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

Solution 6c

Show inequality

${\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 $x\in [x_{i-1},x_{i}]\!\,$, Newton's polynomial interpolation error gives for some $\xi _{i}\in [x_{i-1},x_{i}]\!\,$

${\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

${\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

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

$u\!\,$ needs to be twice differentiable.