DDE-Solver: Diverse Equation Approaches

• DDE-Solver is a polysemous term defining solver families that transform implicit differential or functional equations into numerical trajectories, continuous approximations, or algebraic certificates.
• It encompasses methods for delay differential equations via ODE embeddings, randomized Runge–Kutta, Galerkin projections, and RBF collocation, addressing issues like discontinuities and stiffness.
• It also includes symbolic solvers in combinatorics and learned solvers for diffusion ODE models, showcasing its versatility across numerical and algebraic problem domains.

“DDE-Solver” designates not a single algorithm but several solver families that share a common aim: to turn an implicit differential or functional specification into a computable object such as a numerical trajectory, a continuous approximation, or an annihilating polynomial. In current research usage, the label spans numerical solvers for delay differential equations, symbolic solvers for discrete differential equations in enumerative combinatorics, and, in an extrapolative usage, learned solvers for diffusion ODE sampling (Widmann et al., 2022, Notarantonio, 10 Sep 2025, Wang et al., 27 May 2025).

1. Scope of the term and underlying equation classes

In the delay-differential-equation sense, the basic model is a retarded DDE with constant delay,

$x'(t)=f\bigl(t,x(t),x(t-\tau)\bigr),$

supplemented by a history function on an interval $t\le t_0$. The numerical difficulty is not limited to evaluating $x(t-\tau)$: propagated discontinuities, method-of-steps restrictions, continuous output, floating-point timing of discontinuities, and stiffness all become central solver concerns (Widmann et al., 2022). A more specialized class is the “time-irregular” DDE studied by randomized Runge–Kutta methods, where $f$ is Lipschitz in the current state but only Hölder continuous in time and in the delayed argument (Difonzo et al., 2024).

In the combinatorial sense, a discrete differential equation is a functional equation for a bivariate formal power series

$F(t,u)\in \mathbb{Q}[u][[t]],$

with $u$ a catalytic variable and $\Delta_a$ a divided-difference operator. In that setting, the computational target is not a numerical orbit but an algebraic certificate such as

$R\bigl(t,F(t,a)\bigr)=0,$

which identifies the specialization $F(t,a)$ as an algebraic series (Bostan et al., 2023, Notarantonio, 10 Sep 2025).

This suggests that “DDE-Solver” is best treated as a polysemous technical label. Its meanings are linked by solver design, elimination, and approximation, but the mathematical objects being solved are substantially different.

2. Delay differential equation solvers built from ODE infrastructure

A defining modern architecture is “DelayDiffEq: Generating Delay Differential Equation Solvers via Recursive Embedding of Ordinary Differential Equation Solvers” (Widmann et al., 2022). Its central design is to construct DDE solvers by recursively embedding an ODE solver within itself, reusing OrdinaryDiffEq infrastructure for explicit Runge–Kutta, Rosenbrock, SDIRK, BDF, FIRK, dense output, DAEs, events, AD, and linear solvers. The history function is implemented by querying the dense interpolant of the same integrator, so delayed state evaluation becomes an interpolation call into ODE integrator data structures (Widmann et al., 2022).

This architecture supports the classical method of steps but is not restricted to step sizes below the minimum delay. When a step $[t_{k-1},t_k]$ is larger than the minimum delay, DelayDiffEq treats the DDE step as an implicit equation for the entire step and solves it by a fixed-point iteration on the function space of solutions on $t\le t_0$0. The solver extrapolates a dense solution onto the new step, performs an ODE step using that guess as history, reconstructs a new dense approximation, and repeats until consistency is achieved. For non-FIRK implicit methods, Jacobians are computed and factorized only once per step and then reused across fixed-point corrections, which is especially important for stiff DDEs (Widmann et al., 2022).

The same framework also addresses propagated discontinuities. User-defined discontinuities and those stemming from constant delays are inserted as grid points, while state-dependent delays are handled by a RADAR5-like root-finding procedure for times $t\le t_0$1 satisfying

$t\le t_0$2

where $t\le t_0$3 is a known past discontinuity. The Julia interface reflects the OrdinaryDiffEq style: $F(t,a)$0 or, for stiff problems, $F(t,a)$1 (Widmann et al., 2022).

A distinct solver line addresses low-regularity problems directly. “A Randomized Runge-Kutta Method for time-irregular delay differential equations” (Difonzo et al., 2024) studies deterministic DDEs with constant delay $t\le t_0$4 on multiple delay intervals, under assumptions that include Hölder continuity in time and in the delayed variable. The method is a randomized two-stage Runge–Kutta scheme with i.i.d. $t\le t_0$5, random time

$t\le t_0$6

and random stepsize $t\le t_0$7. The same random variable is used consistently in both the current-state and delayed intermediate stages, and for fixed horizon parameter $t\le t_0$8 the method uses $t\le t_0$9 evaluations of $x(t-\tau)$0 (Difonzo et al., 2024).

3. Alternative numerical architectures for DDE-type solving

One prominent alternative is the Galerkin projection route developed in “Approximating strange attractors and Lyapunov exponents of delay differential equations using Galerkin projections” (Sadath et al., 2018). A scalar DDE

$x(t-\tau)$1

is rewritten as an advection PDE with nonlinear boundary condition by introducing

$x(t-\tau)$2

This yields

$x(t-\tau)$3

with the DDE entering at $x(t-\tau)$4. A Galerkin expansion in shifted Legendre polynomials,

$x(t-\tau)$5

combined with the tau method, produces an ODE system

$x(t-\tau)$6

For smooth solutions, the error decreases exponentially as the number of terms used in the Galerkin approximation increases, and the method was used to approximate strange attractors and Lyapunov exponents with fewer ODEs than a standard method-of-lines discretization (Sadath et al., 2018).

A second alternative is RBF collocation. “Solving delay differential equations through RBF collocation” (Bernal et al., 2017) proposes a method based on Hardy’s multiquadric,

$x(t-\tau)$7

combined with the Residual Subsampling Algorithm of Driscoll and Heryudono for support adaptivity. The core development is for a scalar, first-order, linear DDE with a possibly time-varying delay,

$x(t-\tau)$8

with a history function $x(t-\tau)$9 for $f$0. The method naturally evaluates $f$1 at off-grid points through the RBF interpolant, and the same collocation framework is extended in the paper to neutral DDEs, state-dependent delays, and higher-order equations. For piecewise smooth solutions, the domain is split into smooth subintervals and the solver is applied segment-wise (Bernal et al., 2017).

In a broader differential-equation sense, “D-iteration: application to differential equations” (Hong, 2012) describes a Gauss–Seidel-equivalent but locally adaptive solver for linear systems arising from discretized differential equations. The method replaces row-based “collection” by column-based “diffusion,” maintains a history vector $f$2 and a fluid vector $f$3, and diffuses only where the local residual exceeds a threshold. On the 2D heat equation with a 5-point stencil, the update sends $f$4 to each neighbor and uses an open flag to avoid scanning inactive nodes. This is not a delay-equation solver, but it exemplifies residual-driven solver design for sparse discretizations (Hong, 2012).

4. Convergence, stability, and demanding test regimes

For time-irregular DDEs, the central convergence result is Theorem 4.2 of (Difonzo et al., 2024). With $f$5 and

$f$6

the randomized two-stage Runge–Kutta approximation $f$7 satisfies

$f$8

for each delay interval $f$9. On the first interval the rate is $F(t,u)\in \mathbb{Q}[u][[t]],$0, and when $F(t,u)\in \mathbb{Q}[u][[t]],$1 this gives $F(t,u)\in \mathbb{Q}[u][[t]],$2. The paper emphasizes that the observed negative slopes in experiments are consistently higher than the theoretical lower bound (Difonzo et al., 2024).

A different notion of solver appears in “Stability analysis of delay differential equations via Semidefinite programming” (Su, 2016). There, the DDE is linear, time-invariant, and has multiple constant delays: $F(t,u)\in \mathbb{Q}[u][[t]],$3 After discretization, stability is reduced to the positive definiteness of

$F(t,u)\in \mathbb{Q}[u][[t]],$4

and, in the singular semidefinite case, to the LMI

$F(t,u)\in \mathbb{Q}[u][[t]],$5

The resulting SDP,

$F(t,u)\in \mathbb{Q}[u][[t]],$6

is a stability-analysis solver rather than a trajectory integrator: it certifies asymptotic stability of the zero solution or identifies stabilizing parameters $F(t,u)\in \mathbb{Q}[u][[t]],$7 (Su, 2016).

The nonlinear DDE

$F(t,u)\in \mathbb{Q}[u][[t]],$8

provides a particularly severe qualitative benchmark (Davidow et al., 2017). For $F(t,u)\in \mathbb{Q}[u][[t]],$9 the system is conservative and exhibits no limit cycles; for $u$0, no matter how small, an infinite number of limit cycles exist, their amplitudes going to infinity in the limit as $u$1. Harmonic balance gives

$u$2

with amplitudes growing like $u$3, while a Melnikov-type calculation yields the closely related asymptotic

$u$4

Adding damping,

$u$5

regularizes the singularity through a Hopf threshold and saddle-node curves of limit cycles (Davidow et al., 2017). This makes the equation a demanding benchmark for nonlinear DDE solvers and continuation methods.

5. Symbolic-algebraic DDE-Solver for discrete differential equations

In enumerative combinatorics, DDE-Solver is the proper name of a Maple package for “Discrete Differential Equations” (Notarantonio, 10 Sep 2025). The formal setting uses a bivariate series

$u$6

and the divided-difference operator

$u$7

A DDE of order $u$8 is written as

$u$9

or, after clearing denominators, as a polynomial equation

$\Delta_a$0

(Notarantonio, 10 Sep 2025, Bostan et al., 2023).

The algebraicity theorem of Bousquet-Mélou and Jehanne implies that, for fixed-point type DDEs, the solution $\Delta_a$1 is unique and algebraic. The computational task is therefore to produce a polynomial witness

$\Delta_a$2

The Maple package exposes this through annihilating_polynomial(P, k, var, algorithm, variable), where algorithm can be "duplication", "elimination", "geometry", or "hybrid", and the output is a nonzero polynomial $\Delta_a$3 such that $\Delta_a$4 (Notarantonio, 10 Sep 2025).

The underlying algorithms mirror the structural developments of “Fast Algorithms for Discrete Differential Equations” (Bostan et al., 2023). The duplication algorithm follows the Bousquet‑Mélou–Jehanne system with $\Delta_a$5 copies of the catalytic variable; the elimination algorithm computes a constructible description of fibers with at least $\Delta_a$6 distinct $\Delta_a$7-coordinates; the geometry algorithm uses multiplication matrices and characteristic polynomials and is implemented only for $\Delta_a$8; and the hybrid algorithm combines elimination-based degree bounds with series guessing via Hermite–Padé approximation. The package compares the timings of all these algorithms on DDEs from the literature and shows that there is no universally best algorithm, though elimination is the default and hybrid is particularly effective on difficult instances such as 3-greedy Tamari intervals and 5-constellations (Notarantonio, 10 Sep 2025).

For the 3‑constellation example, the package returns a factorization whose cubic factor annihilates

$\Delta_a$9

thereby turning a catalytic functional equation into a concrete algebraic relation (Notarantonio, 10 Sep 2025).

6. Extrapolative and adjacent uses of the label

An extrapolative use of “DDE-Solver” appears in “Differentiable Solver Search for Fast Diffusion Sampling” (Wang et al., 27 May 2025). There the problem is reverse-diffusion ODE sampling rather than delay equations. The paper argues that Adams-like multistep solvers based only on $R\bigl(t,F(t,a)\bigr)=0,$0-related Lagrange interpolation are suboptimal for diffusion models, defines a compact search space of time steps and solver coefficients, and learns a linear multistep solver by differentiable trajectory matching. The searched solver attains FID $R\bigl(t,F(t,a)\bigr)=0,$1 for SiT-XL/2, FID $R\bigl(t,F(t,a)\bigr)=0,$2 for FlowDCN-XL/2, and FID $R\bigl(t,F(t,a)\bigr)=0,$3 for DiT-XL/2 on ImageNet256 with only 10 steps (Wang et al., 27 May 2025). This is not a delay-differential-equation solver, but it shows the term being repurposed for solver search in diffusion ODEs.

A related adjacent usage is “A polynomial-based constrained solver for fuel-optimal low-thrust trajectory optimization” (Caleb et al., 1 Feb 2025). The proposed DADDy solver combines constrained Differential Dynamic Programming, differential algebra-based polynomial expansions of dynamics and cost, and a polynomial-based Newton solver for high-precision constraint enforcement. The paper reports that the most stable method achieved $R\bigl(t,F(t,a)\bigr)=0,$4 convergence and runtime reductions of $R\bigl(t,F(t,a)\bigr)=0,$5 in the Sun-centered two-body problem, $R\bigl(t,F(t,a)\bigr)=0,$6 to $R\bigl(t,F(t,a)\bigr)=0,$7 in the Earth-Moon CR3BP, and $R\bigl(t,F(t,a)\bigr)=0,$8 to $R\bigl(t,F(t,a)\bigr)=0,$9 in the Earth-centered two-body problem (Caleb et al., 1 Feb 2025). Again, this is not a delay solver in the standard sense. It is, however, another instance where polynomial approximation, elimination of expensive evaluations, and constraint-preserving updates define what a differential-equation solver can mean in contemporary research.

Taken together, these usages indicate that “DDE-Solver” is not a single canonical method. In delay-equation numerics it names architectures such as method-of-steps embeddings, randomized Runge–Kutta, Galerkin reductions, and RBF collocation; in combinatorics it names a Maple package for discrete differential equations; and in adjacent areas it can denote solver-search or polynomial-optimal-control frameworks whose defining feature is systematic solution of differential or functional equations rather than any specific delay formalism (Widmann et al., 2022, Notarantonio, 10 Sep 2025, Wang et al., 27 May 2025).