Source-linked AI summary
An improved fast iterative moment method for the simulation of near-continuum gas flows with Shakhov model
Zhenning Cai, Zhicheng Hu, Guanghan Li, Chunwu Wang, Xiaowen Wang
TL;DR
Near-continuum gas-flow simulation with kinetic models remains computationally challenging, and the original FIM method had accuracy and efficiency limitations. This work extends FIM to the Shakhov model with revised spatial and iterative solvers. The experiments report improved efficiency and robustness, particularly at small Kn.
Problem
Kinetic simulation of near-continuum gas flows is computationally challenging, and the original FIM method was limited in model coverage, spatial accuracy, and hydrodynamic-stage convergence.
Source: §1 · pp. 1–2Method
The solver extends FIM to the Shakhov model, using linear-reconstruction finite volumes, implicit treatment of collisions, and a Newton-multigrid hydrodynamic solver.
Source: §6 · p. 23Results
Numerical experiments show improved efficiency and robustness, particularly at small Kn, with one Newton-multigrid step generally sufficient per alternating iteration.
Source: §6 · p. 23Takeaways & Limitations
The method provides a practical framework for steady-state computation of high-order moment systems in near-continuum gas flows.
Source: §6 · p. 23The Fourier analysis assumes periodic space and neglects boundary effects.
Source: §3.2 · p. 9
Abstract
from arXiv · showhide
The fast iterative moment (FIM) method developed in [G. Li, C. Wang, and Z. Hu, Comput. Fluids, 311 (2026), 107031] achieves fast convergence for near-continuum gas flows modeled by the Boltzmann-BGK equation through alternating iterations between the high-order moment system and the consistent hydrodynamic equations. In this paper, we extend the FIM solver to the Shakhov model and improve it in the following aspects. A finite volume method with linear reconstruction is employed for high-resolution spatial discretization of both systems. Following the original FIM method, an implicit-explicit scheme treating the Shakhov collision term implicitly is adopted for the high-order moment system. A local Fourier analysis of a linearized one-dimensional BGK model reveals the need for more hydrodynamic iterations at smaller Knudsen numbers and their slower convergence on refined grids. Motivated by these findings, a Newton-multigrid method with a new regularization strategy is introduced for the hydrodynamic equations, generally providing an effective hydrodynamic correction with one step per alternating iteration, thereby reducing the need for empirical tuning. Numerical experiments for planar Couette flow, lid-driven cavity flow, and heat transfer in a bottom-heated cavity are carried out to investigate the performance of the improved FIM solver. The results demonstrate improved efficiency and robustness over the original FIM method, particularly at small Knudsen numbers, and favorable performance with both first- and second-order spatial discretizations.
1. Introduction
The original FIM method accelerated near-continuum flow calculations, but its convergence and spatial resolution remained limited. This work extends FIM to the Shakhov model and combines higher-order spatial discretization with revised moment and hydrodynamic solvers.
1. Introduction
- Earlier multigrid solvers became substantially slower as Kn decreased, while integrating velocity-space and spatial coarse corrections remained nontrivial.Source: §1 · p. 2
- The FIM framework alternates between a high-order moment system and consistent hydrodynamic equations, whose lower-cost macroscopic updates accelerate the iteration.Source: §1 · p. 2
- The original FIM was limited to the Boltzmann-BGK equation and first-order spatial discretization, while its hydrodynamic stage could require many iterations at small Kn or on refined grids.Source: §1 · p. 2
- The extended solver uses finite-volume linear reconstruction for second-order spatial accuracy in both the moment and consistent hydrodynamic systems.Source: §1 · p. 2
- An IMEX scheme treats the Shakhov collision term implicitly in the moment system to improve stability at small Kn.Source: §1 · p. 2
- A local Fourier analysis of a linearized one-dimensional BGK model explains why smaller Kn requires more hydrodynamic iterations and why refinement slows their convergence.Source: §1 · p. 2
- One Newton-multigrid step per alternating iteration generally gives an effective hydrodynamic correction and improves efficiency and robustness, particularly at small Kn.The method also reduces the need for parameter tuning.Source: §1 · p. 3
2. Governing equations and FIM framework
This section formulates the Shakhov kinetic model and its high-order moment and hydrodynamic systems, then discretizes them consistently and describes the alternating FIM solution procedure.
2.1. Boltzmann equation with Shakhov collision term
- The Shakhov model modifies the equilibrium distribution to specify the Prandtl number, reducing to BGK when Pr = 1.The collision frequency is inversely proportional to Kn, making the moment system stiff as Kn approaches zero.Source: §2.1 · p. 3, §2.2 · p. 4
2.2. High-order hyperbolic moment system and hydrodynamic equations
- The distribution is approximated by a truncated Hermite expansion of order M, yielding a globally hyperbolic moment system whose conservation-law subsystem gives the equations for density, velocity, and temperature.For the near-continuum regime considered, Grad’s basis uses local mean velocity and temperature.Source: §2.2 · p. 4
2.3. Spatial discretization of the moment system and linear reconstruction
- A unified finite-volume discretization reconstructs cell-interface states, with zero slopes giving first-order accuracy and linear reconstruction providing second-order accuracy.The hydrodynamic discretization is required to match the moment system’s reconstruction, fluxes, and HLL wave speeds.Source: §2.3 · pp. 5–6, §2.4 · p. 7
2.4. FIM framework
- The FIM iteration performs moment updates, extracts macroscopic fields, applies hydrodynamic corrections, then replaces the intermediate solution’s macroscopic fields with the corrected values.The hydrodynamic equations are a lower-cost stage, but small Kn can require many such iterations with the original solver.Source: §2.4 · pp. 6–7
3. Convergence analysis of the FIM method
A local Fourier analysis of a linearized one-dimensional BGK model studies FIM convergence in the near-continuum regime. It links slower convergence at smaller Kn and on refined grids to the need for a more efficient hydrodynamic solver.
3.1. Linearized model and the fast iterative method
- The analysis models FIM as alternating moment-stage and hydrodynamic-stage SGS iterations on a linearized one-dimensional BGK system.It assumes a relaxation parameter independent of Kn and uses a Fourier error-propagation matrix.Source: §3.1 · p. 8, §3.2 · p. 9
3.2. Results of the Fourier analysis and low-frequency behavior
- For fixed relaxation parameter, reducing a fixed low-frequency error mode by a constant factor requires SGS sweeps scaling as Δx^-2.At small Kn, the analysis permits increasing hydrodynamic iterations while reducing moment iterations at fixed total iteration count, lowering cost because the hydrodynamic stage is cheaper.Source: §3.2 · pp. 9–10
- As Kn decreases, the amplification factor approaches its hydrodynamic-stage leading-order value, so more hydrodynamic iterations are needed to maintain the same convergence factor.The analysis attributes the observed need for larger γ2 to the change in amplification behavior.Source: §3.2 · p. 10
3.3. Motivation for an improved hydrodynamic solver
- Because the hydrodynamic iteration count depends on both Kn and grid resolution, the proposed solver uses Newton’s method with multigrid for the resulting linear systems.The multigrid solve is intended to reduce errors across spatial frequencies and mitigate convergence deterioration under refinement.Source: §3.3 · p. 10
4. Improved FIM method
The improved FIM method extends the moment stage to the Shakhov model and replaces the hydrodynamic solver with Newton-multigrid iterations. Its design combines an implicit treatment of the Shakhov collision term with Kn-aware regularization and multigrid correction.
4.1. IMEX scheme for the Shakhov model
- The Shakhov moment iteration treats collision implicitly and convection explicitly, improving stability while keeping per-iteration cost essentially unchanged.The scheme is introduced for the stiffness of the Shakhov collision term in the high-order moment system.Source: §4.1 · pp. 10–11, §4 · p. 10
4.2. Newton-multigrid iteration for hydrodynamic equations
- Kn is incorporated alongside residual-based regularization because moment-supplied σ and q make the hydrodynamic system Kn-dependent.The modified Newton iteration addresses the possibility of a singular Jacobian.Source: §4.2 · p. 12
- Geometric multigrid solves the Newton linear system, using block LU-SGS smoothing and a V-cycle with two pre- and post-smoothing steps.The coarsest grid has 8 cells in one dimension or 8 × 8 cells in two dimensions.Source: §4.2 · p. 12
- One Newton-multigrid step generally corrects primary macroscopic quantities effectively, allowing the default hydrodynamic iteration count to be one without case-dependent tuning.The method combines velocity coarse-level corrections in the outer FIM iteration with spatial corrections in the inner Newton-multigrid iteration.Source: §4.2 · p. 13
4.3. FIM solvers
- Table 1 summarizes FIM variants using different moment-system and hydrodynamic schemes; the study focuses on FIM-3 and FIM-4 with first- and second-order discretizations.With Pr = 1 and first-order discretization, FIM-1, FIM-2, and FIM-3 correspond to solvers studied in, where FIM-3 was reported most efficient.Source: §4.3 · p. 13, Table 1 · p. 13
5. Numerical experiments
The experiments assess first- and second-order FIM solvers on Couette flow, lid-driven cavity flow, and bottom-heated cavity heat transfer. Across these tests, second-order discretization improves spatial resolution, while FIM-4 generally avoids empirical tuning and is more efficient at smaller Knudsen numbers.
5.1.1. Solution validation and spatial accuracy verification
- Second-order FIM solutions agree with DUGKS references for Couette flow, and the first- and second-order schemes achieve their expected convergence orders across the tested Knudsen numbers.The Couette profiles are computed on Nx = 128; grid-refinement errors use the Nx = 4096 solution as the benchmark.Source: §5.1.1 · p. 14, Figure 1 · p. 15, Figure 2 · p. 15
5.1.2. Performance of first-order FIM solvers
- For first-order Couette flow, FIM-4 reduces wall-clock time by factors of approximately 1.4 to 4.0 relative to carefully tuned FIM-3.FIM-4 uses one Newton-multigrid step per alternating iteration and avoids tuning γ2.Source: §5.1.2 · p. 14, Table 2 · p. 16
5.1.3. Performance of second-order FIM solvers
- For second-order Couette flow, FIM-4 is stable in all shown cases and outperforms carefully tuned FIM-3, with speedups rising from about 1.4 at Kn = 0.1 to over 17 at Kn = 0.001.A very large γ2 can cause loss of convergence; FIM-4 uses one Newton-multigrid step per alternating iteration.Source: Table 3 · p. 16, §5.1.3 · pp. 16–17
5.1.4. E
- At second order, the Couette FIM-4 regularization reduces iterations by more than 60% for Kn = 0.1, with the largest benefit in the reported M = 10 and M = 9 cases.The Kn-based strategy also suppresses residual oscillations and has little effect with first-order discretization.Source: Table 4 · p. 17, §5.1.4 · p. 17
5.2.1. Solution validation and spatial DOF comparison
- For the lid-driven cavity at Kn = 0.001, second-order FIM agrees with DUGKS on a 64 × 64 grid, while first-order results still visibly deviate on a 2048 × 2048 grid.This demonstrates substantially improved spatial resolution and reduced spatial degrees of freedom without sacrificing accuracy.Source: §5.2.1 · p. 18, Figure 6 · p. 19, Figure 7 · p. 19
5.2.2. Performance of FIM solvers
- In the first-order lid-driven cavity tests, FIM-4 is slightly slower than tuned FIM-3 at Kn = 0.1, but reduces wall-clock time by about 4× at Kn = 0.01 and maintains about 2× speedup at Kn = 0.001.The one-step Newton-multigrid correction requires no case-dependent γ2 tuning.Source: §5.2.2 · p. 18, Table 5 · p. 20
- In the lid-driven cavity, second-order FIM-4 exceeds tuned FIM-3 speedups of 4 at Kn = 0.01 and about 3 at Kn = 0.001, but is moderately less efficient at Kn = 0.1.The authors note that the full Jacobian system may partly account for the lower efficiency at Kn = 0.1.Source: §5.2.2 · p. 19, Table 6 · p. 20
5.3. Heat transfer in a bottom-heated cavity
- In bottom-heated cavity tests, second-order FIM solutions agree well with DUGKS, and FIM-4 is moderately less efficient at Kn = 0.1 but significantly faster at smaller Kn, particularly at second order.Temperature jumps near the lower corners become more apparent as Kn increases.Source: §5.3 · pp. 21–22, Figure 10 · p. 21, Table 7 · p. 22
6. Conclusion
The paper extends FIM to the Shakhov model, combining higher-resolution discretization with an implicit collision treatment and a Newton-multigrid hydrodynamic solver. Across three flow examples, the improved solver is more efficient and robust, especially at small Knudsen numbers.
6. Conclusion
- The solver extends FIM to the Shakhov model with second-order finite-volume discretization, implicit collision treatment, and a regularized Newton-multigrid hydrodynamic solver.The FIM iteration and Newton-multigrid correction provide acceleration in velocity and spatial spaces, respectively.Source: §6 · p. 23
- A linearized one-dimensional BGK Fourier analysis finds that smaller Knudsen numbers and refined grids require more hydrodynamic iterations and slow their convergence.These findings motivate replacing the time-integration-based hydrodynamic solver with Newton-multigrid.Source: §6 · p. 23
- The improved FIM solver enhances efficiency and robustness across Couette flow and two cavity-flow tests, especially at small Knudsen numbers.Its efficiency advantage over first-order discretization is greater with second-order spatial discretization at small Kn.Source: §6 · p. 23
- Generally, one Newton-multigrid step per alternating iteration maintains fast convergence while reducing case-dependent parameter tuning.The method is presented as a practical framework for steady-state high-order moment computations in near-continuum flows.Source: §6 · p. 23
Appendix A. Expression of the moment system (9)
The appendices specify the Shakhov moment-system terms and hydrodynamic wall-boundary construction. The boundary formulation uses ghost-cell macroscopic quantities and recursively defined auxiliary functions.
Appendix A. Expression of the moment system (9)
- Appendix A presents the moment-system expression and identifies Qα as the Hermite coefficient of the collision term, with negative-index fα set to zero.It then gives the Shakhov-specific collision terms.Source: Appendix A · p. 23
- The displayed moment-system expression includes terms involving θ, velocity components, and moment coefficients indexed by α.Source: Appendix A · p. 23
Appendix B. The hydrodynamic boundary conditions
- For a wall with specified temperature and velocity, ghost-cell macroscopic quantities are obtained from formulas in; the stress and heat flux use interior moment coefficients, with χ=1.The auxiliary functions S and K are defined recursively.Source: Appendix B · p. 24
Appendix C. Derivation of the Fourier convergence estimates
The appendix derives Fourier convergence estimates using the moment system, hydrodynamic variables, and SGS iterations. Its stated parameter conditions include a CFL restriction, while another condition holds automatically for M ≥2.
Appendix C. Derivation of the Fourier convergence estimates
- The derivation uses the moment system (26)–(27), hydrodynamic variables (29), and SGS iterations (31)–(34).Source: Appendix C · p. 24
- The first stated parameter condition is the model's CFL restriction.Source: Appendix C · p. 24
- The second stated condition holds automatically when M ≥2.Source: Appendix C · p. 24
Appendix C.1. Error evolution and the error-propagation matrix
The appendix tracks errors through moment and hydrodynamic SGS stages, then combines them through reconstruction to obtain the outer error propagation matrix. It also partitions moment errors into primary and higher-order components.
Appendix C.1. Error evolution and the error-propagation matrix
- Moment errors are defined before the forward sweep and after the forward and backward sweeps, yielding moment-stage error evolution.Source: Appendix C.1 · pp. 24–25
- The analysis splits each moment error into its first three components and its remaining M −2 components, and defines hydrodynamic-stage errors and initialization.Source: Appendix C.1 · p. 25
- Hydro-SGS updates and reconstruction complete the outer error after one FIM iteration, with the resulting propagation matrix denoted T(ε, k).The appendix shows block structures for matrices arising from the hydrodynamic SGS stage; ∗ marks a potentially nonzero block.Source: Appendix C.1 · p. 25
Appendix C.2. Order analysis of the error-propagation matrix
The error-propagation matrix has an O(ε) lower block row, so bounded errors in the nonequilibrium components become O(ε) after one FIM iteration.
Appendix C.2. Order analysis of the error-propagation matrix
- The lower block-row scaling lets the leading-order amplification analysis focus on hydrodynamic error components.Source: Appendix C.2 · p. 26
- The block-order analysis uses invertibility of the kf block under parameter range (C.1) and an O(ε) inverse Schur complement.The backward-sweep inverse blocks have the same order as those for the forward sweep.Source: Appendix C.2 · p. 26
- T21(ε, k) and T22(ε, k) are O(ε), while T11(ε, k) and T12(ε, k) are O(1).This block scaling implies that the nonequilibrium components of a bounded error become O(ε) after one FIM iteration.Source: Appendix C.2 · p. 26
Appendix C.3. Leading-order amplification factors
The leading-order analysis identifies three hydrodynamic amplification factors; at least M − 2 other eigenvalues are O(ε). The three hydrodynamic factors are then considered for comparison.
Appendix C.3. Leading-order amplification factors
- The leading-order amplification factor associated with an eigenvalue μ of Ah represents one SGS iteration, where Ah is the 3 × 3 upper-left block of A.Source: Appendix C.3 · p. 27
- At least M − 2 eigenvalues of T(ε, k) are O(ε) for sufficiently small ε.The argument uses Gershgorin disks and m ≥ M − 2 zero diagonal entries of D(k).Source: Appendix C.3 · p. 27
- The FIM iteration has three hydrodynamic leading-order amplification factors, with the remaining M − 2 eigenvalues O(ε).The comparison of the three hydrodynamic factors is deferred to the next subsection.Source: Appendix C.3 · p. 27
Appendix C.4. Ordering of the three leading factors
The three leading-factor moduli are even in the wave number and nonincreasing with its absolute value. They are strictly ordered over the stated nonzero range and all equal one at zero.
Appendix C.4. Ordering of the three leading factors
- The modulus |G(0)μ(k)| is even in ζ = k∆x and nonincreasing in |ζ| over 0 ≤ |ζ| ≤ π.Source: Appendix C.4 · p. 27
- For 0 < |k∆x| ≤ 2π/5, the compared factor moduli satisfy a strict inequality; all three moduli equal one at k = 0.Source: Appendix C.4 · p. 28
- For 0 < |k∆x| ≤ 2π/5, positivity of the factors in the subtraction of squared moduli establishes inequality (C.21).Source: Appendix C.4 · pp. 28–29
Appendix C.5. First-order correction for γ1 = 1
For γ1 = 1, the appendix derives a first-order correction to the leading-order amplification factor. The correction can be analyzed using only the first four error components because higher components are not coupled to them through first order in ε.
Appendix C.5. First-order correction for γ1 = 1
- The first-order correction to the amplification factor is derived along an eigenvalue branch with a nonzero limit.The derivation assumes γ1 = 1 and expands the error and intermediate amplitudes along the leading-order branch.Source: Appendix C.5 · p. 29
- Only the first four components of the error vector are needed through first order in ε, since the other components are not coupled to them at that order.The fifth component is O(ε²) under the stated conditions, including M ≥ 4.Source: Appendix C.5 · p. 29
- The first-order correction is obtained by expanding the moment-sweep equations and then incorporating the hydrodynamic stage, which changes only the first three components.The derivation applies the forward and backward sweeps and uses a recurrence for the correction; the branch μ = 0 additionally requires γ2(k∆x)² ≪ 1 for the cited expansion.Source: Appendix C.5 · pp. 29–31