I was looking for some solvers for my example model to find steady state values for states. I found SSRootfind() and DynamicSS(). Upon reading I felt. It is much better to use DynamicSS as for SSRootfind() your initial values assumption for states should be really good.
Here maybe I’m wrong but for my example model as I change initial condition I find different steady state values (I’m using DynamicSS()). In my head I always thought it should be same irrespective of your initial values.
Can someone shed light on this? and also as I’m getting only single value here for states, does this reflect it has only one solution for steadystate?
function F_simple(x, p, t; γ = 1/18, R₀ = 3.0, σ = 1/5.2)
s, e, i, r = x
return [-γ*R₀*s*i; # ds/dt = -γR₀si
γ*R₀*s*i - σ*e;# de/dt = γR₀si -σe
σ*e - γ*i; # di/dt = σe -γi
γ*i; # dr/dt = γi
]
end
i_0 = 0.7
e_0 = 0.1
r_0 = 0.2
x_0 = [s_0, e_0, i_0, r_0] # initial condition
tspan = (0.0, 350.0) # ≈ 350 days
problem = ODEProblem(F_simple, x_0, tspan)
prob = SteadyStateProblem(problem)
sol = solve(prob, DynamicSS(Rodas5P()))
julia> sol = solve(prob, DynamicSS(Rodas5P()))
u: 4-element Vector{Float64}:
0.010268836444676122
1.6803043675976122e-9
1.365844689856922e-7
0.98973102529055
Now if I change initial condition:
i_0 = 0.5
e_0 = 0.4
r_0 = 0.1
output:
julia> sol = solve(prob, DynamicSS(Rodas5P()))
u: 4-element Vector{Float64}:
0.007545374053424994
7.285142390091502e-10
8.02329983116032e-8
1.1924545449850648