OrdinaryDiffEqTaylorSeries

Taylor series methods for ordinary differential equations using automatic differentiation. These methods achieve very high-order accuracy by computing Taylor expansions of the solution using automatic differentiation techniques through TaylorDiff.jl.

Development Status

These methods are still in development and may not be fully optimized or reliable for production use.

Key Properties

Taylor series methods provide:

Very high-order accuracy with arbitrary order capability
Automatic differentiation for derivative computation
Step size control through Taylor series truncation
Natural error estimation from higher-order terms
Excellent accuracy for smooth problems
Single-step methods without requiring history

When to Use Taylor Series Methods

These methods are recommended for:

Ultra-high precision problems where maximum accuracy is needed
Smooth problems with well-behaved derivatives
Scientific computing requiring very low error tolerances
Problems with expensive function evaluations where high-order methods reduce total steps

Mathematical Background

Taylor series methods compute the solution using Taylor expansions: u(t + h) = u(t) + h*u'(t) + h²/2!*u''(t) + h³/3!*u'''(t) + ...

The derivatives are computed automatically using automatic differentiation, allowing arbitrary-order methods without manual derivative computation.

Solver Selection Guide

Available Methods

ExplicitTaylor2: Second-order Taylor series method for moderate accuracy
ExplicitTaylor: Arbitrary-order Taylor series method (specify order with order = Val{p}())

Usage considerations

Smooth problems only: Methods assume the function has many continuous derivatives
Computational cost: Higher orders require more automatic differentiation computations
Memory requirements: Higher orders store more derivative information

Performance Guidelines

When Taylor series methods excel

Very smooth problems where high-order derivatives exist and are well-behaved
High precision requirements beyond standard double precision
Long-time integration where accumulated error matters
Problems where function evaluations dominate computational cost

Problem characteristics

Polynomial and analytic functions work extremely well
Smooth ODEs from physics simulations
Problems requiring very low tolerances (< 1e-12)

Limitations and Considerations

Method limitations

Requires smooth functions - non-smooth problems may cause issues
Memory overhead for storing multiple derivatives
Limited to problems where high-order derivatives are meaningful
Automatic differentiation compatibility - requires functions compatible with TaylorDiff.jl and Symbolics.jl tracing
Long compile times due to automatic differentiation and symbolic processing overhead

When to consider alternatives

Non-smooth problems: Use adaptive Runge-Kutta methods instead
Stiff problems: Taylor methods are explicit and may be inefficient
Large systems: Automatic differentiation cost may become prohibitive
Standard accuracy needs: Lower-order methods are often sufficient

Alternative Approaches

Consider these alternatives:

High-order Runge-Kutta methods (Feagin, Verner) for good accuracy with less overhead
Extrapolation methods for high accuracy with standard function evaluations
Adaptive methods for problems with varying smoothness
Implicit methods for stiff problems requiring high accuracy

Installation and Usage

Taylor series methods require explicit installation of the specialized library:

using Pkg
Pkg.add("OrdinaryDiffEqTaylorSeries")

Then use the methods with:

using OrdinaryDiffEqTaylorSeries: ExplicitTaylor2, ExplicitTaylor,
    ExplicitTaylorAdaptiveOrder
using CommonSolve: solve
using SciMLBase: ODEProblem

# Example: Second-order Taylor method
function f(u, p, t)
    σ, ρ, β = p
    du1 = σ * (u[2] - u[1])
    du2 = u[1] * (ρ - u[3]) - u[2]
    du3 = u[1] * u[2] - β * u[3]
    return [du1, du2, du3]
end

u0 = [1.0, 0.0, 0.0]
tspan = (0.0, 10.0)
p = [10.0, 28.0, 8 / 3]
prob = ODEProblem(f, u0, tspan, p)

# Second-order Taylor method
sol = solve(prob, ExplicitTaylor2(), dt = 0.01)

# Arbitrary-order Taylor method (e.g., 8th order)
sol = solve(prob, ExplicitTaylor(order = Val(8)))

# Adaptive-order Taylor method (orders 4 through 7)
sol = solve(
    prob, ExplicitTaylorAdaptiveOrder(min_order = Val(4), max_order = Val(8)),
    abstol = 1.0e-10, reltol = 1.0e-10,
)

Full list of solvers

OrdinaryDiffEqTaylorSeries.ExplicitTaylor2 — Type

ExplicitTaylor2(;
    stage_limiter! = trivial_limiter!,
    step_limiter! = trivial_limiter!, thread = Serial()
)

Second-order explicit Taylor series method. ExplicitTaylor2 evaluates the first two time derivatives of the solution with TaylorDiff.jl and advances with a fixed step size. Supply dt when calling solve.

The right-hand side must accept TaylorDiff numbers. In particular, avoid branches or external functions that discard or reject Taylor coefficients.

Fields

stage_limiter!: limiter configuration retained by the common explicit-solver interface. The current Taylor step does not invoke this callback.
step_limiter!: deprecated algorithm-level fallback for the accepted-step limiter. When no solve-level step_limiter is supplied, a nontrivial fallback is invoked centrally once per accepted step.
thread: threading configuration retained by the common explicit-solver interface. The second-order implementation currently performs its internal broadcasts serially.

Keywords

stage_limiter! = trivial_limiter!: stage-limiter callback with signature limiter!(u, integrator, p, t).
step_limiter! = trivial_limiter!: deprecated algorithm-level fallback with signature limiter!(u, integrator, p, t). When no solve-level step_limiter is supplied, a nontrivial fallback is invoked centrally once per accepted step.
thread = Serial(): FastBroadcast threading mode stored with the algorithm.

Returns

An ExplicitTaylor2 algorithm object for use with solve or init.

Examples

using OrdinaryDiffEqTaylorSeries: ExplicitTaylor2
using CommonSolve: solve
using SciMLBase: ODEProblem

prob = ODEProblem((u, p, t) -> u, 1.0, (0.0, 1.0))
sol = solve(prob, ExplicitTaylor2(), dt = 0.05)
sol(1.0)

source

OrdinaryDiffEqTaylorSeries.ExplicitTaylor — Type

ExplicitTaylor(;
    order = Val{1}(), stage_limiter! = trivial_limiter!,
    step_limiter! = trivial_limiter!, thread = Serial()
)

Fixed-order explicit Taylor series method with adaptive step-size control. The method symbolically constructs a Taylor jet of the requested order and reuses it while solving. Choose order = Val(N) to use an order-N expansion.

The right-hand side must be traceable by Symbolics.jl and must accept TaylorDiff numbers. Changing the right-hand-side function, parameters, or order constructs a different cached jet, so these should remain stable across repeated solves when compilation cost matters.

Fields

order: compile-time Val containing the Taylor expansion order.
stage_limiter!: limiter configuration retained by the common explicit-solver interface. The current Taylor step does not invoke this callback.
step_limiter!: deprecated algorithm-level fallback for the accepted-step limiter. When no solve-level step_limiter is supplied, a nontrivial fallback is invoked centrally once per accepted step.
thread: FastBroadcast threading mode used for eligible in-place error-estimate broadcasts.

Keywords

order = Val{1}(): Taylor expansion order, represented as Val(N) for a positive integer N.
stage_limiter! = trivial_limiter!: stage-limiter callback with signature limiter!(u, integrator, p, t).
step_limiter! = trivial_limiter!: deprecated algorithm-level fallback with signature limiter!(u, integrator, p, t). When no solve-level step_limiter is supplied, a nontrivial fallback is invoked centrally once per accepted step.
thread = Serial(): FastBroadcast threading mode. Use a threading mode such as Threaded() for eligible broadcasts when Julia has multiple threads.

Returns

An ExplicitTaylor algorithm object for use with solve or init.

Examples

using OrdinaryDiffEqTaylorSeries: ExplicitTaylor
using CommonSolve: solve
using SciMLBase: ODEProblem

prob = ODEProblem((u, p, t) -> u, 1.0, (0.0, 1.0))
sol = solve(
    prob, ExplicitTaylor(order = Val(8)),
    abstol = 1.0e-10, reltol = 1.0e-10,
)
sol(1.0)

source

OrdinaryDiffEqTaylorSeries.ExplicitTaylorAdaptiveOrder — Type

ExplicitTaylorAdaptiveOrder(;
    min_order = Val{1}(), max_order = Val{10}(),
    stage_limiter! = trivial_limiter!, step_limiter! = trivial_limiter!,
    thread = Serial()
)

Explicit Taylor series method that adapts both the expansion order and the step size. At each accepted step, the controller compares nearby orders and selects the estimated lowest-work choice for the next step.

The right-hand side must be traceable by Symbolics.jl and must accept TaylorDiff numbers. The admissible order window is min_order:(max_order - 1) because one additional Taylor coefficient is required to estimate the local error.

Fields

min_order: compile-time Val containing the lowest admissible step order.
max_order: compile-time Val containing the Taylor jet order. The highest admissible step order is one less than this value.
stage_limiter!: limiter configuration retained by the common explicit-solver interface. The current Taylor step does not invoke this callback.
step_limiter!: deprecated algorithm-level fallback for the accepted-step limiter. When no solve-level step_limiter is supplied, a nontrivial fallback is invoked centrally once per accepted step.
thread: FastBroadcast threading mode used for eligible in-place error-estimate broadcasts.

Keywords

min_order = Val{1}(): lowest order considered by the adaptive-order controller, represented as Val(N).
max_order = Val{10}(): Taylor jet order, represented as Val(N). It must be greater than min_order.
stage_limiter! = trivial_limiter!: stage-limiter callback with signature limiter!(u, integrator, p, t).
step_limiter! = trivial_limiter!: deprecated algorithm-level fallback with signature limiter!(u, integrator, p, t). When no solve-level step_limiter is supplied, a nontrivial fallback is invoked centrally once per accepted step.
thread = Serial(): FastBroadcast threading mode. Use a threading mode such as Threaded() for eligible broadcasts when Julia has multiple threads.

Returns

An ExplicitTaylorAdaptiveOrder algorithm object for use with solve or init.

Throws

solve throws ArgumentError when max_order <= min_order.

Examples

using OrdinaryDiffEqTaylorSeries: ExplicitTaylorAdaptiveOrder
using CommonSolve: solve
using SciMLBase: ODEProblem

prob = ODEProblem((u, p, t) -> u, 1.0, (0.0, 1.0))
sol = solve(
    prob, ExplicitTaylorAdaptiveOrder(min_order = Val(4), max_order = Val(9)),
    abstol = 1.0e-10, reltol = 1.0e-10,
)
sol(1.0)

source