Skip to content

New method Rodas6P #2774

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 7 commits into
base: master
Choose a base branch
from
Open

New method Rodas6P #2774

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
7 commits
Select commit Hold shift + click to select a range
bb82cd2
New method Rodas6P
gstein3m Jul 25, 2025
a4f1bff
Update rosenbrock_perform_step.jl
gstein3m Jul 25, 2025
c6c8282
Update lib/OrdinaryDiffEqRosenbrock/src/algorithms.jl
gstein3m Jul 25, 2025
e1c7661
Update lib/OrdinaryDiffEqRosenbrock/src/algorithms.jl
gstein3m Jul 25, 2025
84822c0
Merge remote-tracking branch 'upstream/master'
gstein3m Aug 13, 2025
4f1080d
Update Rodas6P
gstein3m Aug 14, 2025
10be42a
Update ode_rosenbrock_tests.jl
gstein3m Aug 15, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion lib/OrdinaryDiffEqRosenbrock/src/OrdinaryDiffEqRosenbrock.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -267,7 +267,7 @@ end

export Rosenbrock23, Rosenbrock32, RosShamp4, Veldd4, Velds4, GRK4T, GRK4A,
Ros4LStab, ROS3P, Rodas3, Rodas23W, Rodas3P, Rodas4, Rodas42, Rodas4P, Rodas4P2,
Rodas5, Rodas5P, Rodas5Pe, Rodas5Pr,
Rodas5, Rodas5P, Rodas5Pe, Rodas5Pr, Rodas6P,
RosenbrockW6S4OS, ROS34PW1a, ROS34PW1b, ROS34PW2, ROS34PW3, ROS34PRw, ROS3PRL,
ROS3PRL2, ROK4a,
ROS2, ROS2PR, ROS2S, ROS3, ROS3PR, Scholz4_7
Expand Down
4 changes: 3 additions & 1 deletion lib/OrdinaryDiffEqRosenbrock/src/alg_utils.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -33,6 +33,7 @@ alg_order(alg::Rodas5) = 5
alg_order(alg::Rodas5P) = 5
alg_order(alg::Rodas5Pr) = 5
alg_order(alg::Rodas5Pe) = 5
alg_order(alg::Rodas6P) = 6

alg_adaptive_order(alg::Rosenbrock32) = 2
alg_adaptive_order(alg::Rosenbrock23) = 3
Expand All @@ -59,10 +60,11 @@ isfsal(alg::Rodas4) = false
isfsal(alg::Rodas42) = false
isfsal(alg::Rodas4P) = false
isfsal(alg::Rodas4P2) = false
isfsal(alg::Rodas6P) = false

function has_stiff_interpolation(::Union{Rosenbrock23, Rosenbrock32, Rodas23W,
Rodas3P, Rodas4, Rodas4P, Rodas4P2, Rodas5,
Rodas5P, Rodas5Pe, Rodas5Pr})
Rodas5P, Rodas5Pe, Rodas5Pr, Rodas6P})
true
end

Expand Down
3 changes: 2 additions & 1 deletion lib/OrdinaryDiffEqRosenbrock/src/algorithms.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,8 @@ for (Alg, desc, refs, is_W) in [
(:Rodas5, "A 5th order A-stable stiffly stable Rosenbrock method with a stiff-aware 4th order interpolant.", "- Di Marzo G. RODAS5(4) – Méthodes de Rosenbrock d'ordre 5(4) adaptées aux problèmes\n différentiels-algébriques. MSc mathematics thesis, Faculty of Science,\n University of Geneva, Switzerland.", false),
(:Rodas5P, "A 5th order A-stable stiffly stable Rosenbrock method with a stiff-aware 4th order interpolant.\nHas improved stability in the adaptive time stepping embedding.", "- Steinebach G. Construction of Rosenbrock–Wanner method Rodas5P and numerical benchmarks\n within the Julia Differential Equations package.\n In: BIT Numerical Mathematics, 63(2), 2023. doi:10.1007/s10543-023-00967-x", true),
(:Rodas5Pe, "Variant of Rodas5P with modified embedded scheme.", "- Steinebach G. Rosenbrock methods within OrdinaryDiffEq.jl - Overview, recent developments and applications -\n Preprint 2024. Proceedings of the JuliaCon Conferences.\n https://proceedings.juliacon.org/papers/eb04326e1de8fa819a3595b376508a40", true),
(:Rodas5Pr, "Variant of Rodas5P with additional residual control.", "- Steinebach G. Rosenbrock methods within OrdinaryDiffEq.jl - Overview, recent developments and applications -\n Preprint 2024. Proceedings of the JuliaCon Conferences.\n https://proceedings.juliacon.org/papers/eb04326e1de8fa819a3595b376508a40", true)
(:Rodas5Pr, "Variant of Rodas5P with additional residual control.", "- Steinebach G. Rosenbrock methods within OrdinaryDiffEq.jl - Overview, recent developments and applications -\n Preprint 2024. Proceedings of the JuliaCon Conferences.\n https://proceedings.juliacon.org/papers/eb04326e1de8fa819a3595b376508a40", true),
(:Rodas6P, "A 6th order A-stable stiffly stable Rosenbrock method with a stiff-aware 5th order interpolant.", "- Steinebach G. Construction of Rosenbrock–Wanner method Rodas6P.\n to prepare, 2025", true)
]
@eval begin
@doc $(is_W ? rosenbrock_wolfbrandt_docstring(desc, String(Alg), references = refs, with_step_limiter = true) : rosenbrock_docstring(desc, String(Alg), references = refs, with_step_limiter = true)) struct $Alg{CS, AD, F, P, FDT, ST, CJ, StepLimiter, StageLimiter} <:
Expand Down
2 changes: 1 addition & 1 deletion lib/OrdinaryDiffEqRosenbrock/src/interp_func.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,6 @@ function SciMLBase.interp_summary(::Type{cacheType},
cacheType <:
Union{RosenbrockCombinedConstantCache,
RosenbrockCache}}
dense ? "specialized 4rd order \"free\" stiffness-aware interpolation" :
dense ? "specialized 4th (Rodas6P = 5th) order \"free\" stiffness-aware interpolation" :
"1st order linear"
end
5 changes: 3 additions & 2 deletions lib/OrdinaryDiffEqRosenbrock/src/rosenbrock_caches.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -753,9 +753,10 @@ tabtype(::Rodas5) = Rodas5Tableau
tabtype(::Rodas5P) = Rodas5PTableau
tabtype(::Rodas5Pr) = Rodas5PTableau
tabtype(::Rodas5Pe) = Rodas5PTableau
tabtype(::Rodas6P) = Rodas6PTableau

function alg_cache(
alg::Union{Rodas4, Rodas42, Rodas4P, Rodas4P2, Rodas5, Rodas5P, Rodas5Pe, Rodas5Pr},
alg::Union{Rodas4, Rodas42, Rodas4P, Rodas4P2, Rodas5, Rodas5P, Rodas5Pe, Rodas5Pr, Rodas6P},
u, rate_prototype, ::Type{uEltypeNoUnits},
::Type{uBottomEltypeNoUnits}, ::Type{tTypeNoUnits}, uprev, uprev2, f, t,
dt, reltol, p, calck,
Expand All @@ -772,7 +773,7 @@ function alg_cache(
end

function alg_cache(
alg::Union{Rodas4, Rodas42, Rodas4P, Rodas4P2, Rodas5, Rodas5P, Rodas5Pe, Rodas5Pr},
alg::Union{Rodas4, Rodas42, Rodas4P, Rodas4P2, Rodas5, Rodas5P, Rodas5Pe, Rodas5Pr, Rodas6P},
u, rate_prototype, ::Type{uEltypeNoUnits},
::Type{uBottomEltypeNoUnits}, ::Type{tTypeNoUnits}, uprev, uprev2, f, t,
dt, reltol, p, calck,
Expand Down
97 changes: 87 additions & 10 deletions lib/OrdinaryDiffEqRosenbrock/src/rosenbrock_interpolants.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -138,6 +138,8 @@ From MATLAB ODE Suite by Shampine
Θ1 = 1 - Θ
if !hasproperty(cache, :interp_order) || cache.interp_order == 2
@.. Θ1 * y₀ + Θ * (y₁ + Θ1 * (k[1] + Θ * k[2]))
elseif cache.interp_order == 4
@.. Θ1 * y₀ + Θ * (y₁ + Θ1 * (k[1] + Θ * (k[2] + Θ * (k[3] + Θ * k[4]))))
else
@.. Θ1 * y₀ + Θ * (y₁ + Θ1 * (k[1] + Θ * (k[2] + Θ * k[3])))
end
Expand All @@ -151,6 +153,9 @@ end
Θ1 = 1 - Θ
if !hasproperty(cache, :interp_order) || cache.interp_order == 2
@views @.. Θ1 * y₀[idxs] + Θ * (y₁[idxs] + Θ1 * (k[1][idxs] + Θ * k[2][idxs]))
elseif cache.interp_order == 4
@views @.. Θ1 * y₀[idxs] + Θ * (y₁[idxs] + Θ1 * (k[1][idxs] +
Θ * (k[2][idxs] + Θ * (k[3][idxs] + Θ * k[4][idxs]))))
else
@views @.. Θ1 * y₀[idxs] +
Θ * (y₁[idxs] + Θ1 * (k[1][idxs] + Θ * (k[2][idxs] + Θ * k[3][idxs])))
Expand All @@ -165,6 +170,8 @@ end
Θ1 = 1 - Θ
if !hasproperty(cache, :interp_order) || cache.interp_order == 2
@.. out = Θ1 * y₀ + Θ * (y₁ + Θ1 * (k[1] + Θ * k[2]))
elseif cache.interp_order == 4
@.. out = Θ1 * y₀ + Θ * (y₁ + Θ1 * (k[1] + Θ * (k[2] + Θ * (k[3] + Θ * k[4]))))
else
@.. out = Θ1 * y₀ + Θ * (y₁ + Θ1 * (k[1] + Θ * (k[2] + Θ * k[3])))
end
Expand All @@ -179,6 +186,9 @@ end
Θ1 = 1 - Θ
if !hasproperty(cache, :interp_order) || cache.interp_order == 2
@views @.. out = Θ1 * y₀[idxs] + Θ * (y₁[idxs] + Θ1 * (k[1][idxs] + Θ * k[2][idxs]))
elseif cache.interp_order == 4
@views @.. out = Θ1 * y₀[idxs] + Θ * (y₁[idxs] + Θ1 * (k[1][idxs] +
Θ * (k[2][idxs] + Θ * (k[3][idxs] + Θ * k[4][idxs]))))
else
@views @.. out = Θ1 * y₀[idxs] +
Θ * (y₁[idxs] +
Expand All @@ -199,6 +209,11 @@ end
idxs::Nothing, T::Type{Val{1}}, differential_vars)
if !hasproperty(cache, :interp_order) || cache.interp_order == 2
@.. (k[1] + Θ * (-2 * k[1] + 2 * k[2] - 3 * k[2] * Θ) - y₀ + y₁) / dt
elseif cache.interp_order == 4
@.. (k[1] + Θ * (-2 * k[1] + 2 * k[2] +
Θ * (-3 * k[2] + 3 * k[3] +
Θ * (-4 * k[3] + 4 * k[4] - 5 * Θ * k[4]))) -
y₀ + y₁) / dt
else
@.. (k[1] + Θ * (-2 * k[1] + 2 * k[2] +
Θ * (-3 * k[2] + 3 * k[3] - 4 * Θ * k[3])) -
Expand All @@ -214,6 +229,11 @@ end
@views @.. (k[1][idxs] +
Θ * (-2 * k[1][idxs] + 2 * k[2][idxs] - 3 * k[2][idxs] * Θ) -
y₀[idxs] + y₁[idxs]) / dt
elseif cache.interp_order == 4
@views @.. (k[1][idxs] + Θ * (-2 * k[1][idxs] + 2 * k[2][idxs] +
Θ * (-3 * k[2][idxs] + 3 * k[3][idxs] +
Θ * (-4 * k[3][idxs] + 4 * k[4][idxs] - 5 * Θ * k[4][idxs]))) -
y₀[idxs] + y₁[idxs]) / dt
else
@views @.. (k[1][idxs] +
Θ * (-2 * k[1][idxs] + 2 * k[2][idxs] +
Expand All @@ -229,6 +249,11 @@ end
idxs::Nothing, T::Type{Val{1}}, differential_vars)
if !hasproperty(cache, :interp_order) || cache.interp_order == 2
@.. out = (k[1] + Θ * (-2 * k[1] + 2 * k[2] - 3 * k[2] * Θ) - y₀ + y₁) / dt
elseif cache.interp_order == 4
@.. out = (k[1] + Θ * (-2 * k[1] + 2 * k[2] +
Θ * (-3 * k[2] + 3 * k[3] +
Θ * (-4 * k[3] + 4 * k[4] - 5 * Θ * k[4]))) -
y₀ + y₁) / dt
else
@.. out = (k[1] +
Θ * (-2 * k[1] + 2 * k[2] +
Expand All @@ -248,6 +273,11 @@ end
(-2 * k[1][idxs] + 2 * k[2][idxs] -
3 * k[2][idxs] * Θ) -
y₀[idxs] + y₁[idxs]) / dt
elseif cache.interp_order == 4
@views @.. out = (k[1][idxs] + Θ * (-2 * k[1][idxs] + 2 * k[2][idxs] +
Θ * (-3 * k[2][idxs] + 3 * k[3][idxs] +
Θ * (-4 * k[3][idxs] + 4 * k[4][idxs] - 5 * Θ * k[4][idxs]))) -
y₀[idxs] + y₁[idxs]) / dt
else
@views @.. broadcast=false out=(k[1][idxs] +
Θ * (-2 * k[1][idxs] + 2 * k[2][idxs] +
Expand All @@ -262,68 +292,115 @@ end
# Second Derivative
@muladd function _ode_interpolant(Θ, dt, y₀, y₁, k, cache::RosenbrockCombinedConstantCache,
idxs::Nothing, T::Type{Val{2}}, differential_vars)
@inbounds (-2 * k[1] + 2 * k[2] + Θ * (-6 * k[2] + 6 * k[3] - 12 * Θ * k[3])) / dt^2
if cache.interp_order == 4
@inbounds (-2 * k[1] + 2 * k[2] + Θ * (-6 * k[2] + 6 * k[3] +
Θ * (-12 * k[3] + 12 * k[4] - 20 * Θ * k[4]))) / dt^2
else
@inbounds (-2 * k[1] + 2 * k[2] + Θ * (-6 * k[2] + 6 * k[3] - 12 * Θ * k[3])) / dt^2
end
end

@muladd function _ode_interpolant(Θ, dt, y₀, y₁, k, cache::RosenbrockCache, idxs::Nothing,
T::Type{Val{2}}, differential_vars)
@inbounds @.. broadcast=false (-2 * k[1] + 2 * k[2] +
if cache.interp_order == 4
@inbounds @.. broadcast=false (-2 * k[1] + 2 * k[2] + Θ * (-6 * k[2] + 6 * k[3] +
Θ * (-12 * k[3] + 12 * k[4] - 20 * Θ * k[4]))) / dt^2
else
@inbounds @.. broadcast=false (-2 * k[1] + 2 * k[2] +
Θ * (-6 * k[2] + 6 * k[3] - 12 * Θ * k[3]))/dt^2
end
end

@muladd function _ode_interpolant(Θ, dt, y₀, y₁, k,
cache::Union{RosenbrockCombinedConstantCache, RosenbrockCache},
idxs, T::Type{Val{2}}, differential_vars)
@.. broadcast=false (-2 * k[1][idxs] + 2 * k[2][idxs] +
if cache.interp_order == 4
@.. broadcast=false (-2 * k[1][idxs] + 2 * k[2][idxs] + Θ * (-6 * k[2][idxs] + 6 * k[3][idxs] +
Θ * (-12 * k[3][idxs] + 12 * k[4][idxs] - 20 * Θ * k[4][idxs]))) / dt^2
else
@.. broadcast=false (-2 * k[1][idxs] + 2 * k[2][idxs] +
Θ * (-6 * k[2][idxs] + 6 * k[3][idxs] - 12 * Θ * k[3][idxs]))/dt^2
end
end

@muladd function _ode_interpolant!(out, Θ, dt, y₀, y₁, k,
cache::Union{RosenbrockCombinedConstantCache, RosenbrockCache},
idxs::Nothing, T::Type{Val{2}}, differential_vars)
@.. broadcast=false out=(-2 * k[1] + 2 * k[2] +
if cache.interp_order == 4
@.. broadcast=false out=(-2 * k[1] + 2 * k[2] + Θ * (-6 * k[2] + 6 * k[3] +
Θ * (-12 * k[3] + 12 * k[4] - 20 * Θ * k[4]))) / dt^2
else
@.. broadcast=false out=(-2 * k[1] + 2 * k[2] +
Θ * (-6 * k[2] + 6 * k[3] - 12 * Θ * k[3])) / dt^2
end
out
end

@muladd function _ode_interpolant!(out, Θ, dt, y₀, y₁, k,
cache::Union{RosenbrockCombinedConstantCache, RosenbrockCache},
idxs, T::Type{Val{2}}, differential_vars)
@views @.. broadcast=false out=(-2 * k[1][idxs] + 2 * k[2][idxs] +
if cache.interp_order == 4
@views @.. broadcast=false out=(-2 * k[1][idxs] + 2 * k[2][idxs] + Θ * (-6 * k[2][idxs] + 6 * k[3][idxs] +
Θ * (-12 * k[3][idxs] + 12 * k[4][idxs] - 20 * Θ * k[4][idxs]))) / dt^2
else
@views @.. broadcast=false out=(-2 * k[1][idxs] + 2 * k[2][idxs] +
Θ *
(-6 * k[2][idxs] + 6 * k[3][idxs] - 12 * Θ * k[3][idxs])) /
dt^2
end
out
end

# Third Derivative
@muladd function _ode_interpolant(Θ, dt, y₀, y₁, k, cache::RosenbrockCombinedConstantCache,
idxs::Nothing, T::Type{Val{3}}, differential_vars)
@inbounds (-6 * k[2] + 6 * k[3] - 24 * Θ * k[3]) / dt^3
if cache.interp_order == 4
@inbounds (-6 * k[2] + 6 * k[3] + Θ * (-24 * k[3] + 24 * k[4] - 60 * Θ * k[4])) / dt^3
else
@inbounds (-6 * k[2] + 6 * k[3] - 24 * Θ * k[3]) / dt^3
end
end

@muladd function _ode_interpolant(Θ, dt, y₀, y₁, k, cache::RosenbrockCache, idxs::Nothing,
T::Type{Val{3}}, differential_vars)
@inbounds @.. broadcast=false (-6 * k[2] + 6 * k[3] - 24 * Θ * k[3])/dt^3
if cache.interp_order == 4
@inbounds @.. broadcast=false (-6 * k[2] + 6 * k[3] + Θ * (-24 * k[3] + 24 * k[4] - 60 * Θ * k[4])) / dt^3
else
@inbounds @.. broadcast=false (-6 * k[2] + 6 * k[3] - 24 * Θ * k[3])/dt^3
end
end

@muladd function _ode_interpolant(Θ, dt, y₀, y₁, k,
cache::Union{RosenbrockCombinedConstantCache, RosenbrockCache},
idxs, T::Type{Val{3}}, differential_vars)
@.. broadcast=false (-6 * k[2][idxs] + 6 * k[3][idxs] - 24 * Θ * k[3][idxs])/dt^3
if cache.interp_order == 4
@.. broadcast=false (-6 * k[2][idxs] + 6 * k[3][idxs] +
Θ * (-24 * k[3][idxs] + 24 * k[4][idxs] - 60 * Θ * k[4][idxs])) / dt^3
else
@.. broadcast=false (-6 * k[2][idxs] + 6 * k[3][idxs] - 24 * Θ * k[3][idxs])/dt^3
end
end

@muladd function _ode_interpolant!(out, Θ, dt, y₀, y₁, k,
cache::Union{RosenbrockCombinedConstantCache, RosenbrockCache},
idxs::Nothing, T::Type{Val{3}}, differential_vars)
@.. broadcast=false out=(-6 * k[2] + 6 * k[3] - 24 * Θ * k[3]) / dt^3
if cache.interp_order == 4
@.. broadcast=false out=(-6 * k[2] + 6 * k[3] + Θ * (-24 * k[3] + 24 * k[4] - 60 * Θ * k[4])) / dt^3
else
@.. broadcast=false out=(-6 * k[2] + 6 * k[3] - 24 * Θ * k[3]) / dt^3
end
out
end

@muladd function _ode_interpolant!(out, Θ, dt, y₀, y₁, k,
cache::Union{RosenbrockCombinedConstantCache, RosenbrockCache},
idxs, T::Type{Val{3}}, differential_vars)
@views @.. broadcast=false out=(-6 * k[2][idxs] + 6 * k[3][idxs] -
if cache.interp_order == 4
@views @.. broadcast=false out=(-6 * k[2][idxs] + 6 * k[3][idxs] +
Θ * (-24 * k[3][idxs] + 24 * k[4][idxs] - 60 * Θ * k[4][idxs])) / dt^3
else
@views @.. broadcast=false out=(-6 * k[2][idxs] + 6 * k[3][idxs] -
24 * Θ * k[3][idxs]) / dt^3
end
out
end
38 changes: 31 additions & 7 deletions lib/OrdinaryDiffEqRosenbrock/src/rosenbrock_perform_step.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1240,7 +1240,7 @@ end
linsolve_tmp = @.. du + dtd[1] * dT
k1 = _reshape(W \ -_vec(linsolve_tmp), axes(uprev))
# constant number for type stability make sure this is greater than num_stages
ks = ntuple(Returns(k1), 10)
ks = ntuple(Returns(k1), 20)
# Loop for stages
for stage in 2:num_stages
u = uprev
Expand Down Expand Up @@ -1268,11 +1268,26 @@ end
ks = Base.setindex(ks, _reshape(W \ -_vec(linsolve_tmp), axes(uprev)), stage)
integrator.stats.nsolve += 1
end
#@show ks
u = u .+ ks[num_stages]
if (integrator.alg isa Rodas6P)
du = ks[16]
u = uprev
for i in 1:15
u = @.. u + A[16, i] * ks[i]
end
u = u .+ ks[16]
else
du = ks[num_stages]
u = u .+ ks[num_stages]
end

if integrator.opts.adaptive
atmp = calculate_residuals(ks[num_stages], uprev, u, integrator.opts.abstol,
if (integrator.alg isa Rodas5Pe)
du = 0.2606326497975715 * ks[1] - 0.005158627295444251 * ks[2] +
1.3038988631109731 * ks[3] + 1.235000722062074 * ks[4] +
-0.7931985603795049 * ks[5] - 1.005448461135913 * ks[6] -
0.18044626132120234 * ks[7] + 0.17051519239113755 * ks[8]
end
atmp = calculate_residuals(du, uprev, u, integrator.opts.abstol,
integrator.opts.reltol, integrator.opts.internalnorm, t)
integrator.EEst = integrator.opts.internalnorm(atmp, t)
end
Expand Down Expand Up @@ -1381,8 +1396,17 @@ end
@.. $(_vec(ks[stage])) = -linres.u
integrator.stats.nsolve += 1
end
du .= ks[end]
u .+= ks[end]
if (integrator.alg isa Rodas6P)
du .= ks[16]
u .= uprev
for i in 1:15
@.. u += A[16, i] * ks[i]
end
u .+= ks[16]
else
du .= ks[end]
u .+= ks[end]
end

step_limiter!(u, integrator, p, t + dt)

Expand All @@ -1393,7 +1417,7 @@ end
-0.7931985603795049 * ks[5] - 1.005448461135913 * ks[6] -
0.18044626132120234 * ks[7] + 0.17051519239113755 * ks[8]
end
calculate_residuals!(atmp, ks[end], uprev, u, integrator.opts.abstol,
calculate_residuals!(atmp, du, uprev, u, integrator.opts.abstol,
integrator.opts.reltol, integrator.opts.internalnorm, t)
integrator.EEst = integrator.opts.internalnorm(atmp, t)
end
Expand Down
Loading
Loading