Hi ! I wrote a code for a project which concerns differential equations but the problem is I don't have the good curves, even if the code is right (I think). Or it's a problem of syntax but I spend a lot of time searching this, without success. I give you the code accompanied with the formulas and the exit of the code. I can't add hyperlink so I write like this or you can't understand my problem. Sorry for inconveniance !
http:// image.noelshack.com/fichiers/2016/44/1478185128-image-equations-differentielles.png
http:// image.noelshack.com/fichiers/2016/44/1478185610-image-formule-de-lambda-lmd.png
http:// image.noelshack.com/fichiers/2016/44/1478185128-image-graphique-sagemath.png
I'm supposed to get them curves :/ The only problem seems to be black curve (on SageMath) which is the green curve on Excel (third picture) and I don't wanted to put the black curve of Excel (so it's normal if it doesn't appears).
http: //image.noelshack.com/fichiers/2016/44/1478185262-sans-titre.png
y[0] = X (Number of likely ie those not infected) = B - nuX - clmd*X
y[1] = Y (Number of contagious) = cXlmd - Y*(v + nu)
y[2] = Z (Number of non-infectious HIV) = vY(1 - p) - nu*Z
y[3] = A (Number of AIDS patients) = pvY - A*(d + nu)
y[4] = lmd (Proportion of contagious, which made themselves pass HIV on the total population) = (Y*beta)/(X + Y + Z)
params[0] = nu (Natural death rate) == 0.03125 car nu + d ~= d (1 à 1.33/an)
params[1] = d (Death rates by disease) == 1 (1/an)
params[2] = B (Nombre de susceptibles immigrants) == 13333.3/an
params[3] = v (Number of immigrants may) == 0.2/an
params[4] = c (Number of sexual partners) == 48 (2 à 6/mois soit 48/an)
params[5] = p (Proportion of catching HIV) == 0.3 (10 à 30%)
params[6] = beta (Transmission probability) == 0.0014
t (Time)
@interact def interactive_function ( nu = slider (0.001, 0.1, 0.005, default = 0.03125), d = slider (0, 1.5, 0.05, default = 1), B = slider (0, 30000, 500, default = 13333.3), v = slider (0, 1, 0.05, default = 0.2), c = slider (0, 365, 1, default = 48), p = slider (0, 1, 0.05, default = 0.3), beta = slider (0, 1, 0.02, default = 0.014) ) :
def f_1 (t,y,params) :
return [ params[2] - params[0]*y[0] - params[4]*y[4]*y[0],
params[4]*y[0]*y[4] - y[1]*(params[3] + params[0]),
params[3]*y[1]*(1-params[5]) - params[0]*y[2],
params[5]*params[3]*y[1] - y[3]*(params[1] + params[0]),
(y[1]*params[6])/(y[0] + y[1] + y[2]) ]
T = ode_solver() T.function = f_1 T.algorithm = "rk8pd" T.ode_solve (y_0 = [60000,40000,0,0,0.0056], t_span = [0,35], params = [nu,d,B,v,c,p,beta], num_points = 1000)
f = T.solution X = [(x[0],x[1][0]) for x in f] Y = [(x[0],x[1][1]) for x in f] Z = [(x[0],x[1][2]) for x in f] A = [(x[0],x[1][3]) for x in f]
P1 = line (X, rgbcolor='black') P2 = line (Y, rgbcolor='green') P3 = line (Z, rgbcolor='red') P4 = line (A, rgbcolor='blue')
show (P1+P2+P3+P4)
Thanks a lot for the help you can bring to me ! I would really like to understand my problem :/