ERROR: UndefVarError: a71 not defined when using Rodas5P with ContinuousCallback

[vancleve]

Stacktrace starts here:

OrdinaryDiffEq.jl/src/dense/stiff_addsteps.jl

Line 1173 in 7f0eb18

@.. broadcast=false tmp=uprev + a71 * k1 + a72 * k2 + a73 * k3 + a74 * k4 +

Minimish working example

using ComponentArrays
using Parameters: @unpack
using DifferentialEquations

default_params = ComponentVector(
    rb = 1.0,
    μB = 1.0,
    βgi = 0.0025,
    βvi = 0.0025,
    γ = 1.5,
    σ = 1.5,
    ϵ = 0.1,
    Φ = 0.5,
    ω_D = 0.001,
    ωI = 0.0001,
    θ = 0.03,
    D_threshold = 0.85,
    ρ = 1,
    η = -0.005,
    Bx = 1000,
    s = 0.5,
    xe = 0.5,
    Kb = 50000,
    KI = 1000
)

default_u0 = ComponentVector(
    B_g = 20,
    B_v = 1,
    V_D = 1,
    V_I = 1,
    I = 1,
    D = 0.001
)

function virulence_ode!(du, u, p, t)
    @unpack B_g, B_v, V_D, V_I, I, D = u
    @unpack rb, μB, βgi, βvi, γ, σ, ϵ, Φ, ω_D, ωI, θ, ρ, η, Bx, s, xe, Kb, KI = p
    
    du.B_g = rb * B_g * ( 1 - ((B_g + B_v) / Kb)) - (βgi * I * B_g) - μB * B_g * (1 / (1 + exp(η * (B_g - Bx)))) ## dB_g/dt
    du.B_v = rb * exp(-ρ * γ) * B_v * ( 1 - ((B_g + B_v) / Kb)) - (βvi * I * B_v) + μB * B_g * (1 / (1 + exp(η * (B_g - Bx)))) ## dB_v / dt
    du.V_D = (1-s)*γ*B_v - ϵ* V_D
    du.V_I = s*γ*B_v - ϵ*V_I - V_I*I
    du.I = Φ * (B_g + B_v) * (1-(I/KI)) - (σ*V_I * I) - (xe * βgi*I*B_g) - (xe * βvi*I*B_v)  ## dI / dt
    du.D = (ω_D * V_D * (1 - D)) + (ωI * I * (1 - D)) - (θ * D)## dD / dt
end

tspan = (0.0, 50)
prob = ODEProblem(virulence_ode!, default_u0, tspan, ComponentVector(default_params, s=0.999, ρ=0.1, xe=1, Bx=15000))
condition(u, t, integrator) = u.D - default_params.D_threshold
affect!(integrator) = terminate!(integrator)
cb = ContinuousCallback(condition, affect!)
sol = solve(prob, Rodas5P(), callback = cb)

macOS 14.5
Julia 1.10.4
DifferentialEquations v7.13.0

[oscardssmith]

yeah, this looks like a simple typo.

[vancleve]

I also get the error with alg=Rodas5() fwiw.

[gstein3m]

In line 1173 of stiff_addsteps instead of

        @.. broadcast=false tmp=uprev + a71 * k1 + a72 * k2 + a73 * k3 + a74 * k4 +
                                a75 * k5 + a76 * k6 + k7

it must be

    @.. broadcast=false tmp .+=  k7

[oscardssmith]

I believe that is correct. Fix incoming.