unit32.mws

ORDINARY DIFFERENTIAL EQUATIONS POWERTOOL

Unit 32 -- Using Laplace Tranforms to Solve Initial Value Problems

Prof. Douglas B. Meade

Industrial Mathematics Institute

Department of Mathematics

University of South Carolina

Columbia, SC 29208

URL: http://www.math.sc.edu/~meade/

E-mail: meade@math.sc.edu

Copyright © 2001 by Douglas B. Meade

All rights reserved

-------------------------------------------------------------------

>

Outline of Unit 32

>

Initialization

> restart;

> with( DEtools ):

> with( plots ):

> with( linalg ):

> with( inttrans ):

Warning, the name changecoords has been redefined

Warning, the name adjoint has been redefined

Warning, the protected names norm and trace have been redefined and unprotected

Warning, the name hilbert has been redefined

>

32.A Definition of Laplace Transform

The Laplace transform of a function f defined on the interval (0, infinity ) is

F(s) = Int(exp(-s*t)*f(t),t = 0 .. infinity)

for all s for which the (improper) integral exists.

For example, the Laplace transform of the constant function f(t) = 1 is

> L( 1 ) := Int( exp(-s*t)*1, t=0..infinity );

L(1) := Int(exp(-s*t),t = 0 .. infinity)

>

which Maple evaluates to

> value( L( 1 ) );

Definite integration: Can't determine if the integral is convergent.

Need to know the sign of --> s

Will now try indefinite integration and then take limits.

limit(-(exp(-s*t)-1)/s,t = infinity)

>

Maple is unable to evaluate the limit because the result depends on the value of s . In most cases, the transform variable is assumed to be positive.

> assume( s>0 );

> value( L( 1 ) );

1/s

> unassign( 's' ):

>

The built-in Maple command for evaluating Laplace transforms is laplace , found in the inttrans package.

> laplace( 1, t, s );

1/s

>

Note that the laplace command obtains this result without any assumptions on the transform variable.

> about( s );

s:

  nothing known about this object

>

Similarly, the Laplace transform of an exponential function by the definition

> L( exp('a'*t) ) := Int( exp(-s*t)*exp(a*t), t=0..infinity );

L(exp(a*t)) := Int(exp(-s*t)*exp(a*t),t = 0 .. infi...

> value( L( exp('a'*t) ) );

Definite integration: Can't determine if the integral is convergent.

Need to know the sign of --> s-a

Will now try indefinite integration and then take limits.

limit((exp(-s*t+a*t)-1)/(-s+a),t = infinity)

>

requires some knowledge of the value of the growth/decay rate

> assume( s>a ):

> value( L( exp('a'*t) ) );

1/(s-a)

> unassign( 's', 'a' );

>

The built-in command does not require this additional information.

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

1/(s-a)

>

The laplace command applies many of the properties and identities commonly found in tables of Laplace transforms. For example,

> L( Diff( f(t), t ) ) := laplace( diff( f(t), t ), t, s );

L(Diff(f(t),t)) := s*laplace(f(t),t,s)-f(0)

>

The appearance of this result can be improved with the use of the alias command

> alias( F(s) = laplace( f(t), t, s ) ):

> 'L'( Diff( f(t), t ) ) = L( Diff( f(t), t ) );

L(Diff(f(t),t)) = s*F(s)-f(0)

>

The Maple command for the inverse Laplace transform is invlaplace , also in the inttrans package:

> G := 1/(s^2+3*s+2)^2;

G := 1/((s^2+3*s+2)^2)

> g := invlaplace( G, s, t );

g := t*exp(-2*t)+2*exp(-2*t)+t*exp(-t)-2*exp(-t)

>

This result can also be obtained from the partial fraction decomposition of the function in the transform domain

> convert( G, fullparfrac, s );

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

>

32.B Laplace Transform and Initial Value Problems

The Laplace transform of a linear ODE with initial conditions for an unknown function x = x(t) is an algebraic equation for the transform function X = X(s) . The key is to solve this algebraic equation for X, then apply the inverse Laplace transform to obtain the solution to the IVP. An example demonstrates the step-by-step implementation of this procedure:

> ode1 := diff( x(t), t$2 ) + 2*diff( x(t), t ) + 5*x(t) = sin(3*t);

ode1 := diff(x(t),`$`(t,2))+2*diff(x(t),t)+5*x(t) =...

> ic1 := x(0) = 2, D(x)(0) = 0;

ic1 := x(0) = 2, D(x)(0) = 0

>

The Laplace transform of the ODE is

> alias( X(s) = laplace( x(t), t, s ) ):

> L_ode1 := laplace( ode1, t, s );

L_ode1 := s*(s*X(s)-x(0))-D(x)(0)+2*s*X(s)-2*x(0)+5...

>

The initial conditions result in the algebraic equation

> L_ivp1 := eval( L_ode1, {ic1} );

L_ivp1 := s*(s*X(s)-2)-4+2*s*X(s)+5*X(s) = 3*1/(s^2...

>

As a result, it is seen that the Laplace transform of the solution is a rational function

> q1 := isolate( L_ivp1, X(s) );

q1 := X(s) = (3*1/(s^2+9)+4+2*s)/(s^2+2*s+5)

> X_sol1 := simplify( q1 );

X_sol1 := X(s) = (39+4*s^2+2*s^3+18*s)/((s^2+9)*(s^...

>

the inverse Laplace transform of which is

> x_sol1 := invlaplace( X_sol1, s, t );

x_sol1 := x(t) = -3/26*cos(3*t)-1/13*sin(3*t)+55/26...

>

the solution to the IVP.

Note that if initial conditions are not provided, the Laplace transform still can be used to obtain the general solution:

> Xg_sol1 := isolate( L_ode1, X(s) );

Xg_sol1 := X(s) = (3*1/(s^2+9)+D(x)(0)+2*x(0)+s*x(0...

> q2 := invlaplace( Xg_sol1, s, t );

q2 := x(t) = -3/26*cos(3*t)-1/13*sin(3*t)+exp(-t)*x...
q2 := x(t) = -3/26*cos(3*t)-1/13*sin(3*t)+exp(-t)*x...

>

> xg_sol1 := collect( q2, x );

xg_sol1 := x(t) = (1/2*exp(-t)*sin(2*t)+exp(-t)*cos...
xg_sol1 := x(t) = (1/2*exp(-t)*sin(2*t)+exp(-t)*cos...

>

Note that this result is a solution to the ODE

> odetest( xg_sol1, ode1 );

0

>

and the form of the general solution illustrates how the initial conditions enter into the solution of an initial value problem. In particular, when the initial conditions are substituted into the general solution

> eval( xg_sol1, {ic1} );

x(t) = -3/26*cos(3*t)-1/13*sin(3*t)+55/26*exp(-t)*c...

>

the result is idential to the solution obtained previously :

> x_sol1;

x(t) = -3/26*cos(3*t)-1/13*sin(3*t)+55/26*exp(-t)*c...

>

This process can be automated with the use of dsolve and the optional argument method=laplace :

> infolevel[dsolve] := 3:

> dsolve( { ode1, ic1 }, x(t), method=laplace );

> infolevel[dsolve] := 0:

dsolve/inttrans/solveit:   Transform of eqns is   {s*(s*laplace(x(t),t,s)-2)-4+2*s*laplace(x(t),t,s)+5*laplace(x(t),t,s)-3/(s^2+9)}

dsolve/inttrans/solveit:   Algebraic Solution is   {laplace(x(t),t,s) = (2*s^3+18*s+4*s^2+39)/(s^4+14*s^2+2*s^3+18*s+45)}

x(t) = -3/26*cos(3*t)-1/13*sin(3*t)+55/26*exp(-t)*c...

>

Laplace transforms are also useful when the forcing function is not explicitly given. For example,

> ode2 := diff( x(t), t ) + 3*x(t) = f(t);

ode2 := diff(x(t),t)+3*x(t) = f(t)

> dsolve( ode2, x(t), method=laplace );

x(t) = x(0)*exp(-3*t)+int(f(_U1)*exp(-3*t+3*_U1),_U...

>

Careful inspection of this result reminds us that this first-order linear ODE has as its integrating factor

> mu = intfactor( ode2 );

mu = exp(3*t)

>

Note that the result obtained with the Laplace transform provides more detailed information about the dependence on the initial condition and the forcing function than is obtained with the basic dsolve command.

> infolevel[dsolve] := 3:

> dsolve( ode2, x(t) );

> infolevel[dsolve] := 0:

Methods for first order ODEs:

Trying to isolate the derivative dx/dt...

Successful isolation of dx/dt

-> Trying classification methods

trying a quadrature

trying 1st order linear

1st order linear successful

x(t) = (Int(f(t)*exp(3*t),t)+_C1)*exp(-3*t)

>

[Back to ODE Powertool Table of Contents]

>