\[Ax=b\]

Iterative methods

Simple form:

\[x^{(k)}=Bx^{(k-1)}+c\]

B & c independent from k -> stationary

Stationary methods

Jacobi

The Jacobi method is based on solving for every variable locally with respect to the other variables; one iteration of the method corresponds to solving for every variable once. The resulting method is easy to understand and implement, but convergence is slow.

\[x^{(k)}_i = (b_i - \sum_{j\neq i} a_{i,j}x_j)/a_{i,i})\]

Matrix Form:

\[x^{(k)}=D^{-1}(L+U)x^{(k-1)} + D^{-1} b\]

Gauss-Seidel

The Gauss-Seidel method is like the Jacobi method, except that it uses updated values as soon as they are available. In general, it will converge faster than the Jacobi method, though still relatively slowly.

\[x^{(k)}_i = (b_i - \sum_{j < i} a_{i,j}x^{(k)}_j - \sum_{j > i} a_{i,j}x^{(k-1)}_j )/a_{i,i})\]

Matrix form:

\[x^{(k)} = (D-L)^{-1}(Ux^{(k-1)}+b)\]

SOR

Successive Overrelaxation (SOR) can be derived from the Gauss-Seidel method by introducing an extrapolation parameter \(w\). For the optimal choice of \(w\), SOR converges faster than Gauss-Seidel by an order of magni- tude.

\[x^{(k)}_i = w\bar{x}^{(k)}_i + (1-w)x^{(k-1)}_i\]

Matrix Form:

\[x^{(k)} = (D- wL)^{-1}(wU + (1-w)D)x^(k-1)+w(D-wL)^{-1}b\]

Fail to converge when \(w \mbox{ outside of } (0,2)\), not possible to compute in advance the optimal value of \(w\).

• Matrix A is symmetric & positive definite¹ guaranteed to converge with \(w \in (0,2)\)
• Heuristic estimates
• Direct relationship bet spectra of the jacobi and SOR iteration matrix, with spectral radius \(\rho\) of the jacobi matrix

SSOR

Symmetric Successive Overrelaxation (SSOR) has no advantage over SOR as a stand-alone iterative method; however, it is useful as a preconditioner for nonstationary methods.

assume A symmetric => two sweeps

• Forward as in SOR
• Backward SOR

\[x^{(k)} = B_1 B_2 x^{(k-1)} w(2-w)(D-wU)^{-1} D (D-wL)^{-1}b\] where \[B_1=(D-wU)^{-1}(wL+(1-w)D)\] \[B_2=(D-wL)^{-1}(wU+(1-w)D)\]

Nonstationary Methods

Nonstationary methods different from stationary methods in that the computations involve information that changes at each iteration.

Conjugate Gradient (CG)

The conjugate gradient method derives its name from the fact that it gen- erates a sequence of conjugate (or orthogonal) vectors. These vectors are the residuals of the iterates. They are also the gradients of a quadratic functional, the minimization of which is equivalent to solving the linear system. CG is an extremely ective method when the coecientcient matrix is symmetric positive denite, since storage for only a limited number of vectors is required.

iterates are updated bu a multiple of the search direction p \[x^{(i)}=x^{(i-1)}+\alpha_i p^{(i)}\] residuals are updated as \[r^{(i)}=r^{(i-1)} - \alpha q^{(i)}, \mbox{where } q^{(i)}=Ap^{(i)}\]

\[\alpha = \alpha_i = r^T r/ p^T A p \mbox{ minimize } r^T A^{-1} r\]

search direction: \[p^{(i)} = r^{(i)} +\beta_{i-1} p^{(i-1)}\]

\(\beta_i = r^{(i)^T} r^{(i)}/r^{(i-1)^T} r^{(i-1)}\) ensures \(p^{(i)}\) and \(Ap^{(i-1)}\) or \(r^{(i)}\) and \(r^{(i-1)}\) p,r orthogonal to all previous q and r

preconditioner => M²

Minimum Residual (MINRES) and Symmetric LQ (SYMMLQ)

These methods are computational alternatives for CG for coecient matri- ces that are symmetric but possibly indenite. SYMMLQ will generate the same solution iterates as CG if the coecientcient matrix is symmetric positive denite.

Conjugate Gradient on the Normal Equations: CGNE and CGNR

These methods are based on the application of the CG method to one of two forms of the normal equations for Ax=b. CGNE solves the system \((A A^T)y=b\) for y and then computes the solution \(x=A^T y\). CGNR solves \((A^T A)x=\tilde{b}\) for the solution vector x where \(\tilde{b} = A^T b\). When the coefficient matrix A is nonsymmetric and nonsingular, the normal equations matrices \(AA^T\) and \(A^T A\) will be symmetric and positive definite, and hence CG can be applied. The convergence may be slow, since the spectrum of the normal equations matrices will be less favorable than the spectrum of A.

Generalized Minimal Residual (GMRES)

The Generalized Minimal Residual method computes a sequence of orthog- onal vectors (like MINRES), and combines these through a least-squares solve and update. However, unlike MINRES (and CG) it requires storing the whole sequence, so that a large amount of storage is needed. For this reason, restarted versions of this method are used. In restarted versions, computation and storage costs are limited by specifying a xed number of vectors to be generated. This method is useful for general nonsymmetric matrices.

BiConjugate Gradient (BiCG)

The Bicongugate Gradient method generates two CG-like sequences of vec- tors, one based on a system with the original coecientcientcient matrix A, and one on \(A^T\) . Instead of orthogonalizing each sequence, they are made mutually orthogonal, or “bi-orthogonal”. This method, like CG, uses limited storage. It is useful when the matrix is nonsymmetric and nonsingular; however, con- vergence may be irregular, and there is a possibility that the method will break down. BiCG requires a multiplication with the coecientcientcient matrix and with its transpose at each iteration.

Quasi-Minimal Residual (QMR)

The Quasi-Minimal Residual method applies a least-squares solve and up- date to the BiCG residuals, thereby smoothing out the irregular convergence behavior of BiCG. QMR largely avoids the breakdown that can occur in BiCG.

Conjugate Gradient Squared (CGS)

The Conjugate Gradient Squared method is a variant of BiCG that applies the updating operations for the A-sequence and the \(A^T\)-sequences both to the same vectors. Ideally, this would double the convergence rate, but in practice convergence may be much more irregular than for BiCG. An added practical advantage is that the method does not need the multiplications with the transpose of the coecientcientcient matrix.

Biconjugate Gradient Stabilized (Bi-CGSTAB)

The Biconjugate Gradient Stabilized method is a variant of BiCG, like CGS, but using different updates for the \(A^T\)-sequence in order to obtain smoother convergence than CGS.

Chebyshev Iteration

The Chebyshev Iteration recursively determines polynomials with coefficient chosen to minimize the norm of the residual in a min-max sense. The coecientcientcient matrix must be positive denite and knowledge of the extremal eigenvalues is required. This method has the advantage of requiring no inner products.