Section 6.3 & 6.4

Mathematica Basics for Step Functions

Unit Step Function Definition

We need to be able to define step functions like y=(Subscript[u, 3])(t) easily, and here are three ways to do it:

Clear[t]
y1[t_]=HeavisideTheta[t-3]
y2[t_]=UnitStep[t-3]
y3[t_]:=Piecewise[{{0,t<3},{1,t≥3}}]

These are all the same, which plots will show.

Plot[{y1[t]},{t,0,6},PlotStyle->{Red,Thick}]
Plot[{y2[t]},{t,0,6},PlotStyle->{Blue,Thick}]
Plot[{y3[t]},{t,0,6},PlotStyle->{Green,Thick}]

Piecewise is too cumbersome--I will use one of the other two instead.

Unit Step Function To Turn On or Off Functions

Turn cos(t) ON at t=4:

turnon[t_]=HeavisideTheta[t-4]
Plot[turnon[t]*Cos[t],{t,0,12},PlotStyle->{Red,Thick}]

Turn cos(t) OFF at t=4:

turnoff[t_]=1-HeavisideTheta[t-4]
Plot[turnoff[t]*Cos[t],{t,0,12},PlotStyle->{Red,Thick}]

Mathematica Basics for Partial Fractions

Apart[(1/((s^2)((s^2)+1)))]

Apart[((1+s)/((s^2)((s^2)+1)))]

Apart[(1/(s((s^4)+5 (s^2)+4)))]

Solve the IVP Using Tables
y''+y={({{t, 0≤t≤1}, {0,1<t          }}), y(0)=(y^((1)))(0)=0.

Laplace transforms are very good at working with piecewise defined forcing functions, but we need the forcing function in terms of a unit step function. Let's start by doing that.

We want to turn OFF the function t after t=1:

g[t_]=(1-UnitStep[t-1])*t

Let's plot it to see that we got it right:

Plot[g[t],{t,0,2},PlotStyle->Thick]

Now, we need to have the forcing function in a special form if we want to take the Laplace transform of it using the result from Table 6.2.1 (this wouldn't be necessary if we used Mathematica to take the Laplace transform):

L [(Subscript[u, c])(t)f(t-c)]=(e^(-cs))L [f(t)]

So rewrite it as

g(t)=(1-(Subscript[u, 1])(t))t
= t - t (Subscript[u, 1])(t)
=t - t (Subscript[u, 1])(t)+(Subscript[u, 1])(t)-(Subscript[u, 1])(t)
=t -  (Subscript[u, 1])(t)(t-1)-(Subscript[u, 1])(t)

g[t_]=t-UnitStep[t-1]*(t-1)-UnitStep[t-1]
Plot[g[t],{t,0,2},PlotStyle->Thick]

Now we are ready to start by taking the Laplace transform of the differential equation:

L [y'']+L [y]=L [t]-L [  (Subscript[u, 1])(t)(t-1)]-L [  (Subscript[u, 1])(t)]

LaplaceTransform[y''[t],t,s]

Since the initial conditions are all zero, this reduces to L [y'']=(s^2)Y(s). Here are the rest:

LaplaceTransform[t,t,s]

LaplaceTransform[UnitStep[t-1](t-1),t,s]

LaplaceTransform[UnitStep[t-1],t,s]

Notice that Mathemtaica could have done this all at once, we rewrote the g(t) solely to make using Table 6.2.1 easier.

LaplaceTransform[(1-UnitStep[t-1])*t,t,s]

This leads the the following relationship:

(s^2)Y(s)+Y(s)=(1/(s^2))+((e^(-s))/s)-(((e^(-s)) (1+s))/(s^2))-((e^(-s))/s)
(s^2)Y(s)+Y(s)=(1/(s^2))-(((e^(-s)) (1+s))/(s^2))
Y(s)=(1/((s^2)((s^2)+1)))-(((e^(-s)) (1+s))/((s^2)((s^2)+1)))

Now we need to use partial fractions:

Apart[(1/((s^2)((s^2)+1)))]

Apart[(((1+s))/((s^2)((s^2)+1)))]

And our relationship becomes:

Y(s)=(1/(s^2))-(1/((s^2)+1))-(e^(-s))(1/(s^2))-(e^(-s))(1/s)+(e^(-s))((1+s)/(1+(s^2)))

Now, take the inverse Laplace transform:

(L^(-1))[Y(s)]=(L^(-1))[(1/(s^2))]-(L^(-1))[(1/((s^2)+1))]-(L^(-1))[(e^(-s))(1/(s^2))]-(L^(-1))[(e^(-s))(1/s)]+(L^(-1))[(e^(-s))((1+s)/(1+(s^2)))]

Again, you may use Table 6.2.1 or Mathematica to work out the inverse transformations.

InverseLaplaceTransform[(1/(s^2)),s,t]

InverseLaplaceTransform[(1/((s^2)+1)),s,t]

InverseLaplaceTransform[(E^(-s))(1/(s^2)),s,t]

InverseLaplaceTransform[(E^(-s))(1/s),s,t]

InverseLaplaceTransform[(E^(-s))((1+s)/(1+(s^2))),s,t]

Combining all these together, we get the final result for the solution to the initial value problem:

y[t_]=t-Sin[t]-(-1+t) HeavisideTheta[-1+t]-HeavisideTheta[-1+t]+HeavisideTheta[-1+t] (Cos[1-t]-Sin[1-t])//Simplify

Here is a plot of the solution:

plot1=Plot[y[t],{t,0,4Pi},PlotStyle->{Red,Thick}]

Here is Mathematica's DSolve solution, which looks a bit different from ours:

Clear[y]
DSolve[{y''[t]+y[t]==g[t],y[0]==0,y'[0]==0},y[t],t]

plot2=Plot[(Cos[1] Cos[t]-Cos[t] Sin[1]+Cos[1] Sin[t]-((Cos[1])^2) Sin[t]+Sin[1] Sin[t]-((Sin[1])^2) Sin[t])/(((Cos[1])^2)+((Sin[1])^2))+(t-Sin[t]-(Cos[1] Cos[t]-Cos[t] Sin[1]+Cos[1] Sin[t]-((Cos[1])^2) Sin[t]+Sin[1] Sin[t]-((Sin[1])^2) Sin[t])/(((Cos[1])^2)+((Sin[1])^2))) UnitStep[1-t],{t,0,4Pi},PlotStyle->Thick]

But we see that it is actually the same as the one we found:

Show[plot1,plot2]

Solve the IVP Using Mathematica y''+y={({{t, 0≤t≤1}, {0,1<t          }}), y(0)=(y^((1)))(0)=0.

If we are using Mathematica instead of tables, we don’t need to fix up the g(t) to be able to the tables, and we won’t need to use partial fractions to be able to use the tables.

We want to turn OFF the function t after t=1:

g[t_]=(1-UnitStep[t-1])*t

Take Laplace Transform of DE:

algeq=LaplaceTransform[y''[t],t,s]+LaplaceTransform[y[t],t,s]==LaplaceTransform[g[t],t,s]

Solve the algebraic equation for L[y(t)]:

solalgeqn=Solve[algeq,LaplaceTransform[y[t],t,s]]/.y[0]->0/.y'[0]->0

Take the Inverse Laplace transform to get the solution to the IVP:

IVPsol[t_]=InverseLaplaceTransform[solalgeqn[[1,1,2]],s,t]

Let's plot it to see that we got it right:

plot3=Plot[{IVPsol[t],y[t]},{t,0,4Pi},PlotStyle->{Purple,Thick}]

But we see that it is actually the same as the one we found:

Show[plot1,plot2,plot3]

Solve the IVP (y^((4)))+5(y^((2)))+4y=1-(Subscript[u, π])(t), y(0)=(y^((1)))(0)=(y^((2)))(0)=(y^((3)))(0)=0.

Here is the bit we can easily do by hand, which is to take the Laplace transform of the differential equation. You will need to use Table 6.2.1, or Mathematica to determine the Laplace transforms. Table 6.2.1 will be provided on the final exam.

L [(y^((4)))]+5L[(y^((2)))]+4L [y]=L [1]-L [(Subscript[u, π])(t)]

LaplaceTransform[y''''[t],t,s]

Since the initial conditions are all zero, this reduces to L [(y^((4)))]=(s^4)Y(s). Here are the rest:

LaplaceTransform[y''[t],t,s]

LaplaceTransform[1,t,s]

LaplaceTransform[HeavisideTheta[t-Pi],t,s]
LaplaceTransform[UnitStep[t-Pi],t,s]

Notice you have a choice for how you define the unit step function (Subscript[u, π])(t).

This leads the the following relationship:

(s^4)Y(s)+5(s^2)Y(s)+4Y(s)=(1/s)-((e^(-π s))/s)
Y(s)=(1/(s((s^4)+5(s^2)+4)))-((e^(-π s))/(s((s^4)+5(s^2)+4)))

Now we need to use partial fractions:

Apart[1/s/(s^4+5s^2+4)]

And our relationship becomes:

Y(s)=(1/(4 s))-(s/(3 (1+(s^2))))+(s/(12 (4+(s^2))))-((e^(-π s))/(4 s))+(((e^(-π s))s)/(3 (1+(s^2))))-(((e^(-π s))s)/(12 (4+(s^2))))

Now, take the inverse Laplace transform:

(L^(-1))[Y(s)]=(1/4)(L^(-1))[(1/s)]-(1/3)(L^(-1))[(s/(1+(s^2)))]+(1/(12 ))(L^(-1))[(s/(4+(s^2)))]-(1/4)(L^(-1))[((e^(-π s))/s)]+(1/3)(L^(-1))[(e^(-π s))(s/(1+(s^2)))]-(1/12)(L^(-1))[(e^(-π s))(s/(4+(s^2)))]

Again, you may use Table 6.2.1 or Mathematica to work out the inverse transformations.

InverseLaplaceTransform[(1/s),s,t]

InverseLaplaceTransform[(s/(1+(s^2))),s,t]

InverseLaplaceTransform[(s/(4+(s^2))),s,t]

InverseLaplaceTransform[((E^(-π s))/s),s,t]

InverseLaplaceTransform[(((E^(-π s))s)/(1+(s^2))),s,t]

InverseLaplaceTransform[(((E^(-π s))s)/(4+(s^2))),s,t]

Combining all these together, we get the final result for the solution to the initial value problem:

y[t_]=(1/4)-(1/3)Cos[t]+(1/(12 ))Cos[2t]-(1/4)HeavisideTheta[-π+t]-(1/3)Cos[t] HeavisideTheta[-π+t]-(1/12)Cos[2 (-π+t)] HeavisideTheta[-π+t]

Here is a plot of the solution:

plot1=Plot[y[t],{t,0,4Pi},PlotStyle->{Red,Thick}]

Here is Mathematica's DSolve solution, which looks a bit different from ours:

Clear[y]
g[t_]=1-UnitStep[t-Pi]
DSolve[{y''''[t]+5y''[t]+4y[t]==g[t],y[0]==0,y'[0]==0,y''[0]==0,y'''[0]==0},y[t],t]

plot2=Plot[-((2 Cos[t])/3)+(((2 Cos[t])/3)+(1/12) (3-4 Cos[t]+Cos[2 t])) UnitStep[π-t],{t,0,4Pi},PlotStyle->Thick]

But we see that it is actually the same as the one we found:

Show[plot1,plot2]