PICurv 0.1.0
A Parallel Particle-In-Cell Solver for Curvilinear LES
 
Loading...
Searching...
No Matches
Configuration Reference: Solver YAML

For the full commented template, see:

# ==============================================================================
#                 PICurv Master Solver Configuration Template
# ==============================================================================
#
# PURPOSE:
# This file defines the NUMERICAL STRATEGY for the simulation. It is divided
# into user-friendly structured sections for common settings, and a flexible
# 'petsc_passthrough_options' section for advanced users to inject any valid
# PETSc command-line flag.
#
# ==============================================================================

# ==============================================================================
# 1. SOLVER OPERATION MODE
# ==============================================================================
operation_mode:
  # 'solve':      (Default) Advance the flow equations. With case.yml start_step > 0
  #               this is also the restart mode: the saved state is loaded, then solved on.
  # 'load':       Replay stored fields; nothing is solved. Every replayed step needs its
  #               own checkpoint in the --restart-from run (io.data_output_frequency: 1).
  # 'analytical': Evaluate a closed-form field each step (analytical_type below).
  eulerian_field_source: "solve"  # -> -euler_field_source
  # Optional analytical solution selector used when eulerian_field_source = "analytical".
  # Supported values today: "TGV3D", "ZERO_FLOW", "UNIFORM_FLOW".
  # TGV3D currently requires case.yml to use grid.mode: programmatic_c.
  # ZERO_FLOW and UNIFORM_FLOW also support file-based grid ingestion.
  analytical_type: "TGV3D"        # -> -analytical_type
  # Parameters for analytical_type: "UNIFORM_FLOW".
  # Uncomment this block only when analytical_type = "UNIFORM_FLOW".
  # uniform_flow:
  #   u: 0.0                      # [m/s] -> -analytical_uniform_u (divided by velocity_ref)
  #   v: 0.0                      # [m/s] -> -analytical_uniform_v
  #   w: 0.0                      # [m/s] -> -analytical_uniform_w

# Verification-only source overrides.
# Use these only when ordinary end-to-end setup cannot expose the behavior you need to test.
# New verification-pathway source injections must be implemented in
# include/verification_sources.h and src/verification_sources.c.
# verification:
#   sources:
#     diffusivity:
#       mode: "analytical"        # Verification-only selector; currently fixed to analytical
#       profile: "LINEAR_X"       # Supported today: LINEAR_X
#       gamma0: 1.0e-3             # [m^2/s] -> -verification_diffusivity_gamma0
#       slope_x: 2.0e-4            # [m/s], diffusivity per metre -> -verification_diffusivity_slope_x
#
#     scalar:
#       mode: "analytical"        # Verification-only selector; currently fixed to analytical
#       profile: "CONSTANT"       # Supported: CONSTANT, LINEAR_X, SIN_PRODUCT
#       value: 1.0                 # CONSTANT -> -verification_scalar_value
#       # phi0: 0.0                # LINEAR_X -> -verification_scalar_phi0
#       # slope_x: 1.0             # [1/m] LINEAR_X -> -verification_scalar_slope_x
#       # amplitude: 1.0           # SIN_PRODUCT -> -verification_scalar_amplitude
#       # kx: 3.141592653589793    # [1/m] SIN_PRODUCT -> -verification_scalar_kx
#       # ky: 3.141592653589793    # SIN_PRODUCT -> -verification_scalar_ky
#       # kz: 3.141592653589793    # SIN_PRODUCT -> -verification_scalar_kz
#
# Scalar verification is intended for runtime diagnostics such as
# output/analysis/metrics/scatter_metrics.csv. It bypasses only the particle Psi evolution path and
# reuses the production scatter operator.

# --- Scalar Transport Properties ---
# These affect Eulerian diffusivity and particle scalar/Brownian models.
scalar_transport:
  schmidt_number: 1.0              # -> -schmidt_number
  turbulent_schmidt_number: 0.7    # -> -turb_schmidt_number
  iem_constant: 2.0                # -> -iem_constant; C_IEM in the IEM mixing rate (default 2.0)

# ==============================================================================
# 2. PRIMARY SOLVER STRATEGY & TOLERANCES
# ==============================================================================
# --- Main Time-Stepping Scheme ---
strategy:
  # Accepted values:
  #   "Dual Time Picard Jameson RK"  -> implicit dual-time pseudo-stepping (recommended)
  #   "Explicit RK4"                 -> explicit fourth-order Runge-Kutta
  #   "Newton Krylov"                -> matrix-free PETSc SNES/KSP solve
  momentum_solver: "Dual Time Picard Jameson RK" # -> -mom_solver_type
  central_diff: false # [true/false] -> -central

# --- Convergence Criteria for the Momentum Solver ---
tolerances:
  max_iterations: 50   # [Integer] -> -mom_max_pseudo_steps
  # absolute_tol: 1.0e-8   # DEPRECATED -> -mom_atol. Takes no part while a
  #                        # residual tolerance is active (the default).
  relative_tol: 1.0e-5 # [Float] -> -mom_rtol
  residual_absolute_tol: 1.0e-8 # [Float] -> -mom_resid_atol. DIMENSIONLESS:
                                # converges when |R| <= this * a0*|Ucont|_inf/dt. Sufficient alone.
  residual_relative_tol: 1.0e-3 # Use 1.0e-2 for looser exploratory LES
  # step_tol: 0.0       # DEPRECATED -> -imp_stol. Accepted but unused by every
  #                     # momentum solver.

# --- Advanced Momentum Solver Controls (optional) ---
# This section exposes newer dual-time controls directly (no passthrough needed).
momentum_solver:
  # Solver-specific controls for Dual Time Picard Jameson RK.
  # Note: 'dual_time_picard_rk4' is an accepted deprecated alias for this block.
  dual_time_picard_jameson_rk:
    max_pseudo_steps: 50         # -> -mom_max_pseudo_steps
    # absolute_tol: 1.0e-8   # DEPRECATED -> -mom_atol. Takes no part while a
    #                        # residual tolerance is active (the default).
    relative_tol: 1.0e-5         # -> -mom_rtol
    # step_tol: 0.0              # DEPRECATED -> -imp_stol; accepted but unused.
    pseudo_cfl:
      # pseudo_cfl.* values are dimensionless Courant numbers, NOT fractions of dt.
      # The solver computes dtau = pseudo_cfl / lambda_max, where lambda_max is the global
      # maximum spectral radius of the convective operator (sum of |face fluxes| / cell volume,
      # global MPI max). This makes pseudo_cfl flow- and grid-independent.
      # Stability limit for 4-stage Jameson RK: ~2.83 (imaginary-axis CFL). Use 0.5-2.0 in practice.
      initial: 0.5               # -> -pseudo_cfl          (default: 0.5; ~half the stability limit)
      minimum: 0.001             # -> -min_pseudo_cfl      (default: 0.001)
      maximum: 2.0               # -> -max_pseudo_cfl      (default: 2.0; stability limit ~2.83)
      growth_factor: 1.1         # -> -pseudo_cfl_growth_factor   (default: 1.1; must be >= 1)
      reduction_factor: 0.75     # -> -pseudo_cfl_reduction_factor (default: 0.75; must be in (0,1))
    jameson_residual_noise_allowance_factor: 1.1  # -> -mom_dt_jameson_residual_norm_noise_allowance_factor
    # Rejection threshold: a pseudo-time trial is rolled back if the EMA-smoothed step-to-step
    # residual ratio exceeds this value. 1.1 allows 10% residual growth before reducing CFL.
    # Raise toward 1.2-1.5 for flows with non-monotonic residual histories; lower toward 1.05
    # for strict monotonic convergence. Must be >= 1.
    # rk4_residual_noise_allowance_factor is a deprecated alias for the above; use only one.
    ratio_ema_alpha: 0.3          # -> -mom_ratio_ema_alpha (default: 0.3; range [0, 1])
    # Exponential moving average coefficient for the step-to-step residual ratio used in the
    # trial-rejection decision. The smoothed ratio is:
    #   smoothed = alpha * raw_ratio + (1 - alpha) * smoothed_prev
    # alpha = 1.0  : raw ratio (original behavior, most aggressive rejection)
    # alpha = 0.3  : moderate smoothing; ~3-4 consecutive bad trials needed to trigger rejection
    # alpha = 0.0  : ratio never updates (disables ratio-based rejection entirely)
    # Increase alpha if the solver rejects too conservatively on noisy residual histories.

  # Solver-specific controls for Newton Krylov. Uncomment only when
  # strategy.momentum_solver is "Newton Krylov". Omitted fields retain the
  # defaults established by src/momentum_newton_krylov.c and PETSc.
  # newton_krylov:
  #   jacobian:
  #     type: "finite_difference"
  #     finite_difference:
  #       mode: "matrix_free"
  #   preconditioner:
  #     model: "none"
  #     # To enable the frozen-momentum point-block preconditioner, use:
  #     # model: "frozen_momentum_jacobian"
  #     # structure:
  #     #   type: "point_block"
  #   nonlinear_solver:
  #     method: "newtonls"           # -> -mom_nk_snes_type
  #     absolute_tolerance: 1.0e-10  # -> -mom_nk_snes_atol
  #     relative_tolerance: 1.0e-8   # -> -mom_nk_snes_rtol
  #     step_tolerance: 1.0e-12      # -> -mom_nk_snes_stol
  #     max_iterations: 12           # -> -mom_nk_snes_max_it
  #     line_search:
  #       type: "bt"                 # -> -mom_nk_snes_linesearch_type
  #     eisenstat_walker:
  #       enabled: true              # false keeps the KSP tolerance fixed
  #       version: 3
  #       initial_relative_tolerance: 0.3
  #       maximum_relative_tolerance: 0.9
  #       gamma: 1.0
  #       exponent: 1.618033988749895
  #       safeguard_exponent: 1.618033988749895
  #       safeguard_threshold: 0.1
  #   linear_solver:
  #     method: "gmres"              # -> -mom_nk_ksp_type
  #     absolute_tolerance: 1.0e-10  # -> -mom_nk_ksp_atol
  #     relative_tolerance: 1.0e-6   # -> -mom_nk_ksp_rtol
  #     max_iterations: 400          # -> -mom_nk_ksp_max_it
  #     gmres:
  #       restart: 80                # -> -mom_nk_ksp_gmres_restart

# ==============================================================================
# 3. POISSON SOLVER CONFIGURATION
# ==============================================================================
poisson_solver:
  # Solves for pressure correction Phi, then the runtime updates pressure P.
  # The outer linear solver uses PETSc KSP under the hood, but this block uses
  # PICurv-facing names for the common controls.
  # Note: 'pressure_solver' is an accepted deprecated alias for this block name.
  method: "fgmres"             # -> -ps_ksp_type. Options: fgmres (default), cg. gmres, lgmres
                               #    and bcgs are refused: the multigrid preconditioner is not a
                               #    fixed linear operator, which only a flexible method tolerates
  absolute_tolerance: 1.0e-5   # -> -ps_ksp_atol
  relative_tolerance: 1.0e-11  # -> -ps_ksp_rtol
  max_iterations: 50           # -> -ps_ksp_max_it
  gmres:
    # Only valid with method: fgmres (cg has no restart).
    restart: 20                # -> -ps_ksp_gmres_restart

  preconditioner:
    # Only multigrid is supported for the outer Poisson preconditioner today.
    # Other values are rejected until the C runtime grows a non-PCMG path.
    type: "multigrid"          # -> -ps_pc_type mg

  # --- Geometric Multigrid (PCMG) Settings ---
  multigrid:
    levels: 3         # [Integer] -> -mg_level
    # Choose levels so the COARSEST grid lands at roughly 1e3-1e4 unknowns.
    # The coarse solver below replicates an LU factor on every rank, so a
    # coarse grid that is too large costs memory and time on all ranks; a
    # coarse grid that is too small wastes a level. As a rule of thumb a 5M
    # cell grid wants 5 levels, not 4.
    # Coarsenability: each level halves a direction as IM -> (IM+1)/2, so IM
    # must stay odd at every level for the coarsening to be exact. An even IM
    # logs "can't be consistently coarsened further" and keeps going on a
    # slightly misaligned coarse grid.
    # NOTE: levels is bounded by the MPI rank layout, not chosen freely. Every
    # level must leave each rank at least stencil_width nodes per axis, and the
    # stencil width is 3 whenever ANY axis is periodic (2 otherwise). Exceeding
    # it aborts during DM creation with:
    #   Local x-width of domain x 2 is smaller than stencil width s 3
    # Constraint: floor((cells_per_axis / 2^(levels-1) + 1) / ranks_per_axis) >= stencil_width
    # picurv validate cannot catch this; it does not see the rank layout.
    # Worked maxima: docs/pages/25_Pressure_Poisson_GMRES_Multigrid.md
    pre_sweeps: 2     # [Integer] -> -mg_pre_it
    post_sweeps: 2    # [Integer] -> -mg_post_it
    # Sweeps before and after the coarse-grid correction on each smoothed level.
    # When they differ, the post-smoother copies the level solver configured below.
    semi_coarsening:
      i: false        # [true/false] -> -mg_i_semi
      j: false        # [true/false] -> -mg_j_semi
      k: true         # [true/false] -> -mg_k_semi
      
    cycle: "v"        # Currently supported: "v"
    mode: "multiplicative" # Currently supported: "multiplicative"

    # --- Smoother / Coarse Solver Configuration Per MG Level ---
    # PETSc/PICurv level numbering uses level_0 as the coarsest grid; larger
    # level numbers are progressively finer.
    # Friendly aliases: `method` -> ksp_type, `preconditioner` -> pc_type.
    # Other keys are forwarded verbatim under PETSc's level prefix: level_0
    # uses -ps_mg_coarse_<key>; positive levels use -ps_mg_levels_N_<key>.
    # Friendly `max_it`, `rtol`, `atol` map to ksp_max_it, ksp_rtol, ksp_atol.
    # Supported direct keys include ksp_type, pc_type, ksp_max_it, ksp_rtol, ksp_atol.
    level_solvers:
      level_0:
        # level_0 is the COARSE SOLVE at the base of the V-cycle, not a
        # smoother. The smoothers on level_1..N kill high-frequency error;
        # this solve kills the low-frequency error they are blind to.
        #
        # Because it sits inside a preconditioner, it must be a FIXED LINEAR
        # operator: the same input vector must always produce the same output.
        # Krylov methods (gmres, fgmres, cg, bcgs) are not - their subspace
        # adapts to the input, which makes the whole MG preconditioner
        # nonlinear. The outer FGMRES then tolerates it for constructing the
        # solution, but its Arnoldi recurrence decouples from the true
        # residual and the convergence test starts passing on a number that no
        # longer describes b-Ax. Pinning ksp_max_it does NOT fix this.
        # PICurv logs a warning at startup if you configure one here.
        #
        # Recommended by coarse-grid size:
        #   <= ~1e4 unknowns : preonly + redundant   (replicated direct LU)
        #   larger           : preonly + telescope   (LU on a rank subset)
        #   last resort      : a Krylov method with tolerances set against the
        #                      TRUE residual, not the tracked one
        method: "preonly"           # -> -ps_mg_coarse_ksp_type
        preconditioner: "redundant" # -> -ps_mg_coarse_pc_type
        # ksp_max_it: 30           # -> -ps_mg_coarse_ksp_max_it
        # ksp_rtol: 1.0e-3         # -> -ps_mg_coarse_ksp_rtol
        # ksp_atol: 1.0e-8         # -> -ps_mg_coarse_ksp_atol
      # level_1..N are SMOOTHERS. Richardson, Jacobi, Chebyshev and SOR are all
      # fixed linear operators, which is exactly why they belong here.
      level_1:
        method: "richardson"
        preconditioner: "bjacobi"
      level_2:
        method: "richardson"
        preconditioner: "bjacobi"

# ==============================================================================
# 4. INTERPOLATION
#    Controls grid-to-particle interpolation numerics.
# ==============================================================================
interpolation:
  # Method for interpolating Eulerian fields to particle positions.
  # Options:
  #   - "Trilinear"       (default) Direct trilinear from 8 nearest cell centers.
  #                        Measured L2 order 1.97 on TGV3D, uniform Cartesian grid.
  #   - "CornerAveraged"  Legacy two-stage: center->corner average, then trilinear
  #                        from corners. On the same test: 3-4x Trilinear's L2 error,
  #                        L2 order 1.68, maximum-error order 0.96.
  method: "Trilinear"     # -> -interpolation_method

# ==============================================================================
# 5. PETSC PASSTHROUGH OPTIONS (FOR ADVANCED USERS)
# ==============================================================================
petsc_passthrough_options:
  # Use this only for advanced PETSc controls that do not have structured YAML
  # above. Structured values and passthrough flags target the same PETSc options;
  # passthrough wins when the same flag is listed in both places.

  # --- MG level preconditioner tuning examples ---
  # jacobi / sor: simple smoother PCs; no extra nested setup is normally needed.
  # -ps_mg_levels_2_pc_type: "sor"
  # -ps_mg_levels_2_pc_sor_omega: 1.0      # PETSc SOR relaxation; positive real

  # ilu / lu: factor PCs; use PETSc factor controls when pivots are fragile.
  # -ps_mg_coarse_pc_type: "ilu"
  # -ps_mg_coarse_pc_factor_levels: 1    # nonnegative integer
  # -ps_mg_coarse_pc_factor_shift_amount: 1.0e-10 # nonnegative real

  # bjacobi: nested block solves. Inspect exact nested prefixes with -ps_ksp_view
  # when tuning sub-KSP/sub-PC options for your PETSc version.

  # --- Additional Momentum Controls (if you prefer passthrough style) ---
  # Prefer the structured momentum solver blocks above.
  # These passthrough flags override the structured block if both are present.
  # -pseudo_cfl: 0.5          # dimensionless CFL = dtau * lambda_max (spectral-radius-based)
  # -max_pseudo_cfl: 2.0
  # -min_pseudo_cfl: 0.001
  # -pseudo_cfl_growth_factor: 1.1
  # -pseudo_cfl_reduction_factor: 0.75
  # -mom_dt_jameson_residual_norm_noise_allowance_factor: 1.1
  # -mom_ratio_ema_alpha: 0.3

solver.yml controls numerical strategy and solver internals.

1. operation_mode

operation_mode:
eulerian_field_source: "solve"
analytical_type: "TGV3D"
uniform_flow:
u: 0.0
v: 0.0
w: 0.0

Mappings:

  • eulerian_field_source -> -euler_field_source (solve, load, analytical)
  • analytical_type -> -analytical_type
  • uniform_flow.u/v/w -> -analytical_uniform_u/-analytical_uniform_v/-analytical_uniform_w when analytical_type: "UNIFORM_FLOW"; physical velocities, divided by velocity_ref

uniform_flow is only valid when analytical_type: "UNIFORM_FLOW".

2. strategy

strategy:
momentum_solver: "Dual Time Picard Jameson RK"
central_diff: false

Mappings:

  • momentum_solver -> -mom_solver_type (picurv accepts Explicit RK4, Dual Time Picard Jameson RK, or Newton Krylov)
  • central_diff -> -central

central_diff selects the convective flux in Convection() (src/rhs.c). false uses the upwind-biased QUICK scheme, whose interpolation switches on the sign of the face flux; true uses the central average of the two neighbouring cells. Any LES model selects the central flux regardless of this key, because the branch is taken on les || central.

Central convection adds no dissipation of its own. Its transport term cannot see a two-cell oscillation at all - the discrete derivative of such a pattern is zero - so only viscosity, molecular or subgrid, removes one. On a grid whose cell Reynolds number is large along some direction that removal can be slower than the run; the grid generator's [Grid-Scale Damping] report states the time per direction. QUICK damps such a pattern but makes the residual non-differentiable where a face flux changes sign, which a Newton solver's finite-difference Jacobian does not model.

Older boolean toggles are not supported; use strategy.momentum_solver. Only implemented momentum solver values are accepted by picurv and the C runtime. The deprecated Dual Time Picard RK4 display name, dual_time_picard_rk4 solver block, and rk4_residual_noise_allowance_factor key remain readable compatibility aliases; generated controls always use the canonical Jameson names.

3. Momentum Solver Entries

The table below is generated from the canonical mapping in normalize_momentum_solver_type() and regenerated by make docs-inventory.

Value Maps to Status
Dual Time Picard Jameson RKDUALTIME_PICARD_JAMESON_RKsupported
Dual Time Picard RK4DUALTIME_PICARD_JAMESON_RKdeprecated - alias of Dual Time Picard Jameson RK
Explicit RK4EXPLICIT_RKsupported
Newton Krylovnewton_krylovsupported

Each entry follows the capability-entry contract described in Documentation Extension Framework.

Explicit RK4

Identity. strategy.momentum_solver: "Explicit RK4" -> -mom_solver_type EXPLICIT_RK -> MOMENTUM_SOLVER_EXPLICIT_RK -> MomentumSolver_Explicit_RungeKutta4.

What it does. Advances momentum with an explicit four-stage Runge-Kutta step. No pseudo-time iteration and no linear solve: each physical step is a fixed sequence of residual evaluations.

When to choose it. When the timestep is already small for physical reasons and you want the cheapest, most predictable step. It is also the clearest baseline when diagnosing whether a problem lies in the momentum solve or elsewhere, because it has almost no machinery of its own.

Parameters it owns. None. Explicit RK4 takes no solver-specific block; pseudo-CFL and Newton-Krylov controls are rejected for it rather than ignored.

Interactions. Stability is governed by the physical timestep alone, so it is the solver most sensitive to grid refinement. The Poisson solve and boundary treatment are unchanged, and the driven periodic flux controller settles to the same law as under Picard (5.7 Known limitations).

Diagnostics. The startup banner names the solver as Explicit 4 stage Runge-Kutta. A step beyond the stability limit stops the run at the step it occurs, naming the explicit stability limit and the approximate viscous bound dt < 2.8 / (4 nu (1/dx^2 + 1/dy^2 + 1/dz^2)); it no longer surfaces later as a Poisson failure blamed on the multigrid depth.

Evidence. Integration verified - make smoke runs a stable flat-channel step and a step past the limit, and asserts the second stops with the stability message. Analytically verified - explicit-rk4-order-2026-09-21: on the two-dimensional Taylor-Green vortex, velocity converges at second order in time and space and pressure at second order in space; on the driven laminar channel the profile error falls at order 1.97 and the bulk velocity matches the controller law to 1e-4.

Limitations. No step control: the user must keep dt under the viscous and convective limits, and a run past them stops rather than recovers. The projection splitting limits the temporal order to two for velocity and one for pressure, below the four stages' own order. Unsuitable for stiff or strongly convective cases, where the stable dt is far below what the physics needs.

Dual Time Picard Jameson RK

Identity. strategy.momentum_solver: "Dual Time Picard Jameson RK" -> -mom_solver_type DUALTIME_PICARD_JAMESON_RK -> MOMENTUM_SOLVER_DUALTIME_PICARD_JAMESON_RK -> MomentumSolver_DualTime_Picard_JamesonRK.

What it does. Drives each physical step to convergence with a pseudo-time Picard iteration, smoothed by a four-stage Jameson Runge-Kutta scheme, under an adaptive pseudo-CFL controller that accepts or rejects each trial.

When to choose it. The default production solver. It tolerates far larger physical timesteps than Explicit RK4 because the pseudo-time loop absorbs the stiffness, and it does not require the PETSc SNES configuration Newton-Krylov does.

Parameters it owns. The dual_time_picard_jameson_rk block: max_iterations (alias max_pseudo_steps), relative_tol, step_tol, the pseudo_cfl sub-block (initial, minimum, maximum, growth_factor, reduction_factor), and jameson_residual_noise_allowance_factor.

Note
max_iterations caps accepted pseudo-iterations, not attempts. A rejected trial rolls back without consuming the budget; a separate hard cap of 3 x max_iterations bounds total attempts.

Interactions. Pseudo-CFL is exclusive to this solver: it is neither accepted from the Newton-Krylov block nor shown in a Newton-Krylov banner. Acceptance and rollback are global across blocks and MPI ranks. The controller-selected pseudo-CFL carries into the next physical timestep within one process only - it is not written to the checkpoint, so a restarted run begins again from the configured pseudo_cfl.initial.

Diagnostics. The per-trial file log records one row per accepted trial, with the residual ratio and the pseudo-CFL used; rejected trials are not given their own row there, so a quiet log does not mean no rejections occurred. A step that exhausts its budget logs "reached N total attempts without convergence" and continues from the last accepted state - that message counts total attempts, accepted plus rejected.

Evidence. Unit verified - make unit-solver. Integration verified - make smoke. Production exercised in examples/flat_channel and examples/bent_channel. Analytically verified - duct-poiseuille-picard-2026-09-18: steady laminar square-duct Poiseuille flow at second order in space; tgv2d-picard-order-2026-09-18: the two-dimensional Taylor-Green vortex at second order in time and space for velocity and pressure; pipe-poiseuille-curvilinear-2026-09-18: Hagen-Poiseuille flow on a strongly non-orthogonal swept circle at second order; periodic-channel-laminar-picard-2026-09-18: the driven periodic channel at second order, every step converged.

Limitations and full treatment. Control law, cadence, and tuning guidance are at Dual-Time Picard Jameson RK Momentum Solver. The results hold when every step meets its pseudo-time tolerance: the square duct at Re = 100 with dt = 1 hit the iteration cap on every step and finished 18% off in pressure gradient, so read the per-step history before trusting a run. Only central differencing has been measured. The periodic wall-bounded stall once recorded at 5.7 Known limitations was re-characterized on 2026-09-18 and does not reproduce. A wall-modelled LES channel at Re_tau ~ 1000 runs under this solver and reproduces Lee & Moser within the criteria of wmles-channel-retau1000-werner-2026-10-09; all 50,000 steps of its inspected segment were accepted with no rejected pseudo-iteration. Wall-resolved turbulent channels at production resolution remain uncharacterized under it.

Newton Krylov

Identity. strategy.momentum_solver: "Newton Krylov" -> -mom_solver_type newton_krylov -> MOMENTUM_SOLVER_NEWTON_KRYLOV -> MomentumSolver_NewtonKrylov.

What it does. Solves the momentum system with a PETSc SNES Newton-Krylov method using a matrix-free finite-difference Jacobian, rather than iterating in pseudo-time.

When to choose it. When the pseudo-time approach converges slowly or not at all, and you want true Newton convergence on the momentum system. It exposes PETSc's solver machinery directly, which is an advantage when you know what you want from it and a liability when you do not.

Parameters it owns. The newton_krylov block: jacobian (a strict discriminated configuration - explicit type: finite_difference requires finite_difference.mode: matrix_free), preconditioner (model and structure), nonlinear_solver (including line_search), and linear_solver.

Interactions. Accepted combinations are finite-difference/matrix-free with either no preconditioner or a frozen momentum Jacobian / point-block preconditioner. Any other Jacobian type or finite-difference mode is rejected. Raw petsc_passthrough_options are applied last, but an incompatible raw -mom_nk_pc_type override is rejected by the runtime.

Diagnostics. SNES and KSP convergence reasons are reported per step. A DIVERGED_LINEAR_SOLVE or DIVERGED_BREAKDOWN with preconditioner.model: none generally indicates the unpreconditioned matrix-free limitation rather than a configuration error.

Evidence. Unit verified - make unit-newton-krylov. Integration verified - make unit-momentum-newton-boundary-fixedpoint on a production-sized straight duct. Production exercised - turbulent-channel-nk-2026-09-29: 96 x 96 x 256 cells on 144 ranks, all 10,000 inspected steps converged and committed; see 10.2 Retained turbulent-channel campaign (2026-09-29) for the retained measurements and DNS discrepancy. Analytically verified - pipe-poiseuille-nk-three-grid-2026-10-10: Hagen-Poiseuille flow on the swept circle with 9, 17 and 33 cells across, observed order 2.00 then 2.04.

Limitations and full treatment. Newton–Krylov Momentum Solver carries the scope limits, preconditioner findings, and tuning guidance. Supported within that scope. The turbulent-channel campaign used preconditioner.model: none; frozen_momentum_jacobian is supported on the laminar curved-duct campaign humphrey-laminar-bend-nk-pointblock-2026-10-07 (10.1 Retained laminar curved-duct campaign (2026-10-07)), which establishes correct execution, not a speed-up. Quantitative turbulent-channel DNS agreement and parallel scaling remain unverified.

Dual Time Picard RK4 (deprecated)

Identity. strategy.momentum_solver: "Dual Time Picard RK4" - a deprecated alias that normalizes to Dual Time Picard Jameson RK. The C enum keeps MOMENTUM_SOLVER_DUALTIME_PICARD_RK4 as an alias of the Jameson constant.

Status. Deprecated. It remains readable so archived case files continue to load; generated control artifacts always emit the canonical Jameson names.

Migration. Replace the value with Dual Time Picard Jameson RK, rename a dual_time_picard_rk4 solver block to dual_time_picard_jameson_rk, and rename rk4_residual_noise_allowance_factor to jameson_residual_noise_allowance_factor. Do not set canonical keys and their deprecated aliases together - that is a validation error.

4. Source, Interpolation and Convergence Entries

Three further selector families live in solver.yml. Their generated inventories:

Value Maps to
analyticalanalytical
loadload
solvesolve
Value Maps to
CornerAveraged1
Trilinear0
Value Maps to
periodic_deterministicPERIODIC_DETERMINISTIC
statistical_steadySTATISTICAL_STEADY
steady_deterministicSTEADY_DETERMINISTIC
transientTRANSIENT

solve

Identity. operation_mode.eulerian_source: solve -> the run evolves the Eulerian fields numerically through FlowSolver.

What it does. Advances velocity and pressure by integrating the governing equations from the configured initial condition.

When to choose it. The normal mode. Choose load instead when you want to analyse or drive particles through a field somebody already computed, and analytical when you want an exact field with no flow solve at all.

Parameters it owns. None directly; it makes the whole solver.yml numerics block meaningful, which the other two sources largely ignore.

Interactions. Requires a valid initial condition. Restart and continuation apply to this source only.

Diagnostics. Per-step momentum and Poisson convergence logs under <run.runtime_logs>/.

Evidence. Production exercised - examples/flat_channel runs this source. Analytically verified - tgv2d-picard-order-2026-09-18 and duct-poiseuille-picard-2026-09-18 measure the solves it drives; the solver paths carry their full evidence at 3. Momentum Solver Entries.

Limitations. Cost scales with the flow solve; the other two sources are far cheaper when the flow field is not the object of study.

load

Identity. operation_mode.eulerian_source: load -> fields are read from previously written output rather than evolved.

What it does. Replays Eulerian state from an existing run's checkpoints: at each step it reads the committed checkpoint of that same step, so downstream stages - particle transport above all - see the stored sequence rather than a recomputed one.

When to choose it. Post-processing an existing solution, or transporting particles through a stored field without recomputing it.

Parameters it owns. The restart/source directory controls that locate the stored field.

Interactions. The stored field must match the configured grid. Solver tolerances and momentum selection have no effect. Because every step reads its own checkpoint, the source run needs a committed checkpoint at every step the replay serves: write it with io.data_output_frequency: 1 over that span.

Diagnostics. Startup reports the resolved source directory and the step loaded.

Evidence. Regression verified - make smoke restarts the flat-channel particle case with eulerian_field_source: load in its restart-variant sequence. Analytically verified - eulerian-source-domain-modes-2026-09-21: a five-step replay reproduced every stored step bitwise, step 0 to 2e-15.

Limitations. No time evolution of the Eulerian state; the field is what was stored, at the cadence it was stored. A source written at a coarser cadence cannot be replayed step by step.

analytical

Identity. operation_mode.eulerian_source: analytical -> -analytical_type selects a closed-form field.

What it does. Supplies an exact analytical velocity field instead of solving for one.

When to choose it. Verification. An exact field makes an error norm meaningful, which is why the particle verification examples use it.

Parameters it owns. The analytical mode selector and its parameters. See Analytical Solution Modes.

Interactions. File-grid support is limited to the non-custom analytical modes.

Diagnostics. Startup reports the selected analytical mode.

Evidence. Production exercised - examples/drift_uniform_flow runs this source; the verification examples that rely on it are catalogued in Example Catalog. Analytically verified - solution-monitoring-tgv3d-2026-09-18: the TGV3D velocity and pressure match the closed form to 9e-16 at every step; uniform-drift-2026-09-18: a cloud in UNIFORM_FLOW drifts at exactly the carrier velocity.

Limitations. Only the shipped analytical forms are available; there is no user-supplied expression path.

Trilinear

Identity. interpolation.method: Trilinear -> -interpolation_method 0 -> InterpolationMethod -> InterpolateEulerFieldToSwarm.

What it does. Interpolates directly from the eight surrounding cell centres to the particle position.

When to choose it. The default, and the right choice unless you are reproducing legacy behaviour. It interpolates in the cell's own logical coordinates, so the same scheme applies on curvilinear grids; its measured order, 1.97, is from a uniform grid.

Parameters it owns. None.

Interactions. Requires valid DMSwarm_CellID and interpolation weights, which the settle step establishes.

Diagnostics. Interpolation accuracy is measurable with the shipped interpolation_test example.

Evidence. Production exercised - examples/interpolation_test runs this path; see 3. Particle Verification Family for what that case establishes. Analytically verified - tgv-interpolation-2026-09-18: 0.67% relative error against the TGV3D field on the shipped 32^3 grid; interpolation-methods-tgv-2026-09-21: 2.63% and 0.669% on 16^3 and 32^3, order 1.97.

Limitations. Accuracy degrades on strongly distorted cells, as any trilinear scheme does.

CornerAveraged

Identity. interpolation.method: CornerAveraged -> -interpolation_method 1.

What it does. Averages cell-centre values to cell corners first, then interpolates trilinearly from the corners.

When to choose it. Reproducing results from before the direct path existed. For new work prefer Trilinear.

Parameters it owns. None.

Interactions. Same swarm prerequisites as the direct path; the extra averaging stage smooths the field the particle sees.

Diagnostics. As above.

Evidence. Regression verified - make smoke runs a flat-channel particle case with CornerAveraged and asserts the runtime banner reports it. Analytically verified - interpolation-methods-tgv-2026-09-21: relative L2 error 7.99% and 2.49% on 16^3 and 32^3 against the TGV3D field, order 1.68, and maximum error order 0.96 - three to four times the Trilinear error at the same resolution.

Limitations. The additional averaging is diffusive, and it is retained for compatibility rather than accuracy. It is below second order even on a uniform Cartesian grid: the corner average smooths the field before the trilinear step sees it, and the maximum error converges at first order.

steady_deterministic

Identity. monitor.yml -> solution_monitoring.convergence.mode: steady_deterministic -> -solution_convergence_mode STEADY_DETERMINISTIC -> SOLUTION_CONVERGENCE_STEADY_DETERMINISTIC.

What it does. Logs, after every completed step, how much the solution changed since the previous step: absolute and relative L2 changes of velocity and of pressure (with its mean removed), and the domain-mean speed and kinetic energy with their drift. It judges nothing: no tolerance is read and the run is never stopped. Whether the flow has settled is read from <run.runtime_logs>/solution_convergence.log.

When to choose it. Flows that genuinely reach a steady state - laminar channels and ducts below their critical Reynolds number - where a change per step falling toward zero is the convergence signal. It is the default.

Parameters it owns. None; enabled: false turns the writer off.

Interactions. A turbulent case logs changes that never fall; that is the answer, not a malfunction. Physical cells only: ghost layers are excluded from every norm and mean.

Diagnostics. The log's ref column is 0 on the first logged step, which has no previous state to compare with.

Evidence. Analytically verified - solution-monitoring-tgv3d-2026-09-18: on the closed-form TGV3D field every logged drift and mean matched its exact value to 4e-11 in all four modes; the same measurement found and fixed means that counted ghost cells.

Limitations. Reports change per step, not error: a slowly evolving flow shows a small change per step long before it is steady, so read the trend over many steps.

periodic_deterministic

Identity. solution_convergence.mode: periodic_deterministic -> -solution_convergence_mode PERIODIC_DETERMINISTIC -> SOLUTION_CONVERGENCE_PERIODIC_DETERMINISTIC.

What it does. Compares each step with the state one period earlier at the same phase, and logs the same velocity, pressure, speed and energy changes as steady_deterministic together with the phase index. A cycle that repeats drives the logged change to zero.

When to choose it. Flows with a known deterministic period in steps - vortex shedding at low Reynolds number, or another periodically forced flow.

Parameters it owns. periodic_deterministic.period_steps, required and positive: the period as a whole number of steps.

Interactions. Stores one velocity and pressure snapshot per phase, so memory grows with period_steps. The first period has no reference and logs ref 0.

Diagnostics. The ph and per columns give the phase and period of each row.

Evidence. Analytically verified - solution-monitoring-tgv3d-2026-09-18: on the closed-form TGV3D field every logged drift and mean matched its exact value to 4e-11 in all four modes; the same measurement found and fixed means that counted ghost cells.

Limitations. The period must be an integer number of steps and known in advance; a period that is not, or that drifts, never compares like with like.

statistical_steady

Identity. solution_convergence.mode: statistical_steady -> -solution_convergence_mode STATISTICAL_STEADY -> SOLUTION_CONVERGENCE_STATISTICAL_STEADY.

What it does. Records the domain-mean speed and kinetic energy each step, and logs the mean and RMS of each over the latest window_steps samples against the same over the window_steps before them. A statistically steady flow drives those window-to-window drifts toward zero while the instantaneous field keeps changing.

When to choose it. Turbulence, where the instantaneous field never settles but its bulk statistics do - the driven periodic campaigns, for instance.

Parameters it owns. statistical_steady.window_steps, required and positive.

Interactions. Independent of Field Statistics: this mode watches two domain-mean scalars, while field statistics accumulate full fields over their own windows.

Diagnostics. Rows before 2 x window_steps samples exist log ref 0.

Evidence. Analytically verified - solution-monitoring-tgv3d-2026-09-18: on the closed-form TGV3D field every logged drift and mean matched its exact value to 4e-11 in all four modes; the same measurement found and fixed means that counted ghost cells.

Limitations. Watches two domain means only, so a flow whose mean speed and energy have settled can still have an unconverged profile. A window short against the flow's slowest time scale reads as converged when it is not.

transient

Identity. solution_convergence.mode: transient -> -solution_convergence_mode TRANSIENT.

What it does. Logs exactly what steady_deterministic logs; only the mode label differs, so the log states that no steady state is expected.

When to choose it. Decaying turbulence, a startup transient, or any run whose point is the evolution rather than an end state.

Parameters it owns. None.

Interactions. Like every mode, it never stops the run; the run ends on its step count.

Diagnostics. As steady_deterministic.

Evidence. Analytically verified - solution-monitoring-tgv3d-2026-09-18: on the closed-form TGV3D field every logged drift and mean matched its exact value to 4e-11 in all four modes; the same measurement found and fixed means that counted ghost cells.

Limitations. Adds no measure of its own; choose it to label the run, not to change what is computed.

5. tolerances

tolerances:
max_iterations: 50
relative_tol: 1.0e-4
residual_absolute_tol: 1.0e-8 # dimensionless (see below)
residual_relative_tol: 1.0e-3

Mappings:

  • max_iterations -> -mom_max_pseudo_steps
  • relative_tol -> -mom_rtol
  • residual_absolute_tol -> -mom_resid_atol
  • residual_relative_tol -> -mom_resid_rtol

absolute_tol -> -mom_atol is deprecated and no longer appears in the shipped configs. It is still accepted, but the CLI warns that it takes no part in convergence while a residual tolerance is active; it survives only for the legacy update-only branch.

Both residual tolerances now default to enabled (-mom_resid_atol 1e-8, -mom_resid_rtol 1e-3). Setting both non-positive is an explicit opt-out that selects the legacy update-only branch, which can converge falsely when dtau collapses; prefer not to.

When either residual tolerance is positive it decides convergence: the absolute test is sufficient on its own, and the relative test is paired with relative_tol as a guard. absolute_tol then takes no part, because |dU| <= absolute_tol is the disguised, step-size-dependent residual bound |R| <= absolute_tol/dtau. When both residual tolerances are non-positive, the solver preserves its existing update-only convergence criterion (absolute_tol and relative_tol). See 3. Convergence and Adaptive Rollback.

residual_absolute_tol is dimensionless. The residual carries units of volumetric flux per time, so a raw bound on |R| would need retuning for every grid, velocity scale and timestep; across the shipped cases |R| at step 1 spans 1.2e-2 (plane channels) to 2.0 (driven duct). The test applied is

|R| <= residual_absolute_tol * resid_ref,   resid_ref = a0 * |Ucont|_inf / dt

where resid_ref is the magnitude of the residual's own BDF term, recomputed each physical step. Normalising by it removes the dt and velocity-scale factors that dominate that spread, so one value is portable. resid_ref is reported on the per-step Dual-time solver: log line. step_tol/-imp_stol remains accepted as a deprecated compatibility option but is unused by active momentum solvers.

6. momentum_solver (Solver-Specific Block)

momentum_solver:
dual_time_picard_jameson_rk:
max_pseudo_steps: 50
relative_tol: 1.0e-5
pseudo_cfl: # dimensionless Courant number: dtau = pseudo_cfl / lambda_max
initial: 0.5
minimum: 0.001
maximum: 2.0
growth_factor: 1.1
reduction_factor: 0.75
jameson_residual_noise_allowance_factor: 1.1

Mappings include:

  • -pseudo_cfl, -min_pseudo_cfl, -max_pseudo_cfl
  • -pseudo_cfl_growth_factor, -pseudo_cfl_reduction_factor
  • -mom_dt_jameson_residual_norm_noise_allowance_factor

Rule: solver-specific blocks must match selected momentum solver type. Do not set canonical Jameson keys and their deprecated RK4 aliases together. Pseudo-CFL is exclusively a Dual Time Picard–Jameson RK control: it is neither accepted from the structured Newton–Krylov block nor shown in a Newton–Krylov startup banner.

For strategy.momentum_solver: "Newton Krylov", the structured PETSc controls are:

momentum_solver:
newton_krylov:
jacobian:
type: finite_difference
finite_difference:
mode: matrix_free
preconditioner:
model: none
nonlinear_solver:
method: newtonls
absolute_tolerance: 1.0e-10
relative_tolerance: 1.0e-8
step_tolerance: 1.0e-12
max_iterations: 12
line_search:
type: bt
eisenstat_walker:
enabled: true
version: 3
initial_relative_tolerance: 0.3
maximum_relative_tolerance: 0.9
gamma: 1.0
exponent: 1.618033988749895
safeguard_exponent: 1.618033988749895
safeguard_threshold: 0.1
linear_solver:
method: gmres
absolute_tolerance: 1.0e-10
relative_tolerance: 1.0e-6
max_iterations: 400
gmres:
restart: 80

Mappings are nonlinear_solver.method/absolute_tolerance/relative_tolerance/step_tolerance/max_iterations to -mom_nk_snes_type/-mom_nk_snes_atol/-mom_nk_snes_rtol/-mom_nk_snes_stol/-mom_nk_snes_max_it, line_search.type to -mom_nk_snes_linesearch_type, and the corresponding linear_solver fields to -mom_nk_ksp_type/-mom_nk_ksp_atol/-mom_nk_ksp_rtol/-mom_nk_ksp_max_it. gmres.restart maps to -mom_nk_ksp_gmres_restart. The Jacobian fields map to -mom_nk_jacobian_type/-mom_nk_jacobian_fd_mode; the preconditioner fields map to -mom_nk_preconditioner_model/-mom_nk_preconditioner_structure.

eisenstat_walker.enabled: false leaves the configured KSP relative tolerance fixed. enabled: true enables PETSc's inexact-Newton forcing term. Its seven fields map to -mom_nk_snes_ksp_ew_version, -mom_nk_snes_ksp_ew_rtol0, -mom_nk_snes_ksp_ew_rtolmax, -mom_nk_snes_ksp_ew_gamma, -mom_nk_snes_ksp_ew_alpha, -mom_nk_snes_ksp_ew_alpha2, and -mom_nk_snes_ksp_ew_threshold. PETSc supports EW versions 1 through 4; version 3 is therefore selectable directly. When EW is active, linear_solver.relative_tolerance is the initial KSP tolerance and PETSc may replace it before each Newton linear solve.

Newton tolerances are nonnegative, iteration/restart counts are positive integers, and GMRES restart is valid only for gmres, fgmres, or lgmres. Supported combinations are finite difference/matrix free with either no preconditioner or frozen momentum Jacobian/point block. The released linear_solver.preconditioner.type: none remains a deprecated alias for the no-preconditioner model. Raw petsc_passthrough_options are applied last, but an incompatible raw -mom_nk_pc_type override is rejected by the runtime. The Jacobian block is a strict discriminated configuration: an explicit type: finite_difference requires finite_difference.mode: matrix_free. Any other mode or Jacobian type, and irrelevant sibling configuration, is rejected. Advanced SNES, line-search, KSP, GMRES, and PC options that do not have structured PICurv fields remain available through petsc_passthrough_options with the -mom_nk_ prefix. See PETSc's SNES manual, SNES option list, and KSP option list. PETSc startup-only options such as -info belong under monitor.diagnostics.petsc.

7. poisson_solver

Note
**preconditioner.type accepts one value.** multigrid is the only outer preconditioner the Poisson solve supports, because the runtime assumes PETSc PCMG setup. mg and pcmg are accepted spellings for it. There is no choice to make here today, which is why it is documented as a parameter rather than as a capability family; tests/tooling/family_census.json records that classification, and the census fails if a second value is ever added without a family being registered for it.

The outer Krylov method under method is a closed set, fgmres or cg, entered at 7.1 Poisson Method Entries. The per-level methods under multigrid.level_solvers are PETSc KSP tokens passed through rather than a closed PICurv set - any token PETSc accepts is accepted there. PICurv recognises the Krylov subset only to warn when one is configured at level_0, which is the coarse solve rather than a smoother; see 6. Related Pages.

poisson_solver:
method: "fgmres"
absolute_tolerance: 1.0e-5
relative_tolerance: 1.0e-11
max_iterations: 50
gmres:
restart: 20
preconditioner:
type: "multigrid"
multigrid:
levels: 3
pre_sweeps: 2
post_sweeps: 2
cycle: "v"
mode: "multiplicative"
semi_coarsening:
i: false
j: false
k: true
level_solvers:
level_0:
method: "preonly"
preconditioner: "redundant"
level_1:
method: "richardson"
preconditioner: "bjacobi"
level_2:
method: "richardson"
preconditioner: "bjacobi"

Mappings:

  • method -> -ps_ksp_type; fgmres (default) or cg, which also emits -ps_ksp_norm_type unpreconditioned so it stops on the true residual. gmres, lgmres and bcgs are refused: they are not flexible, and the multigrid preconditioner is not a fixed linear operator (see 3. YAML Mapping and PETSc Options)
  • absolute_tolerance -> -ps_ksp_atol
  • relative_tolerance -> -ps_ksp_rtol
  • max_iterations -> -ps_ksp_max_it
  • gmres.restart -> -ps_ksp_gmres_restart; valid only with method: fgmres
  • preconditioner.type -> -ps_pc_type; currently only multigrid is supported
  • multigrid.levels -> -mg_level
  • multigrid.pre_sweeps -> -mg_pre_it
  • multigrid.post_sweeps -> -mg_post_it
  • multigrid.semi_coarsening.i/j/k -> -mg_i_semi/-mg_j_semi/-mg_k_semi
  • multigrid.level_solvers.level_N.method -> -ps_mg_levels_N_ksp_type for N > 0
  • multigrid.level_solvers.level_N.preconditioner -> -ps_mg_levels_N_pc_type for N > 0
  • multigrid.level_solvers.level_N.max_it/rtol/atol -> -ps_mg_levels_N_ksp_max_it/_ksp_rtol/_ksp_atol for N > 0
  • multigrid.level_solvers.level_0.* -> -ps_mg_coarse_*; PETSc names the coarsest solver separately from the positive levels
  • multigrid.cycle and multigrid.mode are validated structured keys; current supported values are v and multiplicative.

Rules:

  • pressure_solver is accepted as a legacy alias, but poisson_solver is preferred because the linear solve computes pressure correction Phi.
  • MG level numbering follows PETSc/PICurv convention: level_0 is the coarsest level and larger numbers are finer.
  • level_0 is the multigrid coarse solve, not a smoother, and the naming hides that. level_1..N are smoothers; level_0 sits at the base of the V-cycle and removes the smooth error the smoothers cannot see. Because multigrid is used here as a preconditioner, level_0 must be a fixed linear operator. A Krylov method there (gmres, fgmres, cg, bcgs, ...) makes the whole preconditioner nonlinear, which decouples the outer FGMRES tracked residual from the true residual b - Ax: the solver then reports convergence on a number that no longer describes the constraint. Pinning ksp_max_it does not fix it. PICurv logs a startup warning if you configure one. Use {method: preonly, preconditioner: redundant} for coarse grids up to roughly 1e4 unknowns, telescope above that. Full discussion, the worked tracked-vs-true residual table, and the level-count sizing rule: Pressure-Poisson, GMRES, and Multigrid.
  • multigrid.levels is bounded by the MPI decomposition, not chosen freely: every level must leave each rank at least stencil_width nodes per axis (3 when any axis is periodic, 2 otherwise). Exceeding it aborts during DM creation with Local x-width of domain ... is smaller than stencil width. picurv run, including --dry-run, refuses such a layout before launch when grid.da_processors_x/y/z is set; a layout left to PETSc is only known at runtime. Formula and worked maxima: Pressure-Poisson, GMRES, and Multigrid.
  • The outer Poisson preconditioner is multigrid-only in the current runtime.
  • pre_sweeps and post_sweeps are applied separately: smoothers run pre_sweeps iterations before the coarse-grid correction and post_sweeps after it. When they differ, the post-smoother copies the configured level solver and reads extra options under -ps_mg_levels_N_up_; see Pressure-Poisson, GMRES, and Multigrid.
  • Advanced PETSc tuning remains available through petsc_passthrough_options; common examples include -ps_mg_levels_N_pc_sor_omega for SOR and -ps_mg_levels_N_pc_factor_shift_amount / -ps_mg_levels_N_pc_factor_levels for factor PCs.

7.1 Poisson Method Entries

Value Maps to
cgcg
fgmresfgmres

fgmres

Identity. poisson_solver.method: fgmres (the default) -> -ps_ksp_type fgmres -> the outer KSP of PoissonSolver_Multigrid.

What it does. Flexible GMRES around the multigrid preconditioner. Flexibility lets the preconditioner change from one iteration to the next, which is what a multigrid cycle with iterative smoothers is, so the residual it tracks is the true residual.

When to choose it. Always, unless measured cost says otherwise. It is the verified default and the only method that also accepts gmres.restart.

Parameters it owns. gmres.restart, the Krylov restart length.

Interactions. Uses the absolute, relative and iteration limits of the block like any method. Stopping at max_iterations is its normal mode and stays silent; a non-finite residual or a failed preconditioner stops the run.

Diagnostics. solver_monitoring.poisson.pic_true_residual logs tracked and true residuals side by side in Poisson_Solver_Convergence_History_Block_*.log.

Evidence. Unit verified - make unit-poisson-rhs solves and projects on a three-level hierarchy. Integration verified - make smoke-driven-periodic asserts the tracked and true residuals agree within 1e-4. Production exercised - every shipped flow case, examples/flat_channel among them. Analytically verified - duct-poiseuille-picard-2026-09-18: the analytic duct pressure gradient at second order; poisson-options-2026-09-21: every multigrid option reproduced this baseline; poisson-persistent-solver-2026-10-10: the solver kept for the run reproduced the previous one and the duct and bend evidence; poisson-option-matrix-2026-10-10: every Poisson setting, varied one at a time around this method, reproduced its solve; hom02-poisson-rewrite-2026-10-10: the HOM02 decaying-turbulence comparison with Wray's DNS reproduced to the printed digits.

Limitations. Stores one vector per iteration up to the restart length, so memory grows with gmres.restart.

cg

Identity. poisson_solver.method: cg -> -ps_ksp_type cg and -ps_ksp_norm_type unpreconditioned.

What it does. Conjugate gradients around the multigrid preconditioner, monitoring the unpreconditioned residual so that it stops on the true residual rather than on the preconditioned one.

When to choose it. When the memory of a Krylov basis matters: CG keeps a fixed handful of vectors whatever the iteration count. CG assumes a symmetric operator; the discrete pressure operator's symmetry on a non-orthogonal grid was not checked, only CG's measured convergence there.

Parameters it owns. None; gmres.restart is refused with it.

Interactions. As fgmres. CG is not flexible in theory; with this preconditioner it held the true residual at 1e-12 or below every step in the measurement below, where gmres, lgmres and bcgs did not.

Diagnostics. As fgmres; watch the true-residual column, which is the norm it stops on.

Evidence. Unit verified - tests/test_config_regressions.py checks that it emits the unpreconditioned norm type. Analytically verified - poisson-options-2026-09-21: on the duct it reproduced the fgmres velocity to 7e-15, and on the curved bent-channel grid it matched fgmres to 4.9e-15 in velocity and 4.9e-14 in pressure with the true residual at 1e-12 or below every step. poisson-option-matrix-2026-10-10: with the rewritten solver it matched the fgmres baseline to 8.7e-10 on the duct, 1.1e-13 on the bent channel and 3.0e-15 on a triply periodic box.

Limitations. Measured on one duct and one curved grid, both small; its cost relative to fgmres was not characterized. On a grid or multigrid configuration the measurement did not cover, confirm with pic_true_residual that the true residual still falls.

8. Physical-solution convergence monitoring

Convergence monitoring is observation policy rather than a numerical solver selection, so it is configured in monitor.yml -> solution_monitoring.convergence rather than here. See Configuration Reference: Monitor YAML for modes, mappings, and defaults. The monitor records every completed timestep; cadence is not a user setting.

solver.yml accepts no solution_convergence key. A file carrying one is rejected by validation with the location to move it to.

Scientific field statistics are likewise a monitor concern, configured at monitor.yml -> field_statistics; see Field Statistics.

9. interpolation

interpolation:
method: "Trilinear" # default; or "CornerAveraged"

Mappings:

  • method -> -interpolation_method (Trilinear = 0, CornerAveraged = 1)

The Trilinear method (default) performs direct trilinear interpolation from the 8 nearest cell centers; on the TGV3D field it converges at order 1.97. The CornerAveraged method is the legacy two-stage path (center-to-corner average, then trilinear from corners); on the same field and a uniform Cartesian grid its L2 error converges at order 1.68 and its maximum error at first order, so it is not second order anywhere measured. See p08_cap_interp_corneraveraged.

See Trilinear Interpolation and Particle-Grid Projection for algorithmic details.

10. scalar_transport

scalar_transport:
schmidt_number: 1.0
turbulent_schmidt_number: 0.7
iem_constant: 2.0

Mappings:

  • schmidt_number -> -schmidt_number
  • turbulent_schmidt_number -> -turb_schmidt_number
  • iem_constant -> -iem_constant, the constant C_IEM in the IEM mixing rate Omega = C_IEM Gamma_eff / Delta^2 (1. IEM Mixing Update In Current Code)

Rules:

  • the Schmidt numbers must be positive; iem_constant must be non-negative, and both configuration validation and the runtime require it to be finite
  • iem_constant: 0 switches micromixing off: the particle update is skipped, so each particle's Psi is carried unchanged and bit-identical - a passive label
  • omitted values use the C runtime defaults: schmidt_number = 1.0, turbulent_schmidt_number = 0.7, and iem_constant = 2.0
  • iem_constant changes nothing while every particle's Psi is equal within each cell; models.physics.particles.fields gives Psi its initial values (1.1 Configured Particle Values (fields)), and the verification scalar source, which also sets Psi, bypasses IEM; see 1. IEM Mixing Update In Current Code
  • use this structured block for ordinary scalar/Brownian transport tuning; reserve petsc_passthrough_options for flags without a YAML schema

11. verification

verification:
sources:
diffusivity:
mode: "analytical"
profile: "LINEAR_X"
gamma0: 1.0e-3
slope_x: 2.0e-4
scalar:
mode: "analytical"
profile: "CONSTANT" # or LINEAR_X, SIN_PRODUCT
value: 1.0
# phi0: 0.0
# slope_x: 1.0
# amplitude: 1.0
# kx: 3.141592653589793
# ky: 3.141592653589793
# kz: 3.141592653589793

Mappings:

  • verification.sources.diffusivity.mode -> -verification_diffusivity_mode
  • verification.sources.diffusivity.profile -> -verification_diffusivity_profile
  • verification.sources.diffusivity.gamma0 -> -verification_diffusivity_gamma0
  • verification.sources.diffusivity.slope_x -> -verification_diffusivity_slope_x
  • verification.sources.scalar.mode -> -verification_scalar_mode
  • verification.sources.scalar.profile -> -verification_scalar_profile
  • verification.sources.scalar.value -> -verification_scalar_value
  • verification.sources.scalar.phi0 -> -verification_scalar_phi0
  • verification.sources.scalar.slope_x -> -verification_scalar_slope_x
  • verification.sources.scalar.amplitude -> -verification_scalar_amplitude
  • verification.sources.scalar.kx/ky/kz -> -verification_scalar_kx/-verification_scalar_ky/-verification_scalar_kz

Units: every value is physical and is converted with the case's scales. gamma0 is a diffusivity (divided by velocity_ref length_ref) and the diffusivity slope_x a diffusivity per length (divided by velocity_ref); the scalar slope_x and kx/ky/kz are per length (multiplied by length_ref); value, phi0, and amplitude are dimensionless. See 3. Input Index.

Rules:

  • this path is verification-only and should be used only when no ordinary end-to-end workflow can expose the behavior under test
  • it is only valid with operation_mode.eulerian_field_source: "analytical"
  • verification.sources.scalar prescribes particle Psi from analytical truth and enables the runtime diagnostic <run.analysis.metrics>/scatter_metrics.csv
  • scalar profiles currently supported are CONSTANT, LINEAR_X, and SIN_PRODUCT
  • new verification source overrides must be implemented in include/verification_sources.h and src/verification_sources.c

Evidence:

  • the diffusivity source, LINEAR_X, drives the gradient drift that diffusivity-gradient-drift-paired-2026-09-18 resolved to 1.1% of a t with paired seeds;
  • the scalar source's CONSTANT profile makes every scatter error term exactly zero, and SIN_PRODUCT gives relative L2 errors of 5.2%, 2.6% and 1.3% on 8, 16 and 32 cells per side, the first order a random cell average must show at fixed particles per cell (scalar-scatter-verification-2026-09-21, examples/scatter_verification);
  • the scalar LINEAR_X profile has not been measured, and neither scatter threshold is gated in CI.

12. petsc_passthrough_options

Advanced escape hatch for raw PETSc flags:

petsc_passthrough_options:
-ps_ksp_type: "gmres"
-ps_pc_type: "ilu"

These are passed into PETSc options DB and consumed by runtime calls like KSPSetFromOptions.

13. Next Steps

Proceed to Configuration Reference: Monitor YAML.

For mapping and extension workflows: