Quasi-Trefftz Methods for PDE Discretization

• Quasi-Trefftz methods are high-order Galerkin discretizations that use local functions tailored to cancel the PDE residual up to a prescribed Taylor order.
• They generalize classical Trefftz methods to variable-coefficient and complex PDE settings, reducing local degrees of freedom while maintaining strong approximation properties.
• Their construction via polynomial bases and generalized plane waves ensures efficient residual cancellation, dimension reduction, and stable discretizations across various operator classes.

Quasi-Trefftz methods are high-order Galerkin discretizations built from local functions that are tailored to the governing partial differential operator but are not, in general, exact local solutions. They generalize classical Trefftz methods—whose elementwise trial and test functions belong to the exact kernel of the PDE—to settings with variable coefficients or other features for which closed-form local solutions are generally not available. The defining idea is to replace exact local solvability by high-order local residual cancellation: on each element, one uses functions whose PDE residual has a vanishing truncated Taylor expansion around a point, equivalently $Lv = O(|x-x_0|^q)$ for a prescribed order $q$. This produces problem-dependent discrete spaces with strong approximation properties and reduced local dimension, and it has been developed for variable-coefficient Helmholtz equations, diffusion–advection–reaction operators, space-time wave equations, first-order Helmholtz systems, anisotropic reaction–diffusion problems, and Oseen linearizations (Fontana et al., 13 Aug 2025, Imbert-Gerard, 24 May 2025, Lederer et al., 2024).

1. Definition and mathematical scope

Classical Trefftz methods use local approximation spaces made of elementwise exact solutions of the target linear homogeneous PDE. For constant-coefficient wave problems, such spaces are readily available: for the homogeneous Helmholtz equation,

$\Delta u + \kappa^2 u = 0,$

classical Trefftz bases are plane waves $e^{i\kappa d\cdot(x-x_0)}$, and formulations such as Trefftz discontinuous Galerkin and ultra-weak variational formulations enforce interface conditions weakly through skeleton terms (Fontana et al., 13 Aug 2025). Quasi-Trefftz methods retain the same medium-aware philosophy, but they relax exact local satisfaction of the PDE when coefficients vary in space.

For a linear differential operator $L$ of order $\gamma$ and a point $x_0$ in an element, a function $v$ is quasi-Trefftz of order $q$ at $x_0$ if the Taylor truncation of $q$0 around $q$1 vanishes up to degree $q$2. In the notation used across the literature,

$q$3

For polynomial constructions of degree $q$4, the natural choice is often $q$5 (Fontana et al., 13 Aug 2025, Imbert-Gerard, 24 May 2025).

This definition has been formulated for general linear scalar differential operators with smooth coefficients, for second-order elliptic diffusion–advection–reaction operators, for the scalar Helmholtz equation with variable refractive index, for the convected Helmholtz equation, for first-order formulations of Helmholtz, and for space-time acoustic wave equations with piecewise-smooth wavespeed (Imbert-Gérard et al., 2024, Imbert-Gérard et al., 10 Sep 2025, Imbert-Gérard et al., 2020). In embedded formulations, the same idea appears in a weak form: the local residual is not required to vanish pointwise, but only to be orthogonal to a local test space. This weak enforcement is also described as quasi-Trefftz in recent DG constructions for Helmholtz and reaction–diffusion problems (Stocker et al., 13 Mar 2026, Gómez et al., 2 Jun 2026).

A recurring consequence is dimension reduction. For a scalar operator of order $q$6, the local polynomial quasi-Trefftz space on $q$7 has dimension

$q$8

so for second-order operators one obtains $q$9 basis functions in two dimensions and $\Delta u + \kappa^2 u = 0,$0 in three dimensions, instead of the full polynomial dimensions $\Delta u + \kappa^2 u = 0,$1 and $\Delta u + \kappa^2 u = 0,$2 (Imbert-Gerard, 24 May 2025, Imbert-Gérard et al., 2024).

2. Local constructions: polynomial, generalized plane wave, and first-order bases

Two principal families dominate the current theory: polynomial quasi-Trefftz spaces and generalized plane wave spaces. In the polynomial case, one seeks $\Delta u + \kappa^2 u = 0,$3 such that

$\Delta u + \kappa^2 u = 0,$4

The resulting operator from $\Delta u + \kappa^2 u = 0,$5 to $\Delta u + \kappa^2 u = 0,$6 is linear, and its graded structure with respect to homogeneous polynomial degree yields a block-triangular system. The coefficients are computed layer by layer from the principal part $\Delta u + \kappa^2 u = 0,$7 of the operator at $\Delta u + \kappa^2 u = 0,$8, provided $\Delta u + \kappa^2 u = 0,$9 is surjective between the relevant homogeneous polynomial layers (Fontana et al., 13 Aug 2025, Imbert-Gerard, 24 May 2025).

The general linear-algebraic formulation rests on the decomposition

$e^{i\kappa d\cdot(x-x_0)}$0

With respect to this grading, the quasi-Trefftz operator takes the form $e^{i\kappa d\cdot(x-x_0)}$1, where the principal operator acts diagonally on degrees and the remainder raises degree. This graded structure implies surjectivity under a nontrivial principal part and gives an explicit constructive right inverse through a finite Neumann series (Imbert-Gerard, 24 May 2025). The same mechanism underlies the basis-construction algorithms for elliptic diffusion–advection–reaction equations and for scalar linear equations in arbitrary dimension and order (Imbert-Gérard et al., 2024).

Generalized plane waves (GPWs) use a phase-based or amplitude-based ansatz. For the scalar Helmholtz operator $e^{i\kappa d\cdot(x-x_0)}$2, the phase-based GPW uses

$e^{i\kappa d\cdot(x-x_0)}$3

with $e^{i\kappa d\cdot(x-x_0)}$4 a polynomial centered at $e^{i\kappa d\cdot(x-x_0)}$5. The residual factors algebraically: $e^{i\kappa d\cdot(x-x_0)}$6 so the quasi-Trefftz condition becomes

$e^{i\kappa d\cdot(x-x_0)}$7

Although the operator is nonlinear in $e^{i\kappa d\cdot(x-x_0)}$8, the same graded, layer-by-layer construction applies after splitting off the principal linear part (Fontana et al., 13 Aug 2025).

Amplitude-based GPWs were introduced to avoid placing high-degree polynomials inside the exponential. Their ansatz keeps a linear phase and moves higher-order terms into a polynomial amplitude,

$e^{i\kappa d\cdot(x-x_0)}$9

with $L$0 chosen so that the residual vanishes to order $L$1. This construction preserves the quasi-Trefftz property and the interpolation order while improving pre-asymptotic behavior relative to phase-based GPWs (Imbert-Gerard, 2020). For the three-dimensional convected Helmholtz equation, amplitude-based GPWs, phase-based GPWs, and polynomial quasi-Trefftz functions were all constructed by explicit layer-wise recursions, using the local second-order coefficient matrix and a direction-selection rule derived from its orthogonal-diagonal decomposition (Imbert-Gerard et al., 2022).

The same Taylor-based principle extends to first-order systems. For a first-order formulation of Helmholtz in pressure–velocity variables, the local space is built from triples $L$2 with pressure in $L$3 and velocity in $L$4, constrained by exact first-order relations with constant coefficients and a truncated Taylor condition on the remaining variable-coefficient residual. The space can be characterized equivalently by a scalar quasi-Trefftz pressure polynomial together with the exact recovery $L$5, and in two dimensions its dimension is again $L$6 (Imbert-Gérard et al., 10 Sep 2025).

3. Discretization frameworks and global coupling

Once local quasi-Trefftz spaces are available, they are inserted into DG-, UWVF-, HDG-, VEM-, or related skeleton-based discretizations. In the standard elementwise construction, one fixes a point $L$7 in each element $L$8, builds a local basis $L$9, and defines

$\gamma$0

The global space is discontinuous across elements, and coupling is enforced weakly on the mesh skeleton through numerical fluxes or interface forms, exactly as in Trefftz DG or UWVF (Fontana et al., 13 Aug 2025).

For variable-coefficient Helmholtz problems, quasi-Trefftz DG schemes use polynomial quasi-Trefftz bases or GPWs inside an interior-penalty and upwind DG framework. In the diffusion–advection–reaction setting, the quasi-Trefftz DG method employs the same SIPG-upwind bilinear forms as a standard DG method, but with local spaces defined by Taylor constraints. Because the trial space for inhomogeneous problems is affine, a local quasi-Trefftz particular solution is constructed elementwise and the global solve is performed for the homogeneous correction (Imbert-Gérard et al., 2024, Perinati, 2023).

Space-time quasi-Trefftz methods apply the same principle to hyperbolic systems. For the acoustic wave equation with piecewise-smooth wavespeed, the local space consists of space-time polynomials $\gamma$1 such that the Taylor polynomial of the wave operator vanishes to degree $\gamma$2 at the element center. The DG method combines space-like upwind traces, time-like averages and jumps, and residual stabilization terms weighted by $\gamma$3 (Imbert-Gérard et al., 2020). Classical space-time Trefftz DG, by contrast, uses exact local wave solutions and eliminates volume terms entirely; its stability and quasi-optimality analysis supplies one of the analytical templates later adapted in quasi-Trefftz settings (Kretzschmar et al., 2015).

A different line of development is the embedded Trefftz or embedded quasi-Trefftz approach. Here one starts from a standard polynomial DG space $\gamma$4 and defines the local constrained space through moment conditions on the residual. For Helmholtz,

$\gamma$5

and the local space is the kernel of these constraints. This avoids explicit construction of plane waves or other analytic Trefftz functions and yields a reduced global system assembled with a simple SIPDG form (Stocker et al., 13 Mar 2026). The same embedded pattern appears for anisotropic reaction–diffusion on quadrilateral meshes, where the residual is enforced orthogonally against a carefully chosen local test space vanishing on the longer edges, and for the Oseen problem, where a local complement space is constructed so that the operator is stably invertible on that complement (Gómez et al., 2 Jun 2026, Stocker et al., 11 Jun 2026).

A unified local–global interpretation has been given for these methods. In that framework, one decomposes the discrete space into a local component and a Trefftz-like component, solves local residual equations elementwise, and couples only the reduced global space. This places exact Trefftz, embedded Trefftz, and Taylor-based quasi-Trefftz methods in a single T-coercive setting (Lederer et al., 2024).

4. Approximation theory, dimension reduction, and conditioning

The main approximation mechanism is Taylor matching. If the coefficients of $\gamma$6 are smooth enough near $\gamma$7 and the exact solution $\gamma$8 is sufficiently regular, then its Taylor polynomial of suitable degree belongs to the local quasi-Trefftz space. For general scalar operators,

$\gamma$9

and the best approximation error in $x_0$0-seminorms has the same order as for full polynomial spaces (Imbert-Gérard et al., 2024). For GPWs, the local space can approximate any smooth solution with order $x_0$1 in the Taylor sense: $x_0$2 for $x_0$3 near $x_0$4, and the same order transfers to DG or UWVF schemes when quadrature and fluxes are chosen consistently (Fontana et al., 13 Aug 2025).

Dimension reduction is one of the central structural advantages. For second-order operators, polynomial quasi-Trefftz spaces scale like $x_0$5 rather than $x_0$6. In two dimensions this gives $x_0$7 local basis functions, and in three dimensions $x_0$8, matching the dimensions of degree-$x_0$9 plane-wave Taylor spans (Imbert-Gerard, 24 May 2025, Fontana et al., 13 Aug 2025). For the first-order Helmholtz formulation, the pressure–velocity quasi-Trefftz space in two dimensions also has dimension $v$0, again reducing the number of local degrees of freedom relative to the ambient polynomial product space (Imbert-Gérard et al., 10 Sep 2025). In anisotropic embedded reaction–diffusion DG on quadrilaterals, the locally constrained space has dimension $v$1, compared to $v$2 for the underlying tensor-product DG space (Gómez et al., 2 Jun 2026).

Conditioning is more nuanced. Wave-like bases, including plane waves and GPWs, can suffer from near-linear dependence as the number of directions or the approximation order increases. The GPW literature states that conditioning depends on the number and distribution of directions, element size $v$3, and coefficient variation; evenly distributed directions and moderate $v$4 improve conditioning, while strong heterogeneities or very high $v$5 may introduce small singular values (Fontana et al., 13 Aug 2025). Amplitude-based GPWs were proposed specifically to tame the poor pre-asymptotic behavior of phase-based GPWs by avoiding high-degree polynomials within an exponential (Imbert-Gerard, 2020).

For the three-dimensional convected Helmholtz equation, numerical condition numbers grow rapidly for plane waves and both GPW families, whereas the polynomial quasi-Trefftz basis remains comparatively well-conditioned. The paper explicitly states that the polynomial basis shows significant promise because it does not suffer from the ill-conditioning issue inherent to wave-like bases (Imbert-Gerard et al., 2022). Related conditioning remedies include orthonormalization of local bases, mild regularization of flux parameters, edgewise orthogonalization-and-filtering in virtual element methods, and spectral filtering or coarse-search-space restriction in VTCR-like solvers (Fontana et al., 13 Aug 2025, Mascotto et al., 2018, Kovalevsky et al., 2016).

5. Applications and operator classes

Quasi-Trefftz methods were introduced, in the case of the Helmholtz equation, to handle wave propagation in inhomogeneous media. The canonical models include variable-index Helmholtz equations,

$v$6

more general second-order operators

$v$7

and convected Helmholtz equations arising in aeroacoustics and plasmas (Fontana et al., 13 Aug 2025). GPWs were originally developed for plasma reflectometry modeled by an O-mode equation with variable $v$8 changing sign, so that propagative and evanescent regions coexist; complex phases or directions naturally capture exponential decay in evanescent regions (Fontana et al., 13 Aug 2025).

The convected Helmholtz equation has been a particularly important testbed. In three dimensions, three families of quasi-Trefftz functions—two generalized-plane-wave families and one polynomial family—were constructed for the convected Helmholtz operator, and all three were shown to satisfy the same local approximation theorem (Imbert-Gerard et al., 2022). The abstract GPW framework also extends to convected Helmholtz by freezing the relevant coefficients at $v$9 in the principal operator and placing the variable-coefficient and nonlinear terms in the graded remainder (Fontana et al., 13 Aug 2025).

Elliptic problems form a second major application area. Polynomial quasi-Trefftz DG has been developed for variable-coefficient diffusion–advection–reaction equations in two and three spatial dimensions, with stability, consistency, quasi-optimality, and high-order convergence established in DG norms (Imbert-Gérard et al., 2024, Perinati, 2023). The same local Taylor-based machinery applies to general scalar linear equations with smooth coefficients and source terms, which broadens the scope beyond wave propagation in the narrow sense (Imbert-Gerard, 24 May 2025).

Embedded quasi-Trefftz methods extend the paradigm to problems for which exact local Trefftz spaces are trivial or difficult to compute. Reaction–diffusion on anisotropic, possibly curved quadrilateral elements is treated by imposing residual orthogonality constraints inside a tensor-product DG space, with anisotropic a priori error estimates and reduced global system size (Gómez et al., 2 Jun 2026). For the Oseen problem, the local quasi-Trefftz constraint is combined with a specially designed complement space and a reduced velocity-only formulation, providing stability and quasi-optimality in standard DG norms (Stocker et al., 11 Jun 2026).

Interface and heterogeneity problems motivate additional variants. For Helmholtz with piecewise constant wave number, a nonconforming Trefftz virtual element method remains exactly Trefftz on elements aligned with material interfaces, while on elements cut by the interface and endowed with an artificial local wave number the method becomes quasi-Trefftz in the sense of a small but nonzero residual (Mascotto et al., 2018). For periodic layered scattering, plane-wave TDG and T-matrix methods are strictly Trefftz when the coefficient is constant per element, but the same source explicitly states that quasi-Trefftz variants would add interior residual stabilizations in variable-$q$0 elements (Monforte et al., 7 Jul 2026).

6. Terminology, neighboring methods, and recurring misconceptions

The term “quasi-Trefftz” is not completely uniform across the literature. In its most precise and now common meaning, it denotes local functions whose PDE residual is zero only up to a prescribed Taylor order or only in a projected weak sense (Fontana et al., 13 Aug 2025, Imbert-Gerard, 24 May 2025). This is the sense used in polynomial quasi-Trefftz spaces, GPWs, first-order Helmholtz constructions, and embedded residual-constrained DG spaces (Imbert-Gérard et al., 10 Sep 2025, Stocker et al., 13 Mar 2026).

A first misconception is to identify every Trefftz-like method with quasi-Trefftz. The nonconforming Trefftz virtual element method for the Helmholtz problem is explicitly “strictly Trefftz”: its local spaces are contained in the kernel of the Helmholtz operator, and the virtual character lies in the projector-and-stabilization machinery rather than in any residual relaxation (Mascotto et al., 2018). Likewise, the plane-wave TDG method for layered periodic scattering is a pure Trefftz method with no interior residual terms (Monforte et al., 7 Jul 2026).

A second misconception concerns embedded methods. In embedded Helmholtz and reaction–diffusion DG, the local PDE is enforced through moment orthogonality rather than pointwise exactness. Those papers call the resulting spaces quasi-Trefftz because the residual is only zero in a weak sense against a local test space (Stocker et al., 13 Mar 2026, Gómez et al., 2 Jun 2026). This usage is consistent with the broader Taylor-based definition, but it emphasizes formulation-level relaxation rather than explicit local approximate solutions.

A third ambiguity appears in indirect wave methods such as the Variational Theory of Complex Rays. There the elementwise fields are exact Herglotz-wave solutions of the homogeneous Helmholtz equation, so the local approximation is Trefftz. The “quasi” aspect enters instead through the construction of a reduced, filtered search subspace and the coarse estimation of an energy norm used for spectral selection (Kovalevsky et al., 2016). In other words, the local PDE satisfaction is exact, but the solver operates in an approximate coarse space.

These distinctions matter because they separate three different mechanisms: exact local kernels, strong approximate local kernels defined by Taylor residual cancellation, and weakly constrained local kernels embedded in larger polynomial spaces. The unifying framework for Trefftz-like discretization methods shows that all three can be analyzed through the same decomposition into local and global problems, but it does not erase the conceptual difference between exact Trefftz and quasi-Trefftz spaces (Lederer et al., 2024).

Taken together, the current literature presents quasi-Trefftz methods as a family of PDE-adapted discretizations rather than a single scheme. Their common invariant is the replacement of generic polynomial approximation by small local spaces tuned to the operator through exact or approximate residual annihilation. The choice between polynomial, GPW, embedded, virtual, or space-time constructions depends on the operator class, coefficient regularity, mesh geometry, desired conditioning, and whether exact local solution spaces are available.