Barry McQuarrie
Aug, 2015

Integral Curves

Example 1

(dy/dx)=((x^2)/(1-(y^2)))

Here is the solution to the DE we saw in class for a variety of choices of c:

integralcurves=Table[(x^3)+(y^3)-3y+c==0,{c,-3,3,0.4}]

Here is the plot of the integral curves:

ContourPlot[Evaluate[integralcurves],{x,-3,3},{y,-3,3},FrameLabel->{"x","y"},BaseStyle->20]

Notice the derivative is infinite when y=±1. This problem was evident in the original ODE. We can exclude those points from the integral curves.

ContourPlot[Evaluate[integralcurves],{x,-3,3},{y,-3,3},FrameLabel->{"x","y"},BaseStyle->20,Exclusions->{y==1,y==-1}]

Example 2

(dy/dx)=(((e^(y^2))sin(x)cos(x))/y)
y (e^(-(y^2)))dy=sin(x)cos(x) dx
∫y (e^(-(y^2)))dy=∫sin(x)cos(x)dx

We can do these integrals by hand (both are u-sub), but let’s use Mathematica instead.

implicitsol=Integrate[Cos[x] Sin[x],x]==Integrate[y Exp[-y^2],y]+c

Here is the solution to the DE we saw in class for a variety of choices of c:

integralcurves=Table[implicitsol,{c,-3,3,0.2}]

Here is the plot of the integral curves:

ContourPlot[Evaluate[integralcurves],{x,-4,4},{y,-4,4},FrameLabel->{"x","y"},BaseStyle->20,PlotPoints->30,Exclusions->{y==0}]

Note that in this problem, there are some vertical asymptotes that we couldn’t identify until we had the solution to the ODE.

Flow In/Out: Mixing Problems

Direction Field

Here is our DE that we built:

Clear[L,Q]
L[Q_]:=  D[Q,t]== (r/4)-(r/100) *Q
L[Q[t]]

Here is the direction field associated with the DE (note we had to set r=1 to get this plot):

field=VectorPlot[{1, (1/4)-(1/100) *Q },{t,0,50},{Q,0,50},VectorScale->{Tiny,Small,None},AspectRatio->1,Frame->True,FrameLabel->{"t","Q"},BaseStyle->20]

It looks like it all solutions approach Q=25lbs. Let's check what happens for different values of r:

field2[r_]:=VectorPlot[{1, (r/4)-(r/100) *Q },{t,0,50},{Q,0,50},VectorScale->{Tiny,Small,None},AspectRatio->1,Frame->True,FrameLabel->{"t","Q"}];

allapproach=Plot[25,{t,0,50},PlotStyle->{Thick,Red}];

Manipulate[Show[field2[r],allapproach,AspectRatio->1,Frame->True,FrameLabel->{"t","Q"},BaseStyle->20],{r,0,10}]

Notice that the solutions approach 25 lbs regardless of the value of r. The solution Q(t)=25 is an equilibrium solution, and it does not depend on the value of r.

Including the full solution in the plot

We solved the DE by hand in class using the integrating factor method and determined that the solution was:

Clear[Q,Q0,r]
Q[t_,r_,Q0_]=(Q0-25)Exp[-r*t/100]+25

Let's include that solution in the direction field, and see what happens to it as r increases.

field2[r_]:=VectorPlot[{1, (r/4)-(r/100) *Q },{t,0,50},{Q,0,50},VectorScale->{Tiny,Small,None}];

sol[r_,Q0_]:=Plot[Q[t,r,Q0],{t,0,50},PlotStyle->{Thick,Blue},BaseStyle->20];

Manipulate[Show[field2[r],sol[r,Q0],allapproach,AspectRatio->1,Frame->True,FrameLabel->{"t","Q"},BaseStyle->20],{r,0,10,0.2},{Q0,0,50}]

Notice that if r=0, all the solutions are equilibrium solutions since we are not adding any brine to the tank, and therefore no extra salt.

So, for this mixing problem, as r increases:

the equilibrium solution remains unchanged at 25 lbs, and

the rate of approach to that equilibrium increases ar r increases.

The graph shows this nicely, and we can see why this happens since the exponential term in the solution goes to zero faster as r increases.

For r>0 the equilibrium solution is called stable, since all solutions that are on either side of the equilibrium will apporach the equilibrium.