This section demonstrates by numerous examples how to apply the Laplace transform for solving the initial value problems for constant coefficient linear differential equations.

Solving IVPs with Laplace transform

The Laplace transform is an alternative approach to the methods for solving initial value problems of linear differential equations with constant coefficients. These were considered in Part IV of this tutorial. The Laplace transform is useful in dealing with discontinuous inputs (closing of a switch) and with periodic functions (sawtooth and rectified waves). Analysis of the effect of such inputs proceeds most smoothly in the frequency domain, that is, in domain of the transform-variable, which we denote by λ.

Application of the Laplace transformation to differential equations is based on the following statements.

Theorem 1: Suppose that f: [0, ∞) → ℝ is a continuous and its derivative f' = df/dt is piecewise continuous on any finite interval 0 ≤ t ≤ b < ∞. Suppose further that there exist constants K, γ, and M such that \( |f(t)| \le K\,e^{\gamma t} \) for t > M. Then the Laplace transform of this function, \( {\cal L}\left[ f(t) \right] (\lambda ) = f^L (\lambda ) \) exists for λ > γ, and more over,
\begin{equation} \label{EqIVP.1} {\cal L}\left[ f'(t) \right] (\lambda ) = \int_0^{\infty} f' (t)\, e^{-\lambda t} {\text d}t = \lambda\,f^L (\lambda ) - f(0^{+}) , \end{equation}
where \( f(0^{+}) = f(0+0) = \lim_{\varepsilon \downarrow 0} f(\varepsilon ) \) is the limit value of f at zero from right.
Theorem 1 states that the Laplace transformation is a spectral representation of the derivative operator \( \texttt{D} = {\text d}/{\text d}t , \) acting in the space of functions defined on semi-line [0, ∞). This means that the Laplace transform maps the unbounded operator \( \texttt{D} \) into multiplication by a spectral parameter λ on a set of functions that vanish at the origin.

Theorem 1 can be extended for arbitrary power of the derivative operator.

Theorem 2: Suppose that f: [0, ∞) → ℝ is continuous along with its derivatives up to the order n−1, and that its n-th derivative \( f^{(n)} \) is piecewise continuous on any finite interval 0 ≤ t ≤ b < ∞. Suppose further that there exist constants K, γ, and M such that \( |f(t)| \le K\,e^{\gamma t} , \quad |f'(t)| \le K\,e^{\gamma t} , \quad \ldots , |f^{(n-1)} (t)| \le K\,e^{\gamma t} \) for t > M. Then the Laplace transform of this function, \( {\cal L}\left[ f(t) \right] (\lambda ) = f^L (\lambda ) \) exists for λ > γ, and moreover,
\begin{equation} \label{EqIVP.2} {\cal L}\left[ f^{(n)}(t) \right] (\lambda ) = \int_0^{\infty} f^{(n)} (t)\, e^{-\lambda t} {\text d}t = \lambda^n f^L (\lambda ) - f^{(n-1)}(0^{+}) - f^{(n-2)}(0^{+}) - \cdots - \lambda^{n-1} f(0^{+}) . \end{equation}

Since the Laplace transform rolls over an unbounded linear differential operator with constant coefficients,

\begin{equation} \label{EqIVP.3} L \left[ \,\texttt{D} \,\right] = a_n \texttt{D}^n + a_{n-1} \texttt{D}^{n-1} + \cdots + a_1 \texttt{D} + a_0 \texttt{I} , \end{equation}
where \( \texttt{I} = \texttt{D}^0 \) is the identity operator, into multiplication by the polynomial \( L \left[ \lambda \right] = a_n \lambda^n + a_{n-1} \lambda^{n-1} + \cdots + a_1 \lambda + a_0 , \) it reduces an initial value problem into an algebraic equation. Upon solving this algebraic equation, we obtain almost immediately the Laplace transform of the unknown function---the solution of the initial value problem. There are no miracles in math, and the price you have to pay for using the beautiful operating method is hidden in the inverse Laplace transform, which is an ill-posed operation. Fortunately, the methods of determination of the inverse Laplace transformation described previously work pretty well in initial value problems for ordinary differential equations.

Schematically, this process can be illustrated starting with, for example, a second order linear differential equation with constant coefficients, as follows.

\begin{eqnarray*} a\,y'' + b\, y' + c\,y &=& f (t) \qquad y(0) = d, \quad y' (0) = v, \\ &\Downarrow & \quad \mbox{apply } {\cal L}\\ {\cal L} \left[ a\,y'' + b\, y' + c\,y \right] &=& {\cal L} \left[ f (t) \right] \\ &\Downarrow & \quad \mbox{solve for } {\cal L}[y] = y^L \\ y^L \equiv {\cal L} \left[ y (t) \right] &=& F(\lambda ) = P(\lambda )/Q(\lambda ), \qquad \mbox{which is a ratio} \\ &\Downarrow & \quad \mbox{apply } {\cal L}^{-1} \\ y(t) &=& {\cal L}^{-1} \left[ F(\lambda ) \right] . \end{eqnarray*}
Of course, some details need to be addressed for this to make sense. In particular, we need to express the left side \( {\cal L} \left[ a\,y'' + b\, y' + c\,y \right] \) in terms of \( {\cal L} \left[ y \right] = y^L . \) We illustrate the method with numerous examples. Here is an example of how you may want to utilize the Laplace transformation to solve the differential equation y'' + 4y = 0:
      We solve and plot with matlab the initial value problem: y'' + 4y = 0, y(0) = 1, y'(0) = 2.
>> y =
cos(2*t) + sin(2*t)
       Graph of the solution.            matlab code

Example 1: Consider the homogeneous differential equation

\[ y'' - y' -2\,y =0 \tag{1.1} \]
and the initial conditions
\[ y(0) = 2, \qquad y' (0) = -1. \tag{1.2} \]
Applying the Laplace transform to Eq.(1.1), we get
\[ {\cal L} \left[ y'' (t) \right] (\lambda ) - {\cal L} \left[ y' (t) \right] (\lambda ) - 2\, {\cal L} \left[ y (t) \right] = 0 . \tag{1.3} \]
Upon using Theorem 2 to express \( {\cal L} \left[ y'' (t) \right] \) and \( {\cal L} \left[ y' (t) \right] \) in terms of \( {\cal L} \left[ y (t) \right] = y^L (\lambda ) , \) we find that the previous equation (1.3) becomes
\[ \lambda^2 y^L - y' (0) - \lambda\, y(0) - \lambda\, y^L + y(0) - 2\,y^L = 0 , \]
or
\[ \left( \lambda^2 - \lambda - 2 \right) y^L = y' (0) + \lambda\, y(0) - y(0) . \]
Solving for yL and using the initial conditions (1.2), we obtain
\[ y^L (\lambda ) = \frac{2\,\lambda -3}{\lambda^2 - \lambda -2} = \frac{2\,\lambda -3}{\left( \lambda -2 \right)\left( \lambda +1 \right)} . \tag{1.4} \]
Since the denominator is a product of simple terms, the residue theorem gives
\[ y(t) = {\cal L}^{-1} \left[ y^L (\lambda ) \right] = \mbox{Res}_{\lambda = 2} \frac{2\,\lambda -3}{\left( \lambda -2 \right)\left( \lambda +1 \right)}\, e^{\lambda t} + \mbox{Res}_{\lambda = -1} \frac{2\,\lambda -3}{\left( \lambda -2 \right)\left( \lambda +1 \right)}\, e^{\lambda t} . \]
Calculations show
\begin{align*} \mbox{Res}_{\lambda = 2} \frac{2\,\lambda -3}{\left( \lambda -2 \right)\left( \lambda +1 \right)}\, e^{\lambda t} &= \lim_{\lambda \to 2} \frac{2\,\lambda -3}{\lambda +1}\, e^{\lambda t} = \frac{1}{3}\, e^{2t} , \\ \mbox{Res}_{\lambda = -1} \frac{2\,\lambda -3}{\left( \lambda -2 \right)\left( \lambda +1 \right)}\, e^{\lambda t} &= \lim_{\lambda \to -1} \frac{2\,\lambda -3}{\lambda -2} \, e^{\lambda t} = \frac{5}{3} \, e^{-t} . \end{align*}
Therefore, the solution of the given initial value problem (1.1), (1.2) becomes
\[ y(t) = \frac{1}{3} \left[ e^{2t} + 5 \, e^{-t} \right] H(t) , \]
where H(t) is the Heaviside function.
Solution of the given homogeneous equation.
   ■
The same procedure can be applied to the general second order linear differential equation with constant coefficients
\begin{equation} \label{EqIVP.4} a\,y'' + b\,y' + c\, y = 0, \qquad y(0) = y_0, \quad y' (0) = y_1 \end{equation}
Assuming that the solution y = ϕ(t) satisfies the conditions of Theorem 2 for n = 2, we can take the Laplace transformation and obtain
\[ a \left[ \lambda^2 y^L - y' (0) - \lambda y(0) \right] + b \left[ \lambda\, y^L - y(0) \right] + c\,y^L = 0. \]
Solving for yL, the Laplace transform of the solution to the IVP \eqref{EqIVP.4}, we get
\begin{equation} \label{EqIVP.5} y^L (\lambda ) = \frac{a\lambda\, y(0) + b\,y(0) + b\,y' (0)}{a \lambda^2 + b\,\lambda + c} . \end{equation}

Example 2: Consider the undamped mechanical oscillator with a forcing function that is a constant f(t) = F0. Recall that it consists of a block of mass m on a table and restrained laterally by an ordinary coil spring. The displacement, denoted as x(t), of the mass (measured as positive to the right) from its equilibrium position: that is, when x = 0 the spring is neither stretched nor compressed. We are now imagining the mass to be disturbed from its equilibrium position by an initial disturbance (position s and velocity v) and applied force f(t).

Since the motion of the mass is constrained to be along a straight line, we need to consider only the forces in the x direction. Neglecting aerodynamic force and the friction force, we are left to consider the only one external (driving) force. Then according to Newton's second law and Hooke's law, the displacement x(t) satisfies the differential equation

\[ m\, x'' + k\, x = f(t) = F_0 , \tag{2.1} \]
where k is a spring constant. Of course, we need to add the initial conditions
\[ x(0) = s, \qquad x' (0) = v , \tag{2.2} \]
In order to apply the Laplace transform, we multiply each term in Eq.(2.1) by the Laplace kernel \( e^{-\lambda t} \) and integrate on t from 0 to ∞. Using the standard notation for the Laplace transform, we get
\[ {\cal L} \left[ m\, x'' + k\, x \right] = {\cal L} \left[ f(t) \right] = {\cal L} \left[ F_0 \right] . \]
Using the linearity of the Laplace transform (which is actually the property of integration), we can rewrite the latter as
\[ {\cal L} \left[ m\, x'' \right] + {\cal L} \left[ k\, x \right] = {\cal L} \left[ F_0 \right] \qquad \mbox{or} \qquad m\, {\cal L} \left[ x'' (t) \right] + k\, {\cal L} \left[ x(t) \right] = {\cal L} \left[ F_0 \right] = \frac{F_0}{\lambda} \]
because the Laplace transform of a constant is just a reciprocal of λ. Recalling Theorem 2, we have
\[ m \left[ \lambda^2 x^L - x'(0) - \lambda\, x(0) \right] - k\, x^L = \frac{F_0}{\lambda} , \tag{2.3} \]
where we denote by xL the Laplace transform of the unknown solution x(t) to the IVP (2.1), (2.2). The point to appreciate is that whereas in the t domain we have the unbounded linear differential operator acting on x(t), in the transform domain we now have the linear algebraic equation (2.3) on xL. By simple algebra, we find its solution almost immediately:
\[ {\cal L} \left[ x(t) \right] = x^L (\lambda ) = \frac{\lambda x(0) + x' (0)}{\lambda^3 + \omega^2} + \frac{1}{m \left( \lambda^2 + \omega^2 \right)}\,{\cal L} \left[ F_0 \right] , \qquad \omega^2 = \frac{k}{m} . \tag{2.4} \]
We now invert the Laplace transform to obtain the required solution
\begin{align*} x(t) &= {\cal L}^{-1} \left[ \frac{\lambda x(0) + x' (0)}{\lambda^2 + \omega^2} \right] + {\cal L}^{-1} \left[ \frac{1}{m \left( \lambda^2 + \omega^2 \right)}\,\frac{F_0}{\lambda} \right] \\ &= x_h (t) + x_p (t) , \end{align*}
where
\[ x_h (t) = {\cal L}^{-1} \left[ \frac{\lambda x(0) + x' (0)}{\lambda^2 + \omega^2} \right] = x(0) {\cal L}^{-1} \left[ \frac{\lambda}{\lambda^2 + \omega^2} \right] + x' (0) {\cal L}^{-1} \left[ \frac{1}{\lambda^2 + \omega^2} \right] \]
is the solution of the IVP for the homogeneous equation (usually is referred to as the complementary function)
\[ m\,x'' + k\, x(t) = 0, \qquad x(0) = s, \quad x' (0) = v . \]
Another component
\[ x_p (t) = {\cal L}^{-1} \left[ \frac{1}{m \left( \lambda^2 + \omega^2 \right)}\,\frac{F_0}{\lambda} \right] \]
is a particular solution of the inhomogeneous equation
\[ m\,x'' + k\, x(t) = F_0, \qquad x(0) = 0, \quad x' (0) = 0 . \]
To determine the explicit expression for xp(t), we have two options. First, we use the explicit formula for the Laplace transform of the forcing function and get
\[ x_p (t) = 2\,\Re\,\mbox{Res}_{\lambda = \omega{\bf j}} \left[ \frac{e^{\lambda t}}{m \left( \lambda^2 + \omega^2 \right)}\,\frac{F_0}{\lambda} \right] + \mbox{Res}_{\lambda = 0} \left[ \frac{e^{\lambda t}}{m \left( \lambda^2 + \omega^2 \right)}\,\frac{F_0}{\lambda} \right] . \]
Here Re z = \( \Re z = \Re \left( a + b{\bf j} \right) = a , \) is the real part of complex number 𝑎 + jb, and j is the unit vector in positive vertical direction on the complex plane ℂ, so j² = −1. Since both singularities λ = ωj and λ = 0 are simple, we just differentiate the denominator and set the value of λ to be the corresponding singularity. So
\begin{align*} ,\mbox{Res}_{\lambda = \omega{\bf j}} \left[ \frac{e^{\lambda t}}{m \left( \lambda^2 + \omega^2 \right)}\,\frac{F_0}{\lambda} \right] &= \lim_{\lambda \to \omega{\bf j}} \, \frac{e^{\lambda t}}{2m\lambda^2} \,F_0 = - \frac{F_0}{2k} \, e^{{\bf j} \omega t} , \\ \mbox{Res}_{\lambda = 0} \left[ \frac{e^{\lambda t}}{m \left( \lambda^2 + \omega^2 \right)}\,\frac{F_0}{\lambda} \right] &= F_0 \lim_{\lambda \to 0} \frac{e^{\lambda t}}{m \left( \lambda^2 + \omega^2 \right)} = \frac{F_0}{k} . \end{align*}
Taking the real part, we obtain
\[ x_p (t) = \frac{F_0}{k} \left( 1 - \cos \omega t \right) H(t) . \]
Note that we always need to multiply the function obtained with the inverse Laplace transform by the Heaviside function H(t).

Another option is to use the convolution rule:

\[ x_p (t) = G(t) \ast f(t) = \int_0^t G\left( t - \tau \right) f(\tau )\,{\text d}\tau = \int_0^t G\left( \tau \right) f(t-\tau )\,{\text d}\tau = f(t) \ast G(t) , \]
where G(t) is the Green function for the harmonic oscillator:
\[ G(t) = {\cal L}^{-1} \left[ \frac{1}{m \left( \lambda^2 + \omega^2 \right)} \right] = {\cal L}^{-1} \left[ \frac{1}{m\lambda^2 + k} \right] . \]
We find the Green function using the residue formula:
\[ G(t) = 2\,\Re \,\mbox{Res}_{\lambda = \omega{\bf j}} \left[ \frac{e^{\lambda t}}{m \left( \lambda^2 + \omega^2 \right)} \right] = 2\,\Re \,\lim_{\lambda \to {\bf j}\omega} \,\frac{e^{\lambda t}}{2m\lambda} = 2\,\Re \,\frac{e^{{\bf j} \omega t}}{2m\,{\bf j}\omega} . \]
Taking the real part, we get
\[ G(t) = \frac{1}{m}\,\frac{\sin \omega t}{\omega} . \]
So
\[ x_p (t) = F_0 \ast \frac{1}{m}\,\frac{\sin \omega t}{\omega} = \frac{F_0}{m}\,\int_0^t \frac{\sin \omega \tau}{\omega} \,{\text d}\tau = \frac{F_0}{m}\,\frac{1 - \cos \omega t}{\omega^2} = \frac{F_0}{k} \left( 1 - \cos \omega t \right) , \]
which we need to multiply by the Heaviside function.

The complementary function xh(t) can be obtained using the table of known Laplace transformations:

\[ x_h (t) = \left[ x(0)\,\cos \omega t + \frac{x' (0)}{\omega}\,\sin \omega t \right] H(t) . \]
When the initial conditions are not specified, the solution formula
\[ x(t) = x_h (t) + x_p (t) = \left[ x(0)\,\cos \omega t + \frac{x' (0)}{\omega}\,\sin \omega t + \frac{F_0}{k} \left( 1 - \cos \omega t \right) \right] H(t) \tag{2.5} \]
gives the general solution because x(0) and x'(0) can be arbitrary. When we choose some numerical values, say s = 2 and v = −1, we get
\[ x(t) = \left[ 2\,\cos \omega t - \frac{\sin \omega t}{\omega} + \frac{F_0}{k} \left( 1 - \cos \omega t \right) \right] H(t) . \]
      We plot the solution with the following parameters:
ω = 3,
k = 1,
F0 = −2.
       Harmonic oscillator.            matlab code
   ■

Example 3: Consider the initial value problem

\[ y'' + 5\, y' +6\,y = \sin 2t , \qquad y(0) =-1, \quad y' (0) =2 . \tag{3.1} \]
Applying the Laplace transform, we get
\[ \lambda^2 y^L - y' (0) - \lambda\,y(0) + 5\,\lambda \, y^L - 5\,y(0) +6\,y^L = \frac{2}{\lambda^2 +4} , \]
>> y =
sin(2*t)/52 - (3*exp(-2*t))/4 - (2*exp(-3*t))/13 - (5*cos(2*t))/52
test =
sin(2*t)
ans =
0

where

\[ y^L ( \lambda ) = {\cal L}\left[ y(t) \right] = \int_0^{\infty} e^{-\lambda\,t} \,y(t)\,{\text d} t \]
is the Laplace transform of the unknown solution y(t). Solving the algebraic equation above, we find
\[ y^L ( \lambda ) = -\frac{\lambda +3}{\lambda^2 +5\lambda +6} + \frac{1}{\lambda^2 +5\lambda +6} \cdot \frac{2}{\lambda^2 +4} = \frac{-1}{\lambda +2} + \frac{1}{\lambda^2 +5\lambda +6} \cdot \frac{2}{\lambda^2 +4} . \tag{3.2} \]
Since it is the sum of two fractions, we apply the inverse Laplace transform to each of them
\begin{eqnarray*} y_h (t) &=& {\cal L}^{-1} \left[ \frac{-1}{\lambda +2} \right] = \mbox{Res}_{\lambda = -2} \, \frac{-1}{\lambda +2} \, e^{\lambda\, t} \\ &=& - e^{-2t} , \end{eqnarray*}
which we have to multiply by the Heaviside function.
>> p(s) = exp(s*t)*(s + 3)
q(s) = s^2 + 5*s + 6
>> ans = (5*exp(2*t))/9
ans = (6*exp(3*t))/11
---->

The denominator in the second term has two real nulls \( \lambda =-2 \) and \( \lambda =-3 \) and two complex conjugate \( \lambda = \pm 2{\bf j} . \) We calculate residues at each pole separately.

\begin{eqnarray*} \mbox{Res}_{\lambda = -2} \, \frac{e^{\lambda\,t}}{\lambda^2 +5\lambda +6} \cdot \frac{2}{\lambda^2 +4} &=& \lim_{\lambda \to -2} \, \frac{e^{\lambda\,t}}{2\lambda +5} \cdot \frac{2}{\lambda^2 +4} = \frac{1}{4}\, e^{-2t} , \\ \mbox{Res}_{\lambda = -3} \, \frac{e^{\lambda\,t}}{\lambda^2 +5\lambda +6} \cdot \frac{2}{\lambda^2 +4} &=& \lim_{\lambda \to -3} \, \frac{e^{\lambda\,t}}{2\lambda +5} \cdot \frac{2}{\lambda^2 +4} = -\frac{2}{13} \, e^{-3t} , \\ \mbox{Res}_{\lambda = 2{\bf j}} \, \frac{e^{\lambda\,t}}{\lambda^2 +5\lambda +6} \cdot \frac{2}{\lambda^2 +4} &=& \lim_{\lambda \to 2{\bf j}} \,\frac{e^{\lambda\,t}}{\lambda^2 +5\lambda +6} \cdot \frac{2}{2\,\lambda} \\ \frac{1}{2(5{\bf j} +1)} \cdot \frac{1}{2{\bf j}} \, e^{2{\bf j}t} = -\frac{1}{4} \left( \frac{5}{52} + \frac{\bf j}{52} \right) e^{2{\bf j}t} . \end{eqnarray*}

Extracting the real part from the latter and multiplying it by 2, we obtain the particular solution

\[ y_p (t) = {\cal L}^{-1} \left[ y_p^L \right] = {\cal L}^{-1} \left[ \frac{1}{\lambda^2 +5\lambda +6} \right] = \left[ \frac{\sin 2t}{52} -\frac{5\,\cos 2t}{52} + \frac{1}{4}\, e^{-2t} - \frac{2}{13} \, e^{-3t} \right] H(t) . \tag{3.3} \]

\[ y_p (t) = {\cal L}^{-1} \left[ y_p^L \right] = {\cal L}^{-1} \left[ \frac{1}{\lambda^2 +5\lambda +6} \right] = \left[ \frac{\sin 2t}{52} -\frac{5\,\cos 2t}{52} + \frac{1}{4}\, e^{-2t} - \frac{2}{13} \, e^{-3t} \right] H(t) . \]
We verify that the function (3.3) satisfies the differential equation (3.1) and homogeneous initial conditions.
Finally, we plot the solution.
Solution of the given IVP (3.1)..

We plot the solution.
   ■

Example 4: Consider the initial value problem

\[ y'' -3\, y' -4\,y = -16\,t , \qquad y(0) =-4, \quad y' (0) =-5 . \tag{4.1} \]
Its solution is the sum \( y(t) = y_h (t) + y_p (t) \) of two functions, one of them yh, is the solution of the homogeneous equation subject to the given initial conditions
\[ y'' -3\, y' -4\,y = 0 , \qquad y(0) =-4, \quad y' (0) =-5 . \]
Using the general formula, we find the Laplace transform of it:
\[ y_h^L = \frac{a \left( \lambda\,y(0) + y' (0) \right) + b\,y(0)}{a\lambda^2 + b\lambda +c} = \frac{7-4\lambda}{\lambda^2 -3\lambda -4} . \]
Another function, yp is the particular solution of the inhomogeneous equation subject to homogeneous initial conditions
\[ y'' -3\, y' -4\,y = -16\,t , \qquad y(0) =0, \quad y' (0) =0 . \]
Its Laplace transform is
\[ y_p^L = \frac{1}{\lambda^2 -3\lambda -4} \cdot {\cal L} \left[ -16\,t \right] = - \frac{16}{\lambda^2 \left( \lambda^2 -3\lambda -4 \right)} . \]
Therefore, to determine the required solution, we need to find the inverse Laplace transforms of these functions and add them. We start with yh:
\begin{eqnarray*} y_h (t) &=& {\cal L}^{-1} \left[ \frac{7-4\lambda}{\lambda^2 -3\lambda -4} \right] = \left( \mbox{Res}_{\lambda = 4} + \mbox{Res}_{\lambda = -1} \right) \frac{7-4\lambda}{\lambda^2 -3\lambda -4} \, e^{\lambda\, t} \\ &=& \left. \frac{7-4\lambda}{2\lambda -3} \, e^{\lambda\, t} \right\vert_{\lambda = 4} + \left. \frac{7-4\lambda}{2\lambda -3} \, e^{\lambda\, t} \right\vert_{\lambda = -1} \\ &=& -\frac{9}{5}\, e^{4t} - \frac{11}{4}\, e^{-t} , \end{eqnarray*}
which we have to multiply by the Heaviside function. Now we find the inverse Laplace transform of the second function:
\begin{eqnarray*} y_p (t) &=& {\cal L}^{-1} \left[ - \frac{16}{\lambda^2 \left( \lambda^2 -3\lambda -4 \right)} \right] = \left( \mbox{Res}_{\lambda = 4} + \mbox{Res}_{\lambda = -1} + \mbox{Res}_{\lambda = 0} \right) \frac{-16}{\lambda^2 \left( \lambda^2 -3\lambda -4 \right)} \, e^{\lambda\, t} \\ &=& \left. \frac{-16}{\lambda^2 \left( 2\lambda -3 \right)} \, e^{\lambda\, t} \right\vert_{\lambda = 4} + \left. \frac{-16}{\lambda^2 \left( 2\lambda -3 \right)} \, e^{\lambda\, t} \right\vert_{\lambda = -1} + \lim_{\lambda \to 0} \,\frac{\text d}{{\text d} \lambda} \frac{-16}{\lambda^2 -3\lambda -4} \, e^{\lambda\, t} \\ &=& \frac{1}{80} \, e^{4t} - \frac{1}{5}\, e^{-t} + \lim_{\lambda \to 0} \, \frac{16 \left( 2\lambda -3 \right)}{\left( \lambda^2 -3\lambda -4 \right)^2} + + \lim_{\lambda \to 0} \, \frac{-16 \,t}{\lambda^2 -3\lambda -4} \, e^{\lambda\, t} \\ &=& \frac{1}{80} \, e^{4t} - \frac{1}{5}\, e^{-t} + \frac{3-4t}{4} . \end{eqnarray*}
Now we add these two functions to obtain the required solution of the given initial value problem:
\[ y(t) = y_h (t) + y_p (t) = \left[ \frac{3-4t}{4} - \frac{143}{80} \, e^{4t} - \frac{59}{20}\, e^{-t} \right] H(t) . \]
Finally, we plot the solution.
Solution of the given IVP (4.1)..
   ■
With examples completed, we can make some observations about the application of the Laplace transform method. Although we considered only second order differential equations, the method can be easily extended to linear differential equations of any order.
  • Upon application of the Laplace transformation, the initial conditions become "build-in."
  • When applying the Laplace transform, we by default assume that the unknown function and all its derivatives are transformable under the Laplace method into holomorphic functions on the half-plane Reλ > γ. However, such restrictions in no way impede the solution steps in finding the Laplace transform of the solution function. If you successfully apply the inverse Laplace transform to obtain the solution, they are no longer relevant.
  • When solving the inhomogeneous equation
    \[ L \left[ \,\texttt{D} \,\right] y (t) = f(t) \]
    for n-th order differential operator \eqref{EqIVP.3} with the aid of the Laplace transform method, we always reduce the corresponding initial value problem to an algebraic equation
    \[ L \left[ \,\lambda \,\right] y^L (\lambda ) = P(\lambda ) + f^L (\lambda ). \]
    Finding y(t) involves determination of the inverse Laplace transform of the expression
    \[ {\cal L}^{-1} \left[ \frac{1}{L \left[ \,\lambda \,\right]}\, f^L (\lambda ) \right] \]
    If the Laplace transform fL of the input function f(t) is not available, the only way to determine the required inverse Laplace transform is to apply the convolution property:
    \begin{equation} \label{EqIVP.6} {\cal L}^{-1} \left[ \frac{1}{L \left[ \,\lambda \,\right]}\, f^L (\lambda ) \right] = G(t) \ast f(t) = \int_0^t G(t-\tau )\,f(\tau )\,{\text d}\tau = f(t) \ast G(t) , \end{equation}
    where G(t) is the Green function. In this case, you don't need to find the Laplace transform of the input function f(t). However, determination of the convolution integral is an ill-posed problem and its practical implementation may lead to some sort of instability.

 

Green Functions

Definition: Given a linear differential operator \eqref{EqIVP.3}, its Green function is
\begin{equation} \label{EqIVP.7} G(t) = {\cal L}^{-1} \left[ \frac{1}{L \left[ \,\lambda \,\right]} \right] . \end{equation}
When the Green function for the linear constant coefficient differential operator \eqref{EqIVP.3} is known, then a particular solution to the inhomogeneous equation
\[ L \left[ \,\texttt{D} \,\right] y(t) = f(t) \]
is obtained by the convolution integral \eqref{EqIVP.6}. The Green function for the linear constant coefficient operator \eqref{EqIVP.3} is the unique solution of the following initial value problem:
\[ L \left[ \,\texttt{D} \,\right] G(t) = 0, \qquad G(0) = 0, \quad G' (0) = 0 , \ldots , \quad \left. \frac{{\text d}^{n-2}G}{{\text d} t^{n-2}} \right\vert_{t=0} = 0, \quad \left. \frac{{\text d}^{n-1}G}{{\text d} t^{n-1}} \right\vert_{t=0} = \frac{1}{a_n} . \]
Theorem 3: Suppose that f(t) is continuous function on [0, ∞) and has a Laplace transform. Then the solution of the initial value problem for the constant coefficient differential equation (𝑎 ≠ 0)
\[ a\,\ddot{y} + b\,\dot{y} + c\,y(t) = f(t) , \qquad y(0) = y_0 , \quad \dot{y}(0) = v, \]
is
\[ y(t) = y(0)\,y_1 (t) + \dot{y}(0)\,g(t) + \frac{1}{a} \int_0^{\infty} g(\tau )\,f(t-\tau )\,{\text d}\tau , \]
where y1(t) and g(t) are solutions of the following IVPs:
\[ a\,\ddot{y}_1 + b\,\dot{y}_1 + c\,y_1(t) = 0, \qquad y_1 (0) = 1, \quad \dot{y}_1 (0) = 0, \]
and
\[ a\,\ddot{g} + b\,\dot{g} + c\,g(t) = 0, \qquad g (0) = 0, \quad \dot{g} (0) =1. \]

Example 5: Consider the linear differential operator

\[ L \left[ \texttt{D} \right] = \texttt{D}^2 + 2\,\texttt{D} -3\,\texttt{I}, \]
where \( \texttt{D} = {\text d}/{\text d}t \) is the derivative operator and \( \texttt{I} \) is the identity operator. Recall that the Green function for the differential operator \( L \left[ \texttt{D} \right] \) is the inverse Laplace transform of the reciprocal to the characteristic polynomial: \( \displaystyle {\cal L}^{-1} \left[ 1/L(\lambda ) \right] . \)

Since the characteristic polynomial L(λ) = (λ-1)(λ+3) has three simple real roots, we find the corresponding Green function using the inverse Laplace transform:

\[ G(t) = {\cal L}^{-1} \left[ \frac{1}{(\lambda -1)(\lambda + 3)} \right] = \mbox{Res}_{\lambda = 1} \, \frac{e^{\lambda\,t}}{(\lambda -1)(\lambda + 3)} + \mbox{Res}_{\lambda = -3} \, \frac{e^{\lambda\,t}}{(\lambda -1)(\lambda + 3)} . \]
For simple roots, we find immediately
\[ G(t) = \frac{1}{4} \left[ e^t - e^{-3t} \right] H(t) . \]
   ■

Example 6: Find the Green function for the following differential operator:

\[ L \left[ \texttt{D} \right] = 4\,\texttt{D}^4 - 4\,\texttt{D}^3 + 33\,\texttt{D}^2 + 38\,\texttt{D} +10\,\texttt{I}, \]
where \( \texttt{D} = {\text d}/{\text d}t \) is the derivative operator and \( \texttt{I} \) is the identity operator. Recall that the Green function for the differential operator \( L \left[ \texttt{D} \right] \) is the inverse Laplace transform of the reciprocal to the characteristic polynomial: \( \displaystyle {\cal L}^{-1} \left[ 1/L(\lambda ) \right] . \)

The Green function for the given differential operator is the solution of the differential equation

\[ L \left[ \texttt{D} \right] G(t) = \delta (t) , \]
where δ(t) is the Dirac delta function. The Green function can be determined by evaluating the inverse Laplace transform:
\[ G(t) = {\cal L}^{-1} \left[ \frac{1}{L(\lambda )} \right] = {\cal L}^{-1} \left[ \frac{1}{4\,\lambda^4 - 4\,\lambda^3 + 33\,\lambda^2 + 38\,\lambda + 10} \right] . \]
Expanding the characteristic polynomial L(λ) with Mathematica, we get
\[ L(\lambda ) = \left( 2\,\lambda +1 \right)^2 \left[ \left( \lambda -1 \right)^2 + 3^2 \right] = 4 \left( \lambda +1/2 \right)^2 \left[ \left( \lambda -1 \right)^2 + 3^2 \right] . \]
Factor[10 + 38 x + 33 x^2 - 4 x^3 + 4 x^4] (1 + 2 x)^2 (10 - 2 x + x^2)
So the characteristic polynomial L(λ) has one double real root λ = -½ and two simple complex conjugate roots λ = 1 ±3j. Recall that j is the unit vector in positive vertical direction on the complex plane ℂ, so j² = -1. To find the explicit expression for the Green function, we need to evaluate two residues:
\begin{align*} \mbox{Res}_{\lambda = -1/2} \frac{e^{\lambda t}}{L(\lambda )} &= \lim_{\lambda \to -1/2} \,\frac{\text d}{{\text d}\lambda} \frac{e^{\lambda t}}{4} \cdot \frac{1}{\left( \lambda -1 \right)^2 + 3^2} \\ &= \frac{t}{4}\, e^{-t/2} \,\lim_{\lambda \to -1/2} \, \frac{1}{\left( \lambda -1 \right)^2 + 3^2} - \frac{1}{4}\, e^{-t/2} \,\lim_{\lambda \to -1/2} \, \frac{2 \left( \lambda -1 \right)}{\left[ \left( \lambda -1 \right)^2 + 3^2 \right]^2} ; \\ \mbox{Res}_{\lambda = 1+ 3{\bf j}} \frac{e^{\lambda t}}{L(\lambda )} &= \lim_{\lambda \to 1+3{\bf j}} \, \frac{e^{\lambda t}}{\left( 2\,\lambda +1 \right)^2 2\,(\lambda -1)} \\ &= \frac{e^{(1+3{\bf j})\, t}}{\left( 3 +6{\bf j} \right)^2 2\times 3{\bf j}} = e^t \, \frac{1}{2\times 3^3} \times \frac{e^{3{\bf j}t}}{{\bf j} \left( 1+ 2{\bf j} \right)^2} = - \frac{e^{t}}{2\times 3^3} \times \frac{e^{3{\bf j}t}}{4+ 3{\bf j}} . \end{align*}
To proceed, we need to perform some arithmetic operations. First, we find the limits
\begin{align*} \lim_{\lambda \to -1/2} \left( \lambda -1 \right) &= -\frac{3}{2} , \\ \lim_{\lambda \to -1/2} \left[ \left( \lambda -1 \right)^2 + 3^2 \right] &= 3^2 \left( \frac{5}{4} \right) , \\ \lim_{\lambda \to -1/2} \frac{2 \left( \lambda -1 \right)}{\left[ \left( \lambda -1 \right)^2 + 3^2 \right]^2} &= - \frac{4^2}{3^3 5^2} = - \frac{16}{675} . \end{align*}
We also need to perform operations with complex numbers:
\[ \frac{e^{3{\bf j}t}}{4+ 3{\bf j}} = \left( \cos 3t + {\bf j}\,\sin 3t \right) \frac{4 - 3{\bf j}}{4^2 + 3^2} . \]
Extracting the real (denoted by ℜ = Re) and imaginary (denoted by ℑ = Im) parts, we get
\begin{align*} \Re\,\frac{e^{3{\bf j}t}}{4+ 3{\bf j}} = \frac{4}{25}\,\cos 3t + \frac{3}{25}\,\sin 3t , \\ \Im \,\frac{e^{3{\bf j}t}}{4+ 3{\bf j}} = - \frac{3}{25}\,\cos 3t + \frac{4}{25}\,\sin 3t . \end{align*}
The imaginary part of the above residue will be canceled out by the residue evaluated at the complex conjugate root 1-3j, so we obtain
\[ \mbox{Res}_{\lambda = 1+ 3{\bf j}} \frac{e^{\lambda t}}{L(\lambda )} + \mbox{Res}_{\lambda = 1-3{\bf j}} \frac{e^{\lambda t}}{L(\lambda )} = - \frac{e^t}{3^3} \left( \frac{4}{25}\,\cos 3t + \frac{3}{25}\,\sin 3t \right) . \]
The residue at the double root becomes
\[ \mbox{Res}_{\lambda = -1/2} \frac{e^{\lambda t}}{L(\lambda )} = \frac{t}{45}\, e^{-t/2} + \frac{4}{675}\,e^{-t/2} . \]
Finally, we get the Green function:
\[ G(t) = \left[ e^{-t/2} \left( \frac{t}{45} + \frac{4}{675} \right) - \frac{e^t}{3^3} \left( \frac{4}{25}\,\cos 3t + \frac{3}{25}\,\sin 3t \right) \right] H(t) , \]
where H(t) is the Heaviside function.

To check the answer, we apply the partial fraction decomposition:

\[ \frac{1}{L(\lambda )} = \frac{a}{(2\lambda +1)^2} + \frac{b}{(2\lambda +1)} + \frac{c+ d\lambda}{\left( \lambda -1 \right)^2 + 9} , \]
where a,b,c,d are constants to be determined. Making a common denominator, we get the equation
\[ 1 = (2\lambda +1)^2 \left( c+ d \lambda \right) + \left( a + b\lambda \right) \left[ \left( \lambda -1 \right)^2 + 9 \right] , \]
Expanding the right-hand side into polynomial form and equating coefficients, we obtain the system of four equations with four unknowns:
\begin{align*} \lambda^0 : \ & 1 = c + 10\,a + 10\,b, \\ \lambda^1 : \ & 0 = 4c+d -2a+18\,b , \\ \lambda^2 : \ & 0 = 4c + 4d -3b +a , \\ \lambda^3 : \ & 0 = 4\,d+ 2\,b . \end{align*}
Solving with Mathematica, we obtain
Solve[{10*a + c + 10*b == 1, 4*c + d - 2*a + 18 b == 0, 4*c + 4*d - 3*b + a == 0, 4*d + 2*b == 0}, {a, b, c, d}]
a -> 4/45, b -> 8/675, c -> -(1/135), d -> -(4/675)
However, Maple can do this job in a blink of an eye:
convert(1/((2*x + 1)^2*(x*x - 2*x + 10)), parfrac, x)
and get
\[ \frac{1}{L(\lambda )} = \frac{4}{45 \left( 2\lambda +1 \right)^2} + \frac{8}{675 \left( 2\lambda +1 \right)} - \frac{4\,\lambda + 5}{675 \left[ \lambda^2 - 2\lambda + 10\right]} . \]
   ■

Example 7: Let us consider the initial value problem (IVP) for the third order differential equation

\[ \frac{{\text d}^3 y}{{\text d} t^3} + 5\,\frac{{\text d}^2 y}{{\text d} t^2} + 2 \,\frac{{\text d} y}{{\text d} t} - 8\, y = e^{-t} \cos (2t) , \qquad y(0) = 1, \quad \dot{y}(0) = -1, \quad \ddot{y}(0) = 2. \tag{7.1} \]
Application of the Laplace transformation yields
\[ \left( \lambda^3 + 5\lambda^2 + 2 \lambda -8 \right) y^L = f^L (\lambda ) + \left( \lambda^2 + 4\lambda -1 \right) , \]
where
\[ y^L = \int_0^{\infty} e^{-\lambda t} y(t)\,{\text d}t \qquad \mbox{and}\qquad f^L = \int_0^{\infty} e^{-\lambda t} e^{-t} \cos (2t) \,{\text d}t = \frac{1+\lambda}{4 + (1+\lambda )^2} . \]
LaplaceTransform[Exp[-t]*Cos[2*t], t, s]
(1 + s)/(4 + (1 + s)^2)
So the Laplace transform yL is expressed via the formula
\[ y^L = \frac{\lambda^2 + 4\lambda -1}{\lambda^3 + 5\lambda^2 + 2 \lambda -8} + \frac{1}{\left( \lambda +4\right) \left( \lambda +2 \right) \left( \lambda -1 \right)} \cdot f^L (\lambda ) = y_h^L + y_p^L , \tag{7.2} \]
where
\[ y_h^L = \frac{\lambda^2 + 4\lambda -1}{\left( \lambda +4\right) \left( \lambda +2 \right) \left( \lambda -1 \right)} , \]
\[ y_p^L = \frac{1}{\left( \lambda +4\right) \left( \lambda +2 \right) \left( \lambda -1 \right)} \cdot f^L (\lambda ) = \frac{1}{\left( \lambda +4\right) \left( \lambda +2 \right) \left( \lambda -1 \right)} \cdot \frac{1+\lambda}{4 + (1+\lambda )^2} . \]
The solution yh(t) of the homogeneous equation subject to the given initial conditions,
\[ \dddot{y} + 5\ddot{y} + 2\dot{y} -8\,y =0, \qquad y(0) = 1, \quad \dot{y}(0) = -1, \quad \ddot{y}(0) = 2. \]
is
\[ y_h (t) = {\cal L}^{-1} \left[ \frac{\lambda^2 + 4\lambda -1}{\left( \lambda +4\right) \left( \lambda +2 \right) \left( \lambda -1 \right)} \right] = \frac{1}{30} \left[ 8\, e^t + 25\, e^{-2t} -3\,e^{-4t} \right] H(t) \tag{7.3} \]
The expression (7.3) is obtained by summing thre residues:
\begin{align*} \mbox{Res}_{\lambda = -4} \frac{\lambda^2 + 4\lambda -1}{\lambda^3 + 5 \lambda^2 +2 \lambda -8} \,e^{\lambda t} &= \lim_{\lambda \to -4} \frac{\lambda^2 + 4\lambda -1}{3\lambda^2 +10\lambda +2} \,e^{\lambda t} = -\frac{1}{10}\, e^{-4t} , \\ \mbox{Res}_{\lambda = -2} \frac{\lambda^2 + 4\lambda -1}{\lambda^3 + 5 \lambda^2 +2 \lambda -8} \,e^{\lambda t} &= \lim_{\lambda \to -2} \frac{\lambda^2 + 4\lambda -1}{3\lambda^2 +10\lambda +2} \,e^{\lambda t} = \frac{5}{6}\, e^{-2t} , \\ \mbox{Res}_{\lambda = 1} \frac{\lambda^2 + 4\lambda -1}{\lambda^3 + 5 \lambda^2 +2 \lambda -8} \,e^{\lambda t} &= \lim_{\lambda \to 1} \frac{\lambda^2 + 4\lambda -1}{3\lambda^2 +10\lambda +2} \,e^{\lambda t} = \frac{4}{15}\, e^{t} . \end{align*}
We check with Mathematica
InverseLaplaceTransform[(s^2 + 4*s - 1)/(s^3 + 5*s^2 + 2*s - 8), s, t];
yh[t_] = 1/30 E^(-4 t) (-3 + 25 E^(2 t) + 8 E^(5 t));
yh[0]
1
D[yh[t], t] /. t -> 0
-1
D[yh[t], t, t] /. t -> 0
2
To determine yp(t), the solution of the nonhomogeneous equation
\[ \dddot{y} + 5\ddot{y} + 2\dot{y} -8\,y = e^{-t} \cos (2t), \qquad y(0) = 0, \quad \dot{y}(0) = 0, \quad \ddot{y}(0) = 0, \]
we have two options. First, we find the Green function:
\[ G(t) = {\cal L}^{-1} \left[ \frac{1}{\left( \lambda +4\right) \left( \lambda +2 \right) \left( \lambda -1 \right)} \right] = \left[ \frac{1}{15}\, e^t - \frac{1}{6}\, e^{-2t} + \frac{1}{10}\, e^{-4t} \right] H(t) . \]
InverseLaplaceTransform[(1)/(s^3 + 5*s^2 + 2*s - 8), s, t]
E^(-4 t)/10 - E^(-2 t)/6 + E^t/15
With this in hand, we obtain the particular solution with the aid of the convolution rule:
\[ y_p (t) = G(t) \ast f(t) = \int_0^t G(t-\tau )\, f(\tau )\,{\text d}\tau = \int_0^t \left[ \frac{1}{15}\, e^{(t-\tau )} - \frac{1}{6}\, e^{-2(t-\tau )} + \frac{1}{10}\, e^{-4(t- \tau )} \right] e^{-\tau} \cos (2\tau )\,{\text d}\tau . \]
This tedious job we delegate to Mathematica:
g[t_] = E^(-4 t)/10 - E^(-2 t)/6 + E^t/15;
f[t_] = Exp[-t]*Cos[2*t];
Integrate[g[t - s]*f[s], {s, 0, t}]
1/780 E^(-4 t) (-18 + 26 E^(2 t) + 13 E^(5 t) - 21 E^(3 t) Cos[2 t] - 27 E^(3 t) Sin[2 t])
\[ y_p (t) = \frac{1}{780}\, e^{-4t} \left[ 26\, e^{2t} -18 + 13\, e^{5t} - 21\,e^{3t}\cos (2t) - 27\, e^{3t} \sin (2t) \right] H(t) . \tag{7.4} \]
Another option to verify this formula is the find the inverse Laplace transform directly with the aid of the residue theorem:
\[ y_p (t) = {\cal L}^{-1} \left[ \frac{1}{\lambda^3 + 5 \lambda^2 + 2\lambda -8} \cdot \frac{1+\lambda}{4 + (1+\lambda )^2} \right] = \mbox{Res}_{\lambda = -4} + \mbox{Res}_{\lambda = -2} + \mbox{Res}_{\lambda = 1} + 2\,\Re \mbox{Res}_{\lambda = -1 + 2{\bf j}} . \]
InverseLaplaceTransform[(1 + s)/(s^3 + 5*s^2 + 2*s - 8)/(4 + (s + 1)^2), s, t]
-(3/130) E^(-4 t) + E^(-2 t)/30 + E^t/60 + (1/520 + I/ 520) E^((-1 - 2 I) t) ((-8 - I) + (1 + 8 I) E^(4 I t))
Since Mathematica provides the expression containing complex numbers, we extract the real part and obtain
Simplify[Assuming[t > 0, ComplexExpand[%]]]
1/780 E^(-4 t) (-18 + 26 E^(2 t) + 13 E^(5 t) - 21 E^(3 t) Cos[2 t] - 27 E^(3 t) Sin[2 t])
This gives us exactly the same answer as in Eq.(7.4). Finally, we verify once more by solving the initial value problem:
DSolve[{u'''[t] + 5*u''[t] + 2*u'[t] - 8*u[t] == Exp[-t]*Cos[2*t], u[0] == 0, u'[0] == 0, u''[0] == 0}, u[t], t]
{{u[t] -> 1/780 E^(-4 t) (-18 + 26 E^(2 t) + 13 E^(5 t) - 21 E^(3 t) Cos[2 t] - 27 E^(3 t) Sin[2 t])}}
      We plot the particular solution yp(t), given in formula (7.4), (in blue) along with the input function f(t) in black. As you can see, there is no relation of the particular solution and the input function.
yp[t_] = 1/ 780 E^(-4 t) (-18 + 26 E^(2 t) + 13 E^(5 t) - 21 E^(3 t) Cos[2 t] - 27 E^(3 t) Sin[2 t]);
f[t_] = Exp[-t]*Cos[2*t];
Plot[{yp[t], f[t]}, {t, 0, 4.5}, PlotStyle -> {{Thickness[0.015], Blue}, {Thickness[0.01], Black}}]
       Plot of the particular solution (7.4) in blue, and f(t) in black.            Mathematica code

   ■

Example 8: Let us consider the initial value problem

\[ y'' + 2\,y' -3\,y = 2x, \qquad y(1) = 3, \quad y' (1) = -1. \tag{8.1} \]
Since the initial conditions are specified at point x = 1, we first make a shift in the independent variable by introducing t = x −1. This yields the IVP:
\[ \ddot{y} + 2\,\dot{y} -3\,y = 2t +2, \qquad y(0) = 3, \quad y' (0) = -1. \tag{8.2} \]
Now we can apply the Laplace transform and obtain
\[ \lambda^2 y^L + 2\lambda\,y^L -3\,y^L = \left( \frac{2}{\lambda^2} + \frac{2}{\lambda} \right) + \dot{y}(0) + \lambda\,y(0) + 2\,y(0) . \]
This linear algebraic equation in the Laplace transform, yL, of the unknown function can be solved. So we get
\[ y^L (\lambda ) = \frac{2 \left( 1 + \lambda \right)}{\lambda^2 \left( \lambda^2 + 2\lambda -3 \right)} + \frac{3\lambda + 5}{\lambda^2 + 2\lambda -3} = y_p^L + y_h^L . \tag{8.3} \]
Application of the inverse Laplace transform yields
InverseLaplaceTransform[ 2*(1 + s)/s^2/(s^2 + 2*s - 3) + (3*s + 5)/(s^2 + 2*s - 3), s, t]
-(10/9) + (10 E^(-3 t))/9 + 3 E^t - (2 t)/3
\[ y(t) = {\cal L}^{-1} \left[ y^L \right] = \left[ \frac{10}{9}\, e^{-3t} - \frac{10}{9} + 3\, e^t - \frac{2t}{3} \right] H(t) . \tag{8.4} \]
We can also use the Green function
\[ G(t) = {\cal L}^{-1} \left[ \frac{1}{\lambda^2 + 2\lambda -3} \right] = \frac{1}{4} \left( e^t - e^{-3t} \right) H(t) \tag{8.5} \]
InverseLaplaceTransform[1/(s^2 + 2*s - 3), s, t]
1/4 E^(-3 t) (-1 + E^(4 t))
Green's function is a solution of the initial value problem
\[ \ddot{G} + 2\dot{G} - 3\,G (t) = 0, \qquad G(0) = 0, \quad \dot{G}(0) = 1. \]
Then the part of the solution that corresponds to the nonhomogeneous equation is
\[ y_p (t) = {\cal L}^{-1} \left[ y_p^L \right] = {\cal L}^{-1} \left[ \frac{2 \left( 1 + \lambda \right)}{\lambda^2 \left( \lambda^2 + 2\lambda -3 \right)} \right] = G \ast \left( 2t + 2 \right) = \int_0^t G(t-\tau )\,2\left( 1 + \tau \right) {\text d}\tau = \frac{1}{9} \left( 3t-1 + e^{-3t} \right) H(t) . \]
Integrate[(1 + x)*(Exp[x - t] - Exp[-3*(t - x)]), {x, 0, t}] /2
1/9 (-1 + E^(-3 t) + 3 t)
   ■