6  Differential Equations

Throughout, we shall continlue to type all Input statements in BasicInput form

6.1  Equations and Solutions

Example 6.1.1  (on page 148)

y[t_] = 1/(1-t^2)

{y'[t], 2t*y[t]^2}

y'[t]- 2t*y[t]^2

y'[t]== 2t*y[t]^2

Plot[ y[t], {t, -3, 3}, PlotRange->{-8, 8} ]

Example 6.1.2  (on page 149)

Clear[y];  y[t_] = 2(t^2 + 1)^(-2) + t^2 + 1

y'[t]+ (4t/(t^2 + 1))y[t]//Simplify

Plot[ y[t], {t, -2, 2}, PlotRange->{0,5}]

Example 6.1.3  (on page 150)

Clear[y];  y[t_] = 1/(c - t^2)

y'[t] - 2t*y[t]^2

solns = y[t] /. c->{-.25,-.1, 0, .5, 2, 4}

Plot[ Evaluate[solns], {t, -2.5, 2.5}, PlotRange->{-11, 10} ]

Example 6.1.4  (on page 150)

Clear[y];  y[t_] = c(t^2 + 1)^(-2) + t^2 + 1

y'[t] + (4t/(t^2 + 1))y[t]  //Simplify

sols = y[t] /. c->{-4, -2, 0, 2, 4}

Plot[ Evaluate[sols], {t, -2, 2}]

Initial-Value Problems

Example 6.1.5   (on page 151)

Clear[y];  y[t_] = 4/(1 - 4t^2)

y'[t]- 2t*y[t]^2

y[0]

Plot[ y[t], {t, -.5, .5}, PlotRange->{0, 20} ]

Examle 6.1.6

Clear[y];   y[t_]  = t^2 + 1 + 8(t^2 + 1)^(-2)

y'[t] +(4t/(t^2 + 1))y[t] //Simplify

y[1]

Plot[ y[t], {t, -4, 4}, PlotRange->{0, 15} ]

Example 6.1.7  (on page 153)

Clear[y];   y[t_, c_] = 10/(1 + c*E^(-t/2))

D[ y[t,c], t] -(1/20)y[t,c](10 - y[t,c])//Simplify

sol = Solve[ y[0,c] == 2, c ]

Flatten[ sol ]

y[t, 4]

Plot[ y[t, 4], {t, -10, 15} ]

6.2  Direction Fields

Needs["VectorFieldPlots`"]

?VectorFieldPlot

VectorFieldPlot[{1,t-y^2},{t,-.5,3},{y,-.7,1.7}]

plotDirectionField[fn_,{x_,a_,b_},{y_,c_,d_},options___]:=(Needs["VectorFieldPlots`"];VectorFieldPlot[{1,Evaluate[fn/.x->t/.y->w]},{t,a,b},{w,c,d},options,ScaleFunction->(1&),PlotPoints->{13,8},ColorFunction->(Hue[0]&)])

plotDirectionField[ t - y^2, {t, -.5,3}, {y,-.7, 1.7}]

plotDirectionField[t-y^2,{t,-.5,3}, {y,-.7,1.7}]

Example 6.2.1  (on page 156)

Clear[y];  y[t_] = c*E^(-t^2/2)

D[ y[t], t] +t*y[t]

dirfield = plotDirectionField[ -t*y, {t, -3, 3}, {y, -2.5, 2.5} ]

curves = Plot[ Evaluate[ Table[ c*E^(-t^2/2), {c, -2.5, 2.5, 1} ] ], {t, -3, 3}, PlotStyle->Thickness[.004]]

Show[ dirfield, curves, Axes->True]

6.3  Euler's Method

eulerStep[ {t_,y_} ] := {t,y} + h{1, f[t,y]}

(Here, delayed set (i.e.   :=  )   is needed because the right-hand member is not known until  f[t, y]  is specified.)

f[t_,y_] = Sin[t*y]

h = 0.1

points = NestList[ eulerStep, {0,2}, 50]

ListPlot[ points ]

ListPlot[points,PlotJoined->True,PlotRange->{0,3}]

Example 6.3.2  (on page 160)

f[t_,y_] = Cos[y] - Sin[t]

h = .1

points = NestList[ eulerStep, {0,1}, 100];

curve=ListPlot[points,PlotJoined->True,AspectRatio->Automatic]

field = plotDirectionField[ f[t,y], {t, 0, 10}, {y, 0, 2.5}]

Show[field, curve, Axes->True]

6.4  Exact Solutions  (starting on page 161)

algEq = Integrate[ y^(-2), y] == Integrate[ Sin[t], t] + c

Solve[ algEq, y ]

Solve[ algEq, y ][[1]]

y[t_] = y /. Solve[ algEq, y ][[1]]

y'[t] - y[t]^2Sin[t]

Example 6.4.2

algEq = Integrate[ y^(-1/2), y ] == Integrate[ -t, t ] + c

eqn = algEq /. {t->0,y->3}

sol = Solve[ eqn, c]

Flatten[ sol ]

sol[[1]]

cValue = c /. sol[[1]]

algEqSubstituted = algEq /. c->cValue

Solve[ algEqSubstituted, y ]

Solve[ algEqSubstituted, y ][[1]]

y[t_] = y /. Solve[ algEqSubstituted, y ][[1]]

Plot[ y[t], {t,0,4}]

Clear[y]

Integration Between Limits  (on page 163)

Example 6.4.3   (on page 163)

algEq = Assuming[ 5<y<100&&t>0,Integrate[ 100/( (100 - u)u ), {u,5,y} ] == Integrate[1/10, {s, 0, t} ] ]

y[t_] = y /. First[Solve[algEq,y]]

Plot[y[t],{t,-20,100}]

As a check, we have

y[0]

y'[t] - (1/10)(1 - y[t]/100)y[t] //Simplify

Clear[y]

First-Order Linear Equations

Example 6.4.4   (on page 164)

algEq = Assuming[ 2 <y < 5 && t > 0, Integrate[ 1/(2-u), {u,5,y}] == Integrate[ 5, {s,0,t}] ]

y[t_] = y /. Solve[ algEq, y][[1]]

Plot[ y[t], {t,0,1}, PlotRange->All ]

As a check, we have

y[0]

y'[t] +5y[t]-10

Clear[y]

Example 6.4.5   (on page 165)

eqn = E^t*y  == Integrate[ E^t*Sin[t], t] + c

E^(-t)*eqn

?Thread

expr2 = Thread[ E^(-t)*eqn, Equal]

expr3 = Simplify[ expr2 ]

y[t_] = y /. Solve[expr3, y][[1]]

As a check, we note that

y'[t] + y[t] - Sin[t] //Simplify

Clear[y]

Example 6.4.6   (on page 166)

phi[t_] =  E^Integrate[ 3/(5 - t), t]

(The integral in the text does not evaluate as desired.  Therefore, proceed as follows.)

y[t_] = phi[t]^(-1)(  phi[0]*y0 + Assuming[ t<5, Integrate[phi[s]*s,{s,0,t}]]) //Simplify

y[t_] = y[t] /. y0->10

Plot[y[t], {t,0,5}]

As a check, we have

y[t]

y[0]

y'[t] + (3/(5-t))y[t] - t  //Simplify

Clear[ y, phi,psi]

Example 6.4.7   (on page 167)

phi[t_] = E^Integrate[k/m, t]

v[t_] = phi[t]^(-1)( phi[0]*v0 + Integrate[ phi[s]*(-g), {s,0,t}]) //Simplify

As a check, we have

v[0]

v'[t] +(k/m)v[t] + g  //Simplify

v[t_] = v[t] /. {k->1, m->1, g->32, v0->0}

Plot[ v[t], {t,0,5}]

Clear[v]

6.5  DSolve andd NDSolve

Using DSolve

Example 6.5.1

Clear[y]

DSolve[ y'[t]+y[t] == Cos[t], y[t],t]

ySoln[t_] = y[t] /. DSolve[  y'[t]+y[t] == Cos[t], y[t],t][[1]]/.C[1]->c

As a check, we have

ySoln'[t] + ySoln[t] - Cos[t]  //Simplify

Example 6.5.2   (on page 168)

Clear[y]

DSolve[ y''[t] + 2y'[t]+2y[t] == t, y[t],t]

y[t_] = y[t]/. DSolve[ y''[t] + 2y'[t]+2y[t] == t, y[t],t][[1]]/.{C[1]->b,C[2]->a}

As a check, we have

y[t]

y'[t]

y''[t]

and

y''[t] +2y'[t]+2y[t]-t //Simplify

Example 6.5.3

Clear[y]

y[t_] = y[t] /. DSolve[{y'[t]+y[t]==E^(-2t), y[0]==c},y[t],t][[1]]

As a check, we have

y[0]

y'[t] + y[t] -E^(-2t) //Simplify

Example 6.5.4   (on page 169)

Do to previous evaluations, it is well to start with an unencombered kernel.

Quit[]

DSolve[ y'[t]+t*y[t]==y[t]^2,y[t],t]

y[t_] = y[t] /. First@DSolve[ {y'[t]+t*y[t]==y[t]^2,y[0]==y0},y[t],t]

?First

As a check, we have

y[0]

y'[t] + t*y[t] - y[t]^2 //Simplify

Example 6.5.5   (on page 169)

Clear[y]

DSolve[ {y'[t] == t^2 - y[t]^3, y[0] == 2}, y[t],t ]

Using NDSolve

Example 6.5.6   (on page 169)

numSol = NDSolve[ {y'[t] == t^2 - y[t]^3, y[0] == 2}, y, {t,0,5}]

yApprox[t_] = y[t] /. numSol[[1]]

Plot[ yApprox[t], {t,0,5}]

Example 6.5.7   (on page 170)

numSolns = Table[  NDSolve[ {y'[t] == 3Cos[t*y[t]],y[0]==c}, y,{t,0,5}], {c,0,5,1}]

yApprox[t_] = y[t] /. numSolns //Flatten

Plot[ Evaluate[ yApprox[t]], {t,0,5}, PlotRange->All]

6.6   Systems of Differential Equations

Example 6.6.1

x[t_] = 2Sin[t] + Cos[t]

y[t_] = Sin[t]

x'[t] -2x[t]+ 5y[t] //Simplify

x[0]

y'[t] - x[t] + 2y[t] //Simplify

y[0]

ellipse = ParametricPlot[{x[t], y[t]}, {t, 0, 2*Pi} ]

Direction Fields

Needs["VectorFieldPlots`"]

field=VectorFieldPlot[{2 x-5 y,x-2 y},{x,-2.5,2.5},{y,-1.5,1.5},PlotPoints->11,ScaleFunction->(1&),Axes->True,ColorFunction->(Blue &)]

Example 6.6.2

field=VectorFieldPlot[{-x-y,2 x-y},{x,-3,3},{y,-2,2},PlotPoints->{13,7},ScaleFunction->(1&),ColorFunction->(Red&)]

Example 6.6.3   (on page 173)

Needs["VectorFieldPlots`"]

field=VectorFieldPlot[{x+y,x-y},{x,-3,3},{y,-2,2},PlotPoints->{13,7},ScaleFunction->(1&),Axes->True,ColorFunction->(Hue[0]&)]

Euler's Method   (page 174)

eulerStep[{x_,y_}] := {x,y} + h{f[x,y],g[x,y]}

Example 6.6.4   (on page 174)

f[x_,y_] = - x - y

g[x_,y_] = 2x - y

h = .05

points = NestList[ eulerStep, {3,-2}, 100];

curve=ListPlot[points,Joined->True]

Show[curve, field]

NDSolve

Example 6.6.5

f[x_,y_] = -y

g[x_,y_] = x

h = .01

points = NestList[eulerStep, {1,0}, 700];

ListPlot[points,Joined->True,AspectRatio->Automatic]

points[[627]]

points[[628]]

points[[629]]

points[[630]]

{1,0} - points[[628]]

{1,0} - points[[629]]

solns = {x[t],y[t]} /. NDSolve[ {x'[t] == -y[t], y'[t] == x[t], x[0] == 1, y[0]==0},{x[t],y[t]},{t,0,7}]//Flatten

ParametricPlot[ Evaluate[solns], {t,0,7}, AspectRatio->Automatic]

{1,0} - solns /. t->2*Pi

Example 6.6.6   (on page 177)

Clear[ solns ]

f[x_,y_] = x(1 - y/9)

g[x_,y_] = -y(1 - x/12)

solns = {x[t],y[t]} /. NDSolve[ {x'[t] == f[x[t],y[t]],y'[t] == g[x[t],y[t]], x[0] == 4, y[0] == 4}, {x[t],y[t]}, {t, 0, 10}]//Flatten

curve = ParametricPlot[ Evaluate[solns], {t, 0, 10}, PlotStyle->Thickness[.006] ]

Needs["VectorFieldPlots`"]

field=VectorFieldPlot[{f[x,y],g[x,y]},{x,0,35},{y,0,25},ScaleFunction->(1&),Axes->True,ColorFunction->(Hue[0]&)]

Show[ curve, field, AspectRatio->Automatic]

The questions for the examination during the class that follows the study of Chapter 6 will be similar to Example 6.5.1 on page 168 and Example 6.5.2 on page 168.

The two most useful MATHEMATICA commands in this chapter are  DSolve  and  NDSolve.  Since NDSolve is a bit more complicated, we shall emphasize  DSolve.

Typical examples for  DSolve   are:

1.  Find the general solution of the differential equation y''''(x) - y(x) = 0.

DSolve[y''''[x] - y[x]==0,y[x],x]

2.  Find the general solution of the differential equation  y''''(x) - y(x) = (13(x^3)+7x+5)(e^(3x+5)) cos (2x+4)

DSolve[ y''''[x]-y[x]==(13x^3+7x+5)Exp[3x+5]*Cos[2x+4], y[x],x]

3.  Find the unique solution of  y'''(x) -3y''(x) + 3y'(x) - y(x) = (x^4)+ 5 that  satisfies  y(0) = 1, y'(0) = 2, and y''(0) = 3.

sol = DSolve[ {y'''[x] -3y''[x]+3y'[x]-y[x]==x^4+5, y[0]==1,y'[0]==2,y''[0]==3},y[x],x]

sol[[1]]

y[x_] = y[x] /. sol[[1]]

4.  Graph the preceding solution over the interval from x = - 1 to x = 1.

Plot[ y[x],{x,-1,1}]