🔢 mathematics

A Discontinuous Galerkin Consistent Splitting Method for the Incompressible Navier-Stokes Equations

This paper introduces a discontinuous Galerkin consistent splitting method for the incompressible Navier-Stokes equations that implicitly enforces divergence-free constraints and consistent boundary conditions to eliminate pressure boundary layers and splitting errors, while achieving optimal convergence and efficient decoupled solves through a semi-implicit treatment of convection and penalty-enhanced mass conservation.

Original authors: Dominik Still, Natalia Nebulishvili, Richard Schussnig, Katharina Kormann, Martin Kronbichler

Published 2026-04-29

📖 4 min read🧠 Deep dive

CC BY 4.0

Original authors: Dominik Still, Natalia Nebulishvili, Richard Schussnig, Katharina Kormann, Martin Kronbichler

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are trying to predict how a crowd of people (the fluid) moves through a complex building (the domain) without anyone bumping into each other or piling up in a corner. This is the challenge of solving the Navier-Stokes equations, the mathematical rules that govern how liquids and gases flow.

For decades, scientists have struggled with a specific problem: how to calculate the "pressure" (the force pushing the crowd) and the "velocity" (how fast they are moving) at the same time without the computer getting stuck in a giant, messy loop.

This paper introduces a new, smarter way to do this calculation using a method called Discontinuous Galerkin (DG) combined with a Consistent Splitting strategy. Here is how it works, broken down into simple concepts:

1. The Old Problem: The "Tangled Knot"

Traditionally, trying to solve for speed and pressure together is like trying to untangle a knot while holding both ends. The computer has to solve one massive, complicated system where every variable depends on every other variable. This is slow, requires huge amounts of computing power, and often forces scientists to use "approximations" at the boundaries (the walls of the room) that create fake ripples or errors in the pressure, known as pressure boundary layers.

2. The New Solution: The "Two-Step Dance"

The authors propose a "splitting" method. Instead of trying to untangle the whole knot at once, they separate the dance into two distinct steps that happen in a specific order every time the clock ticks forward:

Step A (The Speed Dance): First, the computer guesses the new speed of the fluid based on where it was a moment ago and the forces pushing it. Crucially, it uses a "linearly implicit" approach. Think of this as a dancer who plans their next move by looking slightly ahead but doesn't need to solve a complex riddle to do it. It's faster than a fully explicit guess (which is unstable) but doesn't require the heavy lifting of a fully implicit calculation.
Step B (The Pressure Fix): Once the speed is guessed, the computer solves a simpler equation (a Poisson equation) to figure out the pressure needed to make sure the fluid doesn't magically appear or disappear (keeping the "divergence-free" rule).

3. The "Consistent" Secret Sauce

The key innovation here is Consistency. In older splitting methods, the rules for how the fluid hits the walls (boundary conditions) were often "inconsistent" with the math, leading to those fake pressure ripples mentioned earlier.

This paper uses Consistent Boundary Conditions. Imagine a wall that perfectly understands the physics of the fluid hitting it. By carefully designing the math so that the pressure equation "knows" exactly what the fluid is doing at the wall, the method eliminates those fake ripples. This means the simulation is much more accurate, especially over long periods of time.

4. The "Leray Projection" (The Quality Control Check)

Even with a great two-step dance, small errors can accumulate, causing the fluid to slightly "leak" or compress when it shouldn't. To fix this, the authors add a Leray projection.

Analogy: Imagine a bouncer at a club. After the dancers (the fluid) have moved, the bouncer checks the crowd. If anyone is standing too close together (violating the incompressibility rule), the bouncer gently nudges them apart to restore order.
In the paper, this is done mathematically to ensure the fluid remains perfectly incompressible, improving the conservation of mass.

5. Why This Matters (The Results)

The authors tested this method on two famous scenarios:

Flow around a cylinder: Like wind blowing around a flagpole. They showed their method could predict the "drag" (how hard the wind pushes) and "lift" (how it tries to lift the pole) with extreme precision, matching the best existing data.
Taylor-Green Vortex: A complex 3D simulation of swirling turbulence, like a storm inside a box. They proved their method could handle high-speed, chaotic flows without breaking down, even when using larger time steps (making the simulation faster).

The Bottom Line

This paper presents a new way to simulate fluid flow that is:

Faster: It breaks a giant problem into smaller, easier pieces.
More Accurate: It removes the "fake ripples" at the walls that plague older methods.
Robust: It handles high-speed, turbulent flows better than previous splitting techniques.

By using this "two-step dance" with perfect boundary rules, the authors have created a tool that allows scientists to simulate complex fluid dynamics (like air over a plane or blood in a vein) with higher precision and less computational cost, without needing to solve the impossible "tangled knot" of the old methods.

1. Problem Statement

The paper addresses the numerical simulation of transient incompressible Navier–Stokes flows. The primary challenges in this domain are:

Coupling Complexity: Standard monolithic approaches (solving velocity and pressure simultaneously) result in large saddle-point linear systems that are difficult to precondition and solve efficiently, especially for high-order methods.
Splitting Errors: Traditional projection or fractional-step methods decouple velocity and pressure to simplify the system. However, these often introduce "splitting errors," leading to spurious pressure boundary layers and limiting temporal accuracy to first or second order, regardless of the spatial discretization order.
Boundary Conditions: Enforcing consistent boundary conditions for the pressure Poisson equation (PPE) in discontinuous Galerkin (DG) settings, particularly for open and traction boundaries, is non-trivial.
Stability vs. Efficiency: Fully implicit schemes are robust but require solving nonlinear systems. Explicit schemes are cheap but suffer from strict CFL restrictions.

2. Methodology

The authors propose a Discontinuous Galerkin (DG) spatial discretization of the Consistent Splitting Scheme originally proposed by Liu (2009), adapted for high-order accuracy and robustness.

A. Temporal Discretization

The method uses Backward Differentiation Formulae (BDF) for time integration. The scheme decouples the momentum and pressure equations into two sequential steps per time step:

Momentum Step: Solves a linearly implicit convection-diffusion-reaction equation for velocity. The convective term is treated in a linearly implicit manner (using extrapolated velocities), which relaxes the CFL restriction of explicit schemes without requiring the solution of nonlinear systems.
Pressure Step: Solves a Pressure Poisson Equation (PPE) to update the pressure.

B. Consistent Splitting & Boundary Conditions

Unlike classical projection methods, this scheme enforces the divergence-free constraint implicitly via a PPE derived by taking the divergence of the momentum equation.

Consistent BCs: The pressure boundary conditions are derived rigorously from the momentum equation and the traction conditions.
DG Fluxes: By selecting specific fluxes for the divergence of the convective operator and the velocity divergence, the authors show that the consistent pressure boundary condition reduces to contributions solely from the acceleration term and the viscous term. This simplifies the implementation significantly.

C. Mass Conservation & Stabilization

To address the fact that consistent splitting does not strictly enforce point-wise divergence-free velocity:

Leray Projection: A Helmholtz-Leray decomposition is employed. To avoid the computational cost of solving an extra Poisson equation for the projection variable at every step, the authors introduce a Modified Pressure Variable ( $\hat{p}$ ). This variable combines the extrapolated pressure and the Leray correction term, allowing the scheme to maintain high-order accuracy and mass conservation while solving only one PPE per time step.
Penalty Terms: The scheme incorporates divergence penalty and normal continuity penalty terms (similar to $H(\text{div})$ -conforming methods) to improve mass conservation and stability.

D. Spatial Discretization

Method: Symmetric Interior Penalty Galerkin (SIPG) is used for second-order derivatives.
Polynomial Degrees: The method uses inf-sup stable pairings ( $k_u = k_p + 1$ ) to ensure stability, though the consistent splitting nature aims to reduce the strictness of this requirement compared to monolithic methods.
Vorticity Handling: For the viscous term in the boundary conditions (curl-curl form), the authors compute vorticity via $L^2$ projection to avoid reformulating the term into first derivatives, which would degrade accuracy in DG settings.

3. Key Contributions

DG Formulation of Consistent Splitting: This is the first work to present a high-order DG discretization of Liu's consistent splitting scheme, combining the robustness of DG in convection-dominated regimes with the accuracy of consistent splitting.
Modified Pressure Variable: The introduction of the modified pressure variable ( $\hat{p}$ ) eliminates the need for an explicit, separate Leray projection step, reducing computational overhead while preserving optimal convergence rates and mass conservation.
Linearly Implicit Convective Treatment: The scheme achieves a balance between stability and efficiency by treating convection linearly implicitly, allowing for larger time steps than explicit schemes without the cost of nonlinear solvers.
Simplified Boundary Conditions: The derivation shows that for $L^2$ -conforming DG, the complex boundary terms in the PPE cancel out, leaving only acceleration and viscous contributions, which simplifies the implementation of open and traction boundaries.
Higher-Order Stability: The method demonstrates compatibility with high-order BDF schemes (up to BDF-4), overcoming the typical splitting error limitations of lower-order projection methods.

4. Numerical Results

The method was validated through extensive numerical experiments:

Manufactured Solution (2D Taylor-Green Vortex):
- Verified optimal convergence rates in both space ( $O(h^{k+1})$ ) and time ( $O(\Delta t^J)$ ) for BDF-1 to BDF-4.
- Confirmed that the Leray projection (via the modified pressure) is vital; without it, convergence rates degrade.
Flow Around a Cylinder (2D Benchmark):
- Results for drag, lift, and pressure difference matched state-of-the-art reference values (including high-fidelity DG results) with high precision.
- Demonstrated temporal stability for CFL numbers up to 16 with the linearly implicit formulation, significantly outperforming explicit schemes.
- Showed that BDF-3 and BDF-4 provide higher accuracy than BDF-2 for the same time step size.
3D Taylor-Green Vortex (Turbulence Benchmark):
- Simulated the transition to turbulence at $Re=1600$.
- Compared against the Dual Splitting Scheme (Karniadakis et al.) and fully implicit coupled solvers.
- The consistent splitting scheme produced results comparable to the dual splitting scheme and the coupled approach but with the advantage of avoiding splitting errors.
- Demonstrated that using $J_c = J-1$ (lower extrapolation order for the convective term in the PPE) improves stability for high-order BDF schemes without sacrificing accuracy.

5. Significance

This work represents a significant step forward in high-fidelity computational fluid dynamics (CFD):

Efficiency: By decoupling velocity and pressure, the method avoids the expensive saddle-point linear systems, allowing the use of standard, efficient preconditioners (like multigrid) for the pressure Poisson equation.
Accuracy: It eliminates the "splitting error" that plagues traditional projection methods, enabling true high-order temporal accuracy (BDF-3/4) which is crucial for Direct Numerical Simulation (DNS) and Large Eddy Simulation (LES).
Robustness: The linearly implicit treatment of convection makes the method suitable for high-Reynolds number flows where explicit schemes would require prohibitively small time steps.
Scalability: The formulation is well-suited for modern high-performance computing (HPC) architectures, as evidenced by its implementation in the ExaDG framework (based on deal.II), which is designed for exascale computing.

In conclusion, the paper presents a robust, high-order, and efficient numerical framework for incompressible flows that successfully bridges the gap between the accuracy of monolithic solvers and the computational efficiency of splitting methods.