Laplace Transforms, Signals and Plants ====================================== With Laplace transforms, differential equations become algebraic equations. The Frequency Domain -------------------- The :index:`Laplace transform` takes a signal from the :index:`time domain`, in :math:`t`, to the *frequency domain*, using :math:`s` as the symbol in the transform. .. index:: frequency domain .. topic:: Definition of the Laplace transform. Given an integrable function :math:`f(t)` in time :math:`t`, the *Laplace transform* of :math:`f(t)` is .. math:: {\cal L}\{ f \} = \int_0^\infty f(t) e^{-s t} dt = F(s). Because the integral operator is a linear operator, the Laplace transform is linear. In particular, * Let :math:`f(t)` and :math:`g(t)` be two integrable functions. * Let :math:`a` and :math:`b` be two coefficients, independent of :math:`t`. Then .. math:: {\cal L}\{a f(t)+ b g(t)\} = a F(s) + b G(s), where * :math:`F(s)` is the Laplace transform of :math:`f(t)`, and * :math:`G(s)` is the Laplace transform of :math:`g(t)`. Another set of properties involves delay and exponentials. For some constant coefficient :math:`a`, we have: * :math:`f(t - a)` has Laplace transform :math:`e^{-as} F(s)`, * :math:`e^{-at} f(t)` has Laplace transform :math:`F(s+a)`, * :math:`f(a t)` has Laplace transform :math:`F(s/a)/a`. *What is s?* Consider the Laplace transform of :math:`f(t) = 1`: .. math:: \begin{array}{rcl} {\displaystyle \int_0^\infty 1 e^{-s t} dt} & = & {\displaystyle \int_0^\infty \frac{-1}{s} e^{(-s t)} d(-st) } \\ & = & {\displaystyle \left. - \frac{e^{-s t}}{s} \right]_0^\infty} ~ = ~ 0 - \left( -\frac{1}{s} \right) = \frac{1}{s}. \end{array} This leads to .. math:: s \int_0^\infty 1 e^{-s t} dt = 1 \quad \mbox{or} \quad s {\cal L}\{ 1 \} = 1. Multiplying with :math:`s` can be called a *generalized differentiation*. .. topic:: the Laplace transform of the derivative .. math:: {\cal L}\left\{ \frac{d}{d~\!t} f(t) \right\} = s F(s) - f(0), \quad F(s) = {\cal L}\{ f \} Apply integration by parts: :math:`\displaystyle \int u dv = uv - \int v du`. .. math:: \begin{array}{rcl} {\displaystyle {\cal L}\left\{ \frac{d}{d~\!t} f(t) \right\}} & = & {\displaystyle \int_0^\infty e^{-st} \frac{d}{d~\!t} f(t) ~=~ \left. e^{-st} f(t) \vphantom{\frac{1}{2}} \right]_0^\infty - \int_0^\infty f(t) \frac{d}{d~\!t} e^{-st}} \\ & = & {\displaystyle 0 - f(0) - \int_0^\infty f(t) (-s) e^{-st} dt} \\ & = & -f(0) + s F(s) \end{array} Applying the property of the Laplace transform repeatedly, we arrive at the following: .. topic:: the Laplace transform of the :math:`n`-th derivative For :math:`n>0`: .. math:: {\cal L}\left\{ \frac{d^n}{d~\!t^n} f(t) \right\} = s^n F(s) - s^{n-1} f(0) - s^{n-2} f'(0) - \cdots - f^{(n-1)}(0), where :math:`F(s) = {\cal L}\{ f \}`. We next explore symbolic computing, and the Laplace transform in ``SymPy.jl``. Consider: :: @syms s, t one = sympy.laplace_transform(1, t, s) which returns ``(1/s, 0, True)``. :: sympy.inverse_laplace_transform(one[1], s, t) returns :math:`\theta(t)`, defined below: .. topic:: Definition of the Heaviside function The Heaviside function is the piecewise function .. math:: \theta(t) = \left\{ \begin{array}{lcl} 1 & \mbox{if} & t > 0, \\ 1/2 & \mbox{if} & t = 0, \\ 0 & \mbox{if} & t < 0. \end{array} \right. The :index:`Heaviside function` can also be seen as a limit of a smoothed step function, as .. math:: \lim_{k \rightarrow \infty} \left( 1 + e^{-2 k x} \right)^{-1} and shown in :numref:`figHeavisideLimit` for increasing values of :math:`k`. .. _figHeavisideLimit: .. figure:: ./figHeavisideLimit.png :align: center The Heaviside function as the limit of a smoothed step function. The ``sympy.inverse_laplace_transform`` used earlier is formally defined below. .. index:: inverse Laplace transform .. topic:: Definition of the inverse Laplace transform The *inverse Laplace transform* of :math:`F(s)` is :math:`f` with the property .. math:: {\cal L}\{ f \}(s) = F(s). Observe the similarity between the Laplace transform of :math:`f(t)` which is :math:`\displaystyle {\cal L}\{ f \} = \int_0^\infty f(t) e^{-s t} dt = F(s)` with the :index:`continuous Fourier transforms`: .. index:: Fourier transform, inverse Fourier transform .. topic:: Definition of the Fourier transform The *Fourier transform* of :math:`f` is .. math:: {\cal F} \{ f \} = \int_{-\infty}^{\infty} f(t) e^{-i \omega t} dt = F(i \omega). .. topic:: Definition of the inverse Fourier transform The *inverse Fourier transform* of :math:`F` is .. math:: {\cal F}^{-1} \{ F \} = \frac{1}{2 \pi} \int_{-\infty}^{\infty} F(i \omega) e^{i \omega t} d \omega. Operational Calculus -------------------- To solve an Ordinary Differential Equation (ODE), do: 1. Transform the ODE into the frequency domain. 2. Solve the algebraic equation in :math:`s`. 3. Transform the solution into the time domain. Consider the following example: .. math:: \frac{d^2}{d~\!t^2} y(t) + 3 \frac{d}{d~\!t} y(t) + 2 y(t) = 1, \quad y(0) = 0, \quad y'(0) = 0. 1. Denoting :math:`{\cal L}\{y\} = Y(s)`, the transformed equation is .. math:: s^2 Y(s) + 3 s Y(s) + 2 Y(s) = \frac{1}{s}. 2. We solve the algebraic equation: :math:`s ( s+1 ) (s+2) Y(s) = 1`. .. math:: Y(s) = \frac{1}{s (s+1) (s+2)} = \frac{1/2}{s} - \frac{1}{s+1} + \frac{1/2}{s+2} 3. Transformed back into the time domain, the solution is .. math:: y(t) = \frac{1}{2} - e^{-t} + \frac{1}{2} e^{-2t}. Consider the system of first order differential equations: .. math:: \left\{ \begin{array}{rcl} x' & = & -2x - 9y \\ y' & = & \hphantom{-2}x - 2y \end{array} \right. \quad x(0) = -1, \quad y(0) = 1. We apply the 3-step process: 1. Denoting :math:`X(s) = {\cal L}\{x\}` and :math:`Y(s) = {\cal L}\{y\}`, the transformed system is .. math:: \left\{ \begin{array}{rcl} sX + 1 & = & -2X - 9Y \\ sY - 1 & = & \hphantom{-2}X - 2Y \end{array} \right. \quad \mbox{or} \quad \left\{ \begin{array}{rcl} (s + 2)X + 9Y & = & - 1 \\ -X + (s+2)Y & = & 1 \end{array} \right. 2. Use ``SymPy`` to solve the system. :: @syms X, Y eq1 = (s+2)*X + 9*Y + 1 eq2 = -X + (s+2)*Y - 1 sol = solve([eq1, eq2], [X, Y]) has output :: Dict{Any, Any} with 2 entries: Y => (s + 1)/((s + 2)**2 + 9) X => -(s + 11)/(s**2 + 4*s + 13) 3. Transform the solution into the time domain. :: xt = sympy.inverse_laplace_transform(sol[X], s, t) returns :math:`-e^{-2t}(3 \sin(3t) + \cos(3t)) \theta(t)` as :math:`x(t)`. :: yt = sympy.inverse_laplace_transform(sol[Y], s, t) returns :math:`-e^{-2t}(\sin(3t)/3 + \cos(3t)) \theta(t)` as :math:`y(t)`. Note: for :math:`t > 0`, we have :math:`\theta(t) = 1` and we can drop :math:`\theta(t)`. Generalizing the previous example to solving linear systems goes as follows. Let :math:`A` be an *n*-by-*n* matrix and consider the system .. math:: \frac{d}{d~\!t} {\bf z}(t) = A {\bf z}, \quad {\bf z}(0) = {\bf z}_0. 1. Denoting :math:`{\cal L}\{{\bf z}\} = {\bf Z}(s)`, the transformed equation is .. math:: s {\bf Z}(s) - {\bf z}_0 = A. 2. The solution in the frequency domain is .. math:: {\bf Z}(s) = (s I - A)^{-1} {\bf z}_0, where :math:`I` is the identity matrix. If :math:`A` has :math:`k` distinct :index:`eigenvalues`, then .. math:: \det(s I - A) = \prod_{i=1}^k (s - \lambda_k)^{m_k}, where :math:`m_k` is the multiplicity of the *k*-th eigenvalue :math:`\lambda_k`. The matrix :index:`partial fraction expansion` has then the form: .. math:: (s I - A)^{-1} = \sum_{k=1}^r \sum_{j=1}^{m_k} \frac{v_{j,k}}{(s-\lambda_k)^j}. 3. Use the partial fraction expansion to transfer into the time domain. Signals and Plants ------------------ We next define the convolution of continuous functions, make the connection between signals and ordinary differential equations. .. index:: continuous convolution .. topic:: Definition of convolutions of functions Let :math:`f` and :math:`g` be functions of :math:`t`. The *convolution* :math:`f \star g` of :math:`f` and :math:`g` is .. math:: f \star g = \int_0^t f(\tau) g(t - \tau) d \tau. Recall the :index:`discrete convolution` :math:`\displaystyle y_k = \sum_{i=0}^k h_i u_{k-i}`, where * :math:`y_k` is the *k*-th component of the output signal, :math:`y = h \star u`, * :math:`h` is the response of the unit impulse, and * :math:`u` is the input signal. Applying the convolution to continous signals, for continuous input :math:`u(t)`, unit impulse response :math:`h(t)`, the output is then .. math:: y(t) = h(t) \star u(t) is the continuous output of the filter. In the frequency domain: .. math:: Y(s) = H(s) U(s) where the Laplace transforms are .. math:: U(s) = {\cal L}\{ u \}, \quad H(s) = {\cal L}\{ h \}, \quad \mbox{and} \quad Y(s) = {\cal L}\{ y \}. .. topic:: Definition of a plant A *plant* is a linear, time invariant, causal map, which maps signals into signals. Examples of plants are: filters, mechanisms, circuits. By the linearity, time invariance, and causality, the :index:`plant` is determined by its response to the unit impulse. Plants can be defined by ordinary differential equations. For example, consider the plant defined by .. math:: \frac{d}{d~t} y(t) + a y(t) = u(t), where :math:`u(t)` is the input and :math:`y(t)` the output. .. math:: \begin{array}{rcl} {\displaystyle \frac{d}{d~\!t} \left( \vphantom{\frac{1}{2}} e^{a t} y(t) \right)} & = & {\displaystyle a e^{a t} y(t) + e^{a t} \frac{d}{d~\!t} y(t)} \\ & = & {\displaystyle e^{a t} \left( a y(t) + \frac{d}{d~\!t} y(t) \right) ~ = ~ e^{a t} u(t)} \end{array} .. math:: \frac{d}{d~\!t} \left( \vphantom{\frac{1}{2}} e^{a t} y(t) \right) = e^{a t} u(t) \quad \Rightarrow \quad y(t) = e^{-a t} \left( \int_0^t e^{a \tau} u(\tau) d \tau + y(0) \right) Convolutions are connected to the Laplace transform: .. math:: \begin{array}{rcl} y(t) & = & {\displaystyle e^{-a t} \left( \int_0^t e^{a \tau} u(\tau) d \tau + y(0) \right)} \\ & = & {\displaystyle \int_0^t e^{-a(t - \tau)} u( \tau) d \tau + y(0) e^{-a t}} \end{array} Set :math:`y(0) = 0`. The output of the plant is then .. math:: y(t) = e^{-a t} \star u(t). In the frequency domain: .. math:: Y(s) = \left( \frac{1}{s+a} \right) U(s). Exercises --------- 1. Verify the properties, for some constant :math:`a`, and function :math:`f` with Laplace transform :math:`F(s)`: * :math:`f(t - a)` has Laplace transform :math:`e^{-as} F(s)`, * :math:`e^{-at} f(t)` has Laplace transform :math:`F(s+a)`, * :math:`f(a t)` has Laplace transform :math:`F(s/a)/a`. 2. A model of a suspension mechanism, using the equation for damped harmonic motion, written in engineering notation: .. math:: \frac{d^2}{d~\!t^2} y(t) + 2 \zeta \omega \left( \frac{d}{d~\!t} y(t) \right) + \omega^2 y(t) = 0, with initial conditions: .. math:: y(0) = 1 \quad \mbox{and} \quad y'(0) = 0. Apply the 3-step operational calculus method to compute the solution of this model. 3. Consider again .. math:: \left\{ \begin{array}{rcl} x' & = & -2x - 9y \\ y' & = & \hphantom{-2}x - 2y \end{array} \right. \quad x(0) = -1, \quad y(0) = 1. 1. Write the system in matrix format. 2. Compute the eigenvalues of the matrix. 3. Compare the eigenvalues with the previously computed solution, previously computed with the aid of ``SymPy``. Bibliography ------------ 1. Charles R. MacCluer, start of Chapter 10 of *Industrial Mathematics. Modeling in Industry, Science, and Government*, Prentice Hall, 2000. 2. Richard C. Dorf and Robert H. Bishop: *Modern Control Systems.* Nineth Edition, Prentice Hall, 2001.