Numerical Methods Qualification Exam Problems and Solutions (University of Maryland)/August 2006

Problem 1

Consider the definite integral $I(f)=\underset{a}{\overset{b}{\int }}f(x)dx$ Let $Q_n(f)$ denote the approximation of$I(f)$ by the composite midpoint rule with $n$ uniform subintervals. For every $j\in \mathbb{Z}$ set $x_j=a+j\frac{b-a}{n},x_{j+\frac{1}{2}}=\frac{x_j+x_{j+1}}{2}$. Let $K(x)$ be defined by $K(x)=-\frac{1}{2}(x-x_j{)}^2,x_{j-\frac{1}{2}}<x<x_{j+\frac{1}{2}}$. Assume that $f\in C^2([a,b])$.

Problem 1a

Show that the quadrature error $E_n(f)$ satisfies $E_n\equiv Q_n(f)-I(f)=\underset{a}{\overset{b}{\int }}K(x)f^{″}(x)dx$ Hint: Use integration by parts over each subinterval $[x_{j-1},x_j]$.

Solution 1a

Integrating $\int K(x)f^{″}(x)dx$ by parts over arbitrary points $p$ and $q$ gives

$\underset{p}{\overset{q}{\int }}K(x)f^{″}(x)dx={[(x-x_j)f(x)-\frac{1}{2}(x-x_j{)}^2{f}^{\prime }(x)]}_{x=p}^{x=q}-\underset{p}{\overset{q}{\int }}f(x)dx$

Since $K(x)$ is defined on $[x_{j-\frac{1}{2}},x_{j+\frac{1}{2}}]$ we use

$[p,q]=[x_j,x_{j+\frac{1}{2}}]$

and

$[p,q]=[x_{j-\frac{1}{2}},x_j]$

Using the first interval we get

$\frac{1}{2}[-\frac{1}{4}(x_{j+1}-x_j{)}^2{f}^{\prime }(x_{j+\frac{1}{2}})+(x_{j+1}-x_j)f(x_{j+\frac{1}{2}})]\mathbf{(}\mathrm{𝟏}\mathbf{)}$

And for the second we get

$\frac{1}{2}[\frac{1}{4}(x_{j-1}-x_j{)}^2{f}^{\prime }(x_{j-\frac{1}{2}})+(x_j-x_{j-1})f(x_{j-\frac{1}{2}})]\mathbf{(}\mathrm{𝟐}\mathbf{)}$

Since these apply to arbitrary half-subintervals, we can rewrite equation $\mathbf{(}\mathrm{𝟐}\mathbf{)}$ with the its indecies shifted by one unit. The equation for the interval $[x_{j+\frac{1}{2}},x_{j+1}]$ is

$\frac{1}{2}[\frac{1}{4}(x_j-x_{j+1}{)}^2{f}^{\prime }(x_{j+\frac{1}{2}})+(x_{j+1}-x_j)f(x_{j+\frac{1}{2}})]\mathbf{(}{\mathrm{𝟐}}^{\prime }\mathbf{)}$

Combining $\mathbf{(}\mathrm{𝟏}\mathbf{)}$ and $\mathbf{(}{\mathrm{𝟐}}^{\prime }\mathbf{)}$ and writing it in the same form as the integration by parts, we have

$\underset{x_j}{\overset{x_{j+1}}{\int }}K(x)f^{″}(x)dx=(x_{j+1}-x_j)f(x_{j+\frac{1}{2}})-\underset{x_j}{\overset{x_{j+1}}{\int }}f(x)dx$

Then our last step is to sum this over all of our $n$ subintervals, noting that $x_{j+1}-x_j=(b-a)/n$

$\begin{array}{cc}\sum _{j=0}^{n-1}\underset{x_j}{\overset{x_{j+1}}{\int }}K(x)f^{″}(x)dx& =\sum _{j=0}^{n-1}(x_{j+1}-x_j)f(x_{j+\frac{1}{2}})-\sum _{j=0}^{n-1}\underset{x_j}{\overset{x_{j+1}}{\int }}f(x)dx\\ \underset{a}{\overset{b}{\int }}K(x)f^{″}(x)dx& =\frac{b-a}{n}\sum _{j=0}^{n-1}f((x_j+x_{x+1})/2)-\underset{a}{\overset{b}{\int }}f(x)dx\\ \underset{a}{\overset{b}{\int }}K(x)f^{″}(x)dx& =Q_n(f)-I(f)\end{array}$

Problem 1b

Derive a sharp bound on the error of the form $|E_n(f)|\le M_n\Vert f^{″}{\Vert }_{\mathrm{\infty }}$ for every $f\in C^2([a,b])$. Here $\Vert \cdot {\Vert }_{\mathrm{\infty }}$ denotes the maximum norm over $[a,b]$. Recall that the above bound is sharp when the inequality is an equality for some nonzero $f$.

Solution 1b

Applying the result of part (a), the triangle inequality, and pulling out the constant $\Vert f(x){\Vert }_{\mathrm{\infty }}$, we have,

$\begin{array}{cc}|E_n(f)|& =|\underset{a}{\overset{b}{\int }}K(x)f^{\prime \prime }(x)dx|\\ & \le \underset{a}{\overset{b}{\int }}|K(x)||f^{\prime \prime }(x)|dx\\ & \le \Vert f^{\prime \prime }(x){\Vert }_{\mathrm{\infty }}\underset{M_n}{\underbrace{\underset{a}{\overset{b}{\int }}|K(x)|dx}}\end{array}$

$M_n$ is some constant less than infinity since $[a,b]$ is compact and $|K(x)|$ is continuous on each of the finite number of intervals for which it is defined.

The above inequality becomes an equality when

$f^{″}(x)=c,$

where $c$ is any constant.

Problem 2

Consider the (unshifted) $QR-$ method for finding the eigenvalues of an invertible matrix $A\in {ℝ}^{n\times n}$

Problem 2a

Give the algorithm.

Solution 2a

The QR algorithm produces a sequence of similar matrices $\{A_i\}$ whose limit tends towards being upper triangular or nearly upper triangular. This is advantageous since the eigenvalues of an upper triangular matrix lie on its diagonal.

i = 0   
A_1 = A 
while ( error > tolerance  )   
   A_i=Q_i R_i  (QR decomposition/factorization)
   A_{i+1}=R_i Q_i (multiply R and Q, the reverse multiplication)
   i=i+1 (increment)

Problem 2b

Show that each of the matrices $A_n$ generated by this algorithm are orthogonally similar to $A$.

Solution 2b

From the factor step (QR decomposition) of the algorithm, we have

$A_i=Q_i{R}_i$

which implies

$Q_i^{-1}{A}_i=R_i$

Substituting into the reverse multiply step, we have

$\begin{array}{cc}{A}_{i+1}& =R_i{Q}_i\\ & =Q_i^{-1}{A}_i{Q}_i\end{array}$

Problem 2c

Show that if $A$ is upper Hessenberg then so are each of the matrices $A_n$ generated by this algorithm.

Solution 2c

A series of Given's Rotations matrices pre-multiplying $A$, a upper Heisenberg matrix, yield an upper triangular matrix$R$ i.e.

$G_{(n-1,n)}\cdots G_{(2,3)}{G}_{(1,2)}A=R$

Since Givens Rotations matrices are each orthogonal, we can write

$A=\underset{Q}{\underbrace{(G_{(n-1,n)}\cdots G_{(2,3)}{G}_{(1,2)}{)}^T}}R$

i.e.

$A=QR$

If we let $A_0:=A$, we have $A_0=Q_0{R}_0$, or more generally for $k=1,2,3,\dots$

$A_k=Q_k{R}_k$.

In each case, the sequence of Given's Rotations matrices that compose $Q$ have the following structure

$\begin{array}{cc}Q& =G_{(1,2)}^T{G}_{2,3}^T\cdots G_{(n-1,n)}^T\\ & =\left(\begin{array}{cccccc}*& *& 0& 0& \cdots & 0\\ *& *& 0& 0& \cdots & 0\\ 0& 0& 1& 0& \cdots & 0\\ 0& 0& 0& 1& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & \ddots & ⋮\\ 0& 0& \cdots & \cdots & 0& 1\end{array}\right)\left(\begin{array}{cccccc}1& 0& 0& 0& \cdots & 0\\ 0& *& *& 0& \cdots & 0\\ 0& *& *& 0& \cdots & 0\\ 0& 0& 0& 1& \cdots & 0\\ ⋮& ⋮& ⋮& \ddots & \ddots & ⋮\\ 0& 0& \cdots & \cdots & 0& 1\end{array}\right)\cdots \left(\begin{array}{cccccc}1& 0& 0& 0& \cdots & 0\\ 0& 1& 0& 0& \cdots & 0\\ ⋮& \ddots & \ddots & \ddots & \ddots & 0\\ 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& *& *\\ 0& 0& 0& 0& *& *\end{array}\right)\\ & =\left(\begin{array}{ccccc}*& *& *& \cdots & *\\ *& *& *& \cdots & *\\ 0& *& *& \cdots & *\\ ⋮& \ddots & \ddots & \ddots & ⋮\\ 0& 0& 0& *& *\end{array}\right)\end{array}$

So $Q$ is upper Hessenberg.

From the algorithm, we have

$\begin{array}{cc}{A}_{k+1}& =R_k{Q}_k\\ & =\left(\begin{array}{ccccc}*& *& *& \cdots & *\\ 0& *& *& \cdots & *\\ 0& 0& *& \cdots & *\\ ⋮& ⋮& \ddots & \ddots & ⋮\\ 0& 0& 0& 0& *\end{array}\right)\left(\begin{array}{ccccc}*& *& *& \cdots & *\\ *& *& *& \cdots & *\\ 0& *& *& \cdots & *\\ ⋮& \ddots & \ddots & \ddots & ⋮\\ 0& 0& 0& *& *\end{array}\right)\end{array}$

We conclude $A_{k+1}$ is upper Hessenberg because for $j=1,2,\dots ,n-2$ the $j$th column of $A_{k+1}$ is a linear combination of the first $j+1$ columns of $R_k$ since $Q_k$ is also upper Hessenberg.

Problem 2d

Let $\left(\begin{array}{cc}3& 3\\ 1& 5\end{array}\right)$ For this $A$ the sequence $\{A_n\}$ has a limit. Find this limit. Give your reasoning.

Solution 2d

The eigenvalues of $A$ can be computed. They are $6$ and $2$. Furthermore, the result of matrix multiplies in the algorithm shows that the diagonal difference, $(a_{1,2}-a_{2,1})$, is constant for all $i$.

Since the limit of $A$ is an upper triangular matrix with the eigenvalues of $A$ on the diagonal, the limit is

$\left(\begin{array}{cc}6& 2\\ 0& 2\end{array}\right)$

Problem 3

Let $A\in {ℝ}^{N\times N}$ be symmetric and positive definite. Let $b\in {ℝ}^N$. Consider solving $Ax=b$ using the conjugate gradient method. The $n^{th}$ iterate $x^{(n)}$ then satisfies $\Vert x^{(n)}-x{\Vert }_A\le \Vert y-x{\Vert }_A$ for every $y\in x^{(0)}+{𝒦}_n(r^{(0)},A)$, where $\Vert \cdot {\Vert }_A$ denotes the vector A-norm, $r^{(0)}$ is the initial residual, and ${𝒦}_n(r^{(0)},A)=\text{span}\{r^{(0)},A{r}^{(0)},\dots ,A^{n-1}{r}^{(0)}\}$.

Problem 3a

Prove that the error $x^{(n)}-x$ is bounded by $\Vert x^{(n)}-x{\Vert }_A\le \Vert p_n(A){\Vert }_A\Vert x^{(0)}-x{\Vert }_A$, where $p_n(z)$ is any real polynomial of degree $n$ or less that satisfies $p_n(0)=1$. Here $\Vert p_n(A){\Vert }_A$ denotes the matrix norm of $p_n(A)$ that is induced by the vector A-norm.

Solution 3a

We know that $\Vert x^{(n)}-x{\Vert }_A\le \Vert y-x{\Vert }_A$ for every $y\in x^{(0)}+{𝒦}_n(r^{(0)},A)$, so if we can choose a $y\in x^{(0)}+{𝒦}_n(r^{(0)},A)$ such that

$\Vert x^{(n)}-x{\Vert }_A\le \Vert y-x{\Vert }_A\le \Vert p_n(A){\Vert }_A\Vert x^{(0)}-x{\Vert }_A$,

then we can solve part a.

Rewrite r^(0)

First note from definition of $r^{(0)}$

$\begin{array}{cc}{r}^{(0)}& =A{x}^{(0)}-b\\ & =A{x}^{(0)}-Ax\\ & =A(x^{(0)}-x)\end{array}$

Rewrite Krylov space

Therefore, we can rewrite ${𝒦}_n(r^{(0)},A)$ as follows:

$\begin{array}{cc}{𝒦}_n(r^{(0)},A)& ={𝒦}_n(A(x^{(0)}-x),A)\\ & =\text{span}\{A(x^{(0)}-x),A^2(x^{(0)}-x),\dots ,A^n(x^{(0)}-x)\}\end{array}$

Write y explicitly

We can then write $y$ explicitly as follows:

$\begin{array}{cc}y& \in x_0+{𝒦}_n(r^{(0)},A)\\ & \in x_0+{𝒦}_n(A(x^{(0)}-x),A)\\ & =x_0+{\alpha }_{1}A(x^{(0)}-x)+{\alpha }_2{A}^2(x^{(0)}-x)+\dots +{\alpha }_n{A}^n(x^{(0)}-x)\\ & \text{for arbitrary real numbers }{\alpha }_{i}i=1,2,\dots ,n\end{array}$

Substitute and Apply Inequality

We substitute $y$ into the hypothesis inequality and apply a norm inequality of matrix norms to get the desired result.

$\begin{array}{cc}\Vert x^{(n)}-x{\Vert }_A& \le \Vert y-x{\Vert }_A\\ & =\Vert (x_0-x)+{\alpha }_{1}A(x^{(0)}-x)+{\alpha }_2{A}^2(x^{(0)}-x)+\dots +{\alpha }_n{A}^n(x^{(0)}-x){\Vert }_A\\ & =\Vert {\alpha }_n{A}^n(x^{(0)}-x)+{\alpha }_{n-1}{A}^{n-1}(x^{(0)}-x)+\dots +{\alpha }_2{A}^2(x^{(0)}-x)+{\alpha }_{1}A(x^{(0)}-x)+(x_0-x){\Vert }_A\\ & =\Vert P(A)(x^{(0)}-x){\Vert }_A\\ & \le \Vert P(A){\Vert }_A\Vert x^{(0)}-x{\Vert }_A\end{array}$

Problem 3b

Let $T_n(z)$ denote the $n^{th}$ Chebyshev polynomial. Let ${\lambda }_{min}$ and ${\lambda }_{max}$ denote respectively the smallest and largest eigenvalues of $A$. Apply the result of part (a) to $p_n(z)=\frac{T_n(\frac{{\lambda }_{max}+{\lambda }_{min}-2z}{{\lambda }_{max}-{\lambda }_{min}})}{T_n(\frac{{\lambda }_{max}+{\lambda }_{min}}{{\lambda }_{max}-{\lambda }_{min}})}$ to show that $\Vert x^{(n)}-x{\Vert }_A\le \frac{1}{T_n(\frac{{\lambda }_{max}+{\lambda }_{min}}{{\lambda }_{max}-{\lambda }_{min}})}\Vert x^{(0)}-x{\Vert }_A$. You can use without proof the fact that $\Vert p_n(A){\Vert }_A=\max \{|p_n(z)|:z\in Sp\{A\}\}$, where $Sp(A)$ denotes the set of eigenvalues of $A$, and the facts that for every $n\in \mathbb{N}$ the polynomial $T_n(z)$ has degree $n$, is positive for $z>1$, and satisfies $T_n(\cos (\theta ))=\cos (n\theta )\text{for all}\theta \in ℝ$.

Solution 3b

Overview

We want to show that

$\Vert p_n(A){\Vert }_A=\frac{1}{T_n(\frac{{\lambda }_{\max }+{\lambda }_{\min }}{{\lambda }_{\max }-{\lambda }_{\min }})}$

Maximize Numerator of p_n(z)

By hypothesis,

$\Vert p_n(A){\Vert }_A=\max \{|p_n(z)|:z\in Sp\{A\}\}$

Since only the numerator of $p_n(z)$ depends on $z$, we only need to maximize the numerator in order to maximize $|p_n(z)|$. That is find

$z:\underset{z\in [{\lambda }_{\min },{\lambda }_{\max }]}{\max }\left|T_n\left(\frac{{\lambda }_{\max }+{\lambda }_{\min }-2z}{{\lambda }_{\max }-{\lambda }_{\min }}\right)\right|$

Rewrite T_n

Let $\cos (\theta )=x$. Then

$\theta =\arccos (x)$

Hence,

$\begin{array}{cc}{T}_n(\cos \theta )& =T_n(x)\\ & =\cos (n\theta )\\ & =\cos (n\arccos (x))\end{array}$

so

$T_n(x)=\cos (n\arccos (x))$

Max of T_n is 1

Denote the argument of $T_n$ as $x(z)$ since the argument depends on $z$. Hence,

$x(z)=\frac{{\lambda }_{\max }+{\lambda }_{\min }-2z}{{\lambda }_{\max }-{\lambda }_{\min }}$,

Then,

$x({\lambda }_{\min })=\frac{{\lambda }_{\max }-{\lambda }_{\min }}{{\lambda }_{\max }-{\lambda }_{\min }}=1$

$x({\lambda }_{\max })=\frac{{\lambda }_{\min }-{\lambda }_{\max }}{{\lambda }_{\max }-{\lambda }_{\min }}=-1$

Thus $x(z)\in [-1,1]$.

Now, since $\arccos (x)$ is real for $x\in [-1,1]$,

$-1\le \underset{T_n(x)}{\underbrace{\cos (n\stackrel{\in ℝ}{\overbrace{\arccos (x)}})}}\le 1$

Hence,

$\underset{x\in [-1,1]}{\max }{T}_n(x)=\underset{x\in [-1,1]}{\max }|\cos (n\arccos (x))|=1$

Show T_n(1)=1

Let $z={\lambda }_{\min }$, then

$\left|T_n\left(\frac{{\lambda }_{\max }+{\lambda }_{\min }-2({\lambda }_{\min })}{{\lambda }_{\max }-{\lambda }_{\min }}\right)\right|=|T_n(1)|$

Using our formula we have,

$\begin{array}{cc}|T_n(1)|& =|\cos (n\cdot \arccos (1))|\\ & =|\cos (n\cdot 2\pi k)|k\in 𝐙\\ & =1\end{array}$

In other words, if $z={\lambda }_{\min }$, $T_n$ achieves its maximum value of $1$.