SciML / SciML/JumpProcesses.jl

Simulating jump problems require non-integer tspan values

Open
#359 3 comments 0 reactions 0 assignees View on GitHub

Nobody has claimed this yet.

Dominant language
Julia
Stars
150
Forks
41
Avg merge
1d 9h
Merged PRs (30d)
28

Description

Simple example:

using Catalyst, JumpProcesses

rn = @reaction_network begin
    (p,d), 0 <--> X
end
u0 = [:X => 1.0]
tspan = (0.0, 10.0)
p = [:p => 1.0, :d => 0.1]

dprob = DiscreteProblem(rn, u0, tspan, p)
jprob = JumpProblem(rn, dprob, Direct())
sol = solve(jprob, SSAStepper())

If I change to tspan = (0, 10):

using Catalyst, JumpProcesses

rn = @reaction_network begin
    (p,d), 0 <--> X
end
u0 = [:X => 1.0]
tspan = (0, 10)
p = [:p => 1.0, :d => 0.1]

dprob = DiscreteProblem(rn, u0, tspan, p)
jprob = JumpProblem(rn, dprob, Direct())
sol = solve(jprob, SSAStepper())

I instead get an error:

ERROR: InexactError: Int64(1.1)
Stacktrace:

Int64 at [float.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

convert at [number.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

setindex! at [array.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

time_to_next_jump(p::JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, u::Vector{Float64}, params::Vector{Float64}, t::Int64) at [direct.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

generate_jumps!(p::JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, integrator::JumpProcesses.SSAIntegrator{DiscreteFunction{true, true, DiscreteFunction{true, SciMLBase.FullSpecialize, SciMLBase.var"#213#214", Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Nothing, Vector{Symbol}, Symbol, Vector{Symbol}, ModelingToolkit.var"#generated_observed#685"{Bool, JumpSystem{RecursiveArrayTools.ArrayPartition{Any, Tuple{Vector{MassActionJump}, Vector{ConstantRateJump}, Vector{VariableRateJump}}}}, Dict{Any, Any}}, JumpSystem{RecursiveArrayTools.ArrayPartition{Any, Tuple{Vector{MassActionJump}, Vector{ConstantRateJump}, Vector{VariableRateJump}}}}}, Vector{Float64}, Int64, Int64, Vector{Float64}, ODESolution{Float64, 2, Vector{Vector{Float64}}, Nothing, Nothing, Vector{Int64}, Nothing, DiscreteProblem{Vector{Float64}, Tuple{Int64, Int64}, true, Vector{Float64}, DiscreteFunction{true, true, DiscreteFunction{true, SciMLBase.FullSpecialize, SciMLBase.var"#213#214", Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Nothing, Vector{Symbol}, Symbol, Vector{Symbol}, ModelingToolkit.var"#generated_observed#685"{Bool, JumpSystem{RecursiveArrayTools.ArrayPartition{Any, Tuple{Vector{MassActionJump}, Vector{ConstantRateJump}, Vector{VariableRateJump}}}}, Dict{Any, Any}}, JumpSystem{RecursiveArrayTools.ArrayPartition{Any, Tuple{Vector{MassActionJump}, Vector{ConstantRateJump}, Vector{VariableRateJump}}}}}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}}, SSAStepper, SciMLBase.ConstantInterpolation{Vector{Int64}, Vector{Vector{Float64}}}, DiffEqBase.Stats, Nothing}, DiscreteCallback{JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, typeof(SciMLBase.FINALIZE_DEFAULT)}, Nothing, NamedTuple{(:callback,), Tuple{CallbackSet{Tuple{}, Tuple{}}}}, Vector{Int64}}, u::Vector{Float64}, params::Vector{Float64}, t::Int64) at [direct.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

initialize! at [direct.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

AbstractSSAJumpAggregator at [ssajump.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

__init(jump_prob::JumpProblem{true, DiscreteProblem{Vector{Float64}, Tuple{Int64, Int64}, true, Vector{Float64}, DiscreteFunction{true, true, DiscreteFunction{true, SciMLBase.FullSpecialize, SciMLBase.var"#213#214", Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Nothing, Vector{Symbol}, Symbol, Vector{Symbol}, ModelingToolkit.var"#generated_observed#685"{Bool, JumpSystem{RecursiveArrayTools.ArrayPartition{Any, Tuple{Vector{MassActionJump}, Vector{ConstantRateJump}, Vector{VariableRateJump}}}}, Dict{Any, Any}}, JumpSystem{RecursiveArrayTools.ArrayPartition{Any, Tuple{Vector{MassActionJump}, Vector{ConstantRateJump}, Vector{VariableRateJump}}}}}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}}, Direct, CallbackSet{Tuple{}, Tuple{DiscreteCallback{JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, typeof(SciMLBase.FINALIZE_DEFAULT)}}}, JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, Tuple{}, Nothing, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Random.TaskLocalRNG, Base.Pairs{Symbol, Nothing, Tuple{Symbol}, NamedTuple{(:callback,), Tuple{Nothing}}}}, alg::SSAStepper; save_start::Bool, save_end::Bool, seed::Nothing, alias_jump::Bool, saveat::Nothing, callback::Nothing, tstops::Vector{Int64}, numsteps_hint::Int64) at [SSA_stepper.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

__init at [SSA_stepper.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

#init_call#30 at [solve.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

init_call at [solve.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

#init#32 at [solve.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

init at [solve.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

#__solve#221 at [SSA_stepper.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

__solve at [SSA_stepper.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

#solve#47 at [solve.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

solve(prob::JumpProblem{true, DiscreteProblem{Vector{Float64}, Tuple{Int64, Int64}, true, Vector{Float64}, DiscreteFunction{true, true, DiscreteFunction{true, SciMLBase.FullSpecialize, SciMLBase.var"#213#214", Nothing, Nothing, Nothing, Nothing, typeof(SciMLBase.DEFAULT_OBSERVED), Nothing}, Nothing, Vector{Symbol}, Symbol, Vector{Symbol}, ModelingToolkit.var"#generated_observed#685"{Bool, JumpSystem{RecursiveArrayTools.ArrayPartition{Any, Tuple{Vector{MassActionJump}, Vector{ConstantRateJump}, Vector{VariableRateJump}}}}, Dict{Any, Any}}, JumpSystem{RecursiveArrayTools.ArrayPartition{Any, Tuple{Vector{MassActionJump}, Vector{ConstantRateJump}, Vector{VariableRateJump}}}}}, Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}}}, Direct, CallbackSet{Tuple{}, Tuple{DiscreteCallback{JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, typeof(SciMLBase.FINALIZE_DEFAULT)}}}, JumpProcesses.DirectJumpAggregation{Int64, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Tuple{}, Tuple{}, Random.TaskLocalRNG}, Tuple{}, Nothing, MassActionJump{Vector{Float64}, Vector{Vector{Pair{Int64, Int64}}}, Vector{Vector{Pair{Int64, Int64}}}, ModelingToolkit.JumpSysMajParamMapper{Vector{Num}, Vector{SymbolicUtils.BasicSymbolic{Real}}, Float64}}, Random.TaskLocalRNG, Base.Pairs{Symbol, Nothing, Tuple{Symbol}, NamedTuple{(:callback,), Tuple{Nothing}}}}, args::SSAStepper) at [solve.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

top-level scope at [playground.jl](vscode-file://vscode-app/snap/code/140/usr/share/code/resources/app/out/vs/code/electron-sandbox/workbench/workbench.html)

The problem also happens when Catalyst is not used, e.g. this also give a similar error

begin
    using JumpProcesses
    # here we order S = 1, I = 2, and R = 3
    # substrate stoichiometry:
    substoich = [[1 => 1, 2 => 1],    # 1*S + 1*I
        [2 => 1]]            # 1*I
    # net change by each jump type
    netstoich = [[1 => -1, 2 => 1],   # S -> S-1, I -> I+1
        [2 => -1, 3 => 1]]   # I -> I-1, R -> R+1
    # rate constants for each jump
    p = (0.1 / 1000, 0.01)

    # p[1] is rate for S+I --> 2I, p[2] for I --> R
    pidxs = [1, 2]

    maj = MassActionJump(substoich, netstoich; param_idxs = pidxs)

    u₀ = [999, 1, 0]       #[S(0),I(0),R(0)]
    tspan = (0, 250)
    dprob = DiscreteProblem(u₀, tspan, p)

    # use the Direct method to simulate
    jprob = JumpProblem(dprob, Direct(), maj)

    # solve as a pure jump process, i.e. using SSAStepper
    sol = solve(jprob, SSAStepper())
end

This works fine with e.g. ODEs:

rn = @reaction_network begin
    (p,d), 0 <--> X
end
u0 = [:X => 1.0]
tspan = (0, 10)
p = [:p => 1.0, :d => 0.1]

oprob = ODEProblem(rn, u0, tspan, p)
sol = solve(oprob, Tsit5())

Contributor guide

Open the contributing guide

First steps

  1. Read the whole issue, then the project's contributing guide.
  2. Comment on the issue to say you are picking it up — it saves two people doing the same work.
  3. Fork the repository and make your change on a branch.
  4. Open a pull request that references the issue number.

Research direction

Reproduce the reported difference between integer and floating-point tspan values using the JumpProcesses example. Start with time_to_next_jump in direct.jl and the SSA initialization path in SSA_stepper.jl; trace where the integer time type is selected and verify that both examples run without an InexactError.

Written by the indexing model from the issue text.

Assessment

Tech stack
julia
Domain
backend
Issue type
Bug
Difficulty
4/5
Estimated time
3-5 days
Activity status
Stale
Clarity
Needs clarification
Newbie friendliness
30/100

Get new issues in your inbox

A short digest of beginner-friendly GitHub issues.