sm311o_7.mws,wdj,1-22-98
>
with(plots):
with(DEtools):
Example 1 on page 330 of Williamson's book, Introduction to Diff Eqns and Dynamical Systems, McGraw-Hill, 1997.
Two 50 gal tanks connected by flow pipes. Tank 1 has salt soln flowing into it at 1 gal/min with salt concentration of 1 lb/gal. Tank 2 has pure water flowing into it at 1 gal/min. Tank 1 connects to Tank 2 with a pipe in which soln flows at 2 gal/min to Tank 2. Tank 2 connects to Tank 1 with another pipe in which soln flows at 3 gal/min to Tank 1. Soln flows out of Tank 1 at 2 gal/min.
>
sys:=diff(y1(t),t)=-(4/50)*y1(t)+(3/50)*y2(t)+1,diff(y2(t),t)=(2/50)*y1(t)-(3/50)*y2(t);
ic:=y1(0)=0,y2(0)=0;
> phaseportrait([sys],[y1(t),y2(t)],t=0..100,[[ic]],scene=[y1(t),y2(t)],stepsize=.05,y1=0..25,y2=0..18);
The above picture indicates that the limiting amount is around y1 = 23 and y2 = 15. In fact, the limits may be easily read off from the explicit solution to be y1 = 25 and y2 = 50/3.
> dsolve({sys,ic},{y1(t),y2(t)});
Example 2 on page 332 of Williamson's book. Two tanks, Tank 1 of capacity 100 gal, the other 200 gal, are connected by flow pipes. Tank 1 has salt soln flowing into it at 1 gal/min with salt concentration of 1 lb/gal. Tank 1 connects to Tank 2 with a pipe in which soln flows at 1 gal/min to Tank 2. Tank 2 connects to Tank 1 with another pipe in which soln flows at 2 gal/min to Tank 1.
>
sys:=diff(y1(t),t)=-(1/(50+2*t))*y1(t)+(2/(100-t))*y2(t)+1,diff(y2(t),t)=(1/(50+2*t))*y1(t)-(2/(100-t))*y2(t);
ic:=y1(0)=0,y2(0)=0;
> phaseportrait([sys],[y1(t),y2(t)],t=0..25,[[ic]],scene=[y1(t),y2(t)],stepsize=.05,y1=0..25,y2=0..3);
>
dsolve({sys,ic},{y1(t),y2(t)});
No output in this case. This means that MAPLE cannot solve the system of ODEs explicitly.
>
nsoln:=dsolve({sys,ic},{y1(t),y2(t)},type=numeric);
MAPLE can solve the system numerically.
>
odeplot(nsoln,[y1(t),y2(t)],0..25,numpoints=100,view=[0..25,0..3],labels=[y1,y2__]);
The phase portrait
>
odeplot(nsoln,[t,y1(t)],0..25,numpoints=100,view=[0..25,0..25],labels=[t,y1]);
The graph of y1
>
odeplot(nsoln,[t,y2(t)],0..25,numpoints=100,view=[0..25,0..25],labels=[t,y2]);
The graph of y2
>