Solvers

The solver section configures the nonlinear solve, the linear solve nested inside it, and the convergence criteria. It is required for quasi static and newmark, and ignored entirely by central difference, which is explicit.

solver:
  type: newton
  linear solver:
    type: cg
    preconditioner:
      type: jacobi
  termination:
    converge when any:
      - absolute residual: 1.0e-08
      - relative residual: 1.0e-12
    fail when any:
      - maximum iterations: 16

Nonlinear solvers

type is required; there is no default. hessian minimizer is accepted for Norma compatibility and builds a plain NewtonSolver.

Anything else aborts with the list above. In particular lbfgs is not a nonlinear solver type — it is a linear solver, selected with solver.linear solver.type: lbfgs. That mistake used to fall through to Newton silently, producing a run that looked fine but was not the algorithm requested:

ERROR: Unknown solver.type = "lbfgs". Supported: "newton" (aliases ...).
       L-BFGS is a linear solver: set solver.linear solver.type = "lbfgs".

For newton, a linear solver sub-section is required. For the two matrix-free solvers it is optional and defaults to type: none.

Keys common to all nonlinear solvers

Nonlinear CG only

Steepest descent accepts preconditioner with the same meaning. Jacobi is the only preconditioner implemented at this level; asking for anything else (chebyshev, ic, amg) is a hard error rather than silently running Jacobi.

Termination criteria

Carina accepts three syntaxes. Prefer the first.

Preferred: converge when / fail when

  termination:
    converge when any:
      - absolute residual: 1.0e-08
      - relative residual: 1.0e-12
    fail when any:
      - maximum iterations: 16
      - divergence: 1.0e6

Four block keys are recognized, each taking a list:

Each list item is a single-key mapping of test name: value. Groups nest via any: / all: inside a list:

  termination:
    converge when all:
      - minimum iterations: 2
      - any:
          - absolute residual: 1.0e-08
          - relative residual: 1.0e-12

Legacy: typed list

Still supported. Each entry carries an explicit type, with the value under tolerance, value, threshold, or window depending on the test, and combo/tests for nesting:

  termination:
    - type: combo
      combo: and
      tests:
        - type: absolute residual
          tolerance: 1.0e-06
        - type: relative update
          tolerance: 1.0e-12
    - type: maximum iterations
      value: 16

combo must be and or or; the default when it is omitted is or. Any other value is an error — it previously fell through to or, which inverted the meaning of the group without saying so.

Oldest: flat tolerances

If termination is absent entirely, Carina builds OR(absolute residual, relative residual, finite value) from the flat absolute tolerance (default 1e-10) and relative tolerance (default 1e-14) keys.

Available tests

Points worth knowing:

• minimum iterations is a convergence test, not a floor. It reports Converged once the iteration count is reached. Placed in a converge when any block it will force convergence on that iteration. It is only meaningful inside an all group, where it prevents premature convergence.
• Relative tests normalize differently. relative residual divides by the initial residual ‖R₀‖; relative update divides by the current solution norm ‖U‖. Each is inert while its denominator is zero — a relative update test cannot converge on the first iteration from a zero initial guess.
• divergence compares against ‖R₀‖, not the previous iterate, so it detects net growth from the start rather than a single bad step.
• stagnation uses its value twice: as the lookback distance and as the number of consecutive stagnant iterations required to fail. With the default window of 5 it needs at least 11 iterations before it can trigger. Its internal ratio floor is 0.95 — i.e. failure when ‖Rk‖/‖R{k−window}‖ > 0.95 — settable only through the legacy syntax's tolerance key.
• A finite value test is appended automatically to every termination tree, in all three syntaxes. You do not need to add one.
• maximum iterations in the tree takes precedence. If the tree contains a maximum iterations test, its value replaces the flat maximum iterations key as the solver loop bound.
• In an OR group, Converged takes precedence over Failed. If one sub-test converges while another fails on the same iteration, the result is Converged. An AND group short-circuits on the first Failed.

Linear solvers

Configured under solver.linear solver.

Selecting direct on a GPU backend is a hard error:

solver.linear_solver.type = "direct" is CPU-only.

Iterative keys

Inexact Newton forcing terms

By default every Newton iteration solves its linear system to tolerance. That is more accuracy than the early iterations can use: the Newton step they produce is about to be discarded by the next correction. On a 530k-DOF Newmark run, CG spent 181, 205 and 193 iterations on the three solves of a step while ‖R‖ fell 1.6e4 → 1.3e2 → 7.8e-2 → 4.6e-8. Only the last one needed eight digits.

A forcing term replaces the fixed tolerance with a per-iteration ηₖ, chosen from how fast the nonlinear residual is actually falling:

  linear solver:
    type: cg
    tolerance: 1.0e-8
    preconditioner:
      type: jacobi
    forcing term:
      type: eisenstat-walker

Eisenstat–Walker keys

Every bound is enforced. Out-of-range values and unknown types are hard errors rather than silent fallbacks — a typo here would otherwise cost exactly the speedup the block was added to get, with nothing in the log to say the setting never took effect.

How ηₖ is chosen

The base rule is ηₖ = γ (‖Rₖ‖ / ‖Rₖ₋₁‖)^α, guarded three ways:

1. Eisenstat–Walker's own safeguard stops η falling faster than the observed convergence rate justifies: η ← max(η, γ·η_prev^α) whenever that term exceeds 0.1.
2. Kelley's over-solve guard stops the final solve being worked past the point the nonlinear test will reward: η ← max(η, safety·τ/‖Rₖ‖), where τ is the residual norm your termination block is actually aiming at. Carina reads τ from the termination tree itself, so a deck's own converge when thresholds are honored rather than the flat tolerance keys.
3. A clamp to [tolerance, maximum].

When a solve is loosened, the log says so:

        [SOLVE]       inexact Newton: η = 2.00e-01
        [SOLVE]       CG: 16 iters : |r|_CG = 1.07e-01 : [CONV]

What to expect

A forcing term trades linear work for nonlinear work: fewer CG iterations, more Newton iterations. It therefore pays in proportion to how much of your step is the linear solve — most on the GPU matrix-free paths, where the step essentially is the stiffness matvec; least on assembled CPU runs where residual assembly and the line search already dominate.

Measured with the defaults (see benchmark/evidence/inexact_newton.txt):

The quasi-static AMG gain is on top of the ~2× AMG already has over Jacobi.

References

• S. C. Eisenstat and H. F. Walker, Choosing the forcing terms in an inexact Newton method, SIAM Journal on Scientific Computing 17(1):16–32, 1996. The ηₖ rule and safeguard 1.
• C. T. Kelley, Iterative Methods for Linear and Nonlinear Equations, SIAM, 1995, §6.3. The over-solve guard, safeguard 2.

L-BFGS keys

The L-BFGS path always builds a Jacobi preconditioner regardless of any preconditioner sub-section.

Preconditioners

Configured under solver.linear solver.preconditioner.

Chebyshev keys

AMG

  linear solver:
    type: cg
    tolerance: 1.0e-8
    maximum iterations: 500
    preconditioner:
      type: amg

AMG is the only preconditioner whose CG iteration count is nearly independent of the conditioning that defeats Jacobi — on a 530k-DOF torsion problem it holds 5–17 CG iterations across a range of Δt where Jacobi grows from 30 to

It runs on both CPU and GPU. On GPU the hierarchy is built on the host (the sparse pattern lives there anyway) and converted to device CSR; the V(2,2)-cycle then applies entirely on the device through KernelAbstractions kernels, with the fine level smoothed through the matrix-free stiffness action so the fine matrix is never formed on the device. It is the fastest quasi-static option Carina has — see benchmark_report.md.

Two behaviors worth knowing. The near-nullspace (six rigid-body modes) is rebuilt from the current nodal coordinates X + u at every hierarchy build, not frozen at the reference configuration — a frozen reference nullspace degrades badly once the body rotates substantially. And the hierarchy setup is expensive (seconds at 500k DOF), so it is lagged: built once, then rebuilt only when the effective-mass coefficient changes (a Δt change) or when CG iteration growth flags it as stale.

AMG combines well with an Eisenstat–Walker forcing term: the two are independent, and on a 530k-DOF quasi-static torsion run on an A100 the forcing term takes a further 1.59× off the per-step wall on top of what AMG already gives over Jacobi.

AMG targets SPD tangents — quasi-static and moderate-Δt dynamics. At very large Δt on violently dynamic problems the Newmark predictor can overshoot into near-inverted configurations whose tangent is indefinite, which breaks CG regardless of preconditioner; explicit integration is the right tool there.

An unrecognized type is a hard error listing the supported values:

Unknown preconditioner.type = "jacobbi". Supported: "jacobi", "ic" (aliases
"incomplete cholesky", "ildl", "incomplete ldlt"), "chebyshev", "amg"
(aliases "algebraic multigrid", "multigrid"), "none".

Omitting the preconditioner section entirely, or writing type: none, means no preconditioning — that remains a valid, unremarkable choice.

Choosing a combination

Integrator and nonlinear solver

Newton, by device

* L-BFGS is the fastest option for implicit dynamics on GPU, but fails on quasi-static on either device — see the warning above.

Practical guidance:

• CPU, small to medium — direct. No tuning, always converges.
• CPU, large — cg + amg, or cg + ic as a simpler alternative.
• GPU, quasi-static — cg + amg. It is the only GPU option that both converges its linear systems and outperforms the CPU direct solver at scale, and at 1.57M DOF it is the fastest option on either device.
• GPU, implicit dynamics at small Δt — cg + jacobi. The mass term conditions the system, so AMG's iteration reduction does not repay its per-application cost; see benchmark_report.md.
• Never run plain cg with no preconditioner on a real mesh. It is valid and very slow.
• Add a forcing term to any Newton + CG combination, on either device. forcing term: {type: eisenstat-walker} is the least expensive speedup available here: it needs no extra memory, changes no kernel, and cannot make the converged answer less accurate than your tolerance already asked for. It helps most where CG dominates the step — the GPU matrix-free paths — and least where residual assembly and the line search already do.

Combinations that fail or mislead