rkstiff.solveras

Base solver infrastructure for rkstiff adaptive step PDE solvers

Adaptive-step solver infrastructure for stiff time-dependent systems.

This module defines the BaseSolverAS abstract class and its associated SolverConfig, which provide the foundation for adaptive time-stepping solvers such as ETD(3,5) and other embedded Runge–Kutta–type integrators.

Adaptive solvers automatically adjust their time step to maintain a user-specified accuracy tolerance, based on local error estimates computed during integration.

Classes

`BaseSolverAS`(lin_op, nl_func[, config, loglevel])	Abstract base class for adaptive-step stiff solvers.
`SolverConfig`([epsilon, incr_f, decr_f, ...])	Configuration parameters for adaptive-step stiff solvers.

class rkstiff.solveras.SolverConfig(epsilon=1e-4, incr_f=1.25, decr_f=0.85, safety_f=0.8, adapt_cutoff=0.01, minh=1e-16)[source]

Bases: object

Configuration parameters for adaptive-step stiff solvers.

Defines numerical tolerances and scaling factors that control adaptive time step selection. Used by all adaptive solvers derived from BaseSolverAS.

Parameters:

epsilon (float, default=1e-4) – Relative error tolerance for adaptive stepping. Smaller values yield higher accuracy but may reduce efficiency.
incr_f (float, default=1.25) – Step-size increase factor applied after successful steps. Must satisfy incr_f > 1.0.
decr_f (float, default=0.85) – Step-size reduction factor applied after rejected steps. Must satisfy decr_f < 1.0.
safety_f (float, default=0.8) – Safety factor for damping overly aggressive changes in step size. Must satisfy safety_f <= 1.0.
adapt_cutoff (float, default=0.01) – Relative magnitude threshold for excluding small modes when computing adaptive step corrections.
minh (float, default=1e-16) – Minimum allowable step size.

Raises:

ValueError – If any parameter violates its constraint.

Notes

The combination of epsilon, safety_f, and adapt_cutoff primarily determines the balance between stability and adaptivity. Typical recommended settings are:

SolverConfig(epsilon=1e-5, safety_f=0.9, adapt_cutoff=0.001)

__init__(epsilon=1e-4, incr_f=1.25, decr_f=0.85, safety_f=0.8, adapt_cutoff=0.01, minh=1e-16)[source]

Initialize and validate adaptive solver configuration.

property epsilon: float: Relative error tolerance for adaptive stepping.

property incr_f: float: Increment factor for adaptive step sizing (must be > 1.0).

property decr_f: float: Decrement factor for adaptive step sizing (must be < 1.0).

property safety_f: float: Safety factor for adaptive stepping (must be ≤ 1.0).

property adapt_cutoff: float: Cutoff threshold for ignoring small modes during adaptivity.

property minh: float: Minimum allowable step size.

class rkstiff.solveras.BaseSolverAS(lin_op, nl_func, config=SolverConfig(), loglevel='WARNING')[source]

Bases: BaseSolver

Abstract base class for adaptive-step stiff solvers.

Provides an adaptive time-stepping framework for integrating semi-linear systems of the form

\[\frac{\partial \mathbf{U}}{\partial t} = \mathcal{L}\mathbf{U} + \mathcal{N}(\mathbf{U}),\]

where \(\mathcal{L}\) is a (possibly stiff) linear operator and \(\mathcal{N}\) is a nonlinear function.

Subclasses must implement _update_stages() (to perform the embedded Runge–Kutta stage updates), _reset() (to clear solver state), and _q() (to define the order of accuracy for step size control).

Parameters:

lin_op (np.ndarray) – Linear operator \(\mathcal{L}\) defining the stiff part of the system. Can be either 1D (diagonal) or a full 2D matrix.
nl_func (Callable[[np.ndarray], np.ndarray]) – Nonlinear function \(\mathcal{N}(\mathbf{U})\).
config (SolverConfig, optional) – Configuration object specifying tolerances and scaling factors for adaptive step control.
loglevel (str or int, default='WARNING') – Logging level for runtime diagnostics.

config

Adaptive step configuration (tolerances and safety parameters).

Type:: SolverConfig

epsilon

Relative error tolerance (config.epsilon).

Type:: float

minh

Minimum allowable step size (config.minh).

Type:: float

_accept

Flag indicating whether the most recent step was accepted.

Type:: bool

Raises:: ValueError – If lin_op has invalid dimensions or configuration parameters are inconsistent.

Notes

Adaptive time-stepping uses local error control via embedded Runge–Kutta pairs, dynamically adjusting \(h\) such that

\[\| \mathbf{e}_n \| \leq \varepsilon \, \| \mathbf{u}_n \|,\]

where \(\mathbf{e}_n\) is the local truncation error and \(\varepsilon\) is the tolerance.

Tip

Implementations typically define _update_stages() to return both a candidate solution and an error estimate, which are used to adjust the next step size.

exception SolverError[source]

Bases: RuntimeError

Base exception for solver-related runtime errors.

exception MaxLoopsExceeded[source]

Bases: SolverError

Raised when too many attempts are made to find a valid adaptive step.

exception MinimumStepReached[source]

Bases: SolverError

Raised when the adaptive step size falls below the minimum allowed value.

MAX_LOOPS = 50: Maximum retry attempts per adaptive step

MAX_S = 4.0: Maximum allowed step size increase factor

MIN_S = 0.25: Minimum allowed step size reduction factor

__init__(lin_op, nl_func, config=SolverConfig(), loglevel='WARNING')[source]

Initialize an adaptive-step solver with validated configuration.

set_loglevel(loglevel)

Adjust the solver’s logging verbosity at runtime.

Parameters:: loglevel (str or int) – New logging level. Accepts standard string levels or numeric constants from logging.
Return type:: None

Examples

>>> solver.set_loglevel("INFO")
>>> solver.set_loglevel(logging.DEBUG)

property solver_type: SolverType

Return the solver type for adaptive-step solvers.

Returns:: Always returns SolverType.ADAPTIVE_STEP.
Return type:: SolverType

Examples

>>> from rkstiff.etd35 import ETD35
>>> solver = ETD35(lin_op, nl_func)
>>> solver.solver_type == SolverType.ADAPTIVE_STEP
True

reset()[source]

Reset solver and clear adaptive-step state.

Return type:: None

step(u, h_suggest)[source]

Perform one adaptive integration step.

Attempts to advance the solution by time h_suggest and adjusts the step size automatically based on local error estimates.

Parameters:

u (np.ndarray) – Current state vector.
h_suggest (float) – Suggested time step size.

Return type:

Tuple[ndarray, float, float]

Returns:

unew (np.ndarray) – Updated solution vector after one accepted step.
h (float) – Actual step size used.
h_suggest (float) – Suggested step size for the next iteration.

Raises:

MaxLoopsExceeded – If too many attempts are made to find an acceptable step size.
MinimumStepReached – If the step size drops below the configured minimum minh.

Notes

The algorithm follows this pattern:

Try the proposed step.
Estimate local error and compute scaling factor s.
If s < 1 → reject step and reduce h.
If s ≥ 1 → accept step and update h_suggest for next step.

Warning

Very small or divergent s values may indicate instability or excessively tight tolerances.

evolve(u, t0, tf, h_init=None, store_data=True, store_freq=1)[source]

Integrate the system from \(t_0\) to \(t_f\) using adaptive time steps.

Repeatedly applies step() to propagate the solution forward while dynamically adjusting the time step size based on local error estimates.

Parameters:

u (np.ndarray) – Initial solution vector at \(t_0\).
t0 (float) – Initial time.
tf (float) – Final time.
h_init (float, optional) – Initial step size. Defaults to (tf - t0) / 100 if not provided.
store_data (bool, default=True) – Whether to store intermediate results in t and u.
store_freq (int, default=1) – Frequency of data storage; store every store_freq accepted steps.

Returns:

Final solution at \(t = t_f\).

Return type:

np.ndarray

Notes

The evolution proceeds until \(t \geq t_f\), automatically adjusting step sizes as needed.
Stored data is accessible via t and u.

Example

>>> solver = ETD35(lin_op, nl_func)
>>> u_final = solver.evolve(u0, t0=0.0, tf=10.0)
>>> solver.t[-1], np.linalg.norm(solver.u[-1])
(10.0, 0.0134)