1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
//! ETDRK4 exponential time-differencing driver (ADR-0189 §D2).
//!
//! Solves the semilinear initial-value problem
//! ```text
//! ∂_t u = L u + N(u), u(0) = u₀
//! ```
//! where `L` is a linear generator ([`GeneratorAction`]) and `N` is a pointwise
//! nonlinearity ([`Nonlinearity`]).
//!
//! ## Algorithm
//!
//! Cox–Matthews (2002) fourth-order exponential time differencing (ETDRK4).
//! The linear part is handled exactly via φ-functions; the nonlinear part is
//! evaluated at four explicit stages. No implicit solve is required.
//!
//! ## Usage
//!
//! ```rust,ignore
//! use semiflow::{Etdrk4, AllenCahn, scratch::ScratchPool};
//!
//! // Build driver (validates dim > 0)
//! let driver = Etdrk4::new(my_op, AllenCahn::<f64>::new(), h)?;
//!
//! // Single step
//! driver.step(&u, &mut u_next, &mut scratch)?;
//!
//! // Multiple steps; final state written to `out`
//! driver.integrate(&u0, 100, &mut out, &mut scratch)?;
//! ```
//!
//! ## References
//!
//! - Cox & Matthews (2002) *J. Comput. Phys.* **176**, 430–455.
//! - ADR-0189; `contracts/semiflow-core.math.md` §58.
use crate::{
error::SemiflowError, etdrk4_helpers::etdrk4_step as step_impl, float::SemiflowFloat,
generator_action::GeneratorAction, nonlinearity::Nonlinearity, phi_krylov::PhiSolver,
scratch::ScratchPool,
};
// ---------------------------------------------------------------------------
// Etdrk4 driver struct
// ---------------------------------------------------------------------------
/// ETDRK4 semilinear time-stepping driver.
///
/// Wraps a linear generator `Op` and a nonlinearity `Nl` together with a
/// fixed step size `h`. After construction via [`Etdrk4::new`], call
/// [`step`](Etdrk4::step) for single steps or [`integrate`](Etdrk4::integrate)
/// for multiple steps.
///
/// Generic parameters:
/// - `F`: scalar float type (`f32` or `f64`).
/// - `Op`: linear generator implementing [`GeneratorAction<F>`].
/// - `Nl`: nonlinearity implementing [`Nonlinearity<F>`].
pub struct Etdrk4<F: SemiflowFloat, Op: GeneratorAction<F>, Nl: Nonlinearity<F>> {
op: Op,
nl: Nl,
h: F,
n: usize,
solver: PhiSolver,
}
impl<F, Op, Nl> Etdrk4<F, Op, Nl>
where
F: SemiflowFloat,
Op: GeneratorAction<F>,
Nl: Nonlinearity<F>,
{
/// Construct the ETDRK4 driver.
///
/// # Errors
/// Returns [`SemiflowError::DomainViolation`] if `op.dim() == 0`.
pub fn new(op: Op, nl: Nl, h: F) -> Result<Self, SemiflowError> {
let n = op.dim();
if n == 0 {
return Err(SemiflowError::DomainViolation {
what: "Etdrk4::new: op.dim() must be > 0",
value: 0.0,
});
}
Ok(Self {
op,
nl,
h,
n,
solver: PhiSolver::Taylor,
})
}
/// Choose the φ evaluator (ADR-0211): [`PhiSolver::Taylor`] (default, exact
/// sweep counts) or [`PhiSolver::Krylov`] (adaptive Arnoldi, for stiff
/// non-symmetric generators where `h‖G‖ ≫ 10²`).
#[must_use]
pub fn with_phi_solver(mut self, solver: PhiSolver) -> Self {
self.solver = solver;
self
}
/// Advance one step: `u_next ← ETDRK4(u, h)`.
///
/// Borrows temporary buffers from `scratch`; all are returned before
/// this call exits.
///
/// # Errors
/// Propagates any error from the generator or nonlinearity.
pub fn step(
&self,
u: &[F],
u_next: &mut [F],
scratch: &mut ScratchPool<F>,
) -> Result<(), SemiflowError> {
step_impl(&self.solver, &self.op, &self.nl, self.h, u, u_next, scratch)
}
/// Integrate `n_steps` steps from `u0`; final state written to `out[..n]`.
///
/// Uses one additional scratch buffer of size `n` for the ping-pong swap.
///
/// # Errors
/// Propagates any error from [`step`](Etdrk4::step).
pub fn integrate(
&self,
u0: &[F],
n_steps: usize,
out: &mut [F],
scratch: &mut ScratchPool<F>,
) -> Result<(), SemiflowError> {
let n = self.n;
out[..n].copy_from_slice(&u0[..n]);
let mut tmp = scratch.take_vec(n);
for _ in 0..n_steps {
self.step(out, &mut tmp, scratch)?;
out[..n].copy_from_slice(&tmp[..n]);
}
scratch.return_vec(tmp);
Ok(())
}
}
// ---------------------------------------------------------------------------
// Unit tests (fast)
// ---------------------------------------------------------------------------
#[cfg(test)]
include!("etdrk4_tests_mod.rs");