ODE-system solving with initial conditions (Maple)

Question

This is my problem:

I have to solve a differential-equation-system, and then to plot the results. The system is :

{ S'=-r*S*J,

J'=r*S*J-a*J,

R'=a*J }.

It's about epidemiology, S, J, and R are the Susceptible, Infected, and Retired people, while a and r are two parameters that I determined before.

The thing is, when I give as initial conditions the 3 values at t=0, everything is OK, but I don't know these 3 values ! So I'd like to give Maple other values, as the derivate of one of the functions at some point, or this sort of thing... But Maple doesn't accept, whatever my syntax is. I've tried many ways, and already checked the online help that they provide, but didn't find anything helpful. Does anybody know how to fix that issue?

Thank you :)

Answer 1

Using $t=100$ as an effective value of $\infty$ and using the other values that you gave in the comment, here's how to code this in Maple:

sys:= {
     diff(S(t),t) = -r*S(t)*J(t),
     diff(J(t),t) = r*S(t)*J(t) - a*J(t),
     diff(R(t),t) = a*J(t),
     # Boundary conditions:
     S(100) = 50000, R(0) = 0, R(100) = 40000
}: 

Sol:= dsolve(eval(sys, [r=4.9e-6, a=1/3]), numeric):

plots:-odeplot(
     Sol, [[t,S(t)], [t,J(t)], [t,R(t)]], t= 0..100, 
     legend= ["susceptible","infected","retired"]
);

Answer 2

Here's what I tried, following your advice :

sys := {J(0) = 1, R(0) = 0, S(0) = 89999, diff(J(t), t) = r*S(t)*J(t)-a*J(t), diff(R(t), t) = a*J(t), diff(S(t), t) = -r*S(t)*J(t)};

Sol:= dsolve(eval(sys, [r=4.9e-6, a=1/3]), numeric);

plots:-odeplot( Sol, [[t,S(t)], [t,J(t)], [t,R(t)]], t= 50..150, legend= ["Susceptible","Infected","Retired"], color=["Blue", "Red", "Green"]);

As you see, the plot is the same, just shifted to the right for approximately 50 days. And no Newton iteration convergence problem anymore ! I think I'll use this solution : by superimposing it over the previous one, I can easily detremine the "true" initial conditions that I need.

EDIT: I searched for the common maximum for J on the two plots, hence I deduced : OldSol(t)=NewSol(t+40)