Here is the code I have,
using Pkg
using PyCall
using CSV; using DataFrames; using Dierckx;using DifferentialEquations;using LSODA;using MATLABDiffEq; using ParameterizedFunctions;using IfElse;
p = [
1.35e-2 #k1
6.33e-1 #k2
5.00e-5 #k3
1.00e-3 #k4
3.80e-3 #k5
5.82e-1 #k6
2.20e-2 #k7
4.71 #k8
1.08e-2 #k9
#2.6 REMOVED
1.35 #sigma /index 10
0.63 #Km /index 11
5.0 #Gb /index 12
8.58 #Ib /index 13
]
input = [
75000, # D
70, # Mb
0 # t_meal_start
]
c = [0.00551,
17/70,
0.043,
1,
31,
3,
13/70,
9,
30,
0.1,
0.043*(p[11]+p[12])/p[12]-p[5]*1*p[13],
p[7]*p[12]/(1*31*p[13])*30
]
function differentialequations_1(du,u,p,t)
Mg,Gpl,Ipl,Int=u
mg = IfElse.ifelse(t> input[3], p[10].*p[1].^p[10].*(t-input[3]).^(p[10]-1) .* exp(-(p[1].*(t-input[3])).^p[10]) .* input[1],0)
mgpl = p[2] .* Mg
du[1] = mg - mgpl
gli = c[3] - p[3] .* (Gpl-p[12])-p[4] .* c[4] .* (Ipl-p[13])
ggu = p[2] .* (c[1]/(c[2]*input[2])) .* Mg;
gnoni = c[11] .* (Gpl./(p[11]+Gpl));
gi = p[5] .* c[4] .* Ipl .* (Gpl./(p[11]+Gpl))
gre= IfElse.ifelse(Gpl > c[8], c[10] ./ (c[2]*input[2]) .* (Gpl - c[8]), 0 )
du[2] = gli + ggu - gnoni - gi - gre
t_lowerbound= t- c[9];
if (time > c[9]) && (length(t_saved)>1) && (length(t_saved) == length(Gpl_saved))
spl=Spline1D(t_saved,Gpl_saved)
Gpl_lowerbound= spl(t_lowerbound);
else
Gpl_lowerbound = Gpl_saved[1];
end
du[4] = (Gpl-p[12]) - (Gpl_lowerbound-p[12]);
ipnc = (c[4].^-1) .* (p[6].*(Gpl-p[12]) +(p[7]/c[5]).*Int + (p[7]/c[5]) .*p[12]* c[9] + (p[8].*c[6]) .* du[2])
iliv = c[12].*Ipl;
iif = p[9].*(Ipl-p[13])
du[3] = ipnc - iliv - iif
end
function compute(u0, local_G_saved,local_I_saved)
global time=0
global t_saved=[0]
convert(Any, local_G_saved)
convert(Any, local_I_saved)
while time < 120
tspan = (time, time+1)
prob = ODEProblem(differentialequations_1,u0,tspan,p)
sol=solve(prob, MATLABDiffEq.ode15s())
time= time+1
convert(Any, sol)
push!(t_saved, time)
push!(local_G_saved, last(sol[2,:]))
push!(local_I_saved,last(sol[3,:]))
u0=last(sol)
end
return local_G_saved, local_I_saved, u0
end
global Gpl_saved= 5.0
global Ipl_saved=11.5
u0=[0.0,5.0,11.5,0.0]
p[1]=0.00933
p[5]=0.01
p[6]=0.3
p[8]=3.3
simulated_G_saved, simulated_I_saved, u= compute(u0, Gpl_saved,Ipl_saved);
which when run provides the error
Unable to resolve the name mxsol.stats.
Now, I have tried reproducing the error is a smaller version but unfortunately I could not. With every attempt there were some other issues popping up. In some cases, the 4th u variable entirely disappeared in the 2nd iteration. Please run the code and see if you can reproduce the error and see where the flaw is.
I tried to follow the julia recommended coding practices, but please feel free to provide additional pointers if necessary.
Thanks in advance