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)} >