Preliminary Examination, Numerical Analysis, Spring 2006

Sensitivity Analysis

(a) Let \(X\in \mathbb {R}^{n\times n}\) such that \(\lVert X\rVert <1\) for some induced norm. Show that \(I-X\) is nonsingular and \[ \lVert (I-X)^{-1}\rVert \le \frac {1}{1-\lVert X\rVert }. \]

(b) Assume that \(A\in \mathbb {R}^{n\times n}\) is nonsingular, and the linear system \[ Ax=b \] is perturbed to \[ (A+\delta A)(x+\delta x)=b+\delta b, \] and the perturbations are small enough so that \[ \lVert \delta A\rVert \le \epsilon \lVert A\rVert ,\qquad \lVert \delta b\rVert \le \epsilon \lVert b\rVert , \] for some small \(\epsilon >0\), and \(\epsilon \kappa (A)=r<1\), where \(\kappa (A)\) is the condition number of \(A\). Show that \(A+\delta A\) is nonsingular.

(c) Let \(\tilde {x}=x+\delta x\), show that \[ \lVert \tilde {x}\rVert \le \frac {1+r}{1-r}\lVert x\rVert ; \] and \[ \frac {\lVert \delta x\rVert }{\lVert x\rVert } \le \frac {2\epsilon }{1-r}\kappa (A). \]

SVD Let \(A\) be a real \(m\times n\) matrix \((m\ge n)\) of full rank. We are interested in the singular value decomposition of \(A=U\Sigma V^T\). It is known that algorithms for SVD can be derived by turning the SVD problem into an eigenvalue problem. One such approach is to find the eigenvalue decomposition of the \((m+n)\times (m+n)\) symmetric matrix \[ \begin {bmatrix} 0 & A^T\\ A & 0 \end {bmatrix}. \] Establish the connections between the singular values, left and right singular vectors of \(A\) and the eigenvalues, eigenvectors of the above matrix. You should pay particular attention to the case \(m>n\). Comment on the advantages and disadvantages of this approach compared with an approach based on the eigendecomposition of \(A^T A\).

Numerical Quadrature Consider the quadrature formula \[ \int _a^b w(x)f(x)\,dx=\sum _{i=1}^n \alpha _i f(\xi _i)+E, \] for the integral of \(f\) with a given positive weight function \(w(x)>0\). We would like to use \(n\) knot values \((f(\xi _i), i=1,\ldots ,n)\) to achieve the exact integral if \(f\) is a polynomial of degree up to \(2n-1\).

(a) Start from the special case \(a=-1\), \(b=1\) and \(w(x)=1\), define the orthogonal polynomials over \([-1,1]\).

(b) Use those orthogonal polynomials and polynomial division to simplify the integral. How would you pick the knots \(\xi _i\) and the weights \(\alpha _i\) such that the formula is exact (i.e., \(E=0\)) for all polynomials of degree up to \(2n-1\)? Verify that your choice works.

(c) Show that you can extend the approach to more general integrals, with arbitrary \(a<b\) and positive weight function \(w\).

ODE Consider the initial value problem for a system of ODEs \[ y'=Ay \] where \(A\) has eigenvalues ranging from \(-10^4\) to \(-1\). Which of the following schemes would be the best choice for solving this problem? Justify your answer in terms of stability, accuracy, and efficiency. \[ y^{n+1}=y^n+hAy^n, \tag {1} \] \[ y^{n+1}=y^n+hAy^{n+1}, \tag {2} \] \[ y^{n+1}=y^n+\frac {h}{2}(Ay^n+Ay^{n+1}). \tag {3} \]

Elliptic Equations Consider the standard five-point difference approximation (centered difference for both the gradient and divergence operators) for the variable coefficient Poisson’s equation \[ -\nabla \cdot (a\nabla u)=f \] with Dirichlet boundary conditions, in a two-dimensional rectangular region. We assume that \(a(\mathbf {x})\ge \delta >0\). The approximate solution \(\{u_{i,j}\}\) satisfies a linear system \(Au=b\).

(a) State and prove the maximum principle for the numerical solution \(u_{i,j}\).

(b) Use the maximum principle to show that the coefficient matrix \(A\) is nonsingular.

(c) Derive the matrix \(A\) in the one-dimensional case and show that it is symmetric and positive definite.

(d) With this linear system, do we get convergence with Jacobi iterations? What about Gauss-Seidel iterations? Give reasons for your answers.

Wave Equation For the PDE \(v_t+cv_x=0\), derive the Lax-Wendroff scheme \[ u_j^{n+1}=u_j^n -\frac {1}{2}\left (\frac {c\Delta t}{\Delta x}\right ) (u_{j+1}^n-u_{j-1}^n) +\frac {1}{2}\left (\frac {c\Delta t}{\Delta x}\right )^2 (u_{j-1}^n-2u_j^n+u_{j+1}^n). \] Analyze this scheme in terms of accuracy and stability. Does this scheme satisfy a maximum principle? If it does, explain why. If it does not, give an example to show this.

Convection-Diffusion Equation For the convection-diffusion equation \[ v_t+av_x=bv_{xx}, \] where \(a>0\) and \(b>0\), a natural Crank-Nicolson scheme is \[ \begin {aligned} \frac {u_j^{n+1}-u_j^n}{\Delta t} &+\frac {a}{2}\left ( \frac {u_{j+1}^n-u_{j-1}^n}{2\Delta x} +\frac {u_{j+1}^{n+1}-u_{j-1}^{n+1}}{2\Delta x} \right ) \\ &= \frac {b}{2}\left ( \frac {u_{j+1}^n-2u_j^n+u_{j-1}^n}{(\Delta x)^2} +\frac {u_{j+1}^{n+1}-2u_j^{n+1}+u_{j-1}^{n+1}}{(\Delta x)^2} \right ). \end {aligned} \] Let \(\nu =b\Delta t/(\Delta x)^2\) and \(\mu =a\Delta t/(2\Delta x)\). Show that the scheme is stable in the maximum norm if \(\mu \le \nu \le 1\). Also use Fourier analysis to derive the stability condition in the \(\ell _2\) norm and show that the scheme is always stable in the \(\ell _2\) norm. Suppose \(u\) stands for some volume fraction of a certain species, which norm is more appropriate? Explain why this scheme is so restrictive that it is not a good choice for many applications.