|
| 1 | +""" |
| 2 | + GlobalRichardson(alg) |
| 3 | +
|
| 4 | +Wrap an ODE algorithm with global Richardson extrapolation. |
| 5 | +""" |
| 6 | +struct GlobalRichardson{A <: SciMLBase.AbstractODEAlgorithm} <: GlobalDiffEqAlgorithm |
| 7 | + alg::A |
| 8 | +end |
| 9 | + |
| 10 | +# Forward algorithm traits to the wrapped algorithm |
| 11 | +# This allows GlobalRichardson to inherit capabilities from the inner algorithm |
| 12 | +SciMLBase.allows_arbitrary_number_types(alg::GlobalRichardson) = |
| 13 | + SciMLBase.allows_arbitrary_number_types(alg.alg) |
| 14 | +SciMLBase.allowscomplex(alg::GlobalRichardson) = |
| 15 | + SciMLBase.allowscomplex(alg.alg) |
| 16 | +SciMLBase.isautodifferentiable(alg::GlobalRichardson) = |
| 17 | + SciMLBase.isautodifferentiable(alg.alg) |
| 18 | + |
| 19 | +function SciMLBase.__solve( |
| 20 | + prob::Union{SciMLBase.AbstractODEProblem, SciMLBase.AbstractDAEProblem}, |
| 21 | + alg::GlobalRichardson, args...; |
| 22 | + dt, kwargs... |
| 23 | + ) |
| 24 | + opt = Dict(kwargs) |
| 25 | + otheropts = delete!(copy(opt), :dt) |
| 26 | + tstops = get(opt, :tstops, range(prob.tspan[1], stop = prob.tspan[2], step = dt)) |
| 27 | + local sol |
| 28 | + val, |
| 29 | + err = Richardson.extrapolate( |
| 30 | + dt, rtol = get(opt, :reltol, 1.0e-3), |
| 31 | + atol = get(opt, :abstol, 1.0e-6), contract = 0.5 |
| 32 | + ) do _dt |
| 33 | + sol = solve(prob, alg.alg, args...; dt = _dt, adaptive = false, otheropts...) |
| 34 | + # Convert Vector{Vector{T}} to Matrix{T} for Richardson.jl compatibility |
| 35 | + reduce(hcat, sol.(tstops)) |
| 36 | + end |
| 37 | + return sol |
| 38 | +end |
0 commit comments