High-Order Strong Formulation

• High-Order Strong Formulation is a method that directly discretizes differential operators in PDEs using high-order polynomial reconstructions, enabling precise derivative matching.
• It employs finite difference, finite element, and particle-based schemes with techniques like SBP-SAT and WENO to ensure stability, conservation, and accurate shock capturing.
• Advanced temporal integration and adaptive boundary treatments make these methods essential for robust simulations of conservation laws and nonlinear, discontinuous phenomena.

A high-order strong formulation refers to the discretization of partial differential equations (PDEs) using high-order approximations directly to their strong (differential) form rather than their weak (integral or variational) form. Such schemes are characterized by high-order spatial and/or temporal accuracy, rigorous representation of derivatives, and in many cases, careful consideration of stability, conservation, and boundary treatment. High-order strong form methods span finite difference, finite element, particle, and reconstruction-based approaches, and are central to the accurate and robust solution of conservation laws, balance laws, and other physically relevant PDEs, especially in regimes requiring the precise capture of shocks, interfaces, or complex boundary dynamics.

1. Defining High-Order Strong Formulations

High-order strong formulations approximate the spatial and/or temporal derivatives in the governing equations using polynomial or other high-order local reconstructions, typically adhering to the following schematic PDE: $$\frac{\partial u}{\partial t} + \mathcal{L}[u] = S(u),$$ where $\mathcal{L}[u]$ is a spatial differential operator of high order, and $S(u)$ is a source term. Instead of using an integration-by-parts approach as in weak formulations, strong form methods directly discretize the derivatives, potentially leading to more compact stencils, higher accuracy, and opportunities for enforcing specific types of conservation or boundary conditions at a discrete level. These features enable the construction of schemes with spatial accuracy at $\mathcal{O}(h^k)$ for $k\ge3$ in appropriately smooth regions.

Prominent families within this framework include:

• High-order finite difference schemes with compact stencils and strong-form WENO reconstructions, crucial for non-conservative systems such as the Z4 formulation of general relativity (Balsara et al., 8 Jun 2024).
• High-order polynomial reconstruction schemes for finite element or flux reconstruction methods, allowing multi-moment or Hermite-constrained strong-form updates (Xiao et al., 2012).
• High-order Lagrangian particle methods where derivatives are reconstructed using high-order kernel-weighted fits, as in Consistent Particle Hydrodynamics in Strong Form (CPHSF) (Yamamoto et al., 2017).
• High-order finite difference conservative operators, such as the IDO-CF for plasma transport (Ludvig-Osipov et al., 13 Jun 2024).
• GRP (Generalized Riemann Problem) solvers for balance laws with high-order interface time-derivative matching (Qian et al., 2013).

2. Core Mathematical Principles and Discretization Strategies

The construction of high-order strong form schemes emphasizes the following analytic and computational elements:

1. Polynomial Reconstruction and Derivative Matching:

Local (often cell- or particle-centered) polynomials approximate the solution and its derivatives up to order $n$, typically by solving moment-matching or minimization problems. For example, in MMC-FR (Multi-Moment Constrained Flux Reconstruction), the flux polynomial on each cell is uniquely reconstructed from pointwise (Lagrange) values and boundary/derivative (Hermite) constraints, forming a closed linear system (Xiao et al., 2012). In CPHSF, high-order kernel-weighted least-squares fits of intensive variables enable direct estimation of derivatives at particle locations (Yamamoto et al., 2017).

2. Summation-by-Parts (SBP) and Energy Stability:

In finite element and finite difference contexts, SBP operators mimic integration by parts properties at the discrete level. When coupled with Simultaneous-Approximation-Terms (SAT) for weak boundary conditions, provable energy stability results can be obtained for both single and multiple elements, with the mass/conservation matrix and differentiation operator designed to satisfy $Q + Q^T = B$ (Malan et al., 2023).

3. Conservative and Non-Conservative Formulations:

Both conservative (divergence) and non-conservative (non-symmetric) strong forms are addressed. IDO-CF in (Ludvig-Osipov et al., 13 Jun 2024) combines a 4th-order interpolated central-difference-like operator with explicit inclusion of cell-integrated balance equations to secure exact (machine-precision) conservation, while the Z4 and CCZ4 formulations for the Einstein field equations employ non-conservative strong forms with strong hyperbolicity and auxiliary variables (Balsara et al., 8 Jun 2024, Dumbser et al., 2017).

4. High-Order Temporal Integration:

Strong-form spatial discretization is often matched with strong-stability-preserving (SSP) or implicit Runge-Kutta time integration, with careful consideration of consistency and stability for stiff or nonlinear systems (Ludvig-Osipov et al., 13 Jun 2024, Balsara et al., 8 Jun 2024).

3. Advanced Boundary Condition and Shock Treatment

Boundary and discontinuity features are critical in high-order strong formulations:

• Boundary Condition Imposition:

Direct enforcement of physical boundary conditions is accomplished via special polynomial constraints (e.g., modified fit for ghost or boundary particles in CPHSF), SAT penalty terms in SBP-SAT CG methods, or flux/derivative matching at element boundaries in MMC-FR (Yamamoto et al., 2017, Malan et al., 2023, Xiao et al., 2012).

• Artificial Viscosity and Dissipation:

To suppress oscillations and handle shocks, high-order strong form methods employ targeted dissipation mechanisms. In SBP-SAT CG, Galerkin-weighted artificial dissipation operators—built from higher derivatives—are activated based on smoothness indicators to localize diffusion near discontinuities (Malan et al., 2023). CPHSF applies a directionally adaptive tensor bulk viscosity that acts along the most compressive eigenvector of the local strain-rate tensor (Yamamoto et al., 2017). FD-WENO schemes for Z4 systems use fluctuation-based dissipation compatible with non-conservative principal terms (Balsara et al., 8 Jun 2024).

• Limiter and Subcell Approaches:

For discontinuous Galerkin (DG) and related schemes, subcell WENO limiters can be incorporated to stabilize the evolution near singularities or mergers (e.g., moving black holes in CCZ4) without degrading high-order accuracy elsewhere (Dumbser et al., 2017).

4. Conservation, Consistency, and Convergence Properties

Rigorous analysis is conducted for high-order strong form discretizations:

• Conservation:

Schemes such as IDO-CF enforce conservation in the spatial variable up to machine precision by directly incorporating cell-integral constraints, whereas particle-based CPHSF achieves small conservation error due to its Lagrangian character, despite non-conservation at the particle level (Ludvig-Osipov et al., 13 Jun 2024, Yamamoto et al., 2017).

• Consistency and Accuracy:

Formally, spatial truncation errors can be established as $\mathcal O(h^{k+1})$ locally for degree-$k$ polynomial reconstructions, with global errors in cell averages, function values, or derivatives scaling as $\mathcal O(h^k)$ or higher, depending on super-convergence effects (e.g., $O(h^{p+2})$ in SBP-SAT CG for $p$-th order Lagrange polynomials) (Malan et al., 2023). Numerical experiments in both linear and nonlinear contexts confirm machine-precision conservation, fourth-order spatial convergence, and maintained high-order rates even in the presence of strong discontinuities (Ludvig-Osipov et al., 13 Jun 2024, Balsara et al., 8 Jun 2024, Qian et al., 2013).

• Stability:

Energy stability is obtained by construction in SBP-SAT CG, with global stability carrying through to multi-element assemblies. Implicit time-stepping (L-stable IRK) robustly handles stiff nonlinearities and avoids restrictive stability limits (Malan et al., 2023, Ludvig-Osipov et al., 13 Jun 2024).

5. Representative Algorithms and Formulations

Scheme	Core Strategy	Key Properties
CPHSF (Yamamoto et al., 2017)	Kernel-weighted high-order fit on particles	Lagrangian, high-order, boundary adaptable
MMC-FR (Xiao et al., 2012)	Multi-moment polynomial flux reconstruction	Arbitrary order, Hermite/Lagrange blend
SBP-SAT CG (Malan et al., 2023)	Strong form CG via SBP, SAT, artificial dissipation	Energy stable, WENO-like at shock, super-convergent
IDO-CF (Ludvig-Osipov et al., 13 Jun 2024)	4th-order interpolated difference, cell-integral coupling	Exact conservation, robust for stiff PDEs
FD–WENO (Balsara et al., 8 Jun 2024)	Non-conservative strong-form with WENO, well-balanced	Arbitrary high-order, robust for Einstein equations
GRP (Qian et al., 2013)	High-order Taylor expansion at interfaces via strong-form characteristic analysis	Second/third-order interface solvers, robust through discontinuities

Each approach addresses the strong form through context-specific discretization and constraint design, tailored for the physical PDE structure.

6. Applications and Numerical Evidence

High-order strong formulation methods have demonstrated substantial accuracy and robustness across application domains:

• Astrophysical and Cosmological Simulations:

CPHSF has been successfully employed for extremely long-term disk and star/disk formation simulations with large deformation and negligible angular momentum loss, achieving order-of-magnitude improvements in lifetime over SPH (Yamamoto et al., 2017). High-order SBP-SAT CG and ADER-DG methods have been critical for robust, high-accuracy gravitational wave and moving black hole simulations in fully general-relativistic settings (Dumbser et al., 2017).

• Plasma and Transport Modeling:

IDO-CF achieves machine-level spatial conservation and maintains fourth-order spatial error up to round-off for stiff nonlinear plasma transport equations, enabling accurate transient simulation on coarse grids (Ludvig-Osipov et al., 13 Jun 2024).

• General Balance Laws and Hyperbolic Systems:

High-order GRP solvers ensure temporal and spatial accuracy for nonlinear conservation laws, with the third-order quadratic solver providing essential accuracy in strong discontinuity regimes (Qian et al., 2013).

• Engineering and Free-Surface Flow:

CPHSF and high-order WENO strong-form schemes outperform classical SPH and low-order methods for dam break, sloshing, and multiphase interface capturing, particularly when accurate enforcement of physical boundaries is needed (Yamamoto et al., 2017).

7. Limitations and Open Challenges

Despite their advantages, high-order strong formulation methods entail several challenges:

• Computational Complexity:

High-order polynomial reconstruction and the inversion of local moment matrices incur higher per-step computational cost relative to low-order or weak form alternatives. Nonlinearities in flux or diffusion operators add algebraic complexity to the coefficient computation (as in IDO-CF (Ludvig-Osipov et al., 13 Jun 2024)).

• Conservation and Re-gridding in Lagrangian Methods:

Particle-based strong-formulation methods, such as CPHSF, are not strictly conservative and require periodic particle rearrangement (re-gridding) based on local ordering criteria (Yamamoto et al., 2017).

• Shock and Topology Management:

Approaches for capturing shocks and handling topological changes (e.g., merging droplets) require problem-specific limiters, surface-bulk switching, or hybrid DG/WENO strategies, and implicit time stepping is often necessary near boundaries or interfaces (Yamamoto et al., 2017, Dumbser et al., 2017).

• Extension to Multidimensional and Nonlinear Systems:

While the above techniques generalize naturally to multidimensional and strongly nonlinear systems, detailed stability and accuracy analysis remains nontrivial, particularly in the presence of complex geometries or non-conservative terms (Balsara et al., 8 Jun 2024).

8. Outlook and Research Directions

Current research focuses on:

• Integration of high-order strong-form discretization with adaptive mesh refinement (AMR), local time stepping (LTS), and subcell limiting to maximize computational efficiency for multi-scale problems (Dumbser et al., 2017).
• Unification of strong-form and weak-form discretizations, particularly in multi-moment and hybrid frameworks that allow problem-specific tuning between conservation, stability, and accuracy as in MMC-FR (Xiao et al., 2012).
• Extension to mixed non-conservative and conservative systems, with rigorous well-balancing and algebraic constraint management, as required in general-relativistic simulations or stiff plasma models (Balsara et al., 8 Jun 2024, Ludvig-Osipov et al., 13 Jun 2024).
• Advanced boundary treatment and interface capturing through high-order polynomial constraints and adaptive dissipation for complex, evolving domains (Yamamoto et al., 2017, Malan et al., 2023).

The continuing development of high-order strong formulation schemes plays a pivotal role in advancing the fidelity and robustness of computational models for complex physical phenomena.