laplace[1].mws

>

MA 507

Chapter 5: Laplace Transforms

16 October 2000

Richard Hitt

Introduction

This worksheet gives a few of the basics on using Maple with Laplace transforms.

Getting Started

We have to load the integral transforms package if we want to use the built-in Maple function for the Laplace transform and its inverse.

> restart;

> with(inttrans);

[addtable, fourier, fouriercos, fouriersin, hankel,...

>

The Basics

Recall the definition of the Laplace transform:

L(f)(s) = int(e^(-s*t)*f(t),t = 0 .. infinity) .

We can use Maple to calculate the Laplace transforms of functions. The syntax is illustrated below where we let Maple calculate the transform of f(t) = t .

> laplace(t,t,s);

1/(s^2)

Of course, Maple can calculate transforms of many functions.

> laplace(sin(t)*exp(t),t,s);

1/((s-1)^2+1)

> laplace(BesselJ(0,t),t,s);

1/(sqrt(s^2+1))

> laplace(diff(y(t),t)+2*y(t)=sin(t),t,s);

s*laplace(y(t),t,s)-y(0)+2*laplace(y(t),t,s) = 1/(s...

Maple can compute inverse Laplace transforms as well.

> invlaplace(s/(s^2-s+1),s,t);

exp(1/2*t)*cos(1/2*sqrt(3)*t)+1/3*sqrt(3)*exp(1/2*t...

Solving an IVP with Laplace Transforms

Maple can be used to solve differential equations via Laplace transforms. For the first example, we use Maple to perform each step along the way.

Example 1

> de1 := diff(y(t),t)+2*y(t)=sin(t);

de1 := diff(y(t),t)+2*y(t) = sin(t)

We can use laplace( , , ) to transform the entire DE.

> lode1:=laplace(de1,t,s);

lode1 := s*laplace(y(t),t,s)-y(0)+2*laplace(y(t),t,...

We can also impose an initial condition like y(0) = 1 .

> livp1:=subs(y(0)=1,lode1);

livp1 := s*laplace(y(t),t,s)-1+2*laplace(y(t),t,s) ...

If we let Y1 denote the laplace transform of y , the above equation is s*Y1-1+2*Y1 = 1/(s^2+1) . In order to solve the equation for Y1 , we can use

> Y1 :=solve(livp1,laplace(y(t),t,s));

Y1 := (s^2+2)/(s^3+s+2*s^2+2)

We can find the solution to the original IVP by taking the inverse Laplace transform of Y1 :

> invlaplace(Y1,s,t);

6/5*exp(-2*t)-1/5*cos(t)+2/5*sin(t)

Maple can automate the entire process using the method=laplace option in dsolve.

> dsolve({de1,y(0)=1},y(t), method=laplace);

y(t) = 6/5*exp(-2*t)-1/5*cos(t)+2/5*sin(t)

Note: You will get the same result in this case if you don't specify "method=laplace" in the dsolve command. This is not always the case, though.

>

Other Examples

Example 2

> ode1:=diff(y(t),t$2)+y(t)=Heaviside[t-Pi]*cos(t);

ode1 := diff(y(t),`$`(t,2))+y(t) = Heaviside[t-Pi]*...

> ic1:=y(0)=0,D(y)(0)=0;

ic1 := y(0) = 0, D(y)(0) = 0

> soln1:=rhs(dsolve({ode1,ic1},y(t),method=laplace));

soln1 := (-1/2*sin(t)*Pi+1/2*sin(t)*t)*Heaviside(t-...

> plot(soln1,t=0..30);

[Maple Plot]

>

Example 3

> ode2:=diff(y(t),t$2)+0*diff(y(t),t)+4*y(t)=2*Heaviside[t-Pi]-Heaviside(t-2*Pi);

ode2 := diff(y(t),`$`(t,2))+4*y(t) = 2*Heaviside[t-...

> ic2:=y(0)=1,D(y)(0)=0;

ic2 := y(0) = 1, D(y)(0) = 0

> soln2:=rhs(dsolve({ode2,ic2},y(t),method=laplace));

soln2 := -1/2*sin(1/2*sqrt(4)*(-t+2*Pi))^2*Heavisid...

> plot(soln2,t=0..15);

[Maple Plot]

>

Square Waves

We define a square wave function in Maple as follows

> sqw := proc(t,d,T) piecewise(t-floor(t/T)*T<(d/100)*T, 1) end;

sqw := proc (t, d, T) piecewise(t-floor(t/T)*T < 1/...

so that sqw[t,d,T] gives a function having value 1 d% of the time with period T. For example,

> plot(sqw(t,50,4),t=0..10);

[Maple Plot]

> dsolve({diff(y(t),t$2) + 2*diff(y(t),t) + 4*y(t) = sqw(t,50,4), y(0)=1, D(y)(0)=2},y(t),method=laplace);

y(t) = Int(1/3*cos(sqrt(3)*u)*Heaviside(2-u+4*floor...
y(t) = Int(1/3*cos(sqrt(3)*u)*Heaviside(2-u+4*floor...

But Maple is not able to handle the integration directly. So we can define a square wave over a finite interval using sums of Heaviside functions. For example

> E := sum(Heaviside(t-2*k)-Heaviside(t-2*k-1),k=0..10);

E := Heaviside(t)-Heaviside(t-1)+Heaviside(t-2)-Hea...
E := Heaviside(t)-Heaviside(t-1)+Heaviside(t-2)-Hea...
E := Heaviside(t)-Heaviside(t-1)+Heaviside(t-2)-Hea...
E := Heaviside(t)-Heaviside(t-1)+Heaviside(t-2)-Hea...

> plot(E,t=0..30);

[Maple Plot]

Then Maple can solve a DE built with such a finite square wave.

> ode3:=diff(y(t),t$2) +0*diff(y(t),t) + 4*y(t) = E;

ode3 := diff(y(t),`$`(t,2))+4*y(t) = Heaviside(t)-H...
ode3 := diff(y(t),`$`(t,2))+4*y(t) = Heaviside(t)-H...
ode3 := diff(y(t),`$`(t,2))+4*y(t) = Heaviside(t)-H...
ode3 := diff(y(t),`$`(t,2))+4*y(t) = Heaviside(t)-H...

> sol3:=dsolve({ode3, y(0)=1, D(y)(0)=1},y(t), method=laplace);

sol3 := y(t) = 1/2*Heaviside(t-4)*sin(1/2*sqrt(4)*(...
sol3 := y(t) = 1/2*Heaviside(t-4)*sin(1/2*sqrt(4)*(...
sol3 := y(t) = 1/2*Heaviside(t-4)*sin(1/2*sqrt(4)*(...
sol3 := y(t) = 1/2*Heaviside(t-4)*sin(1/2*sqrt(4)*(...
sol3 := y(t) = 1/2*Heaviside(t-4)*sin(1/2*sqrt(4)*(...
sol3 := y(t) = 1/2*Heaviside(t-4)*sin(1/2*sqrt(4)*(...
sol3 := y(t) = 1/2*Heaviside(t-4)*sin(1/2*sqrt(4)*(...
sol3 := y(t) = 1/2*Heaviside(t-4)*sin(1/2*sqrt(4)*(...

> plot(rhs(sol3),t=0..20);

[Maple Plot]

>

An RLC Circuit

We wish to analyze a series RLC circuit having an applied voltage which is a square wave. Specifically, we would like to find the initial conditions which completely eliminate any transient in the circuit. The relevant differential equation is L*diff(q(t),`$`(t,2))+R*diff(q(t),t)+q(t)/C = E(t) with initial conditions q(0) = q[0] and D(q)(0) = q[1] where q = q(t) is the charge on the capacitor. We will assume L=1 henry, R = 3 ohms, C = 1/2 farad, and E(t) is a square wave function starting with value 1 and alternating from 1 to 0 to 1 across each unit interval.

Exercise: By trial and error, find a set of initial conditions that give a periodic solution, i.e., that produce no transient.

Dirac Delta

For convenience, we define a shorthand name for the Heaviside function:

> H := t -> Heaviside(t);

H := Heaviside

This pulse function corresponds to what we will refer to as I[h](t-t[0]) in class if we use t[0] = 0 :

> pulse := (t,h) -> (1/h)*H(t)*H(h-t);

pulse := proc (t, h) options operator, arrow; H(t)*...

Let's graph the pulse for several values of h:

> plot( [pulse(t,1),pulse(t,1/2),pulse(t,1/4),pulse(t,1/5),pulse(t,1/6)], t=-0.5..1.5, color=[red,blue,green,black,brown], title=`Pulse functions for h = 1, 1/2, 1/4, 1/5, 1/6`);

[Maple Plot]

Maple can compute the Laplace transforms of Dirac delta functions.

> laplace(Dirac(t),t,s);

1

> laplace(Dirac(t-2),t,s);

exp(-2*s)

Here is an example.

> ode10:= diff(y(t),t$2)+2*diff(y(t),t)+4*y(t)=Dirac(t-1);

ode10 := diff(y(t),`$`(t,2))+2*diff(y(t),t)+4*y(t) ...

> sol10:=dsolve([ode10,y(0)=0,D(y)(0)=0],y(t));

sol10 := y(t) = -1/3*sqrt(3)*Heaviside(t-1)*exp(1-t...

> plot(rhs(sol10),t=0..10);

[Maple Plot]

>

Additional Exercises

In addition to the exercise in the RLC circuit section, work the following.

Exercises 6(g,j), 7(g,j) on page 275, and 1(c,f) on page 280. Graph all the solutions.