Past exam of the mathematics course of the University of Cambridge / 2003 / iii / Paper 67 / 7

The multigrid method addresses a characteristic failure of single-grid relaxation: its cheapest updates rapidly remove oscillatory error but barely change smooth error. For example, the one-dimensional Poisson equation matrix has diagonal $2/h^2$. The weighted Jacobi method with weight $\omega$ has Fourier error factorAt low physical frequency, $\theta =O(h)$ and this factor is $1-O(h^2)$, requiring $O(h^{-2})$ sweeps for a fixed error reduction. With $\omega =2/3$, high grid frequencies $\pi /2\le |\theta |\le \pi$ are multiplied by at most $1/3$ in modulus. Smooth error on a fine mesh is well represented on a coarser mesh, where it is relatively more oscillatory and hence cheaper to remove. These two complementary mechanisms motivate smoothing and coarse-grid correction.

Write the discrete problem as $A_h{u}_h=f_h$. For an approximation $v_h$, the residual is $r_h=f_h-A_h{v}_h$ and the error $e_h=u_h-v_h$ satisfies $A_h{e}_h=r_h$. It is this error equation that is transferred to a coarse grid; restricting the approximate solution alone would not solve the missing error. With prolongation $P$ and restriction $R$, solveFor nested Cartesian grids with $H=2h$, linear interpolation copies a coarse value at a coincident fine node and averages the two nearest coarse values at an intermediate node. Its two-dimensional counterpart is bilinear interpolation. Full-weighting restriction uses weights $(1,2,1)/4$ in one dimension andin two dimensions. With these choices $R=2^{-d}{P}^T$ in ordinary coordinate inner products. A Galerkin coarse-grid operator is $A_H=R{A}_{h}P$; alternatively one can rediscretize the PDE on the coarse mesh when that choice remains compatible with the transfers.

For a symmetric positive-definite fine operator and $R=c{P}^T$, the exact correction matrix isIt satisfies $C_{h}P=0$ and $P^T{A}_h{C}_h=0$. Thus it removes the coarse-space component of the error and is an orthogonal projection in the energy inner product. With pre- and post-smoothing matrices $S_{\mathrm{pre}},S_{\mathrm{post}}$, the exact two-grid error matrix is $S_{\mathrm{post}}{C}_h{S}_{\mathrm{pre}}$. Effective smoothing on the complementary space and accurate coarse representation of smooth error are both needed for contraction.

A concrete harmonic calculation makes the mechanism quantitative. On the one-dimensional constant-coefficient problem, couple the two harmonics $\theta$ and $\theta +\pi$, put $s={\sin }^2(\theta /2)$ and $c=1-s$, and use the above transfers. Linear interpolation has harmonic weights $(c,s{)}^T$, while the fine operator has diagonal symbol $4\mathrm{diag}(s,c)/h^2$. Therefore its Galerkin correction on this pair isOne pre- and one post-sweep of weighted Jacobi with $\omega =2/3$ gives $E=SCS$, $S=\mathrm{diag}(1-4s/3,1-4c/3)$. This rank-one matrix has eigenvalues zero andThe two-grid Poisson factor with weighted Jacobi is therefore mesh-independent in this model. A periodic constant nullspace must first be removed; a Dirichlet problem has no such nullspace. This is evidence for the design, not a claim that the same factor applies to every PDE or every smoother.

The geometric multigrid V-cycle replaces the exact coarse solve by one recursively defined coarse cycle. A complete implementation on level $\ell$ proceeds as follows:

Boundary values of the solution are imposed on the fine level; the error equation has homogeneous boundary conditions when those values were already imposed exactly. Residuals must be computed after pre-smoothing, not before. The zero initialization on each recursive error solve makes the cycle an error correction rather than an accidental addition of an old coarse solution. For an energy-symmetric cycle usable in conjugate gradient method, choose compatible adjoint transfers and symmetric smoothing, such as forward Gauss-Seidel before correction and its backward sweep afterward.

With bounded stencil size and fixed sweep counts, a sweep, residual calculation or transfer costs $O(N_{\ell })$. Uniform coarsening gives $N_{\ell -1}\approx 2^{-d}{N}_{\ell }$, soIf the full V-cycle has an error reduction factor $q<1$ independent of mesh size, reducing an arbitrary initial algebraic error by a factor $\epsilon$ takes $O(\log (1/\epsilon ))$ cycles and $O(N\log (1/\epsilon ))$ work. Linear work per cycle alone does not prove the mesh-independent factor; it rests on the smoothing and approximation properties just illustrated. Full multigrid starts on a coarse mesh, interpolates its solution to each successive finer mesh, and uses a few cycles per level. When the interpolation accuracy and cycle contraction match the discretization error, this reaches discretization-level accuracy in $O(N)$ work.

Performance depends on the operator. Strong anisotropy can leave error oscillatory in one direction but weakly penalized by the differential operator, defeating point smoothing; line or plane relaxation and selective coarsening address this. Large coefficient jumps demand compatible coarse spaces. Nonrectangular meshes motivate algebraic rather than purely geometric coarsening. Nonlinear problems use a full-approximation formulation instead of the linear residual correction verbatim. Singular Neumann or periodic operators need right-hand-side compatibility and consistent nullspace treatment on every level. Multigrid is fast because its components address different error scales, not because coarsening alone guarantees a good solver.