Suggested Solution to Homework 4

1

To find an upper bound of the eigenvalues of \(B\), we apply the Gerschgorin theorem, which is also known as the Gerschgorin circle theorem.

\[ B = \begin{pmatrix} 2 & 1 & 2 \\ 1 & 4 & 1 \\ 2 & 3 & 6 \end{pmatrix} \] \[ |\lambda - b_{ii}| \leq \sum_{j\not=i} |b_{ij}| \]

More concretely, for every eigenvalue of \(B\), we obtain an inequality according to each row of \(B\).

\begin{align} |\lambda-2| & \leq |1| + |2| = 3 \\ |\lambda-4| & \leq |1| + |1| = 2 \\ |\lambda-6| & \leq |2| + |3| = 5 \end{align}

Hence, we can determine an upper bound of the eigenvalues of \(B\).

\[|\lambda| \leq 11 \quad \forall \text{ eigenvalues \(\lambda\) of $B$}\]

2

2.1

\[ A = \begin{pmatrix} 1 & \alpha & \alpha \\ \alpha & 1 & \alpha \\ \alpha & \alpha & 1 \end{pmatrix} \] Note that \(A\) is symmetric, which means \(A\) is diagonalizable and all eigenvalues are real. For \(A\) to be positive definite, we need all the eigenvalues of \(A\) to be positive.

From the characteristic polynomial of \(A\), which is \((1-\lambda-\alpha)^2(2\alpha+1-\lambda)\), we obtain three eigenvalues.

\begin{align} \lambda_1 & = 1-\alpha \\ \lambda_2 & = 1-\alpha \\ \lambda_3 & = 2\alpha+1 \end{align}

Now, we need all \(\lambda_i\) to be positive. Essentially, we have \(-\frac{1}{2} < \alpha < 1\).

2.2

In the Gauss-Seidel method, we take \(N\) to be the lower triangular part of \(A\).

\[ N = \begin{pmatrix} 1 & 0 & 0 \\ \alpha & 1 & 0 \\ \alpha & \alpha & 1 \end{pmatrix} \]

Given \(A\) is positive definite, we see that both \(A\) and \(N^*+N-A\) are self-adjoint and positive definite matrices.

\[ N^*+N-A = \begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix} \]

Hence, we can apply the Household-John theorem, which implies that the method will converge.

3

3.1

First, we decompose \(T\) as a sum of three matrices: a diagonal \(D\), a lower triangular matrix \(L\), and an upper triangular matrix \(U\).

\[ T = D + L + U \]

Let's look at the eigenvalues of the matrix \(\alpha D^{-1} L + \frac{1}{\alpha} D^{-1} U\), where \(\alpha \not = 0\).

\[ \begin{pmatrix} 0 & 0 & -\frac{1}{2\alpha} \\ 0 & 0 & 0 & -\frac{1}{2\alpha} \\ -\frac{\alpha}{2} & 0 & 0 & 0 & -\frac{1}{2\alpha} \\ & -\frac{\alpha}{2} & 0 & 0 & 0 \\ & & -\frac{\alpha}{2} & 0 & 0 \end{pmatrix} \]

A careful computation should yield the following characteristic polynomial.

\[ -\lambda^5 +\frac{3}{4}\lambda^3 -\frac{1}{8}\lambda \]

This shows that the eigenvalues of \(\alpha D^{-1} L + \frac{1}{\alpha} D^{-1} U\) are independent of \(\alpha\). Hence, we have shown that \(T\) is consistenly ordered.

3.2

Now that \(T\) is consistently ordered. Let \(M_J\) be the matrix in Jacobi method, i.e. \(M_J = D^{-1}(D-A)\). Note that \(M_J\) is symmetric, so \(M_J\) has only real eigenvalues according to spectral theorem. In particular, from the above calculation, we see that the eigenvalues of \(M_J\) are \(0\), \(\pm \frac{1}{2}\), \(\pm \frac{1}{\sqrt{2}}\), so we have \(\rho(M_J) = \frac{1}{\sqrt{2}} < 1\).

Using the Theorem by D. Young (Lecture 17), the spectral radius \(\rho(M_{\text{SOR}}, \omega) < 1\) for \(\omega \in (0, 2)\). Hence, the condition \(\omega \in (0, 2)\) is not only the necessary condition, but also the sufficient condition for the SOR method to converge.

Meanwhile, the optimal \(\omega\) for the SOR method is given by the following.

\[\omega = \frac{2}{1+\sqrt{1-\rho(M_j)^2}}\]

Hence, we have \(\omega_{\text{opt}} = \frac{2}{1+\frac{1}{\sqrt{2}}}\).

3.3

From the above, we see that \(\rho(M_J) < 1\), so the Jacobi method converges.

On the other hand, the Gauss-Seidel method is a special case of the SOR method with \(\omega = 1\), which, from the above analysis, is within the range for the method to converge. Hence, it also converges.

4

4.1

For the first method to converge, we need the spectral radius \(\rho(G) < 1\). But \(\rho(G) = c\). So, we need \(0 \leq c < 1\).

4.2

We can express the second method in another form.

\[ u^{n+1} = \left[ \gamma G + (1-\gamma)I \right] u^{n} + \gamma v \]

Let \(H = \gamma G + (1-\gamma)I\). Note that the eigenvalues of \(G\) lie in the interval \([-c, 0]\).

\[ -c \leq \lambda_i \leq 0 \]

While the eigenvalues of \(H\) are given by \(\sigma_i := \gamma \lambda_i + (1-\gamma)\). Consider \(\gamma > 0\), we have the following inequality.

\[ -c \gamma + 1 - \gamma \leq \sigma_i \leq 1 - \gamma \]

(If \(\gamma < 0\), then one could easily get \(\rho(H) > 1\). So, we only consider \(\gamma > 0\).)

Obviously, we have \(\sigma_i < 1\). So, for the spectral radius \(\rho(H) < 1\), we need \(-1 < -c \gamma + 1 - \gamma\), that is, \(0 < \gamma < \frac{2}{1+c}\).

4.3

Note that \(\rho(H)\) is smallest when we have the following equality.

\[ -(-c \gamma + 1-\gamma) = 1-\gamma\]

Solving the equation gives the optimal value of \(\gamma\), namely \(\gamma^*\).

\[\gamma^* = \frac{2}{2+c}\]

4.4

Note that the second scheme always converges when we take \(\gamma = \gamma^*\). But for both methods to converge, we consider \(c \in [0, 1)\) and compare the spectral radii of the two methods.

\[ \left\{ \begin{array}{l} \rho(G) = c \\ \rho(H) = \frac{c}{2+c} \end{array} \right. \]

This shows that \(\rho(H)\) is always smaller than \(\rho(G)\). So, the rate of convergence of the second method is much better than the first one. For example, when \(\rho(G)\) is close to \(1\), \(\rho(H)\) is less than \(0.334\).

(Remark: in a more rigorous sense, one can replace the estimate \(-c \leq \lambda_i \leq 0\) by \(\lambda_1 \leq \lambda_i \leq \lambda_n\) and derive similar inequalities; nevertheless, we have \(\rho(H) \leq \frac{c}{2+c}\) with \(\gamma^* = \frac{2}{2-\lambda_1-\lambda_n}\).)

Created: 2018-04-03 Tue 11:57

Emacs 25.3.1 (Org mode 8.2.10)

Validate