Skip to content

Commit

Permalink
Merge pull request #2058 from SciML/dw/rodas5
Browse files Browse the repository at this point in the history
Fix Rodas5 and Rodas5P
  • Loading branch information
ChrisRackauckas authored Nov 14, 2023
2 parents 11a1a20 + 16764b1 commit 6b1437f
Show file tree
Hide file tree
Showing 2 changed files with 16 additions and 0 deletions.
6 changes: 6 additions & 0 deletions src/dense/stiff_addsteps.jl
Original file line number Diff line number Diff line change
Expand Up @@ -838,6 +838,9 @@ function _ode_addsteps!(k, t, uprev, u, dt, f, p, cache::Rosenbrock5Cache,
veck8 = _vec(k8)
@.. broadcast=false veck8=-vecu

# https://github.com/SciML/OrdinaryDiffEq.jl/issues/2055
tmp = linsolve_tmp

@unpack h21, h22, h23, h24, h25, h26, h27, h28, h31, h32, h33, h34, h35, h36, h37, h38, h41, h42, h43, h44, h45, h46, h47, h48 = cache.tab
@.. broadcast=false tmp=h21 * k1 + h22 * k2 + h23 * k3 + h24 * k4 + h25 * k5 +
h26 * k6 + h27 * k7 + h28 * k8
Expand Down Expand Up @@ -1092,6 +1095,9 @@ function _ode_addsteps!(k, t, uprev, u, dt, f, p, cache::Rosenbrock5Cache{<:Arra

@unpack h21, h22, h23, h24, h25, h26, h27, h28, h31, h32, h33, h34, h35, h36, h37, h38, h41, h42, h43, h44, h45, h46, h47, h48 = cache.tab

# https://github.com/SciML/OrdinaryDiffEq.jl/issues/2055
tmp = linsolve_tmp

@inbounds @simd ivdep for i in eachindex(u)
tmp[i] = h21 * k1[i] + h22 * k2[i] + h23 * k3[i] + h24 * k4[i] + h25 * k5[i] +
h26 * k6[i] + h27 * k7[i] + h28 * k8[i]
Expand Down
10 changes: 10 additions & 0 deletions test/integrators/event_detection_tests.jl
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
using StaticArrays
using OrdinaryDiffEq
using DiffEqDevTools
using Test

@inbounds @inline function ż(z, p, t)
Expand Down Expand Up @@ -60,6 +61,15 @@ prob = ODEProblem(f, u0, tspan, p)
sol = solve(prob, Tsit5(), callback = cb2)
@test minimum(Array(sol)) > -40

# https://github.com/SciML/OrdinaryDiffEq.jl/issues/2055
for alg in (Rodas5(), Rodas5P())
sol2 = solve(prob, alg; callback = cb2)
sol3 = appxtrue(sol, sol2)
@test sol3.errors[:L2] < 1e-5
@test sol3.errors[:L∞] < 5e-5
@test sol3.errors[:final] < 1e-5
end

function fball(du, u, p, t)
du[1] = u[2]
du[2] = -p
Expand Down

0 comments on commit 6b1437f

Please sign in to comment.