 1
 0
Hello,
This is my first post and hopefully my question has not been answered elsewhere already as I realize it is annoying to answer the same type of posts over and over again.
I am working on a system of PDEs with one ODE, coupled. It is an SEIR model with one extra class for the group of people who have received treatment.
I believe I am having a problem because of the multiple dependencies in each of the equations and conditions.
Ideally I would like the system to be solved and plotted however I have no idea how to go about this with what I have here. I have searched the internet for other example but I have found nothing helpful. I have tried everything in the Maple help and used all options available to pdsolve. If somebody could take a look at the worksheet and guide me in the right direction that would be very helpful.
I will paste the contents of my worksheet below as well, I will attach it to this post if that is easier to read.
The version of Maple that I am using is Maple 10 student edition.
Thanks,
Shannon
> restart
> with(PDEtools); with(plots);
Warning, the name changecoords has been redefined
> with(DEtools)
> #Ideally I would like a system solved analytically and the results plotted. I would like the results for S,E,I,R,D plotted on the same graph.I have no idea what I have done wrong but I believe my issue lies in the fact that the system I am trying to solve is recursive.
>
>
> #alpha := 0.3257;
> #d := 0.00761;
> #mu := 0.0015;
> #omega := 1/30;
> #lambda := 0.01028;
> #gam := 0.49089;
> #R0 := 1.42;
> #p := 0.1;
> #beta := (R0 * d * (d + mu + gam) * (d + alpha))/(lambda*alpha);
> "#v:=0.8;"
> int(Ip(a, t), a = 0 .. infinity) := Ip(a, t);
> int(Del(a, t), a = 0 .. infinity) := Del(a, t);
> int(R(a, t), a = 0 .. infinity) := R(a, t);
int(Ip(a, t), a = 0 .. infinity) := Ip(a, t)
int(Del(a, t), a = 0 .. infinity) := Del(a, t)
int(R(a, t), a = 0 .. infinity) := R(a, t)
> "#ode:=diff(S(a,t),t)  lambda+d*S(a,t) + S(a,t)*beta*int(Ip(a,t),a=0..infinity)+ v*p*S(a,t)  omega*int(Del(a,t),a=0..infinity)f*int(R(a,t),a=0..infinity)=0; "
> "#dsolve(ode, S(a,t));"
> [/ d \ / d \
> pdesys := [ E(a, t) +  E(a, t) + d E(a, t) + alpha E(a, t) = 0,
> [\ dt / \ da /
>
> / d \ / d \
>  Ip(a, t) +  Ip(a, t) + d Ip(a, t) + gam Ip(a, t) + mu Ip(a, t) = 0,
> \ dt / \ da /
>
> / d \ / d \
>  R(a, t) +  R(a, t) + f R(a, t) + d R(a, t) = 0,
> \ dt / \ da /
>
> / d \ / d \ ]
>  Del(a, t) +  Del(a, t) + d Del(a, t) = 0]
> \ dt / \ da / ]
> print();
[/ d \ / d \
[ E(a, t) +  E(a, t) + d E(a, t) + alpha E(a, t) = 0,
[\ dt / \ da /
/ d \ / d \
 Ip(a, t) +  Ip(a, t) + d Ip(a, t) + gam Ip(a, t) + mu Ip(a, t) = 0,
\ dt / \ da /
/ d \ / d \
 R(a, t) +  R(a, t) + f R(a, t) + d R(a, t) = 0,
\ dt / \ da /
/ d \ / d \ ]
 Del(a, t) +  Del(a, t) + d Del(a, t) = 0]
\ dt / \ da / ]
> IBCs := {E(a, 0) = 0, Ip(a, 0) = 0.05, R(a, 0) = 0, Del(a, 0) = 0,
>
> E(0, t) = beta Ip(0, t), Ip(0, t) = alpha E(0, t), R(0, t) = gam Ip(0, t),
>
> Del(0, t) = v p}
> print();
{E(0, t) = beta Ip(0, t), Ip(0, t) = alpha E(0, t), R(0, t) = gam Ip(0, t),
Del(0, t) = v p, E(a, 0) = 0, Ip(a, 0) = 0.05, R(a, 0) = 0, Del(a, 0) = 0}
> pdesol := pdsolve(pdesys)
> print();
{Del(a, t) = _F1(t  a) exp(d a), E(a, t) = _F2(t  a) exp((d + alpha) a),
R(a, t) = _F3(t  a) exp((f + d) a),
Ip(a, t) = _F4(t  a) exp((d + gam + mu) a)}
>
This is my first post and hopefully my question has not been answered elsewhere already as I realize it is annoying to answer the same type of posts over and over again.
I am working on a system of PDEs with one ODE, coupled. It is an SEIR model with one extra class for the group of people who have received treatment.
I believe I am having a problem because of the multiple dependencies in each of the equations and conditions.
Ideally I would like the system to be solved and plotted however I have no idea how to go about this with what I have here. I have searched the internet for other example but I have found nothing helpful. I have tried everything in the Maple help and used all options available to pdsolve. If somebody could take a look at the worksheet and guide me in the right direction that would be very helpful.
I will paste the contents of my worksheet below as well, I will attach it to this post if that is easier to read.
The version of Maple that I am using is Maple 10 student edition.
Thanks,
Shannon
> restart
> with(PDEtools); with(plots);
Warning, the name changecoords has been redefined
> with(DEtools)
> #Ideally I would like a system solved analytically and the results plotted. I would like the results for S,E,I,R,D plotted on the same graph.I have no idea what I have done wrong but I believe my issue lies in the fact that the system I am trying to solve is recursive.
>
>
> #alpha := 0.3257;
> #d := 0.00761;
> #mu := 0.0015;
> #omega := 1/30;
> #lambda := 0.01028;
> #gam := 0.49089;
> #R0 := 1.42;
> #p := 0.1;
> #beta := (R0 * d * (d + mu + gam) * (d + alpha))/(lambda*alpha);
> "#v:=0.8;"
> int(Ip(a, t), a = 0 .. infinity) := Ip(a, t);
> int(Del(a, t), a = 0 .. infinity) := Del(a, t);
> int(R(a, t), a = 0 .. infinity) := R(a, t);
int(Ip(a, t), a = 0 .. infinity) := Ip(a, t)
int(Del(a, t), a = 0 .. infinity) := Del(a, t)
int(R(a, t), a = 0 .. infinity) := R(a, t)
> "#ode:=diff(S(a,t),t)  lambda+d*S(a,t) + S(a,t)*beta*int(Ip(a,t),a=0..infinity)+ v*p*S(a,t)  omega*int(Del(a,t),a=0..infinity)f*int(R(a,t),a=0..infinity)=0; "
> "#dsolve(ode, S(a,t));"
> [/ d \ / d \
> pdesys := [ E(a, t) +  E(a, t) + d E(a, t) + alpha E(a, t) = 0,
> [\ dt / \ da /
>
> / d \ / d \
>  Ip(a, t) +  Ip(a, t) + d Ip(a, t) + gam Ip(a, t) + mu Ip(a, t) = 0,
> \ dt / \ da /
>
> / d \ / d \
>  R(a, t) +  R(a, t) + f R(a, t) + d R(a, t) = 0,
> \ dt / \ da /
>
> / d \ / d \ ]
>  Del(a, t) +  Del(a, t) + d Del(a, t) = 0]
> \ dt / \ da / ]
> print();
[/ d \ / d \
[ E(a, t) +  E(a, t) + d E(a, t) + alpha E(a, t) = 0,
[\ dt / \ da /
/ d \ / d \
 Ip(a, t) +  Ip(a, t) + d Ip(a, t) + gam Ip(a, t) + mu Ip(a, t) = 0,
\ dt / \ da /
/ d \ / d \
 R(a, t) +  R(a, t) + f R(a, t) + d R(a, t) = 0,
\ dt / \ da /
/ d \ / d \ ]
 Del(a, t) +  Del(a, t) + d Del(a, t) = 0]
\ dt / \ da / ]
> IBCs := {E(a, 0) = 0, Ip(a, 0) = 0.05, R(a, 0) = 0, Del(a, 0) = 0,
>
> E(0, t) = beta Ip(0, t), Ip(0, t) = alpha E(0, t), R(0, t) = gam Ip(0, t),
>
> Del(0, t) = v p}
> print();
{E(0, t) = beta Ip(0, t), Ip(0, t) = alpha E(0, t), R(0, t) = gam Ip(0, t),
Del(0, t) = v p, E(a, 0) = 0, Ip(a, 0) = 0.05, R(a, 0) = 0, Del(a, 0) = 0}
> pdesol := pdsolve(pdesys)
> print();
{Del(a, t) = _F1(t  a) exp(d a), E(a, t) = _F2(t  a) exp((d + alpha) a),
R(a, t) = _F3(t  a) exp((f + d) a),
Ip(a, t) = _F4(t  a) exp((d + gam + mu) a)}
>
Attachments

105.9 KB Views: 345