Splitting Scalar Auxiliary Variable (SSAV)

• SSAV scheme is a structure-preserving numerical method that reformulates nonlinear energies using a scalar auxiliary variable to simplify integration.
• It splits the evolution into subproblems by handling stiff linear dynamics and nonlinear energy tracking separately, ensuring a modified discrete energy law.
• The approach is versatile, extending to kinetic Langevin dynamics and stochastic PDEs while delivering stability, convergence, and effective long-time behavior.

The Splitting Scalar Auxiliary Variable (SSAV) scheme is a class of structure-preserving numerical methods built around the scalar auxiliary variable methodology and a splitting of the underlying evolution into subproblems that are easier to integrate or analyze. In its most concrete current formulation, SSAV denotes an explicit, energy-structured integrator for kinetic Langevin dynamics with possibly non-globally Lipschitz, superlinear potentials, where a modified scalar auxiliary variable is combined with a Hamiltonian/Ornstein–Uhlenbeck splitting (Dai et al., 4 Sep 2025). More broadly, the surrounding literature shows the same design pattern across gradient flows and stochastic PDEs: one rewrites the nonlinear energy contribution through one or more scalar variables, advances the stiff linear or exactly solvable part separately, and enforces a discrete modified-energy law or a close surrogate of the original energy law (Lam et al., 15 Jan 2025, Liu et al., 2022, Cui et al., 4 Mar 2026).

1. Conceptual setting and relation to the SAV methodology

The starting point is the classical SAV idea: split an energy into a quadratic part and a nonlinear part, then introduce a scalar variable such as

$r(t)\approx \sqrt{E_1(\phi(t))+C},$

so that the difficult nonlinear contribution is carried by a scalar factor rather than by a fully nonlinear implicit solve. In gradient-flow settings this yields linear or linearly implicit schemes with a modified energy of the form

$$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$

and the modified energy can often be shown to be unconditionally dissipative (Lam et al., 15 Jan 2025, Wang et al., 2020).

This suggests that “splitting” in SSAV is best understood structurally rather than narrowly. In some works it means a direct operator splitting between linear and nonlinear parts, or between different physical fields; in others it appears as a decomposition into a PDE step and an auxiliary-scalar correction, or as a stochastic reformulation in which the stiff linear part is propagated by a semigroup while the nonlinear energy is tracked by a scalar variable (Liu et al., 2022, Kemmochi et al., 2021, Cui et al., 4 Mar 2026). In that sense, SSAV is not a single formula but a design principle: isolate the nonquadratic energy in scalar form, preserve a modified energy identity, and combine this with a split treatment of the dominant operators.

A common misconception is that SAV-based methods automatically preserve the original energy law. The literature does not support that conclusion. Standard SAV schemes preserve a modified energy, and discrete equivalence between the auxiliary scalar and its defining energy functional is generally lost after time discretization. This is exactly why relaxed SAV, improved SAV, and related variants were introduced (Jiang et al., 2021, Chen et al., 2024).

2. Core construction principles

Across the available SSAV-type constructions, several ingredients recur. First, the nonlinear energy is rewritten through a scalar surrogate. Second, the evolution is split so that the stiff or linear part is treated implicitly, exactly, or by a semigroup, while the nonlinear contribution is frozen, extrapolated, or otherwise reduced to a low-dimensional correction. Third, the scalar update is designed so that the combined step satisfies a discrete energy identity or inequality.

For gradient flows, this often starts from

$E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$

with $\mathcal L$ handled implicitly and $E_1$ handled through a scalar. For more complicated energies, the nonlinear part may itself be decomposed into components, producing multiple auxiliary variables or higher-rank corrections (Lam et al., 15 Jan 2025, Zhang et al., 17 Jun 2026). For stochastic problems, Itô correction terms may have to be incorporated into the scalar update to retain consistency with the stochastic energy evolution law (Cui et al., 4 Mar 2026).

A second recurring principle is that the scalar update should remain algebraic or low-dimensional. Relaxed SAV schemes choose the new scalar as a convex combination of a predicted scalar and the exact energy surrogate, with the relaxation parameter constrained by a discrete energy inequality. Exponential SAV and exponential semi-implicit SAV schemes update an exponential energy variable through a one-step scalar formula, avoiding bounded-below assumptions on the energy. Constant SAV and improved SAV schemes modify the scalar dynamics to improve consistency with the original energy or to recover original-energy stability (Liu et al., 2022, Liu et al., 2020, Zhang et al., 2024, Chen et al., 2024).

A plausible implication is that SSAV methods can be organized along two almost independent axes: the energy tracker used for the nonlinear part, and the split operator used for the state update. The pullback-corrected SAV framework makes this explicit by separating scalar energy tracking from the rank of the correction applied to the state equation, while retaining a single scalar variable (Zhang et al., 17 Jun 2026).

3. Explicit SSAV for kinetic Langevin dynamics

The most explicit and technically developed SSAV construction in the current literature is the scheme for the kinetic Langevin dynamics

$$\left\{ \begin{aligned} \mathrm{d}v(t) &= -\kappa \nabla \Phi(u(t))\,\mathrm{d}t - \gamma v(t)\,\mathrm{d}t + \Gamma \,\mathrm{d}W_t,\ \mathrm{d}u(t) &= v(t)\,\mathrm{d}t, \end{aligned} \right.$$

with random initial data in $\mathbb R^m\times\mathbb R^m$, where $\nabla\Phi$ may be non-globally Lipschitz and have superlinear growth (Dai et al., 4 Sep 2025).

For $\Gamma=\sqrt{2\kappa\gamma}I$, the invariant measure has density

$$\mu_\infty(x,y)\propto \exp\Big(-\tfrac{|y|^2}{2\kappa}\Big)\exp\big(-\Phi(x)\big).$$

The numerical difficulty is that explicit Euler-type methods may explode with positive probability, while the full kinetic system fails global monotonicity even when $$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$0 is monotone. The SSAV construction addresses this by rewriting the Hamiltonian

$$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$1

through a modified auxiliary variable

$$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$2

where $$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$3 and $$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$4 are chosen so that

$$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$5

The corresponding modified energy is

$$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$6

and satisfies $$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$7.

The splitting is then into two subproblems. The first is a deterministic Hamiltonian/SAV part, chosen so that $$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$8 is preserved at the continuous level. The second is an Ornstein–Uhlenbeck part containing friction and noise, which is exactly solvable. This decomposition is central: the deterministic substep preserves the rewritten Hamiltonian, while the stochastic dissipative substep reproduces the energy growth law of the original dynamics.

The nominally implicit deterministic SAV step can be solved in closed form, which makes the overall method fully explicit. With

$$\widetilde E(\phi,r)=\tfrac12(\phi,\mathcal L\phi)+r^2,$$9

the practical SSAV iteration uses

$E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$0

followed by the exact Ornstein–Uhlenbeck velocity update. No nonlinear solve is required.

The Hamiltonian/SAV substep satisfies the exact discrete conservation law

$E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$1

and the full method reproduces the exact linear-in-time growth bound

$E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$2

4. Stability, convergence, and long-time behavior

The analytical novelty of the kinetic-Langevin SSAV scheme lies in its combination of explicitness with exponential integrability and convergence theory under superlinear, non-globally Lipschitz drift (Dai et al., 4 Sep 2025). The key estimate is an exponential integrability bound for both the exact solution and the numerical approximation: $E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$3 This supplies control of superlinear terms without assuming a global monotonicity condition for the full kinetic system.

The convergence results separate several notions. On the grid, the method first satisfies an order-$E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$4 moment estimate

$E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$5

and in continuous time the associated extension satisfies the analogous order-$E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$6 bound. More strongly, the grid-point error in maximal norm obeys

$E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$7

which is the order-one strong convergence statement emphasized by the paper.

For long-time analysis, exponential-in-$E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$8 bounds are replaced by polynomial growth estimates: $E(\phi)=\tfrac12(\phi,\mathcal L\phi)+E_1(\phi),$9 This is what enables weak error estimates with only polynomial dependence on the final time. For test functions $\mathcal L$0 with polynomially growing derivatives, the invariant-measure approximation satisfies

$\mathcal L$1

The weak order in $\mathcal L$2 is therefore $\mathcal L$3, while the constant grows polynomially rather than exponentially in $\mathcal L$4.

A second misconception addressed by this analysis is that polynomial moment growth would make an explicit scheme ineffective for invariant-measure sampling. The paper shows the opposite in this setting: despite only polynomial growth in the error constants, the method remains computationally effective for approximating the invariant distribution of an exponentially ergodic kinetic Langevin process.

5. Variants and neighboring formulations

The broader SAV literature has produced a family of methods that illuminate what can be changed inside an SSAV design without losing linearity or structure. The main variations concern how the scalar energy tracker is defined, how closely it is forced to follow the original energy, and whether the state equation receives a rank-one or higher-rank correction.

Variant	Defining device	Representative paper
RSAV	Relax $\mathcal L$5 toward $\mathcal L$6 by a convex combination constrained by a discrete energy inequality	(Jiang et al., 2021)
E-SAV / RE-SAV	Exponentiate the total or nonlinear energy and optionally relax $\mathcal L$7 toward $\mathcal L$8	(Liu et al., 2022)
ESI-SAV	Exponential semi-implicit scalar variable with one constant-coefficient linear solve per step	(Liu et al., 2020)
CSAV	Replace the square-root variable by a constant scalar auxiliary variable $\mathcal L$9 at the continuous level, with a stabilization parameter $E_1$0	(Zhang et al., 2024)
iSAV	Keep full linearity while ensuring stability of the original energy rather than only a modified energy	(Chen et al., 2024)
PB-SAV	Keep a single scalar tracker but replace the rank-one SAV correction by a pullback correction of rank up to the number of components	(Zhang et al., 17 Jun 2026)
Stochastic SSAV	Reformulate the stochastic system with Itô-corrected scalar updates so the modified SAV energy asymptotically preserves the stochastic energy evolution law	(Cui et al., 4 Mar 2026)

Relaxed SAV for Cahn–Hilliard systems with bounded mass source is particularly relevant to SSAV because it shows how a scalar relaxation step can be inserted after a linearized PDE solve while preserving a discrete energy inequality in a non-dissipative setting (Lam et al., 15 Jan 2025). Exponential SAV and exponential semi-implicit SAV show that the scalar variable can be based on $E_1$1 rather than $E_1$2, eliminating bounded-from-below restrictions and still requiring only one constant-coefficient linear system per step (Liu et al., 2022, Liu et al., 2020). Improved SAV and constant SAV target the discrepancy between modified and original energies from two different directions: one by original-energy-stable stabilization, the other by making the continuous scalar dynamics trivial and forcing the discrete scalar to stay close to $E_1$3 (Chen et al., 2024, Zhang et al., 2024).

The pullback-corrected SAV framework adds a different perspective. A first-order SAV step can be interpreted as a semi-implicit base step plus a rank-one positive semidefinite correction; PB-SAV generalizes that correction to rank up to the number of energy components while keeping a single scalar variable. The paper shows that this correction is the Gauss–Newton matrix of a least-squares representation of the nonlinear energy, and that refinement over component decompositions increases the correction by an explicit weighted variance term (Zhang et al., 17 Jun 2026). This suggests a route to SSAV schemes with adaptive-rank or componentwise corrections without proliferating scalar bookkeeping.

6. Applications, numerical evidence, and outlook

The direct application area of the explicit kinetic-Langevin SSAV scheme is molecular dynamics and Hamiltonian or underdamped Langevin Monte Carlo. The reported numerical tests use a Gaussian mixture potential, a double-well potential in dimensions $E_1$4, and a two-dimensional bimodal “banana-shaped” target. For $E_1$5, the strong errors in $E_1$6 and the energy errors follow slopes approximately $E_1$7 before the optimal estimate and approximately $E_1$8 after it; weak errors show slope approximately $E_1$9. In long-time tests on $$\left\{ \begin{aligned} \mathrm{d}v(t) &= -\kappa \nabla \Phi(u(t))\,\mathrm{d}t - \gamma v(t)\,\mathrm{d}t + \Gamma \,\mathrm{d}W_t,\ \mathrm{d}u(t) &= v(t)\,\mathrm{d}t, \end{aligned} \right.$$0, the approximation to $$\left\{ \begin{aligned} \mathrm{d}v(t) &= -\kappa \nabla \Phi(u(t))\,\mathrm{d}t - \gamma v(t)\,\mathrm{d}t + \Gamma \,\mathrm{d}W_t,\ \mathrm{d}u(t) &= v(t)\,\mathrm{d}t, \end{aligned} \right.$$1 becomes very good after about $$\left\{ \begin{aligned} \mathrm{d}v(t) &= -\kappa \nabla \Phi(u(t))\,\mathrm{d}t - \gamma v(t)\,\mathrm{d}t + \Gamma \,\mathrm{d}W_t,\ \mathrm{d}u(t) &= v(t)\,\mathrm{d}t, \end{aligned} \right.$$2, and the method remains accurate and stable for larger time steps where Euler–Maruyama distorts the stationary distribution or becomes unstable (Dai et al., 4 Sep 2025).

Outside kinetic Langevin dynamics, SSAV-type constructions appear in several nearby settings. A stochastic SSAV system for the stochastic Cahn–Hilliard equation uses an exponential Euler update for the stiff linear part, an SSAV linearization of the polynomial nonlinearity, and explicit Itô correction terms in the scalar update; it achieves strong order $$\left\{ \begin{aligned} \mathrm{d}v(t) &= -\kappa \nabla \Phi(u(t))\,\mathrm{d}t - \gamma v(t)\,\mathrm{d}t + \Gamma \,\mathrm{d}W_t,\ \mathrm{d}u(t) &= v(t)\,\mathrm{d}t, \end{aligned} \right.$$3 and asymptotic preservation of the averaged energy evolution law (Cui et al., 4 Mar 2026). In Cahn–Hilliard systems with additional fields, the tumor-growth discretization described in the relaxed-SAV literature already has a concrete splitting form: first update the nutrient equation, then solve the SAV phase-field step, while proving a joint discrete energy estimate (Lam et al., 15 Jan 2025). For higher-order gradient flows such as the square phase field crystal equation, a second-order BDF2–SAV scheme with a carefully chosen energy decomposition yields unconditional energy stability and FFT-based constant-coefficient solves; the paper explicitly identifies this as a blueprint for future splitting SAV designs (Wang et al., 2020).

The current literature also delineates the main open directions. The kinetic-Langevin SSAV work remarks that the framework could be extended to other stochastic Hamiltonian systems and to higher-order schemes (Dai et al., 4 Sep 2025). Exponential SAV and ESI-SAV papers point toward Lie or Strang-type splittings between linear stiff operators and nonlinear SAV substeps (Liu et al., 2022, Liu et al., 2020). PB-SAV suggests adaptive decompositions and low-rank corrections guided by weighted-variance identities and Gauss–Newton structure (Zhang et al., 17 Jun 2026). Taken together, these results suggest that future SSAV schemes will likely be judged along three axes simultaneously: preservation of the relevant energy structure, separation of scalar tracking from correction rank or splitting architecture, and long-time accuracy for invariant measures or metastable dynamics.