Appendix 4.  The Laplace Transform

你4.1 Introduction

    The Laplace transform is one of the most useful tools we'll see in this course for solving differential equations with initial conditions.  Applying a Laplace transform to a differential equation with initial conditions effectively hides the fact that we're differentiating and turns everything into expressions without any derivatives in sight.  We can then manipulate those expressions using little more than high school algebra.  Finally, an inverse Laplace transform takes us back to our original differential equation, except that the algebraic manipulations we did cause the solution of our equation to pop out with little extra work on our part.  The chart below illustrates this process:

 

(1)
Start with differential equation with initial conditions:

y'' + y = et
y(0) = 0
y'(0) = 0

(2)
Take Laplace transform and obtain subsidiary equation:

(s2 + 1)Y(s) = 1/(s - 1)

   
(4)
Use inverse Laplace transform to get back to solution:

y(t) = (1/2)et - (1/2)cos(t) - (1/2)sin(t)

(3)
Solve subsidiary equation for Y:

Y(s) = (1/2)/(s - 1) + (-1/2)s/(s2 + 1)

+ (-1/2)/(s2 + 1)

 

你4.2 The Basics

    In this lab we will be using symbolic variables Y2, s, and t.  So the first thing to do after running MATLAB is to enter

>> syms Y2 s t

You need only do this once (unless you overwrite s or t later, or restart MATLAB).

Example A4.1  We show how to use MATLAB to compute L{cos t}, the Laplace transform of cos t.  The following command tells MATLAB to calculate the Laplace transform of cos t with respect to t, and give the result in terms of s:

>> Y = laplace(cos(t),t,s)

The result is:

Y  =  s/(s^2+1)

MATLAB can also compute inverse Laplace transforms.  (Notice that the order of s and t is reversed in this command!)

>> y = ilaplace(Y,s,t)

Depending on what kind of mood MATLAB is in, you may either get what you expect:

y = cos(t)

or you might get something weird like:

y = 1/2*i*(-i*exp(i*t)-i*exp(-i*t))

If the latter occurs, the command

>> simplify(y)

should tell you that y is in fact equal to cos(t).  (You can check this by using Euler's formula, eit = cos t + i sin t.)

Exercise A4.1
Include input and output. (a)  Use MATLAB to compute the Laplace transforms L{1} and L{sin2(t)}.  (If MATLAB gives you something like a/b/c this means a/(b*c).)
Include full-sentence response. Do work by hand. (b)  Explain how, using the results from part (a), we can at this point compute L{(cos2(t)) without using MATLAB or doing an integral.  What do you get?  (Hint: use the trig identity sin2(t) + cos2(t) = 1).
Include input and output. Do work by hand. (c)  Now use MATLAB to calculate L{cos2(t)}, and verify your answer to part (b).  You may have to do some additional scratch work to see that your answers match.
Include input and output. (d)  Use MATLAB to compute the inverse Laplace transform of 4/(9+s2).  If your answer looks complicated, try using the simplify command on it.

你4.3 Heaviside

    The Laplace transform is a powerful tool, but we do have to acknowledge that it applies only in limited situations.  Namely, we can't use the Laplace transform to solve differential equations unless such equations are linear equations with constant coefficients.  Furthermore, the solution method is most beneficial when we have initial conditions.

    As explained above, the Laplace transform, in a sense, "takes the calculus out" of solving differential equations, and thus may be a good technique for solving differential equations by hand.  Of course, we're interested in ways MATLAB can be of use.  If all we deal with are linear equations with constant coefficients, then we know ode45 can already do the job.  So why go to the trouble of using MATLAB to accomplish the four steps of the Laplace transform outlined in the beginning of this assignment?

    The answer is that ode45 is often not able to solve non-homogeneous equations where the non-homogeneous part is a discontinuous function.  In fact, such equations are even more resistant to the methods you learned in the lectures for solving them by hand.  One of the greatest accomplishments of the Laplace transform is that it provides us a relatively easy way to solve the types of equations that come up in real-world applications where functions aren't always continuous.

Example A4.2  Imagine we have an electric circuit to which we attach a power source (like a battery) and a switch that either lets current through or prevents it from flowing depending on whether the switch is closed or open.  We'll start with the switch open so that no current can pass through it.

    Now when we close the switch, current generated by the battery immediately begins to flow.  The electromotive force in such a situation might look like:

Of course, this function is discontinuous, but nevertheless it models a situation which occurs all the time in the study of electric circuits.

    The function illustrated above is called the Heaviside function, after the English engineer and mathematician Oliver Heaviside.  Due to inconsistencies between versions of MATLAB, we're going to create an m-file to define the Heaviside function instead of using the built-in Heaviside function.

    As in the last assignment, make sure that you are saving your files into MATLAB's current directory.  Open up a new M-file and type

1 function y=Heaviside(x)
2 y=(x>=0);

Remark A4.1  True statements have the value 1 and false statements have the value 0.  Therefore, this function just says that if x is nonnegative then the function evaluates to one, and not, it evaluates to zero, exactly as it should.

    Save the M-file as "Heaviside.m" .

Exercise A4.2 
Include input and output. Include plots and graphs. (a)  Test your M-file by typing

>> f='Heaviside(t)'
>> ezplot(f,[-5,5])

Include your input and the resulting graph in your Word document.  Note that even though the function is discontinuous, ezplot tries to draw it in a continuous manner by connecting the two pieces together.  This won't have any effect on our calculations though.

Include input and output. Include plots and graphs. Include full-sentence response. (b)  Note that the Heaviside function models a situation where the force is applied at time zero.  In what follows below we will want to be able to apply such a force at any time.  Try the commands

>> f='Heaviside(t+1)'
>> ezplot(f,[-5,5])

Include your input and the resulting graph.  At what value of t does the function jump from 0 to 1 in this case?

Include input and output. Include plots and graphs. Include full-sentence response. (c)  What variation of the Heaviside function should we use if we wanted the function to jump at t = 4 instead?  Test your answer with MATLAB by graphing it as in the above two exercises.

你4.4 Electrical Circuit

In this section we will work with a differential equation that models the following electrical circuit.

 

A differential equation that describes charge in a capacitor is given by the following formula

Ly" + Ry' + (1/C)y = E(t)

where L is the inductance, R the resistance, C the capacitance, and E(t) the electromotive force.  (This is described in detail in the textbook in section 3.8.)  The example that we will explore will use the following characteristics:

L = 3

R = 1

C = 1/5

E(t)   =

{

1,   40t50

0,   otherwise

Remark A4.2  The electromotive force function tells us that the battery is supplying electricity to the circuit only for the time period 40t50.  The battery is disconnected at all other times (i.e., 0t < 40 and 50 < t).

    After we substitute the constants in we get the following equation:

3y'' + y' + 5y = E(t)

where the capacitor has an initial charge of 2 and the initial current is 0 (i.e., y(0) = 2, y'(0) = 0).

    Before we can solve this equation, we have to tweak the Heaviside function so that MATLAB will be able to recognize our E(t).  Here is what we have to do.

Example A4.3 

    You already know how to shift the Heaviside function.  What if we try to subtract two Heaviside functions?  Let's type the following into MATLAB

>> d = 'Heaviside(t)-Heaviside(t-1)'
>> ezplot(d,[-5,5,-1,2])

Your graph should look like this

Remark A4.3  All we have done in the above example is subtract two functions.  We can also add and multiply these functions (dividing is a bit tricky since we can't divide by 0), which will result in a variety of piecewise defined functions.

Exercise A4.3
Include input and output. Include plots and graphs. (a)  Using Example A3.3 as a template, define a function h in terms of the Heaviside function that will describe the following behavior:

Now use the following command to graph your function.

>> ezplot(h,[-3,3,-0.5,1.5])

 

Include input and output. Include plots and graphs.

(b)  Now define a function g which when graphed looks like this:

 

(Hint: sometimes it helps to go one step at a time adding the results at the end. Remember to use parentheses to denote the order of operations.)

Use the following command to graph your function.

>> ezplot(g, [-3,3,-0.5,1.5])

 

    We are now ready to solve our differential equation. Follow the procedure outlined below.

Step (1).  We need to input our differential equation.  To simplify the commands, we'll start off by assigning our differential equation to a variable:

>> ode='3*D(D(y))(t)+D(y)(t)+5*y(t)=Heaviside(t-40)-Heaviside(t-50)'

Step (2).  Now take the Laplace transform.  Again, for convenience, we assign the result to a new variable:

>> LTode=laplace(ode,t,s)

LTode =

3*s*(s*laplace(y(t),t,s)-y(0))-3*D(y)(0)+s*laplace(y(t),t,s)-y(0)+5*laplace(y(t),t,s) = exp(-40*s)/s-exp(-50*s)/s
 

    Now we need to substitute our initial conditions.  While we're at it, we'll also substitute the variable Y2 for MATLAB's awkward output laplace(y(t),t,s).  All of this can be done in one slick step.

>> LTode=subs(LTode,{'laplace(y(t),t,s)','y(0)','D(y)(0)'},{Y2,2,0})

LTode =

3*s*(s*Y2-2)-2+s*Y2+5*Y2=exp(-40*s)/s-exp(-50*s)/s
 

Step (3).  The above is our subsidiary equation which we now solve for Y2.

>> Y2=solve(LTode,Y2)

Y2 =

(6*s^2+2*s+exp(-40*s)-exp(-50*s))/s/(3*s^2+s+5)
 

Step (4).  Now we use the inverse Laplace transform to get back from Y(s) to y(t).

>> y2=ilaplace(Y2,s,t)

y2 =

Heaviside(t-40)*(-1/5*exp(-1/6*t+20/3)*cos(1/6*59^(1/2)*(t-40))-1/295*59^(1/2)*exp(-1/6*t+20/3)*sin(1/6*59^(1/2)*(t-40))+1/5)+(1/5*exp(-1/6*t+25/3)*cos(1/6*59^(1/2)*(t-50))+1/295*59^(1/2)*exp(-1/6*t+25/3)*sin(1/6*59^(1/2)*(t-50))-1/5)*Heaviside(t-50)+2*exp(-1/6*t)*cos(1/6*59^(1/2)*t)+2/59*59^(1/2)*exp(-1/6*t)*sin(1/6*59^(1/2)*t)
 

(It's a bit messy.  You can see the advantage of doing this on a computer.)

Include input and output. Include plots and graphs. Include full-sentence response. Exercise A4.4  Plot the function y2 with the commands

>> y2=char(y2)
>> fplot(y2,[0,70,-3,3])

Looking at our initial value problem, try to explain the shape that results from graphing the function y2.  Why does the function tend toward zero, then jump up a little, and then flatten again?  (Hint: look at the initial conditions and remember that capacitors hold charge for some time after the electricity has been disconnected.)  Include all work done for this problem in your Word document.

你4.5 Conclusion

    The Laplace transform is an effective computational tool for solving linear equations with constant coefficients in the presence of initial conditions.  MATLAB helps to sort out all the tedious algebra involved in using this method.  When solving non-homogeneous equations involving discontinuous functions, the Laplace transforms is the only tool we have at our disposal and, again, MATLAB makes the job much simpler.


Last Modified:

04/15/2008