What is the advantage of multigrid over domain decomposition preconditioners, and vice versa?

Question

This is mostly aimed for elliptic PDEs over convex domains, so that I can get a good overview of the two methods.

Answer 1

Multigrid and multilevel domain decomposition methods have so much in common that each can usually be written as a special case of the other. The analysis frameworks are somewhat different, as a consequence of the different philosophies of each field. Generally speaking, multigrid methods use moderate coarsening rates and simple smoothers while domain decomposition methods use extremely rapid coarsening and strong smoothers.

Multigrid (MG)

Multigrid uses moderate coarsening rates and achieves robustness through modification of interpolation and smoothers. For elliptic problems, the interpolation operators should be "low energy", such that they preserve the near-null space of the operator (e.g. rigid body modes). An example geometric approach to these low energy interpolants is Wan, Chan, Smith (2000), compare to the algebraic construction of smoothed aggregation Vaněk, Mandel, Brezina (1996) (parallel implementations in ML and PETSc via PCGAMG, the replacement for Prometheus). Trottenberg, Oosterlee, and Schüller's book is a good general reference on Multigrid methods.

Most multigrid smoothers involve pointwise relaxation, either additively (Jacobi) or multiplicatively (Gauss Seidel). These correspond to tiny (single node or single element) Dirichlet problems. Some spectral adaptivity, robustness, and vectorizability can be achieved using Chebyshev smoothers, see Adams, Brezina, Hu, Tuminaro (2003). For non-symmetric (e.g. transport) problems, multiplicative smoothers like Gauss-Seidel are generally necessary, and upwinded interpolants may be used. Alternatively, smoothers for saddle point and stiff wave problems can be constructed by transforming via Schur-complement-inspired "block preconditioners" or by the related "distributed relaxation", into systems in which simple smoothers are effective.

Textbook multigrid efficiency refers to solving to discretization error in a small multiple of the cost of a few residual evaluations, as few as four, on the fine grid. This implies that the number of iterations to a fixed algebraic tolerance goes down as the number of levels in increased. In parallel, the time estimate involves a logarithmic term arising due to synchronization implied by the multigrid hierarchy.

Domain Decomposition (DD)

The first domain decomposition methods had only one level. With no coarse level, the condition number of the preconditioned operator cannot be less than $\mathcal O\big( \frac{L^2}{H^2}\big)$ where $L$ is the diameter of the domain and $H$ is the nominal subdomain size. In practice, condition numbers for one-level DD falls between this bound and $\mathcal O\big( \frac{L^2}{hH} \big)$ where $h$ is the element size. Note that the number of iterations needed by a Krylov method scales as the square root of the condition number. Optimized Schwarz methods (Gander 2006) improve the constants and dependence on $H/h$ relative to Dirichlet and Neumann methods, but generally do not include coarse levels and thus degrade in the case of many subdomains. See the books by Smith, Bjørstad, and Gropp (1996) or Toselli and Widlund (2005) for a general reference to domain decomposition methods.

For optimal or quasi-optimal convergence rates, multiple levels are necessary. Most DD methods are posed as two-level methods and some are very difficult to extend to more levels. DD methods can be classified as overlapping or non-overlapping.

Overlapping

These Schwarz methods use overlap and are generally based on solving Dirichlet problems. The strength of the methods can be increased by increasing the overlap. This class of methods is usually robust, do not require local null-space identification or technical modifications for problems with local constraints (common in engineering solid mechanics), but involve extra work (especially in 3D) due the overlap. Additionally, for constrained problems like incompressible, the inf-sup constant of the overlapping strip usually appears, leading to suboptimal convergence rates. Modern overlapping methods using similar coarse spaces to BDDC/FETI-DP (discussed below) are developed by Dorhmann, Klawonn, and Widlund (2008) and Dohrmann and Widlund (2010).

Non-overlapping

These methods usually solve Neumann problems of some sort, which means that unlike Dirichlet methods, they cannot work with a globally assembled matrix, and instead require unassembled or partially assembled matrices. The most popular Neumann methods either enforce continuity between subdomains by balancing at every iteration or by Lagrange multipliers that will enforce continuity only once convergence is reached. The early methods of this sort (Balancing Neumann-Neumann and FETI) require precise characterization of the null space of each subdomain, both to construct the coarse level and to make the subdomain problems non-singular. Later methods (BDDC and FETI-DP) select subdomain corners and/or edge/face moments as coarse level degrees of freedom. See Klawonn and Rheinbach (2007) for an in-depth discussion of coarse space selection for 3D elasticity. Mandel, Dohrmann, and Tazaur (2005) showed that BDDC and FETI-DP have all the same eigenvalues, except for possible 0 and 1.

More than two levels

Most DD methods are only posed as two-level methods, and some select coarse spaces that are inconvenient for use with more than two levels. Unfortunately, especially in 3D, the coarse level problems quickly become a bottleneck, limiting the problem sizes that can be solved. Additionally, the condition numbers of the preconditioned operators, especially for the DD methods based on Neumann problems, tend to scale as

$$ \kappa( P_{\text{DD}}^{-1} A ) = C \Big(1 + log \frac H h \Big)^{2(L-1)} $$

where $L$ is the number of levels. As long as aggressive coarsening is used, this is perhaps not so critical because $L \le 4$ should be capable of solving problems with more than $10^{12}$ degrees of freedom, but is certainly a concern. See this question for further discussion of this limitation.

Answer 2

This is an excellent writeup but I think saying that (multilevel) DD and MG have a lot in common is not accurate, or at least not useful. The methods are very different and I don't think that expertise in one is very useful in the other.

First, the two communities use different definitions of complexity: DD optimizes the condition number of the preconditioned systems and MG optimizes work/memory complexity. This is a big fundamental difference -- "optimality" has a totally different meaning in these two contexts. Things don't change when you add in parallel complexity (although you get a log term added in MG). The two communities are almost speaking different languages.

Second, MG has multilevel built into it and multilevel DD methods have all been developed with two level theory and implementations. This limits the space of coarse grid spaces that you can use in MG -- they must be recursive. For instance, You can not implement FETI in an MG framework. People do some multilevel DD methods as Jed mentioned but at least some of the current popular DD methods do not seem to be implementable recursively.

Third, I see the algorithms themselves, as practiced, as very different. Qualitatively speaking I would say that DD methods project onto domain boundaries and solve this interface problem. MG works directly with the native equations. Avoiding this projection allows MG to be applied to nonlinear and unsymmetrical problems easily. Although the theory all but goes away for nonlinear and unsymmetric problems they have worked for a lot of people. MG also explicitly decouples the problem into two parts: the coarse grid space for scaling and an iterative solver (the smoother) to solve the physics. This is critical in understanding and working with MG and is an attractive property to me.

Although theoretically the smoothers and coarse grid spaces are tightly coupled, in practice you can often swap different smoother in and out as an optimization parameter. As Jed mentioned point or vertex smoothers are popular and usually faster, but for challenging problems heavier smoothers can be useful. This plot is from my dissertation showing the solve time as a function of Poisson ratio for Jacobi, block Jacobi and "additive Schwarz" (overlapped). Its a little hard to read but at the highest Poisson ratio (0.499) overlapping Schwarz is about 2x faster than (vertex) Jocobi whereas it is about 3x slower at pedestrian Poisson ratios.

Answer 3

According to Jed's answer, MG uses moderate coarsening while DD uses rapid coarsening. I think this makes a difference when they are parallelized. There will be multiples of communications and synchronizations for MG to go through many levels of coarsening that are equivalent to a single coarsening of DD. Another point from Jed's answer is MG uses cheap smoother and DD uses strong smoother. Considering the two points, it has been reported that MG at coarse levels will have bad communication/computation ratios. So according to Amdahl's law, the parallel speedup is not good. A remedy of this is parallel coarse grids correction such as BPX preconditioner. Besides, MG can use DD as smoother as Adams pointed out, and MG can also be used within subdomains of DD. Based on the considerations Barker have pointed out, I guess using MG within DD is better, which exploits both parallelsim of DD and optimal complexity of MG.

Answer 4

I want to make one small addition to Jed's excellent answer, namely that the motivations behind the two approaches are (or at least were) different.

Domain decomposition is motivated as a technique for parallel computing. Especially for one-level methods DD is very natural to implement on a parallel machine - you divide the domain into pieces and give each piece to a different processor. In some sense the motivation behind DD is to divide arithmetic operations between processors.

Good parallel multigrid implementations exist, but it is often less natural to do in parallel. Instead, the motivation behind multigrid is to do less arithmetic operations in the first place.

Answer 5

@Jed Brown About the condition number of one-level DD, I think the size of the subdomain $H$ should be in the denominator to multiply with the mesh size $h$. It is common for the one-level method that using more subdomains needs more iterations. One can also refer to the book of Toselli and Widlund. For example, the condition number of the Schur complement is $\mathcal{O}(\frac{1}{Hh})$ from page 98 of the book.

Another point is about the condition number of optimized Schwarz methods. Actually, it improves the exponent of $h$ for one-level method (non-overlapping and overlapping) and in addition the exponent of $H$ for two-level method.