rkstiff.models

Common benchmark models and initial conditions for stiff PDE solvers

This module provides utilities for constructing canonical test equations used to validate time-integration schemes such as exponential, integrating-factor, and Runge–Kutta methods.

Each model supplies:

• The governing PDE

• The corresponding linear and nonlinear operators

• Representative initial conditions

Models included

• Korteweg–de Vries (KdV) equation — soliton dynamics

• Viscous Burgers equation — nonlinear advection–diffusion

• Allen–Cahn equation — bistable phase-field dynamics

Functions

allen_cahn_ops(x, d_cheb_matrix[, epsilon])	Return the linear and nonlinear operators for the Allen–Cahn equation.
burgers_ops(kx, mu)	Return the linear and nonlinear operators for the viscous Burgers equation.
kdv_multi_soliton(x, ampl, x0[, t])	Construct a multi-soliton superposition for the KdV equation.
kdv_ops(kx)	Return the linear and nonlinear operators for the KdV equation.
kdv_soliton(x[, ampl, x0, t])	Return the analytic single-soliton solution of the Korteweg–de Vries (KdV) equation.

rkstiff.models.kdv_soliton(x, ampl=0.5, x0=0.0, t=0.0)

Return the analytic single-soliton solution of the Korteweg–de Vries (KdV) equation.

The KdV equation is:

\[\frac{\partial u}{\partial t} + 6 u \frac{\partial u}{\partial x} + \frac{\partial^3 u}{\partial x^3} = 0.\]

Its single-soliton solution is given by:

\[u(x,t) = \tfrac{1}{2} a^2 \operatorname{sech}^2\! \left[\tfrac{1}{2} a (x - x_0 - a^2 t)\right],\]

where \(a\) is the amplitude parameter.

Parameters:

• x (np.ndarray) – Spatial grid points.

• ampl (float, optional) – Soliton amplitude \(a\). Default is 0.5.

• x0 (float, optional) – Initial position \(x_0\). Default is 0.0.

• t (float, optional) – Time. Default is 0.0.

Returns:

Soliton profile \(u(x,t)\).

Return type:

np.ndarray

rkstiff.models.kdv_multi_soliton(x, ampl, x0, t=0.0)

Construct a multi-soliton superposition for the KdV equation.

Parameters:

• x (np.ndarray) – Spatial grid points.

• ampl (Sequence[float]) – Sequence of soliton amplitudes \((a_1, a_2, \dots, a_m)\).

• x0 (Sequence[float]) – Initial positions \((x_{0,1}, x_{0,2}, \dots, x_{0,m})\).

• t (float, optional) – Time. Default is 0.0.

Returns:

Sum of all soliton profiles at time t.

Return type:

np.ndarray

Raises:

ValueError – If ampl and x0 have mismatched lengths.

rkstiff.models.kdv_ops(kx)

Return the linear and nonlinear operators for the KdV equation.

Governing PDE:

\[\frac{\partial u}{\partial t} = -6u \frac{\partial u}{\partial x} - \frac{\partial^3 u}{\partial x^3}.\]

Parameters:

kx (np.ndarray) – Wavenumber array in Fourier space.

Return type:

Tuple[ndarray, Callable[[ndarray], ndarray]]

Returns:

• lin_op (np.ndarray) – Linear operator in spectral space \(L = i k_x^3\).

• nl_func (Callable[[np.ndarray], np.ndarray]) – Nonlinear term \(N(u) = -6\,\mathcal{F}[u u_x]\), where \(u_x = \mathcal{F}^{-1}[i k_x \hat{u}]\).

rkstiff.models.burgers_ops(kx, mu)

Return the linear and nonlinear operators for the viscous Burgers equation.

Governing PDE:

\[\frac{\partial u}{\partial t} + u \frac{\partial u}{\partial x} = \mu \frac{\partial^2 u}{\partial x^2},\]

where \(\mu > 0\) is the kinematic viscosity.

Parameters:

• kx (np.ndarray) – Wavenumber array in Fourier space.

• mu (float) – Viscosity coefficient.

Return type:

Tuple[ndarray, Callable[[ndarray], ndarray]]

Returns:

• lin_op (np.ndarray) – Linear operator \(L = -\mu k_x^2\).

• nl_func (Callable[[np.ndarray], np.ndarray]) – Nonlinear term \(N(u) = -\mathcal{F}[u u_x]\).

Notes

In Fourier space, differentiation corresponds to multiplication by \(i k_x\), allowing both linear diffusion and nonlinear advection terms to be computed spectrally.

rkstiff.models.allen_cahn_ops(x, d_cheb_matrix, epsilon=0.01)

Return the linear and nonlinear operators for the Allen–Cahn equation.

Governing PDE:

\[\frac{\partial u}{\partial t} = \epsilon \frac{\partial^2 u}{\partial x^2} + u - u^3,\]

where \(\epsilon\) is a small positive diffusion coefficient.

Parameters:

• x (np.ndarray) – Spatial grid (Chebyshev–Gauss–Lobatto points).

• d_cheb_matrix (np.ndarray) – Chebyshev differentiation matrix \(D\).

• epsilon (float, optional) – Diffusion parameter \(\epsilon\). Default is 0.01.

Return type:

Tuple[ndarray, Callable[[ndarray], ndarray]]

Returns:

• lin_op (np.ndarray) – Linear operator \(L = \epsilon D^2 + I\) with boundary rows removed.

• nl_func (Callable[[np.ndarray], np.ndarray]) – Nonlinear term \(N(u) = u - u^3\).

Notes

Boundary points are excluded from lin_op since Dirichlet or Neumann boundary conditions are typically enforced separately.