Why Multiscale States Need Deliberate Tolerances
A mechanistic model can include states with very different typical magnitudes: a drug concentration might be tens of nanomoles per litre while a signaling intermediate is near \(10^{-6}\) nanomoles per litre. An adaptive ODE solver uses error estimates to decide whether a proposed step is accurate enough. Its decisions depend on the tolerances supplied by the modeler. If those tolerances do not reflect each state’s scale and scientific relevance, a solver may accept consequential error in a small state or spend excessive effort resolving changes that do not matter.
A tolerance is a numerical error-control setting, not a statement about measurement precision, biological variability, or parameter uncertainty. It also does not guarantee that the error in the final solution is bounded by the tolerance. It sets the scale against which a solver judges its estimated local integration error.
How the Error Scale Works
Suppose an adaptive solver estimates the local error for component \(i\) on a proposed step as \(e_i\). A common normalized error ratio is
Here \(y_{i,n}\) and \(y_{i,n+1}\) are the component values at the start and proposed end of the step. The solver combines the ratios across the state vector using its chosen norm, then typically accepts the step if the combined ratio is no greater than one. Details differ among solvers: some use a root-mean-square norm, and some use a maximum norm. Thus, a component ratio above one does not necessarily mean rejection when a solver uses an aggregate norm.
When \(rtol\,|y_i|\) is much larger than \(atol_i\), the error scale is approximately relative: the allowed estimated error grows with the state magnitude. When \(y_i\) is near zero, the absolute term dominates and prevents the error scale from shrinking to zero. This makes \(atol_i\) particularly important for states that approach zero or are naturally small.
The units must be consistent. Since \(rtol\) is dimensionless, both \(atol_i\) and \(rtol\,|y_i|\) have the units of \(y_i\). If states use different units, a single absolute tolerance generally cannot have the right units for every component; use component-specific tolerances or appropriately scaled state variables.
Choosing Tolerances From State Scales
Begin by listing each state’s units, typical magnitude over the simulation, and the smallest value or change that could affect a scientific conclusion. Use that information to choose \(atol_i\) and \(rtol\). Relative tolerance sets a broadly shared fractional target; absolute tolerance sets a component-specific floor in the state’s units. Neither should be chosen solely from the initial state if the trajectory later enters a very different range.
A useful first check is to compare \(atol_i\) with the smallest scientifically meaningful value of \(y_i\). If \(atol_i\) is larger than that value, the solver may accept estimated local errors of a size that obscures it, especially near zero. If \(atol_i\) is much smaller than any relevant scale, it may force unnecessary work. These comparisons guide choices but do not replace checking how predictions change when tolerances are tightened.
Worked Example: Setting Scales for Two Concentration States
Consider two hypothetical concentration states, both measured in nM. During the simulation, \(y_1\) is about \(800\) nM and \(y_2\) is about \(0.003\) nM. Choose \(rtol=10^{-4}\), \(atol_1=0.01\) nM, and \(atol_2=10^{-8}\) nM. Estimate each component’s error scale at those representative values.
For the larger state, the relative contribution is \(10^{-4}\times800=0.08\) nM. Therefore,
For the smaller state, the relative contribution is \(10^{-4}\times0.003=3\times10^{-7}\) nM. Therefore,
These scales show how the same relative tolerance gives a larger absolute error scale for the larger concentration, while the component-specific absolute tolerance helps set the scale for the small state. They are targets for judging estimated local error, not guarantees that the numerical solution differs from the exact solution by less than these amounts.
Why Absolute Tolerance Matters Near Zero
For a state close to zero, relative error alone is not a useful control target: the fractional error is undefined at zero and can become very large for a tiny absolute discrepancy. The absolute tolerance supplies a finite scale there. Choose it in the state’s units and in relation to the smallest magnitude or change that matters for the model’s intended use.
Worked Example: A State Approaching Zero
Suppose a hypothetical mediator concentration is \(2\times10^{-5}\) nM during part of a simulation. Its values near zero may still matter down to \(10^{-7}\) nM. Set \(rtol=10^{-3}\) and consider \(atol=10^{-8}\) nM. At \(2\times10^{-5}\) nM, the error scale is
At \(y=0\), the scale is \(10^{-8}\) nM. This absolute floor is smaller than the stated \(10^{-7}\) nM scale of scientific interest, so it is a plausible starting choice. By contrast, \(atol=10^{-4}\) nM would make the near-zero error scale about one thousand times larger than \(10^{-7}\) nM. Such a choice could fail to resolve changes that matter. This comparison identifies a risk; it does not by itself establish the actual error in the computed trajectory.
Deriving the Role of Scaling
The normalized ratio follows directly from comparing an estimated error with a permitted error scale. If the target scale for component \(i\) is \(s_i=atol_i+rtol\,|y_i|\), a solver asks whether its error estimate is small relative to that scale:
A ratio of one means that the estimated error equals the scale; a ratio below one means it is smaller. Substituting \(s_i\) gives the componentwise expression used above. This also explains the failure of a large shared absolute tolerance for states of different magnitudes: it makes \(s_i\) large even for a small state, reducing \(r_i\) and making that state’s estimated error less influential in the solver’s decision.
For example, if \(y_i=0.003\) nM and \(atol=0.1\) nM with \(rtol=10^{-4}\), the scale is \(0.1+0.0000003=0.1000003\) nM. An estimated error of \(0.001\) nM gives
The small ratio may make this component appear unimportant to error control, even though the error is one third of the state’s current magnitude. A more suitable component-specific \(atol_i\) can prevent that mismatch.
Checking Whether the Choice Is Working
Tolerance selection is an iterative modeling decision. After choosing plausible values, rerun the simulation with tighter tolerances and compare quantities relevant to the model’s use, such as state trajectories, peak values, or timing of a response. Inspect the small states as well as the large ones; a visually smooth trajectory is not evidence that error control is appropriate for every component.
Tightening tolerances can increase computational cost. In a smooth region where a method’s local error estimate behaves approximately as \(C h^{p+1}\), with step size \(h\), constant \(C\), and method order \(p\), reducing the target error scale by a factor \(q\) can require a step-size factor near \(q^{1/(p+1)}\). The number of steps may then rise substantially. This is an approximation, not a universal cost law: stiffness, changing dynamics, error norms, and solver implementation all affect the result.
Worked Example: Deciding Whether Tighter Tolerances Are Useful
In a hypothetical test, a model is run with tolerances \((rtol, atol)=(10^{-4}, 10^{-8})\), then with both values reduced by a factor of 100. The predicted peak concentration changes from \(4.250\) nM to \(4.251\) nM, while the scientific comparison requires resolving differences of at least \(0.01\) nM. The tighter run takes 2.4 times as long.
The change in the peak is \(4.251-4.250=0.001\) nM. Relative to the comparison scale, this is \(0.001/0.01=0.1\), or one tenth of the smallest difference the analysis needs to resolve. The result suggests that the original setting may already be adequate for this particular peak-based use, while the tighter setting costs more. It does not establish adequacy for other outputs, states, or time points. Check those quantities too, and investigate any unstable differences rather than relying on this one comparison.
Alternatives, Diagnostics, and Limitations
Component-specific absolute tolerances are often the most direct choice for states with different magnitudes but compatible units. If states have different units, or their scales are awkward to manage, nondimensionalizing or rescaling the state variables can provide a common numerical scale. Such transformations must be applied consistently to the equations and interpreted outputs; they do not remove the need to check predictions under tighter settings. A single scalar absolute tolerance is simpler, but is appropriate only when one error scale makes sense for all components.
Warning signs include a small state changing substantially under modest tolerance tightening, a near-zero state becoming negative when nonnegative values are expected, or a large increase in run time with negligible changes in quantities of interest. These signs motivate investigation rather than identify a unique cause: model behavior, solver choice, and other numerical difficulties can also matter. Tolerances alone do not repair an incorrect model, guarantee conservation, or ensure accurate results at every time point. The local error estimate itself is approximate, and the norm used to combine components may allow some component ratios to exceed one.
Check Your Understanding
Answer each question, then check the short solutions.
- A state has units of nM, \(rtol=10^{-3}\), \(atol=10^{-7}\) nM, and value \(0.002\) nM. What is its error scale?
- Why does relative tolerance alone not provide a useful finite error scale at a state value of zero?
- A solver uses a shared \(atol=0.05\) nM for states near \(100\) nM and \(0.001\) nM. Give one potential problem.
- Does a solver tolerance guarantee that the final solution’s error is less than the tolerance? Explain briefly.
- A tighter run changes an important prediction by more than the smallest difference the analysis aims to resolve. What should you do?
Solutions. (1) \(10^{-7}+10^{-3}(0.002)=2.1\times10^{-6}\) nM. (2) At zero, the relative contribution \(rtol|y|\) is zero, so it supplies no positive scale; the absolute tolerance provides the floor. (3) The absolute scale may be too large relative to the \(0.001\) nM state, allowing consequential estimated error in that component. (4) No. It scales the solver’s estimated local error; it is not a guaranteed bound on final solution error. (5) Investigate the discrepancy, assess the relevant states and outputs, and choose settings that make the numerical change acceptably small for the intended analysis.