Skip to content

Commit ac13d59

Browse files
committed
Another test
1 parent f46663d commit ac13d59

File tree

3 files changed

+53
-3
lines changed

3 files changed

+53
-3
lines changed

src/systems/diffeqs/sdesystem.jl

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -558,7 +558,7 @@ function DiffEqBase.SDEProblem{iip}(sys::SDESystem, u0map = [], tspan = get_tspa
558558
sparsenoise === nothing && (sparsenoise = get(kwargs, :sparse, false))
559559

560560
noiseeqs = get_noiseeqs(sys)
561-
if noiseeqs isa AbstractVector
561+
if noiseeqs isa AbstractVector || size(noiseeqs, 2) == 1
562562
noise_rate_prototype = nothing
563563
elseif sparsenoise
564564
I, J, V = findnz(SparseArrays.sparse(noiseeqs))

src/systems/systems.jl

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -77,11 +77,12 @@ function structural_simplify(sys::AbstractSystem, io = nothing; simplify = false
7777
_iszero(dvar) && break
7878
g_row = get(dvar2eq, dvar, 0)
7979
iszero(g_row) && error("$dvar isn't handled.")
80+
g_row > size(g, 1) && continue
8081
@views copyto!(sorted_g_rows[i, :], g[g_row, :])
8182
end
8283

83-
return SDESystem(equations(ode_sys), sorted_g_rows, get_iv(ode_sys),
84-
states(ode_sys), parameters(ode_sys);
84+
return SDESystem(full_equations(ode_sys), sorted_g_rows,
85+
get_iv(ode_sys), states(ode_sys), parameters(ode_sys);
8586
name = nameof(ode_sys))
8687
end
8788
end

test/sdesystem.jl

Lines changed: 49 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -614,3 +614,52 @@ sys2 = SDESystem(drift_eqs, diffusion_eqs, t, sts, ps, name = :sys1)
614614
prob = SDEProblem(sys1, sts .=> [1.0, 0.0, 0.0],
615615
(0.0, 100.0), ps .=> (10.0, 26.0))
616616
@test_nowarn solve(prob, LambaEulerHeun(), seed = 1)
617+
618+
@variables t
619+
D = Differential(t)
620+
sts = @variables S(t) E(t) I(t) R(t) D1(t) D2(t) β(t) logβ(t)
621+
ps = @parameters σ ρ γ λ N
622+
@brownian α
623+
eqs = Equation[
624+
β ~ logβ #exp(logβ)
625+
D(logβ) ~ sqrt(0.2) * α
626+
D(S) ~ -β * S * I / N
627+
E ~ N - (S + I + R + D1 + D2)
628+
D(I) ~ σ * E - γ * I
629+
D(R) ~ (1 - ρ) * γ * I
630+
D(D1) ~ ρ * γ * I - λ * D1
631+
D(D2) ~ λ * D1]
632+
@named mecha_bayes = System(eqs, t)
633+
mecha_bayes = structural_simplify(mecha_bayes)
634+
dE = 0.4
635+
dI = 2.0
636+
Rh = 3.0
637+
NN = 1
638+
using Random
639+
Random.seed!(1)
640+
641+
E0 = rand(Uniform(0, 0.02NN))
642+
I0 = rand(Uniform(0, 0.02NN))
643+
R0 = rand(Uniform(0, 0.02NN))
644+
D10 = rand(Uniform(0, 0.02NN))
645+
D20 = rand(Uniform(0, 0.02NN))
646+
S0 = 1 - (E0 + I0 + R0 + D10 + D20)
647+
u0 = [
648+
S => S0
649+
E => E0
650+
I => I0
651+
R => R0
652+
D1 => D10
653+
D2 => D20
654+
β => 0.0
655+
logβ => log(rand(Gamma(1, dI / Rh)))
656+
]
657+
658+
p = [N => NN,
659+
σ => rand(Gamma(100, 100dE)),
660+
ρ => rand(Beta(1, 99)),
661+
γ => rand(Gamma(100, 100dI)),
662+
λ => rand(Gamma(10, 10*25))]
663+
prob = SDEProblem(mecha_bayes, u0,
664+
(0.0, 10.0), p)
665+
@test_nowarn sol = solve(prob, SRIW1(), reltol=1e-3, abstol=1e-3, seed = 1);

0 commit comments

Comments
 (0)