Extended-Domain Methods

• Extended-domain methods are reformulation strategies that embed physical problems into larger computational or spectral spaces to simplify complex geometries and interactions.
• They employ techniques like fictitious-domain, diffuse-boundary, spectral, and augmented Lagrangian approaches to handle irregular interfaces and improve numerical stability.
• These methods enhance simulation accuracy and reduce computational challenges in applications such as full-waveform inversion, fluid–structure interaction, and multirate system modeling.

An extended-domain method is a family of formulations in which a problem posed on a physical domain, a physical source space, or a physical set of active degrees of freedom is recast on a larger, simpler, or otherwise augmented domain, while the original physics is recovered by restriction, weak enforcement, stabilization, projection, or perturbative splitting. In the literature sampled here, the term covers several non-identical constructions: fictitious-domain and unfitted-mesh methods for PDEs on complex geometries, diffuse-domain and smoothed-boundary formulations, spectral extensions on enlarged boxes, wavefield- or source-extended formulations for full-waveform inversion, spatial extension of perturbative–numerical splitting in cosmological $N$-body simulation, and space–time domain decomposition with enhanced velocity spaces (Court et al., 2014, Yang et al., 2022, Tassev et al., 2015, 1803.02075).

1. Terminological scope and unifying principle

The phrase “extended-domain method” does not denote a single algorithm. In unfitted finite element and fictitious-domain formulations for Stokes or fluid–structure interaction, the actual fluid region is embedded in a larger simple computational domain $\mathcal{O}$, and the mesh is built on $\mathcal{O}$ rather than on the physical fluid domain $\mathcal{F}$; interface conditions are imposed weakly, often with Lagrange multipliers and stabilization (Court et al., 2013, Court et al., 2014). In spectral methods for elliptic PDEs, the physical domain $\Omega$ is embedded in an enlarged box $R$ or $\tilde I$, the basis is defined on that larger domain, and the PDE is enforced only on the physical region by collocation (1803.02075). In diffuse-domain formulations, the physical domain is represented by a smooth domain parameter $\psi$, and the PDE is solved on the whole computational box, with boundary conditions encoded by terms involving $\nabla\psi$ and $|\nabla\psi|$ (Yu et al., 2011).

In inverse problems, especially full-waveform inversion, the extension is not primarily geometric. Yang and Zhou formulate wavefield reconstruction inversion and extended source inversion by enlarging the optimization space from the model $\mathcal{O}$0 alone to $\mathcal{O}$1 or $\mathcal{O}$2, so that the wave equation is no longer enforced as a hard constraint at each iteration (Yang et al., 2022). Aghamiry and coauthors recast the same idea through an augmented Lagrangian and then move the extended variables from full wavefield space into data space, yielding a time-domain algorithm with data-sized dual variables rather than wavefield-sized multipliers (Gholami et al., 2020). Guo et al. retain the extended-source viewpoint and focus on practical approximation of the inverse data-domain Hessian used to deblur data residuals before back-propagation (Guo et al., 2023).

A plausible implication is that “extended domain” is best understood structurally rather than lexically: the physical equations remain central, but they are embedded into a larger representational setting where geometry, nonlinearity, interface motion, or multi-rate coupling become easier to handle.

2. Geometric embedding: fictitious domains, diffuse domains, and unfitted discretizations

A canonical geometric formulation appears in fictitious-domain Stokes methods. Let $\mathcal{O}$3 be a fixed simple domain, $\mathcal{O}$4 a rigid solid with boundary $\mathcal{O}$5, and $\mathcal{O}$6 the fluid domain. The mesh is built on $\mathcal{O}$7, not on $\mathcal{O}$8, and $\mathcal{O}$9 cuts the background mesh arbitrarily. In the XFEM interpretation, basis functions are multiplied by a Heaviside function $\mathcal{O}$0 that equals $\mathcal{O}$1 in $\mathcal{O}$2 and $\mathcal{O}$3 in $\mathcal{O}$4; degrees of freedom fully outside $\mathcal{O}$5 are removed, whereas cut-element degrees of freedom are retained as “virtual” DOFs (Court et al., 2013, Court et al., 2014).

In these Stokes formulations, the interface $\mathcal{O}$6 is localized by a level-set function, and the Dirichlet condition on $\mathcal{O}$7 is imposed by a Lagrange multiplier $\mathcal{O}$8. The multiplier is identified with the interface traction

$\mathcal{O}$9

so the method does not merely impose kinematic compatibility; it also provides the force per unit area exerted by the fluid on the structure. Because arbitrary cut configurations make the velocity–multiplier coupling delicate, the method adds an augmented-Lagrangian stabilization term

$\mathcal{F}$0

a Barbosa–Hughes-type device that improves the approximation of the normal trace of the Cauchy stress tensor and yields a discrete inf–sup condition under explicit assumptions on interface trace norms (Court et al., 2014).

The same extended-domain logic appears in other embedded formulations. Yu, Chen, and Thornton introduce a domain parameter $\mathcal{F}$1 that is approximately $\mathcal{F}$2 inside the physical domain, $\mathcal{F}$3 outside, and smoothly varying across a thin transition layer. For diffusion,

$\mathcal{F}$4

multiplication by $\mathcal{F}$5 and repeated use of product rules produce volumetric terms localized in the diffuse layer, with Neumann conditions represented through $\mathcal{F}$6 and Dirichlet conditions through terms proportional to $\mathcal{F}$7 (Yu et al., 2011). The unit normal is recovered as $\mathcal{F}$8, so geometry and boundary conditions are both encoded by $\mathcal{F}$9.

Optimization-based embedding takes yet another route. The Smooth Selection Embedding Method embeds $\Omega$0 in a simple periodic box $\Omega$1, defines the unknown on all of $\Omega$2, imposes the original PDE only on $\Omega$3 and $\Omega$4, and then selects from the resulting affine family of discrete extensions the one minimizing a Sobolev-type norm induced by $\Omega$5 (Agress et al., 2018). This avoids explicit extensions of the data and leaves the PDE unmodified on the physical domain.

Unfitted hybridizable DG methods for elasticity fit the same pattern. The X-HDG method uses a background triangulation not fitted to the interface or boundary and replaces standard element spaces on cut elements by split spaces such as

$\Omega$6

with separate polynomial fields on each side of the interface. The method is HDG because stress and displacement are represented inside elements while numerical traces on faces and on the interface or boundary provide the globally coupled unknowns (Han et al., 2020).

3. Spectral extended-domain constructions and their stability theory

A different branch of the subject uses the extended domain primarily to support global spectral bases. In the Laplacian-eigenbasis collocation method, the physical domain $\Omega$7 is embedded in a rectangle $\Omega$8, which is itself enlarged to

$\Omega$9

The basis is defined by the Dirichlet Laplacian eigenproblem on $R$0,

$R$1

with separated sine eigenfunctions on the extended box. The approximate solution is expanded in those eigenfunctions and the PDE is enforced by collocation on the physical region (1803.02075).

The rationale is regularity preservation. Embedded-boundary strategies often lose smoothness when a solution is extended from a curved $R$2 to a simple box by ad hoc continuation. By instead constructing the approximation in a smooth basis on the enlarged box and restricting it back to $R$3, the method seeks to retain the exponential decay of Fourier-type coefficients for analytic solutions. The 2018 analysis ties the collocation error to the Fourier best approximation through a Lebesgue constant and derives, for the interpolation setting under consideration, an upper bound

$R$4

showing directly how domain stretching influences stability and accuracy (1803.02075).

Subsequent analysis reveals that this picture is operator-dependent. “On the Lebesgue Constant of Extended-Domain Spectral Methods for Elliptic PDEs” proves a stability dichotomy: for the self-adjoint Poisson equation, the method is unstable, with a Lebesgue constant that grows super-polynomially, whereas for the convection–diffusion operator $R$5 the Lebesgue constant is polynomially bounded, with a proved $R$6 bound in one dimension under an explicit coercivity condition on $R$7 or $R$8 (Wu, 10 Sep 2025). This is not a generic statement about all spectral embeddings; it is a result for the particular extended-domain collocation construction analyzed there. Still, it materially changes the interpretation of earlier empirical success, because it separates smooth-solution accuracy from worst-case $R$9-stability.

The Projection Extension method places the spectral extension idea into an explicitly least-squares framework. One embeds $\tilde I$0 in a simple $\tilde I$1, chooses a global spectral basis $\tilde I$2 on $\tilde I$3, and then minimizes a residual functional that combines PDE mismatch in $\tilde I$4 with boundary mismatch on $\tilde I$5. In its forcing-extension form, the coefficients $\tilde I$6 satisfy normal equations built from

$\tilde I$7

where $\tilde I$8 (Qadeer et al., 2021). The resulting method is stable under iterative application and is used for elliptic, parabolic, Newtonian fluid, and viscoelastic problems on arbitrary domains, including disconnected or non-smooth ones.

4. Extended domains in full-waveform inversion

In full-waveform inversion, the “domain” being extended is the optimization space. Classical FWI solves

$\tilde I$9

with the wave equation enforced exactly. Extended-domain variants enlarge the admissible set by allowing non-physical wavefields or sources and penalizing those deviations rather than eliminating them outright (Yang et al., 2022).

Wavefield Reconstruction Inversion replaces the hard constraint by

$\psi$0

whereas Extended Source Inversion enforces $\psi$1 exactly but treats the source $\psi$2 as unknown and penalizes $\psi$3. Yang and Zhou show that both can be written as weighted least-squares FWI in data space: $\psi$4 with model-dependent weightings built from the modeling operator $\psi$5 (Yang et al., 2022). They emphasize that this enlarges the basin of attraction and mitigates cycle skipping, but does not fundamentally remove local minima.

The main algorithmic obstacle in time-domain WRI and ESI is storage. The extended source $\psi$6, wavefield $\psi$7, and adjoint field $\psi$8 live in full space–time, and matrix-free conjugate-gradient solutions of the normal equations still require repeated forward and adjoint propagations together with access to large time histories. Yang and Zhou identify the huge storage demand to store time-domain wavefields through iterations as the dominant practical bottleneck, and discuss two workarounds—extracting sparse frequential wavefields and solving in a data-space formulation—while concluding that these options should be explored more intensively for tractable workflows (Yang et al., 2022).

Aghamiry and collaborators develop an augmented-Lagrangian variant in the time domain. Starting from the penalty formulation

$\psi$9

they derive a data-space system

$\nabla\psi$0

and use it to formulate Data Reconstruction Inversion, in which the only stored dual variable is data sized (Gholami et al., 2020). The practical DRI iteration requires two forward and two backward propagations per outer iteration, a cost comparable to conventional adjoint-state FWI, while avoiding wavefield-sized multipliers.

Guo et al. remain within time-domain extended-source FWI but focus on the inverse data-domain Hessian

$\nabla\psi$1

which maps FWI residuals to the deblurred residuals whose back-propagation yields source extensions (Guo et al., 2023). Their implementation approximates $\nabla\psi$2 by matching filters in Fourier and short-time Fourier domains and refines the result with conjugate-gradient iterations when necessary. They combine this with the augmented Lagrangian method, multiscale frequency continuation, grid coarsening, and total-variation regularization for large-contrast reconstruction (Guo et al., 2023).

5. Spatial and space–time extensions of operator splitting and domain decomposition

The term also appears when an existing decomposition is generalized from one coordinate domain to another. In sCOLA, the original COLA method already split temporal evolution into an LPT part and an $\nabla\psi$3-body residual part. sCOLA extends that split to the spatial domain by restricting the $\nabla\psi$4-body solver to a small spatial COLA box while using Lagrangian perturbation theory to encode the far field generated by the rest of the universe (Tassev et al., 2015). The central modification is the replacement of the full inverse Laplacian by a COLA-box inverse Laplacian,

$\nabla\psi$5

so that the near field is computed numerically and the far field is injected perturbatively. This lets a small simulation volume reproduce the results of a standard $\nabla\psi$6-body run for the same small volume embedded inside a much larger simulation, provided a sufficiently large buffer surrounds the region of interest (Tassev et al., 2015).

In porous-media flow, the enhanced velocity mixed finite element method extends from spatial domain decomposition to space–time domain decomposition. The global domain $\nabla\psi$7 is partitioned into non-overlapping space–time blocks $\nabla\psi$8, each with its own spatial and temporal resolution. A naive direct sum of local $\nabla\psi$9-type velocity spaces is not $|\nabla\psi|$0-conforming across non-matching interfaces, so the method adds enhanced velocity basis functions on interface subelements created by intersecting the traces of the neighboring space–time meshes (Singh et al., 2018). The resulting space $|\nabla\psi|$1 is $|\nabla\psi|$2-conforming, preserves local mass conservation, and supports a monolithic fully coupled system that does not require subdomain iterations (Singh et al., 2018).

A similar structural theme appears in unfitted HDG elasticity. There the extension is not temporal but geometric and algebraic: polynomial spaces on cut elements are split by characteristic functions $|\nabla\psi|$3, and the interface or boundary is treated by HDG numerical traces and stabilization terms of order $|\nabla\psi|$4 rather than by mesh fitting (Han et al., 2020). This suggests that “extended-domain” can denote not just a larger box, but any systematic enlargement of the representational space needed to decouple mesh generation from physical interfaces.

6. Advantages, limitations, and recurrent misconceptions

Across these literatures, the principal advantage is geometric or algorithmic decoupling. Embedded-domain PDE solvers avoid body-fitted meshing, remeshing, or interface-aligned quadrature rules on moving geometries; Stokes fictitious-domain methods on $|\nabla\psi|$5 allow moving rigid solids to be handled without remeshing and provide direct access to the hydrodynamic force through the stabilized multiplier $|\nabla\psi|$6 (Court et al., 2014). Smoothed-boundary methods replace sharp surface conditions by volumetric terms on simple grids, which is particularly useful for image-based microstructures and evolving boundaries (Yu et al., 2011). Spectral embeddings make global high-order bases available on irregular domains (1803.02075), while PE avoids domain decomposition and remains stable under iterative use (Qadeer et al., 2021). Inverse-problem extensions enlarge the linear regime and can reduce starting-model sensitivity (Yang et al., 2022, Gholami et al., 2020). Space–time EVMFE preserves local mass conservation on non-matching grids and enables multirate, monolithic simulation (Singh et al., 2018). sCOLA localizes gravity for the $|\nabla\psi|$7-body solver, supports zoom-in calculations, and can make the code embarrassingly parallel (Tassev et al., 2015).

The limitations are equally recurrent. Unfitted and fictitious-domain finite element methods require cut-element integration, interface quadrature, and stabilization choices; in the Stokes multiplier formulation, $|\nabla\psi|$8 must be “large enough for stabilizing, but not too much in order to keep the system coercive” (Court et al., 2014). Diffuse-domain methods introduce interface-thickness error and terms involving $|\nabla\psi|$9 or $\mathcal{O}$00, which require regularization near $\mathcal{O}$01 (Yu et al., 2011). Optimization-based embeddings such as SSEM can become ill-conditioned for large smoothing orders and need preconditioning or explicit factorization strategies (Agress et al., 2018). In extended-domain FWI, storage rather than raw floating-point cost is the dominant obstacle in the time domain (Yang et al., 2022). In the spectral collocation setting, the 2025 Lebesgue-constant analysis shows that excellent $\mathcal{O}$02 convergence for smooth Poisson problems does not imply uniform $\mathcal{O}$03-stability (Wu, 10 Sep 2025).

A common misconception is that extended-domain methods are defined by one enforcement mechanism. The literature here includes Lagrange multipliers, augmented Lagrangians, least-squares projections, diffuse-interface source terms, perturbative–numerical operator splitting, and enhanced interface basis functions. Another misconception is that extension automatically resolves nonconvexity: in FWI, WRI and ESI can enlarge the basin of attraction, but they remain weighted least-squares data-misfit formulations and can still cycle-skip (Yang et al., 2022). A further misconception is that a larger computational domain is always the central ingredient; in some branches, especially inverse problems, the decisive extension is of the variable space or source space rather than of physical geometry.

Taken together, these works support a broad but technically precise characterization: an extended-domain method is a reformulation strategy that embeds the original problem into a larger computational, geometric, temporal, spectral, or optimization setting so that discretization, coupling, or inversion becomes tractable, while the physical solution is recovered by restriction and by explicit mechanisms that control the discrepancy introduced by the extension.