Using Tau-leaping
bug
It seems that tau-leaping methods are currently broken? I might be defining something incorrectly, though the errors are a bit cryptic... Using Catalyst v16.1.1 and JumpProcesses v9.26.0 in Julia v1.12.6.
```julia
using JumpProcesses
using Catalyst
using StochasticDiffEq
sir_model = @reaction_network begin
β, S + I --> 2I
ν, I --> R
end
p = (:β => 0.1 / 1000, :ν => 0.01)
u₀ = [:S => 990, :I => 10, :R => 0]
tspan = (0.0, 250.0)
jump_prob = JumpProblem(sir_model, u₀, tspan, p, aggregator=PureLeaping())
sol = solve(jump_prob, TauLeaping())
```
```julia
ERROR: MethodError: no method matching zero(::Nothing)
The function `zero` exists, but no method is defined for this combination of argument types.
Closest candidates are:
zero(::Type{Union{}}, Any...)
@ Base number.jl:315
zero(::Type{Missing})
@ Base missing.jl:104
zero(::Type{Dates.Time})
@ Dates ~/.julia/juliaup/julia-1.12.6+0.aarch64.apple.darwin14/Julia-1.12.app/Contents/Resources/julia/share/julia/stdlib/v1.12/Dates/src/types.jl:460
...
Stacktrace:
[1] alg_cache(alg::TauLeaping, prob::DiscreteProblem{…}, u::Vector{…}, ΔW::Nothing, ΔZ::Nothing, p::MTKParameters{…}, rate_prototype::Vector{…}, noise_rate_prototype::Nothing, jump_rate_prototype::Nothing, ::Type{…}, ::Type{…}, ::Type{…}, uprev::Vector{…}, f::DiscreteFunction{…}, t::Float64, dt::Float64, ::Type{…}, verbose::DEVerbosity{…})
@ StochasticDiffEqLeaping ~/.julia/packages/StochasticDiffEqLeaping/YRM93/src/caches/tau_caches.jl:27
[2] _sde_init(_prob::JumpProblem{…}, alg::TauLeaping; saveat::Tuple{}, tstops::Tuple{}, d_discontinuities::Tuple{}, save_idxs::Nothing, save_everystep::Bool, save_noise::Bool, save_on::Bool, save_start::Bool, save_end::Nothing, callback::Nothing, dense::Bool, calck::Bool, dt::Float64, adaptive::Bool, gamma::Rational{…}, abstol::Nothing, reltol::Nothing, qmin::Rational{…}, qmax::Rational{…}, qsteady_min::Int64, qsteady_max::Int64, beta2::Nothing, beta1::Nothing, qoldinit::Rational{…}, controller::Nothing, fullnormalize::Bool, failfactor::Int64, delta::Rational{…}, maxiters::Int64, dtmax::Float64, dtmin::Float64, internalnorm::typeof(DiffEqBase.ODE_DEFAULT_NORM), isoutofdomain::typeof(DiffEqBase.ODE_DEFAULT_ISOUTOFDOMAIN), unstable_check::typeof(DiffEqBase.ODE_DEFAULT_UNSTABLE_CHECK), verbose::Bool, force_dtmin::Bool, timeseries_errors::Bool, dense_errors::Bool, advance_to_tstop::Bool, stop_at_next_tstop::Bool, initialize_save::Bool, progress::Bool, progress_steps::Int64, progress_name::String, progress_message::typeof(DiffEqBase.ODE_DEFAULT_PROG_MESSAGE), progress_id::Symbol, userdata::Nothing, save_discretes::Bool, initialize_integrator::Bool, seed::UInt64, rng::Nothing, alias::Nothing, initializealg::DefaultInit, kwargs::@Kwargs{})
@ StochasticDiffEqCore ~/.julia/packages/StochasticDiffEqCore/qCsWr/src/solve.jl:521
[3] __init(_prob::JumpProblem{…}, alg::TauLeaping; kwargs::@Kwargs{…})
@ StochasticDiffEqCore ~/.julia/packages/StochasticDiffEqCore/qCsWr/src/solve.jl:97
[4] __solve(prob::JumpProblem{…}, alg::TauLeaping; merge_callbacks::Bool, kwargs::@Kwargs{})
@ StochasticDiffEqCore ~/.julia/packages/StochasticDiffEqCore/qCsWr/src/solve.jl:33
[5] __solve(prob::JumpProblem{…}, alg::TauLeaping)
@ StochasticDiffEqCore ~/.julia/packages/StochasticDiffEqCore/qCsWr/src/solve.jl:27
[6] solve(prob::JumpProblem{…}, args::TauLeaping; kwargs::@Kwargs{})
@ DiffEqBase ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:656
[7] solve(prob::JumpProblem{…}, args::TauLeaping)
@ DiffEqBase ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:655
[8] top-level scope
@ ~/Projects/draft/src/tau_leaping_test.jl:19
```
I get a different error trying to use `CaoTauLeaping` (which should still be compatible?):
```julia
sol = solve(jump_prob, CaoTauLeaping())
```
```julia
ERROR: FieldError: type Nothing has no field `dW`; Nothing has no fields at all.
Stacktrace:
[1] getproperty
@ ./Base_compiler.jl:54 [inlined]
[2] perform_step!(integrator::OrdinaryDiffEqCore.ODEIntegrator{…}, cache::StochasticDiffEqLeaping.TauLeapingCache{…})
@ StochasticDiffEqLeaping ~/.julia/packages/StochasticDiffEqLeaping/YRM93/src/perform_step/tau_leaping.jl:25
[3] solve!(integrator::OrdinaryDiffEqCore.ODEIntegrator{…})
@ OrdinaryDiffEqCore ~/.julia/packages/OrdinaryDiffEqCore/rnOL4/src/solve.jl:889
[4] __solve(prob::JumpProblem{…}, alg::CaoTauLeaping; merge_callbacks::Bool, kwargs::@Kwargs{})
@ StochasticDiffEqCore ~/.julia/packages/StochasticDiffEqCore/qCsWr/src/solve.jl:34
[5] __solve
@ ~/.julia/packages/StochasticDiffEqCore/qCsWr/src/solve.jl:27 [inlined]
[6] #solve#31
@ ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:656 [inlined]
[7] solve(prob::JumpProblem{…}, args::CaoTauLeaping)
@ DiffEqBase ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:655
[8] top-level scope
@ ~/Projects/draft/src/tau_leaping_test.jl:20
```
On a different note, trying to use `SSAStepper` also leads to an error, which makes sense given the `PureLeaping` aggregator, but isn't clear from the error message:
```julia
sol = solve(jump_prob, SSAStepper())
```
```julia
ERROR: BoundsError: attempt to access Tuple{} at index [0]
Stacktrace:
[1] getindex(t::Tuple, i::Int64)
@ Base ./tuple.jl:33
[2] __init(jump_prob::JumpProblem{…}, alg::SSAStepper; save_start::Bool, save_end::Bool, seed::Nothing, alias_jump::Bool, saveat::Nothing, callback::Nothing, tstops::Nothing, numsteps_hint::Int64)
@ JumpProcesses ~/.julia/packages/JumpProcesses/sLhTC/src/SSA_stepper.jl:0
[3] init_call(_prob::JumpProblem{…}, args::SSAStepper; merge_callbacks::Bool, kwargshandle::Nothing, kwargs::@Kwargs{})
@ DiffEqBase ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:87
[4] init_call(_prob::JumpProblem{…}, args::SSAStepper)
@ DiffEqBase ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:67
[5] init(prob::JumpProblem{…}, args::SSAStepper; kwargs::@Kwargs{})
@ DiffEqBase ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:106
[6] init(prob::JumpProblem{…}, args::SSAStepper)
@ DiffEqBase ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:105
[7] __solve(jump_prob::JumpProblem{…}, alg::SSAStepper; kwargs::@Kwargs{})
@ JumpProcesses ~/.julia/packages/JumpProcesses/sLhTC/src/SSA_stepper.jl:131
[8] __solve(jump_prob::JumpProblem{…}, alg::SSAStepper)
@ JumpProcesses ~/.julia/packages/JumpProcesses/sLhTC/src/SSA_stepper.jl:129
[9] solve(prob::JumpProblem{…}, args::SSAStepper; kwargs::@Kwargs{})
@ DiffEqBase ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:656
[10] solve(prob::JumpProblem{…}, args::SSAStepper)
@ DiffEqBase ~/.julia/packages/DiffEqBase/CJTzr/src/solve.jl:655
[11] top-level scope
@ ~/Projects/draft/src/tau_leaping_test.jl:21
```
4 条评论