Innumerical analysis, amultigrid method (MG method) is analgorithm for solvingdifferential equations using ahierarchy ofdiscretizations. They are an example of a class of techniques calledmultiresolution methods, very useful in problems exhibitingmultiple scales of behavior. For example, many basicrelaxation methods exhibit different rates of convergence for short- and long-wavelength components, suggesting these different scales be treated differently, as in aFourier analysis approach to multigrid.^[1] MG methods can be used as solvers as well aspreconditioners.

The main idea of multigrid is to accelerate the convergence of a basic iterative method (known as relaxation, which generally reduces short-wavelength error) by aglobal correction of the fine grid solution approximation from time to time, accomplished by solving acoarse problem. The coarse problem, while cheaper to solve, is similar to the fine grid problem in that it also has short- and long-wavelength errors. It can also be solved by a combination of relaxation and appeal to still coarser grids. This recursive process is repeated until a grid is reached where the cost of direct solution there is negligible compared to the cost of one relaxation sweep on the fine grid. This multigrid cycle typically reduces all error components by a fixed amount bounded well below one, independent of the fine grid mesh size. The typical application for multigrid is in the numerical solution ofelliptic partial differential equations in two or more dimensions.^[2]

Multigrid methods can be applied in combination with any of the common discretization techniques. For example, thefinite element method may be recast as a multigrid method.^[3] In these cases, multigrid methods are among the fastest solution techniques known today. In contrast to other methods, multigrid methods are general in that they can treat arbitrary regions andboundary conditions. They do not depend on theseparability of the equations or other special properties of the equation. They have also been widely used for more-complicated non-symmetric and nonlinear systems of equations, like theLamé equations ofelasticity or theNavier-Stokes equations.^[4]

There are many variations of multigrid algorithms, but the common features are that a hierarchy ofdiscretizations (grids) is considered. The important steps are:^[5][6]

• Smoothing – reducing high frequency errors, for example using a fewiterations of theGauss–Seidel method.
• Residual Computation – computingresidual error after the smoothing operation(s).
• Restriction – downsampling theresidual error to a coarser grid.
• Interpolation orprolongation – interpolating a correction computed on a coarser grid into a finer grid.
• Correction – Adding prolongated coarser grid solution onto the finer grid.

There are many choices of multigrid methods with varying trade-offs between speed of solving a single iteration and the rate of convergence with said iteration. The 3 main types are V-Cycle, F-Cycle, and W-Cycle. These differ in which and how many coarse-grain cycles are performed per fine iteration. The V-Cycle algorithm executes one coarse-grain V-Cycle. F-Cycle does a coarse-grain V-Cycle followed by a coarse-grain F-Cycle, while each W-Cycle performs two coarse-grain W-Cycles per iteration. For adiscrete 2D problem, F-Cycle takes 83% more time to compute than a V-Cycle iteration while a W-Cycle iteration takes 125% more. If the problem is set up in a 3D domain, then a F-Cycle iteration and a W-Cycle iteration take about 64% and 75% more time respectively than a V-Cycle iteration ignoringoverheads. Typically, W-Cycle produces similar convergence to F-Cycle. However, in cases ofconvection-diffusion problems with highPéclet numbers, W-Cycle can show superiority in its rate of convergence per iteration over F-Cycle. The choice of smoothing operators are extremely diverse as they includeKrylov subspace methods and can bepreconditioned.

Any geometric multigrid cycle iteration is performed on a hierarchy of grids and hence it can be coded using recursion. Since the function calls itself with smaller sized (coarser) parameters, the coarsest grid is where the recursion stops. In cases where the system has a highcondition number, the correction procedure is modified such that only a fraction of the prolongated coarser grid solution is added onto the finer grid.

functionphi=V_Cycle(phi,f,h)% Recursive V-Cycle Multigrid for solving the Poisson equation (\nabla^2 phi = f) on a uniform grid of spacing h% Pre-Smoothingphi=smoothing(phi,f,h);% Compute Residual Errorsr=residual(phi,f,h);% Restrictionrhs=restriction(r);eps=zeros(size(rhs));% stop recursion at smallest grid size, otherwise continue recursionifsmallest_grid_size_is_achievedeps=coarse_level_solve(eps,rhs,2*h);elseeps=V_Cycle(eps,rhs,2*h);end% Prolongation and Correctionphi=phi+prolongation(eps);% Post-Smoothingphi=smoothing(phi,f,h);end

functionphi=F_Cycle(phi,f,h)% Recursive F-cycle multigrid for solving the Poisson equation (\nabla^2 phi = f) on a uniform grid of spacing h% Pre-smoothingphi=smoothing(phi,f,h);% Compute Residual Errorsr=residual(phi,f,h);% Restrictionrhs=restriction(r);eps=zeros(size(rhs));% stop recursion at smallest grid size, otherwise continue recursionifsmallest_grid_size_is_achievedeps=coarse_level_solve(eps,rhs,2*h);elseeps=F_Cycle(eps,rhs,2*h);end% Prolongation and Correctionphi=phi+prolongation(eps);% Re-smoothingphi=smoothing(phi,f,h);% Compute residual errorsr=residual(phi,f,h);% Restrictionrhs=restriction(r);% stop recursion at smallest grid size, otherwise continue recursionifsmallest_grid_size_is_achievedeps=coarse_level_solve(eps,rhs,2*h);elseeps=V_Cycle(eps,rhs,2*h);end% Prolongation and Correctionphi=phi+prolongation(eps);% Post-smoothingphi=smoothing(phi,f,h);end

functionphi=W_cycle(phi,f,h)% Recursive W-cycle multigrid for solving the Poisson equation (\nabla^2 phi = f) on a uniform grid of spacing h% Pre-smoothingphi=smoothing(phi,f,h);% Compute Residual Errorsr=residual(phi,f,h);% Restrictionrhs=restriction(r);eps=zeros(size(rhs));% stop recursion at smallest grid size, otherwise continue recursionifsmallest_grid_size_is_achievedeps=coarse_level_solve(eps,rhs,2*h);elseeps=W_cycle(eps,rhs,2*h);end% Prolongation and correctionphi=phi+prolongation(eps);% Re-smoothingphi=smoothing(phi,f,h);% Compute residual errorsr=residual(phi,f,h);% Restrictionrhs=restriction(r);% stop recursion at smallest grid size, otherwise continue recursionifsmallest_grid_size_is_achievedeps=coarse_level_solve(eps,rhs,2*h);elseeps=W_cycle(eps,rhs,2*h);end% Prolongation and correctionphi=phi+prolongation(eps);% Post-smoothingphi=smoothing(phi,f,h);end

This sectiondoes notcite anysources. Please helpimprove this section byadding citations to reliable sources. Unsourced material may be challenged andremoved.(September 2021) (Learn how and when to remove this message)

This approach has the advantage over other methods that it often scales linearly with the number of discrete nodes used. In other words, it can solve these problems to a given accuracy in a number of operations that is proportional to the number of unknowns.

Assume that one has a differential equation which can be solved approximately (with a given accuracy) on a grid${\displaystyle i}$ with a given grid point density${\displaystyle N_i}$. Assume furthermore that a solution on any grid${\displaystyle N_i}$ may be obtained with a given effort${\displaystyle W_i=\rho K{N}_i}$ from a solution on a coarser grid${\displaystyle i+1}$. Here, is the ratio of grid points on "neighboring" grids and is assumed to be constant throughout the grid hierarchy, and${\displaystyle K}$ is some constant modeling the effort of computing the result for one grid point.

The following recurrence relation is then obtained for the effort of obtaining the solution on grid${\displaystyle k}$:$${\displaystyle W_k=W_{k+1}+\rho K{N}_k}$$

And in particular, we find for the finest grid${\displaystyle N_1}$ that$${\displaystyle W_1=W_2+\rho K{N}_1}$$Combining these two expressions (and using${\displaystyle N_k={\rho }^{k-1}{N}_1}$) gives$${\displaystyle W_1=K{N}_1\sum _{p=0}^n{\rho }^p}$$

Using thegeometric series, we then find (for finite${\displaystyle n}$)

that is, a solution may be obtained in${\displaystyle O(N)}$ time. It should be mentioned that there is one exception to the${\displaystyle O(N)}$ i.e. W-cycle multigrid used on a 1D problem; it would result in${\displaystyle O(N\log N)}$ complexity.^{[citation needed]}

A multigrid method with an intentionally reduced tolerance can be used as an efficientpreconditioner for an external iterative solver, e.g.,^[7] The solution may still be obtained in${\displaystyle O(N)}$ time as well as in the case where the multigrid method is used as a solver. Multigrid preconditioning is used in practice even for linear systems, typically with one cycle per iteration, e.g., inHypre. Its main advantage versus a purely multigrid solver is particularly clear for nonlinear problems, e.g.,eigenvalue problems.

If the matrix of the original equation or an eigenvalue problem is symmetric positive definite (SPD), the preconditioner is commonly constructed to be SPD as well, so that the standardconjugate gradient (CG)iterative methods can still be used. Such imposed SPD constraints may complicate the construction of the preconditioner, e.g., requiring coordinated pre- and post-smoothing. However,preconditionedsteepest descent andflexible CG methods for SPD linear systems andLOBPCG for symmetric eigenvalue problems are all shown^[8] to be robust if the preconditioner is not SPD.

Originally described in Xu’s Ph.D. thesis^[9] and later published in Bramble-Pasciak-Xu,^[10] the BPX-preconditioner is one of the two major multigrid approaches (the other is the classic multigrid algorithm such as V-cycle) for solving large-scale algebraic systems that arise from the discretization of models in science and engineering described by partial differential equations. In view of the subspace correction framework,^[11] BPX preconditioner is a parallel subspace correction method whereas the classic V-cycle is a successive subspace correction method. The BPX-preconditioner is known to be naturally more parallel and in some applications more robust than the classic V-cycle multigrid method. The method has been widely used by researchers and practitioners since 1990.

Multigrid methods can be generalized in many different ways. They can be applied naturally in a time-stepping solution ofparabolic partial differential equations, or they can be applied directly to time-dependentpartial differential equations.^[12] Research on multilevel techniques forhyperbolic partial differential equations is underway.^[13] Multigrid methods can also be applied tointegral equations, or for problems instatistical physics.^[14]

Another set of multiresolution methods is based uponwavelets. These wavelet methods can be combined with multigrid methods.^[15][16] For example, one use of wavelets is to reformulate the finite element approach in terms of a multilevel method.^[17]

Adaptive multigrid exhibitsadaptive mesh refinement, that is, it adjusts the grid as the computation proceeds, in a manner dependent upon the computation itself.^[18] The idea is to increase resolution of the grid only in regions of the solution where it is needed.

Practically important extensions of multigrid methods include techniques where no partial differential equation nor geometrical problem background is used to construct the multilevel hierarchy.^[19] Suchalgebraic multigrid methods (AMG) construct their hierarchy of operators directly from the system matrix. In classical AMG, the levels of the hierarchy are simply subsets of unknowns without any geometric interpretation. (More generally, coarse grid unknowns can be particular linear combinations of fine grid unknowns.) Thus, AMG methods become black-box solvers for certain classes ofsparse matrices. AMG is regarded as advantageous mainly where geometric multigrid is too difficult to apply,^[20] but is often used simply because it avoids the coding necessary for a true multigrid implementation. While classical AMG was developed first, a related algebraic method is known as smoothed aggregation (SA).

In an overview paper^[21] by Jinchao Xu and Ludmil Zikatanov, the "algebraic multigrid" methods are understood from an abstract point of view. They developed a unified framework and existing algebraic multigrid methods can be derived coherently. Abstract theory about how to construct optimal coarse space as well as quasi-optimal spaces was derived. We note that this result appeared first in a note on Algebraic Multigrid by Brannick and Zikatanov and was just rewritten in the overview paper. Also, they proved that, under appropriate assumptions, the abstract two-level AMG method converges uniformly with respect to the size of the linear system, the coefficient variation, and the anisotropy. Their abstract framework covers most existing AMG methods, such as classical AMG, energy-minimization AMG, unsmoothed and smoothed aggregation AMG, and spectral AMGe.

Multigrid methods have also been adopted for the solution ofinitial value problems.^[22]Of particular interest here are parallel-in-time multigrid methods:^[23]in contrast to classicalRunge–Kutta orlinear multistep methods, they can offerconcurrency in temporal direction.The well knownParareal parallel-in-time integration method can also be reformulated as a two-level multigrid in time.

Nearly singular problems arise in a number of important physical and engineering applications. Simple, but important example of nearly singular problems can be found at the displacement formulation oflinear elasticity for nearly incompressible materials. Typically, the major problem to solve such nearly singular systems boils down to treat the nearly singular operator given by${\displaystyle A+\epsilon M}$ robustly with respect to the positive, but small parameter${\displaystyle \epsilon }$. Here${\displaystyle A}$ is symmetricsemidefinite operator with largenull space, while${\displaystyle M}$ is a symmetricpositive definite operator. There were many works to attempt to design a robust and fast multigrid method for such nearly singular problems. A general guide has been provided as a design principle to achieve parameters (e.g., mesh size and physical parameters such asPoisson's ratio that appear in the nearly singular operator) independent convergence rate of the multigrid method applied to such nearly singular systems,^[24] i.e., in each grid, a space decomposition based on which the smoothing is applied, has to be constructed so that the null space of the singular part of the nearly singular operator has to be included in the sum of the local null spaces, the intersection of the null space and the local spaces resulting from the space decompositions.

• Links to AMG presentations

• Numerical analysis
• Partial differential equations
• Wavelets

• CS1: long volume value
• Articles with short description
• Short description is different from Wikidata
• Articles needing additional references from September 2021
• All articles needing additional references
• All articles with unsourced statements
• Articles with unsourced statements from September 2021
• CS1 maint: bot: original URL status unknown