Skip to content

Simulation

The simulation block defines how the reaction network is solved.

Basic Form

simulation {
  base @ transient_solve(
    model=thermal,
    solver={method=cvode_bdf, atol=1e-10, rtol=1e-10, t_end=1e6},
    conditions=[{T=600}, {T=650}]
  )
}

Simulation Entries

Each entry has:

  • a name
  • a simulation model
  • an argument list

Simulation Models

  • transient_solve
  • steady_state_solve

What Each Simulation Model Does

transient_solve

transient_solve integrates the kinetic model forward in time.

In this mode, MKMCXX computes the time evolution of the system starting from the initial values defined in the components block. The result is a transient trajectory that shows how the system approaches its later-time behavior.

This model is useful when:

  • the time evolution itself is physically important
  • the system may take a long time to approach steady state
  • a robust route to a stationary solution is needed
  • sensitivity analyses such as orders, DRC, TDRC, Eapp, or selectivity are to be computed

The run stops when either:

  • the final simulation time t_end is reached, or
  • the optional steady_state_threshold criterion is satisfied

In practical terms, transient_solve answers the question:

"How does the system evolve over time under these conditions?"

Example

simulation {
  baseline @ transient_solve(
    model=thermal,
    solver={
      method=cvode_bdf,
      atol=1e-10,
      rtol=1e-10,
      t_end=1e6,
      steady_state_threshold=1e-8
    },
    analyses=[
      {type=reaction_orders, enabled=true},
      {type=drc, enabled=true}
    ],
    conditions=[
      {T=500},
      {T=550, t_end=5e6},
      {T=600, atol=1e-12, rtol=1e-12}
    ]
  )
}

steady_state_solve

steady_state_solve solves the stationary nonlinear problem directly.

In this mode, MKMCXX does not integrate the system in time. Instead, it seeks a solution for which all time derivatives vanish, that is, a state satisfying:

\[ \frac{dy}{dt} = 0 \]

This model is useful when:

  • only the stationary solution is needed
  • transient behavior is not of interest
  • a direct steady-state solve is faster than long time integration
  • a good initial guess is already available

In practical terms, steady_state_solve answers the question:

"What stationary state satisfies the model under these conditions?"

Fallback Mechanism

In the current implementation, steady_state_solve does not rely on a single nonlinear solve attempt.

MKMCXX first tries a direct KINSOL solve using the configured strategy. If that fails in a recoverable way, it tries the alternate globalization strategy. If the direct nonlinear solve still does not converge, MKMCXX falls back to a long transient preconditioning run and then retries the steady-state solve from that warmed-up state.

Conceptually, the ladder is:

  1. direct solve with the configured strategy
  2. direct solve with the alternate strategy
  3. transient warm-start followed by the configured strategy
  4. transient warm-start followed by the alternate strategy

This means that a reported steady-state solution may still be valid even if the purely direct algebraic route did not converge on its own.

Path In The Convergence Report

The convergence report contains a Path column for steady_state_solve. This column tells you how the final steady-state solution was obtained.

  • direct The solve converged with the configured nonlinear strategy and no fallback was needed.
  • alternate The configured strategy did not converge, but the alternate globalization strategy did.
  • warm_start The direct solve did not converge cleanly, so MKMCXX first performed a transient warm-up and then converged with the configured strategy.
  • warm_start+alt The direct solve and the configured warm-start retry were insufficient, but the warm-started solve converged with the alternate strategy.

For transient_solve, the Path entry is reported as - because no steady-state fallback ladder is involved.

In practice:

  • direct is the cleanest outcome for a steady-state solve
  • alternate means the problem was still solved directly, but needed a more robust globalization strategy
  • warm_start and warm_start+alt indicate that the stationary problem was numerically more difficult and needed transient preconditioning

When interpreting results, the Path column is therefore a useful diagnostic. It tells you not only that a steady state was found, but also how difficult the nonlinear solve was.

Example

simulation {
  baseline @ steady_state_solve(
    model=thermal,
    solver={
      method=kinsol_newton,
      fnormtol=1e-10,
      scsteptol=1e-12,
      max_iter=50,
      strategy=line_search
    },
    conditions=[
      {T=500},
      {T=550, max_iter=75},
      {T=600, fnormtol=1e-12, scsteptol=1e-14}
    ]
  )
}

If a run finishes with output such as:

Job ID    SimID       End t         Stop          Path              Threshold     Max|dy/dt|
------------------------------------------------------------------------------------------------
1         T=550       n/a           steady_state  direct            1.000e-10     8.252e-11
2         T=600       n/a           steady_state  warm_start        1.000e-10     6.818e-11

then both conditions reached a stationary solution, but the 600 K case needed the transient fallback path whereas the 550 K case converged directly.

Choosing Between Them

Use transient_solve when you want the physical time-dependent behavior or a robust path toward a stationary state.

Use steady_state_solve when you want the final stationary solution directly.

For difficult systems, a common workflow is:

  1. run transient_solve
  2. use the final state as the starting point for steady_state_solve

Runtime Model

The current production workflow documented here uses:

  • model=thermal

model

  • model Selects the runtime mode for the simulation entry.

At present, the production workflow documented in this manual uses model=thermal.

vib_cutoff

  • vib_cutoff Optional vibrational cutoff in cm^-1 used by ArrheniusFreq.

When provided, MKMCXX replaces every vibrational frequency \tilde{\nu} used by ArrheniusFreq with

\[ \tilde{\nu}^\mathrm{eff} = \max(\tilde{\nu}, \tilde{\nu}_\mathrm{min}) \]

where \tilde{\nu}_\mathrm{min} is the configured vib_cutoff. The same effective frequency is used in both the vibrational partition-function term and the ZPE correction.

This parameter is useful when low-frequency modes would otherwise produce an unphysically large harmonic vibrational contribution.

Solver Configuration

transient_solve

The solver object for transient_solve controls the time-integration settings.

  • method
  • atol
  • rtol
  • t_end
  • steady_state_threshold
  • jacobian
  • linear_algebra

Variable Meanings

  • method Selects the transient solver method. The current supported method is cvode_bdf.
  • atol Absolute tolerance for the ODE solver. This controls the allowed absolute numerical error.
  • rtol Relative tolerance for the ODE solver. This controls the allowed error relative to the magnitude of the solution.
  • t_end Final simulation time for the transient integration.
  • steady_state_threshold Optional early-stop threshold. If the maximum absolute derivative in the system falls below this value, the run may stop before reaching t_end.
  • jacobian Controls how the transient solver obtains Jacobian information. Supported values are none and analytic. The default is none.
  • linear_algebra Controls the matrix/storage backend used by the transient solver. Supported values are dense and sparse. The default is dense. Sparse mode requires MKMCXX to be linked against a KLU-enabled Sundials installation.
  • max_order Optional upper bound on the multistep order used by the transient solver.
  • max_num_steps Optional upper bound on the number of internal solver steps.
  • initial_step Optional initial time step suggested to the solver.
  • min_step Optional lower bound for the time step.
  • max_step Optional upper bound for the time step.

steady_state_solve

The solver object for steady_state_solve controls the nonlinear stationary solve.

  • method
  • fnormtol
  • scsteptol
  • max_iter
  • strategy
  • jacobian
  • linear_algebra

Variable Meanings

  • method Selects the steady-state solver method. The current supported method is kinsol_newton.
  • fnormtol Tolerance on the residual norm of the nonlinear system.
  • scsteptol Tolerance used for the nonlinear step-size convergence criterion.
  • max_iter Maximum number of nonlinear iterations.
  • strategy Nonlinear solution strategy. Supported values are none and line_search.
  • jacobian Controls how the steady-state solver obtains Jacobian information. Supported values are none and analytic. The default is none.
  • linear_algebra Controls the matrix/storage backend used by the steady-state solver. Supported values are dense and sparse. The default is dense. Sparse mode requires MKMCXX to be linked against a KLU-enabled Sundials installation.

Conditions

The conditions list defines per-run settings.

Each condition must define at least:

  • T

Conditions may also override solver settings for individual runs.

Condition Variables

  • T Temperature of the run condition.

For transient_solve, a condition may also override:

  • t_end Final integration time for that specific run.
  • atol Absolute tolerance for that specific run.
  • rtol Relative tolerance for that specific run.
  • max_order Maximum solver order for that specific run.
  • max_num_steps Maximum number of internal steps for that specific run.
  • initial_step Initial time step for that specific run.
  • min_step Minimum time step for that specific run.
  • max_step Maximum time step for that specific run.
  • steady_state_threshold Early-stop steady-state threshold for that specific run.
  • jacobian Jacobian mode override for that specific transient run.
  • linear_algebra Linear algebra backend override for that specific transient run.

For steady_state_solve, a condition may also override:

  • fnormtol Residual-norm tolerance for that specific run.
  • scsteptol Step-size convergence tolerance for that specific run.
  • max_iter Maximum number of nonlinear iterations for that specific run.
  • jacobian Jacobian mode override for that specific steady-state run.
  • linear_algebra Linear algebra backend override for that specific steady-state run.

Other Simulation Keys

analyses

  • analyses Defines the list of additional analyses to compute together with the base simulation.

This key is supported for transient_solve. The current steady-state workflow does not support analyses.