Ensemble algorithm for parametrized flow problems with energy stable open boundary conditions

Technical Summary: Ensemble Algorithm for Parametrized Flow Problems with Energy Stable Open Boundary Conditions

Problem Statement
The paper addresses the computational challenge of solving ensembles of incompressible Navier-Stokes equations subject to varying initial conditions, body forces, and viscosity coefficients ($\nu_j$). In applications such as sensitivity analysis, reduced-order modeling, and ensemble forecasting, the number of ensemble members ($J$) can be large. Standard semi-implicit methods require assembling and solving a distinct linear system for each member $j$ at every time step due to the dependence of the system matrix on $\nu_j$ and the nonlinear terms. This results in prohibitive computational costs. While fully explicit approaches allow for a single system assembly, they impose restrictive time-step conditions, particularly for high Reynolds number flows. Previous ensemble methods for constant viscosity or specific viscosity variations often required restrictive assumptions on the viscosity coefficients (e.g., bounded relative differences) to ensure stability.

Methodology
The authors propose a novel ensemble calculation method that decouples the linear system matrix from the ensemble member index $j$, allowing for the solution of a single linear system with multiple right-hand sides at each time step. The core of the methodology involves:

1. Splitting Strategy: The nonlinear convective term is split into an ensemble mean component and a fluctuating component. The mean component is treated semi-implicitly, while the fluctuating component is treated explicitly.
2. Viscosity Treatment: Unlike previous works (e.g., [8]) that required specific bounds on the variation of viscosity coefficients to ensure stability, this paper introduces a different treatment of the diffusive term. By defining a maximum viscosity $\nu_\infty = \max_j \nu_j$ and splitting the diffusion term, the authors derive a scheme where the linear system matrix is independent of $j$ without imposing restrictions on the magnitude of the viscosity variations.
3. Open Boundary Conditions (OBC): The method is extended to domains with open boundaries ($\Gamma_N$). The boundary is partitioned into outflow and backflow regions. The authors utilize an energy-stable boundary condition involving a characteristic length scale $L$ and a Heaviside-like function to handle the convective nature of the outflow.
4. Discretization: The paper presents both first-order (Algorithms 1–3) and second-order (Algorithms 4–6) time-stepping schemes. The spatial discretization uses conforming finite element spaces (specifically Taylor-Hood $P_2/P_1$ pairs) with skew-symmetrized trilinear forms to ensure energy stability.
5. Stability Analysis: Theoretical proofs establish energy stability for the first-order schemes under specific CFL conditions involving velocity fluctuations. For the second-order schemes, while a rigorous stability proof without viscosity restrictions is not provided (due to the complexity of the $\nabla \partial_{tt} u$ term), numerical experiments suggest stability under the same conditions as the first-order schemes.

Key Contributions

• Unrestricted Viscosity Stability: The primary theoretical contribution is an ensemble scheme for variable viscosity coefficients that guarantees stability without the restrictive assumption $|\nu_j - \nu|/\nu \leq \sqrt{\mu}$ found in prior literature.
• Open Boundary Extension: The authors extend ensemble methods to problems with open boundaries, deriving weak formulations and proving energy stability for first-order schemes under specific time-step constraints.
• Algorithmic Efficiency: The proposed algorithms result in a single coefficient matrix per time step, solvable efficiently using multi-right-hand-side solvers, significantly reducing the computational cost compared to solving $J$ independent systems.

Results
Numerical experiments were conducted using the FreeFem++ package:

1. Convergence Study: A convergence test on a domain with exact solutions confirmed the expected second-order convergence rates in both velocity and pressure for the proposed second-order schemes. The accuracy was comparable to independent sequential simulations.
2. Flow Around a Cylinder: The algorithms were tested on a benchmark problem involving flow past a cylinder with three different viscosity values. The results (drag, lift, and pressure drop) matched reference values from literature using Dirichlet and "do-nothing" outflow conditions. The open boundary condition successfully allowed vortices to exit the domain naturally, avoiding the artificial reflection issues associated with prescribed parabolic outflow profiles.
3. Channel Flow with Contraction: A complex flow scenario with two outlets and varying viscosities was simulated. The ensemble results qualitatively matched independent Direct Numerical Simulation (DNS) runs, despite the ensemble simulations being performed on a coarser mesh stabilized by an adaptive nonlinear filter.

Significance and Claims
The paper claims to offer a robust and efficient framework for parametrized flow simulations. By removing the restrictive assumptions on viscosity variations required by previous ensemble methods, the proposed algorithm broadens the applicability of ensemble simulations to problems with highly variable physical parameters. The extension to open boundary conditions with provable energy stability is a significant step toward realistic simulations of unbounded flows. The authors conclude that while the ensemble schemes perform well qualitatively and quantitatively against independent simulations, further research is needed to fully quantify computational advantages, particularly regarding CPU time comparisons on adaptively refined meshes for high Reynolds number flows.