Analysis gaps in multi-dimensional high-order weighted compact central schemes: mixed-derivative averaging, limiter parameter selection, order-dependent CFL restriction, and pointwise convergence under rotated characteristic decomposition

The weighted compact central schemes (WCCS) of Shen, Al Jahdali and Parsani (J. Comput. Phys. 466:111370, 2022) evolve, on each cell, a polynomial whose degrees of freedom are the cell average and all spatial derivatives, updated on space-time staggered meshes; they bypass Riemann solvers and achieve arbitrarily high, uniform space-time order on a super-compact stencil in a single explicit time step. In multiple dimensions and at high order, however, the methodology rests on empirical heuristics whose consequences the authors explicitly left unanalyzed: the arithmetic averaging of two discrepant discrete mixed derivatives ( $u_{xy}$ versus $u_{yx}$ ), a WENO-type limiter whose hand-tuned parameter $\alpha$ leaves unestablished whether its non-oscillatory property holds for arbitrarily strong shocks, a maximum stable CFL number that shrinks as the order grows, and an unexplained loss of pointwise ( $L_\infty$ ) convergence for the second-order limited scheme in two dimensions, whose reason the authors state is "not completely clear because of the complexity of the rotated characteristic decomposition". A 2025-2026 follow-up line by the same group removed the Cauchy-Kovalewski tensor-algebra bottleneck (Jacobian-free and recursively Jacobian-free variants) and replaced the limiter for orders above four, but left the four analysis questions open. The problem posed here is to replace the surviving heuristics with quantified theory: characterize (or remove by design) the effect of mixed-derivative averaging on accuracy, conservation, and stability, derive a principled limiter-parameter selection strategy, quantify the order-dependent timestep restriction, and explain the observed pointwise convergence behavior under the rotated characteristic projection.

Background

Systems of hyperbolic conservation laws, $\partial_t U + \nabla\cdot F(U) = 0$ with $U(x,t)$ vector-valued and the flux Jacobians of $F$ having real eigenvalues (the compressible Euler equations of gas dynamics are the canonical example), develop shock waves and contact discontinuities from smooth initial data. Numerical methods must therefore be simultaneously high-order accurate in smooth regions and free of spurious oscillations near discontinuities. Classical Godunov-type finite-volume schemes compute the flux across each cell interface from the solution of a Riemann problem — the initial-value problem for two constant states — typically through an approximate Riemann solver. Central schemes avoid this by staggering the control volumes in space-time so that the flux is continuous across interfaces: the Lax-Friedrichs scheme is the first-order ancestor, and space-time conservation element and solution element (CESE) schemes are a modern staggered central family.

The weighted compact central schemes (WCCS) of Shen, Al Jahdali and Parsani (2022) combine this staggered central flux with a compact polynomial representation: on each cell the solution is approximated by a polynomial of degree $P$ whose degrees of freedom — the cell average and every spatial derivative up to order $P$ — are stored and updated separately. The cell average is updated by a finite-volume flux on the staggered mesh, the $k$ -th spatial derivatives are updated by central differences of the $(k-1)$ -th derivatives at the cell vertices, and all time derivatives are obtained from spatial ones by the Cauchy-Kovalewski (Lax-Wendroff) procedure, which repeatedly substitutes the PDE to convert time derivatives into spatial ones. The result is arbitrarily high, uniform order in space and time on a stencil of only the neighboring cells, in a single explicit time step and without Riemann solvers. For nonlinear systems the closure is automated with a generalized Leibniz rule technique (Dyson, 2002), also used in ADER schemes (one-step schemes built on the same Lax-Wendroff idea). To capture discontinuities, a WENO-type limiter — a nonlinear convex combination of the high-order polynomial with lower-degree polynomials, weighted by smoothness indicators — is applied component-wise to characteristic variables (the variables that diagonalize the flux Jacobian, in which each component obeys a simple advection equation); in two dimensions a rotated characteristic decomposition performs one projection along the local gradient direction instead of one per mesh direction, saving roughly 40% of the characteristic-decomposition cost. Linear von Neumann stability analysis gives maximum CFL numbers (the bound $\nu = a\Delta t/\Delta x$ on the timestep) that shrink as the order grows: $\nu \le 0.5$ , $0.384$ , and $0.304$ for the 2nd-, 3rd-, and 4th-order members respectively, obtained numerically because of the complexity of the amplification matrix; the paper's numerical examples use the more conservative CFL numbers 0.4, 0.3, and 0.25.

The authors explicitly flag, or leave open, what the construction does not establish. In two dimensions the vertex-based update produces two different discrete approximations of the same mixed derivative, $u_{xy}$ and $u_{yx}$ , and the method simply averages them (Remark 5 of the paper), with no analysis of the effect on accuracy or conservation; whether the limiter preserves its non-oscillatory property for arbitrarily strong shocks is left open, and its parameter $\alpha$ (which controls the strength of the fallback to first order and must satisfy $\alpha\ge 1/2$ to preserve order at smooth extrema) is set by hand to $\alpha=2$ ; the timestep restriction tightens with order, with the maximum stable CFL numbers known only from numerical von Neumann computations for the low-order members; and in the two-dimensional isentropic-vortex test the pointwise ( $L_\infty$ ) convergence of the second-order limited scheme fell below its design order without the characteristic decomposition technique — an observation whose reason the authors state is "not completely clear because of the complexity of the rotated characteristic decomposition". A 2025-2026 follow-up line by the same group (Jacobian-free and recursively Jacobian-free Cauchy-Kovalewski variants) has since removed the tensor-algebra bottleneck of the closure for general systems, so the questions that survive concern the multi-dimensional accuracy, stability, and limiter theory of the family.

Problem Statement

The problem is to replace the surviving empirical heuristics of the multi-dimensional, arbitrarily high-order weighted compact central scheme (WCCS) family — the arithmetic averaging of discrepant mixed derivatives, the hand-tuned WENO-type limiter parameter, the empirically characterized order-dependent timestep restriction, and the unexplained pointwise convergence behavior under the rotated characteristic decomposition — with quantified theory or analysis-backed constructions, so that the accuracy, conservation, stability, and shock-capturing behavior of the family in multiple dimensions are predicted rather than tuned. The objective comprises the following four components.

(1) Mixed-derivative consistency. Determine the effect of the arithmetic averaging of $u_{xy}$ and $u_{yx}$ in the two-dimensional vertex-based derivative update on the formal order of accuracy, on the conservative character of the staggered cell-average update, and on linear stability, at arbitrary order — or construct a consistent mixed-derivative update that removes the $u_{xy}$ / $u_{yx}$ ambiguity, and quantify the accuracy, conservation, and stability properties of the corrected scheme.

(2) Limiter parameters and strong shocks. Provide a principled, automated strategy for selecting the WENO-type limiter's parameters (in particular $\alpha$ , hand-set to 2 in the source) — or a redesigned limiter for the same scheme family — together with an analysis of the conditions under which spurious oscillations near very strong discontinuities can and cannot be eliminated, while preserving the design order in smooth regions including at smooth extrema.

(3) Order-dependent timestep restriction. Explain and quantify how the maximum stable CFL number of the family depends on the order $P$ as $P$ grows, extending beyond the numerically computed von Neumann limits reported for the 2nd-, 3rd-, and 4th-order members.

(4) Pointwise convergence and the rotated characteristic decomposition. Explain the observation that in the two-dimensional isentropic-vortex test the $L_\infty$ convergence of the second-order limited scheme WCCS-2 falls below second order without the characteristic decomposition technique — the interaction between the limiter, the rotated characteristic projection, and pointwise error at smooth extrema — and give conditions under which the limited multi-dimensional scheme attains, or loses, its design pointwise order.

Each component is to be addressed for the WCCS family itself, or for a correction of it that preserves the family's defining properties (staggered Riemann-solver-free flux, super-compact stencil, single explicit time step, uniform space-time order); substituting a different scheme family, or confining the analysis to a fixed low order or to linear scalar problems only, does not answer the question — in component (4) the fixed second-order instance is the observation to be explained, and the objective there is the general conditions on the limited multi-dimensional scheme, not an analysis of that single instance alone.

Each component is posed at the following generality. (1) The effect of the averaging on formal order, conservation, and linear stability at general order $P$ , or a consistent mixed-derivative update with those properties established. (2) A parameter-selection rule that is derived rather than fitted, or a redesigned limiter, together with conditions characterizing when oscillations near strong discontinuities are and are not eliminated while the design order is preserved at smooth extrema. (3) An order law for the maximum stable CFL number valid for all $P$ . (4) The mechanism of the observed pointwise order loss together with the general conditions on the limited multi-dimensional scheme. The objective in each case is the general characterization rather than a numerical demonstration at the orders examined in the source.

References

Hua Shen, Rasha Al Jahdali, Matteo Parsani (2022). A class of high-order weighted compact central schemes for solving hyperbolic conservation laws. Journal of Computational Physics 466:111370. DOI: 10.1016/j.jcp.2022.111370 (arXiv:2104.04347). https://doi.org/10.1016/j.jcp.2022.111370
Run Li, Hua Shen (2025). Jacobian-free high-order weighted compact central schemes for hyperbolic conservation laws. Journal of Scientific Computing 105:38. DOI: 10.1007/s10915-025-03064-z. https://doi.org/10.1007/s10915-025-03064-z
Run Li, Yuxuan Zhou, Hua Shen (2026). Increasingly high-order compact central schemes with a recursively Jacobian-free Cauchy-Kovalevskaya method. Journal of Computational Physics 567:115300. DOI: 10.1016/j.jcp.2026.115300. https://doi.org/10.1016/j.jcp.2026.115300
Yuxuan Zhou, Hua Shen (2026). Arbitrarily High-Order Compact Central Schemes for Solving a Class of Conservative Partial Differential Equations. Journal of Scientific Computing 107:63. DOI: 10.1007/s10915-026-03260-5. https://doi.org/10.1007/s10915-026-03260-5
Hua Shen, Matteo Parsani (2021). A Rotated Characteristic Decomposition Technique for High-Order Reconstructions in Multi-dimensions. Journal of Scientific Computing 88:87. DOI: 10.1007/s10915-021-01602-z. https://doi.org/10.1007/s10915-021-01602-z
D. Zorío, A. Baeza, P. Mulet (2017). An Approximate Lax-Wendroff-Type Procedure for High Order Accurate Schemes for Hyperbolic Conservation Laws. Journal of Scientific Computing 71:246-273. DOI: 10.1007/s10915-016-0298-2. https://doi.org/10.1007/s10915-016-0298-2
R. W. Dyson (2002). Technique for very High Order Nonlinear Simulation and Validation. Journal of Computational Acoustics 10(2):211-229. DOI: 10.1142/S0218396X02001553. https://doi.org/10.1142/S0218396X02001553
Yang Cao, Yong Luo, Hu Li, Ruixuan Ma, Hua Shen, Conghai Wu (2026). An Upwind-Type Rotated Characteristic Decomposition Method for Steady-State Simulations of the Euler Equations. Journal of Computational Physics 562:115025. DOI: 10.1016/j.jcp.2026.115025. https://doi.org/10.1016/j.jcp.2026.115025