Source-linked AI summary
A rapidly converging domain decomposition method for the Helmholtz equation
Christiaan C. Stolk
TL;DR
The paper addresses efficient solution of heterogeneous 2-D and 3-D Helmholtz systems, especially with many thin subdomains. It introduces PML-based transmission conditions and multiplicative sweeps, achieving rapid convergence and near-linear multifrontal cost in the reported setting. The approach remains subject to PML, right-hand-side, cavity, and parallelism constraints.
Problem
The paper targets Helmholtz domain decomposition with many thin subdomains, where iteration counts can otherwise grow with the number of subdomains.
Method
The method uses non-reflecting PML transmission conditions, ordered upward and downward subdomain solves, and GMRES preconditioning.
Results
Fewer than 10 iterations generally reduced the residual by 10^-6, while multifrontal subdomain solves produced near-linear behavior with iteration counts essentially independent of problem size and subdomain count.
Takeaways & Limitations
The method offers a rapidly converging Helmholtz solver for thin-layer decompositions and is one of the few reported methods with near-linear behavior in this setting.
Takeaways & Limitations
Near-linear performance requires sufficiently many right-hand sides, increasing PML thickness, and does not hold for cavity problems; the multiplicative method is inherently sequential.
Abstract
from arXiv · showhide
A new domain decomposition method is introduced for the heterogeneous 2-D and 3-D Helmholtz equations. Transmission conditions based on the perfectly matched layer (PML) are derived that avoid artificial reflections and match incoming and outgoing waves at the subdomain interfaces. We focus on a subdivision of the rectangular domain into many thin subdomains along one of the axes, in combination with a certain ordering for solving the subdomain problems and a GMRES outer iteration. When combined with multifrontal methods, the solver has near-linear cost in examples, due to very small iteration numbers that are essentially independent of problem size and number of subdomains. It is to our knowledge only the second method with this property next to the moving PML sweeping method.
1. Introduction
The paper introduces PML-based transmission conditions and a multiplicative double-sweep domain decomposition method for heterogeneous 2-D and 3-D Helmholtz problems. Numerical examples report rapid convergence largely independent of domain size and subdomain count, with near-linear cost when paired with multifrontal subdomain solves.
- 1. Introduction: The method constructs PML-based transmission conditions that are non-reflecting and match outgoing waves from one subdomain to incoming waves in its neighbor.Single-layer potentials and PML boundary layers provide the transmission mechanism.
- 1. Introduction: The domain is divided into many thin quasi-2-D layers, while consecutive upward and downward subdomain solves let information cross the full domain in two solves per subdomain.The resulting multiplicative procedure is used as a GMRES preconditioner.
- 1. Introduction: Fewer than 10 iterations generally reduced the residual by 10^-6, with iteration counts essentially independent of domain size and subdomain number.The experiments varied grid distance, frequency with constant ωh, and the number of subdomains.
- 1. Introduction: Multifrontal solves for quasi-2-D subdomains yield near-linear behavior, with 3-D factorization cost O(N^4/3) instead of O(N^2) for direct multifrontal treatment.The solve cost is O(N log N), ignoring additional logarithmic factors noted elsewhere.
- 1. Introduction: The near-linear claim requires increasing PML thickness, sufficiently many right-hand sides, and excludes cavity problems; the multiplicative method is also inherently sequential.The estimated PML thickness growth is O(log N), adding logarithmic cost factors.
- 1.2. Results: The method is theoretically exact for a constant-coefficient strip with perfect PML, solving it in one upward and one downward sequence of subdomain solves.The upward and downward sequences can be performed simultaneously with careful handling at their crossing point.
2. The method
The method partitions a rectangular Helmholtz domain into subdomains with PML-extended operators and uses ordered forward and backward solves to construct a GMRES preconditioner. A reduced interface-space formulation can restrict right-preconditioned GMRES to boundary-adjacent grid layers.
- PML-modified subdomain operators extend the core subdomains at internal boundaries to provide the transmission framework.The domain is divided along the x-axis, and each core subdomain is padded with PML layers.
- Forward subdomain solves are performed consecutively, with each solve using information from the previously solved neighbor.The subdomain solutions are assembled into a global function piecewise over the core subdomains.
- A backward sequence solves subdomain restrictions of the residual, producing a correction that is combined with the forward solution.The correction is computed in reverse subdomain order and added to the first solution.
- The map f 7→v + w defines the approximate solution and can serve as either a left or right preconditioner for GMRES.This map is the central preconditioning operator for the iterative solution method.
- Right-preconditioned GMRES can be restricted to 2(J −1) interface-adjacent grid layers because the residual is nonzero only near internal boundaries.The original solution is recovered from the reduced-space solution using the right-preconditioned formulation.
3. Theoretical results
For a constant-coefficient strip with perfectly non-reflecting PML boundaries, upward and downward subdomain sweeps recover the full Helmholtz solution in one iteration. The analysis derives this result through Fourier-domain transmission conditions and induction.
- 3.1. Multiplicative domain decomposition with upward and downward sweeps in 1-D: For the upward sweep, interface conditions propagate the neighboring subdomain data through successive solves, with the subdomain Helmholtz equation enforced on each interval.The construction is established by induction from the interface condition at the preceding boundary.
- 3.1. Multiplicative domain decomposition with upward and downward sweeps in 1-D: The concurrently executed upward and downward sweeps also produce the solution of the original problem.The proof uses the solution formula and verifies the result interval by interval.
- 3.2. PML based transmission on the strip: In the Fourier domain, perfect non-reflecting boundaries are represented mode by mode, with incoming-wave data encoded separately from homogeneous boundary conditions.The spatial-domain operator is non-local after inverse Fourier transformation, which distinguishes these conditions from Robin conditions in higher dimensions.
- 3.1. Multiplicative domain decomposition with upward and downward sweeps in 1-D: The downward sweep solves residual-driven subdomain problems and imposes matching interface conditions on the corrected field.The proof combines the forward solution with the backward correction and completes the argument by induction.
4. Numerical results
Numerical tests in 2-D and 3-D show rapid convergence across heterogeneous and constant media, with iteration counts largely insensitive to domain size and subdomain count. Increasing PML thickness can maintain this behavior as grids grow, whereas Robin transmission conditions lose the same scaling.
- Experimental setup: The experiments use sparse multifrontal solvers, with MUMPS for 3-D subdomain systems and UMFPACK for the 2-D examples.The large computations were implemented under MPI, while the Robin comparison used Matlab.
- Example 1: Marmousi: 4–9 iterations reduce the residual by 10^-6 in the 2-D Marmousi tests, with only mild dependence on grid size and subdomain count.The experiments vary h from 1 to 16 m and use 3 to 300 subdomains while keeping hω constant.
- Example 2: random medium in 3-D: 6–8 iterations suffice for the 3-D random-medium examples, and wpml = 4 is sufficient in these tests.The examples range from 100^3 to 400^3 and use layers of thickness 10 grid points excluding PML layers.
- Example 3: constant medium in 3-D: Increasing wpml somewhat reduces iterations in larger constant-medium examples, while fixed wpml yields increasing iteration counts as grid size grows.The reported data indicate that wpml must increase with grid size to preserve approximately constant iteration counts.
- Comparison between Robin and PML-based transmission conditions: Robin transmission conditions lose the PML method’s scaling: iteration counts grow roughly linearly with Nx in 2-D, despite subdomain solves being about four times cheaper.The comparison uses constant and random media with layers 10 grid points thick.
5. Discussion
The method combines PML-based transmission conditions with multiplicative double sweeps and achieves fast convergence even for thin-layered subdomains. Its sequential structure limits parallelism, while Robin alternatives do not reproduce the near-linear scaling.
- Main discussion: The proposed method has remarkably fast convergence, including for thin-layered subdomains, when paired with sparse direct subdomain solvers.The paper focuses on this combination in its discussion of the method’s practical behavior.
- Comparison with Robin conditions: PML-based transmission conditions retain near-linear cost in the reported setting, unlike the compared Robin double-sweep method.The Robin method can be attractive for small J but does not have near-linear cost like the PML-based method.
- Parallelism and limitations: Multiplicative domain decomposition is intrinsically sequential, so parallel performance requires parallel factorization or pipelining multiple right-hand sides.The pipelined strategy is especially relevant when factorization setup costs can be amortized across multiple right-hand sides.
- Broader applicability: The techniques are expected to apply to time-harmonic Maxwell and linear-elastic wave equations because their solutions behave similarly to Helmholtz solutions.This is stated as an expectation rather than demonstrated numerically in the paper.