Engineering Mathematics 1: Ordinary Differential Equations (ODE)

Engineering mathematics is an unavoidable subject for any engineering student. In particular, Ordinary Differential Equations (ODEs) form the mathematical foundation of virtually every engineering discipline — from circuit analysis and control system design to mechanical vibration analysis. This article systematically explains ODEs from the basics to Laplace transforms, centered on engineering application examples.

1. Fundamentals of Differential Equations

What is a Differential Equation?

A differential equation is an equation that expresses a relationship between an unknown function and its derivatives. For example:

$\frac{dy}{dx} = 2x, \quad y'' + 3y' + 2y = 0, \quad \frac{\partial u}{\partial t} = k\frac{\partial^2 u}{\partial x^2}$

Differential equations are a core tool for mathematically modeling natural phenomena and engineering systems.

Classification: ODE vs PDE

Ordinary Differential Equation (ODE): Contains derivatives with respect to only one independent variable.

$\frac{d^2x}{dt^2} + 2\frac{dx}{dt} + x = 0$

Partial Differential Equation (PDE): Contains partial derivatives with respect to two or more independent variables.

$\frac{\partial u}{\partial t} = c^2\frac{\partial^2 u}{\partial x^2}$

This article focuses on ODEs.

Order and Degree

The order is the highest derivative appearing in the equation.

• First-order ODE: $y' + y = x$ (highest derivative is $y'$)
• Second-order ODE: $y'' + 2y' + y = 0$ (highest derivative is $y''$)

The degree is the exponent of the highest-order derivative. In most engineering problems, the degree is 1.

Linear vs Nonlinear

Linear ODE: The unknown function and its derivatives appear only as first-degree terms.

$a_n(x)y^{(n)} + a_{n-1}(x)y^{(n-1)} + \cdots + a_1(x)y' + a_0(x)y = f(x)$

Nonlinear ODE: Contains nonlinear terms such as $y^2$, $y \cdot y'$, $\sin(y)$, etc.

$y'' + \sin(y) = 0 \quad \text{(pendulum equation)}$

General Solution, Particular Solution, Singular Solution

• General solution: A complete solution containing arbitrary constants
• Particular solution: A solution obtained by applying initial conditions to determine the constants
• Singular solution: A special solution that cannot be obtained from the general solution

ODEs in Engineering

Most engineering systems are described by ODEs.

RC charging circuit: $RC\frac{dv}{dt} + v = V_s$

Spring-mass-damper system: $m\ddot{x} + c\dot{x} + kx = F(t)$

Population growth model: $\frac{dP}{dt} = kP$

RL circuit: $L\frac{di}{dt} + Ri = E(t)$

All of these are ODEs. Learning how to solve them allows you to fully analyze the dynamic behavior of engineering systems.

2. First-Order Ordinary Differential Equations

2.1 Separation of Variables

Form: $\frac{dy}{dx} = f(x)g(y)$

Solution procedure:

$\frac{dy}{g(y)} = f(x)dx$

Integrate both sides:

$\int \frac{dy}{g(y)} = \int f(x)dx + C$

Example: RC Charging Circuit

Circuit equation: $\frac{dv}{dt} = \frac{V_s - v}{RC}$

Separation of variables: $\frac{dv}{V_s - v} = \frac{dt}{RC}$

Integration: $-\ln|V_s - v| = \frac{t}{RC} + C_1$

Applying initial condition $v(0) = 0$:

$V_s - v = V_s e^{-t/(RC)}$

$\boxed{v(t) = V_s\left(1 - e^{-t/(RC)}\right)}$

This result represents the transient response of an RC circuit. The voltage converges to $V_s$ at a rate determined by the time constant $\tau = RC$.

Physical interpretation:

• At $t = \tau$: $v \approx 0.632 V_s$ (63.2% charged)
• At $t = 5\tau$: $v \approx 0.993 V_s$ (almost fully charged)

2.2 Integrating Factor Method

Standard form of a first-order linear ODE: $\frac{dy}{dx} + P(x)y = Q(x)$

Integrating factor: $\mu(x) = e^{\int P(x)dx}$

Derivation:

Multiply both sides by $\mu(x)$: $\mu(x)\frac{dy}{dx} + \mu(x)P(x)y = \mu(x)Q(x)$

The left side, by the product rule: $\frac{d}{dx}[\mu(x)y] = \mu(x)Q(x)$

Integrating: $\mu(x)y = \int \mu(x)Q(x)dx + C$

$\boxed{y = \frac{1}{\mu(x)}\left[\int \mu(x)Q(x)dx + C\right]}$

Example: RL Circuit

$L\frac{di}{dt} + Ri = E$

Converting to standard form: $\frac{di}{dt} + \frac{R}{L}i = \frac{E}{L}$

Integrating factor: $\mu(t) = e^{(R/L)t}$

$\frac{d}{dt}\left[e^{(R/L)t} \cdot i\right] = \frac{E}{L}e^{(R/L)t}$

Integrating: $e^{(R/L)t} \cdot i = \frac{E}{R}e^{(R/L)t} + C$

Applying initial condition $i(0) = 0$:

$\boxed{i(t) = \frac{E}{R}\left(1 - e^{-(R/L)t}\right)}$

The current converges to steady state $E/R$ with time constant $\tau = L/R$.

2.3 Exact ODEs

Form: $M(x,y)dx + N(x,y)dy = 0$

Exactness condition: $\frac{\partial M}{\partial y} = \frac{\partial N}{\partial x}$

When this condition is satisfied, there exists a potential function $F$ such that $F(x,y) = C$.

Solution procedure:

1. Integrate $\frac{\partial F}{\partial x} = M$ with respect to $x$ to find $F$
2. Use the condition $\frac{\partial F}{\partial y} = N$ to determine the integration constant
3. General solution: $F(x,y) = C$

Example:

$(2xy + y^2)dx + (x^2 + 2xy)dy = 0$

Checking exactness: $\frac{\partial M}{\partial y} = 2x + 2y = \frac{\partial N}{\partial x}$

Since it is exact: $F = \int M\,dx = \int (2xy + y^2)dx = x^2y + xy^2 + g(y)$

Since $\frac{\partial F}{\partial y} = x^2 + 2xy + g'(y) = N = x^2 + 2xy$, we have $g'(y) = 0$, $g(y) = C$

General solution: $x^2y + xy^2 = C$

2.4 Bernoulli Equation

Form: $\frac{dy}{dx} + P(x)y = Q(x)y^n \quad (n \neq 0, 1)$

Substitution: Setting $v = y^{1-n}$ transforms this into a linear ODE.

Since $v' = (1-n)y^{-n}y'$, dividing the original equation by $y^n$:

$y^{-n}\frac{dy}{dx} + P(x)y^{1-n} = Q(x)$

$\frac{1}{1-n}\frac{dv}{dx} + P(x)v = Q(x)$

$\frac{dv}{dx} + (1-n)P(x)v = (1-n)Q(x)$

This is a first-order linear ODE in $v$, solvable by the integrating factor method.

Application: The logistic growth model $\frac{dP}{dt} = rP\left(1 - \frac{P}{K}\right)$ is a Bernoulli equation with $n=2$.

3. Second-Order Linear Ordinary Differential Equations

3.1 Homogeneous Equations

Standard form: $ay'' + by' + cy = 0$

Characteristic equation derivation: Substitute $y = e^{rx}$.

$a r^2 e^{rx} + b r e^{rx} + c e^{rx} = 0$

$\boxed{ar^2 + br + c = 0}$

Three cases arise depending on the sign of the discriminant $\Delta = b^2 - 4ac$.

Case 1: Two distinct real roots ($\Delta > 0$)

For $r_1 \neq r_2$: $y = C_1 e^{r_1 x} + C_2 e^{r_2 x}$

Case 2: Repeated real root ($\Delta = 0$)

For $r_1 = r_2 = r = -b/(2a)$: $y = (C_1 + C_2 x)e^{rx}$

The term $C_2 x e^{rx}$ must be added so that the Wronskian is nonzero, constructing an independent second solution.

Case 3: Complex conjugate roots ($\Delta < 0$)

For $r = \alpha \pm i\beta$ (where $\alpha = -b/(2a)$, $\beta = \sqrt{4ac - b^2}/(2a)$): $y = e^{\alpha x}(C_1\cos\beta x + C_2\sin\beta x)$

This form is converted to real solutions using Euler's formula $e^{i\beta x} = \cos\beta x + i\sin\beta x$.

3.2 Complete Analysis of RLC Circuits

Circuit equation (in terms of charge $q(t)$): $L\frac{d^2q}{dt^2} + R\frac{dq}{dt} + \frac{q}{C} = 0$

In standard form: $q'' + \frac{R}{L}q' + \frac{1}{LC}q = 0$

Characteristic equation: $r^2 + \frac{R}{L}r + \frac{1}{LC} = 0$

$r = -\frac{R}{2L} \pm \sqrt{\left(\frac{R}{2L}\right)^2 - \frac{1}{LC}}$

Defining the damping coefficient $\alpha = R/(2L)$ and resonant frequency $\omega_0 = 1/\sqrt{LC}$:

$r = -\alpha \pm \sqrt{\alpha^2 - \omega_0^2}$

Overdamped: $\alpha > \omega_0$ ($R > 2\sqrt{L/C}$)

$q(t) = C_1 e^{r_1 t} + C_2 e^{r_2 t} \quad (r_1, r_2 < 0)$

A sum of two exponential terms — decays monotonically without oscillation.

Critically Damped: $\alpha = \omega_0$ ($R = 2\sqrt{L/C}$)

$q(t) = (C_1 + C_2 t)e^{-\alpha t}$

Returns to equilibrium most rapidly without oscillation. Very important in control system design.

Underdamped: $\alpha < \omega_0$ ($R < 2\sqrt{L/C}$)

Defining the damped oscillation frequency $\omega_d = \sqrt{\omega_0^2 - \alpha^2}$:

$q(t) = e^{-\alpha t}(C_1\cos\omega_d t + C_2\sin\omega_d t)$

Oscillates while decaying. Commonly seen in communication circuits and vibration systems.

Summary of physical interpretation:

3.3 Non-homogeneous Equations

Form: $ay'' + by' + cy = f(x)$

Complete solution: $y = y_h + y_p$

• $y_h$: General solution of the homogeneous equation (complementary function)
• $y_p$: Particular solution

Method of Undetermined Coefficients

Guess the form of $y_p$ based on the form of $f(x)$.

Note: If the guessed form of $y_p$ overlaps with a term in $y_h$, multiply by $x$.

Example: $y'' + 4y = 3\cos 2x$

Homogeneous solution: $y_h = C_1\cos 2x + C_2\sin 2x$

Since the guessed form of $y_p$ overlaps with $y_h$, modify: $y_p = x(A\cos 2x + B\sin 2x)$

After substitution: $A = 0$, $B = 3/4$

$y_p = \frac{3}{4}x\sin 2x$

This is resonance, where the amplitude increases proportionally with time.

Forced vibration response:

$m\ddot{x} + c\dot{x} + kx = F_0\cos\omega t$

Steady-state particular solution: $x_p(t) = \frac{F_0/k}{\sqrt{(1-r^2)^2 + (2\zeta r)^2}}\cos(\omega t - \phi)$

Where frequency ratio $r = \omega/\omega_n$ and damping ratio $\zeta = c/(2\sqrt{mk})$.

The amplitude is maximum at $r = 1$ (resonance condition).

Variation of Parameters

Used for general $f(x)$ where the method of undetermined coefficients cannot be applied.

Given $y_h = C_1 y_1 + C_2 y_2$:

$y_p = -y_1 \int \frac{y_2 f}{W}dx + y_2 \int \frac{y_1 f}{W}dx$

Wronskian: $W = y_1 y_2' - y_1' y_2$

This formula applies to any continuous function $f(x)$.

4. Systems of Differential Equations

4.1 Matrix Representation

A system of $n$ first-order ODEs:

$\frac{d\mathbf{X}}{dt} = A\mathbf{X}$

Where $\mathbf{X} = [x_1, x_2, \ldots, x_n]^T$ and $A$ is an $n \times n$ coefficient matrix.

4.2 Eigenvalue/Eigenvector Solution

Form of solution: $\mathbf{X} = \mathbf{v}e^{\lambda t}$

Substituting: $\lambda \mathbf{v} = A\mathbf{v}$

So $\lambda$ is an eigenvalue of $A$, and $\mathbf{v}$ is the corresponding eigenvector.

Characteristic equation: $\det(A - \lambda I) = 0$

Example: Two connected tanks system

$\frac{d}{dt}\begin{pmatrix} x_1 \\ x_2 \end{pmatrix} = \begin{pmatrix} -2 & 1 \\ 1 & -2 \end{pmatrix}\begin{pmatrix} x_1 \\ x_2 \end{pmatrix}$

Characteristic equation: $(\lambda+2)^2 - 1 = 0$

$\lambda_1 = -1, \quad \lambda_2 = -3$

Eigenvectors: $\mathbf{v}_1 = [1, 1]^T$, $\mathbf{v}_2 = [1, -1]^T$

General solution: $\mathbf{X} = C_1 \begin{pmatrix} 1 \\ 1 \end{pmatrix} e^{-t} + C_2 \begin{pmatrix} 1 \\ -1 \end{pmatrix} e^{-3t}$

4.3 Application to Circuit Network Analysis

Circuit with two loops (Kirchhoff's voltage law):

$L_1\frac{di_1}{dt} + R_1 i_1 - M\frac{di_2}{dt} = E_1(t)$ $L_2\frac{di_2}{dt} + R_2 i_2 - M\frac{di_1}{dt} = E_2(t)$

Where $M$ is the mutual inductance. This can be written in matrix form and solved using the eigenvalue method or Laplace transform.

5. Laplace Transform

5.1 Definition and Convergence

Definition: $\mathcal{L}\{f(t)\} = F(s) = \int_0^\infty f(t)e^{-st}dt$

$s = \sigma + j\omega$ is a complex variable. For convergence, $\text{Re}(s) > \sigma_c$ (half-plane of convergence) must be satisfied.

Why use the Laplace transform?

• Differential operations are converted to algebraic operations
• Initial conditions are automatically incorporated
• Complex ODEs can be solved as simple algebraic equations
• It forms the basis of the transfer function concept

5.2 Key Transform Pairs

Derivation example: $\mathcal{L}\{e^{at}\}$

$\mathcal{L}\{e^{at}\} = \int_0^\infty e^{at} e^{-st} dt = \int_0^\infty e^{-(s-a)t} dt = \left[\frac{e^{-(s-a)t}}{-(s-a)}\right]_0^\infty = \frac{1}{s-a}$

Converges when $s > a$.

5.3 Key Properties

Linearity: $\mathcal{L}\{af(t) + bg(t)\} = aF(s) + bG(s)$

Time shift (First shifting theorem): $\mathcal{L}\{e^{at}f(t)\} = F(s-a)$

Frequency shift (Second shifting theorem): $\mathcal{L}\{f(t-a)u(t-a)\} = e^{-as}F(s)$

Differentiation transform (most important property): $\mathcal{L}\{f'(t)\} = sF(s) - f(0)$ $\mathcal{L}\{f''(t)\} = s^2F(s) - sf(0) - f'(0)$ $\mathcal{L}\{f^{(n)}(t)\} = s^nF(s) - s^{n-1}f(0) - s^{n-2}f'(0) - \cdots - f^{(n-1)}(0)$

Initial conditions are automatically included in the s-domain.

Integration transform: $\mathcal{L}\left\{\int_0^t f(\tau)d\tau\right\} = \frac{F(s)}{s}$

Convolution theorem: $\mathcal{L}\{(f*g)(t)\} = F(s)G(s)$

Where $(f*g)(t) = \int_0^t f(\tau)g(t-\tau)d\tau$

This theorem mathematically proves that the output of a system is the convolution of the impulse response and the input signal.

Final value theorem: $\lim_{t\to\infty} f(t) = \lim_{s\to 0} sF(s)$

Allows determination of the steady-state value without solving in the t-domain.

Initial value theorem: $\lim_{t\to 0^+} f(t) = \lim_{s\to\infty} sF(s)$

5.4 Inverse Laplace Transform

Partial Fraction Decomposition:

For $F(s) = P(s)/Q(s)$, factor the denominator and decompose into a sum of table transform pairs.

Case 1: Distinct real poles

$F(s) = \frac{N(s)}{(s-p_1)(s-p_2)\cdots(s-p_n)} = \frac{A_1}{s-p_1} + \frac{A_2}{s-p_2} + \cdots$

$A_k = \lim_{s\to p_k}(s-p_k)F(s)$

Case 2: Repeated poles ($s = p$ is an $m$-fold root)

$F(s) = \cdots + \frac{A_m}{(s-p)^m} + \frac{A_{m-1}}{(s-p)^{m-1}} + \cdots + \frac{A_1}{s-p}$

$A_k = \frac{1}{(m-k)!}\left[\frac{d^{m-k}}{ds^{m-k}}(s-p)^m F(s)\right]_{s=p}$

Case 3: Complex conjugate poles

$\frac{As + B}{s^2 + 2\alpha s + (\alpha^2 + \beta^2)} \Rightarrow e^{-\alpha t}(C_1\cos\beta t + C_2\sin\beta t)$

Heaviside Expansion Theorem:

When the denominator is a product of distinct linear factors: $F(s) = \sum_{k=1}^n \frac{N(p_k)}{Q'(p_k)} \cdot \frac{1}{s - p_k}$

5.5 Solving ODEs with the Laplace Transform

Example: RLC Series Circuit Transient Response

$L\frac{d^2i}{dt^2} + R\frac{di}{dt} + \frac{i}{C} = E_0\delta(t)$

Initial conditions: $i(0) = 0$, $i'(0) = 0$ (no charge on capacitor)

Solution:

Applying the Laplace transform to both sides:

$I(s) = \frac{E_0/L}{s^2 + (R/L)s + 1/(LC)}$

Defining $\alpha = R/(2L)$, $\omega_0 = 1/\sqrt{LC}$, $\omega_d = \sqrt{\omega_0^2 - \alpha^2}$ (underdamped case):

$I(s) = \frac{E_0/L}{(s+\alpha)^2 + \omega_d^2}$

Inverse Laplace transform (from table: $\mathcal{L}\{e^{-\alpha t}\sin\omega_d t\} = \omega_d/[(s+\alpha)^2+\omega_d^2]$):

$i(t) = \frac{E_0}{L\omega_d}e^{-\alpha t}\sin\omega_d t$

This is the impulse response of the RLC circuit, which is also its Green's function.

6. Numerical ODE Solutions with Python

6.1 Using scipy.integrate.solve_ivp

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt

# RLC circuit simulation (2nd order ODE -> system of 1st order ODEs)
# L * q'' + R * q' + q/C = Vs (step input)
# State variables: y[0] = q (charge), y[1] = i = dq/dt (current)

def rlc_circuit(t, y, R, L, C, Vs):
    q, dq_dt = y
    # q'' = (Vs - R*i - q/C) / L
    d2q_dt2 = (Vs - R * dq_dt - q / C) / L
    return [dq_dt, d2q_dt2]

# Circuit parameters
R = 100      # Resistance (Ohm)
L = 0.1      # Inductance (H)
C = 1e-6     # Capacitance (F)
Vs = 10.0    # Supply voltage (V)

# Calculate resonant frequency and damping ratio
omega_0 = 1.0 / np.sqrt(L * C)
alpha = R / (2 * L)
zeta = alpha / omega_0
print(f"Resonant frequency: {omega_0/(2*np.pi):.1f} Hz")
print(f"Damping ratio zeta: {zeta:.3f}")
if zeta > 1:
    print("-> Overdamped")
elif zeta == 1:
    print("-> Critically damped")
else:
    print("-> Underdamped")

# Numerical solution
y0 = [0.0, 0.0]  # Initial conditions: q(0) = 0, i(0) = 0
t_span = (0, 0.01)
t_eval = np.linspace(0, 0.01, 2000)

sol = solve_ivp(
    rlc_circuit, t_span, y0,
    t_eval=t_eval,
    args=(R, L, C, Vs),
    method='RK45',
    rtol=1e-8,
    atol=1e-10
)

# Voltage calculation: v_C = q/C
v_C = sol.y[0] / C
i = sol.y[1]

# Visualization
fig, axes = plt.subplots(2, 1, figsize=(10, 8))

axes[0].plot(sol.t * 1e3, v_C, 'b-', linewidth=2, label='Capacitor voltage')
axes[0].axhline(y=Vs, color='r', linestyle='--', label=f'Steady state {Vs}V')
axes[0].set_xlabel('Time (ms)')
axes[0].set_ylabel('Voltage (V)')
axes[0].set_title('RLC Series Circuit - Step Response')
axes[0].legend()
axes[0].grid(True, alpha=0.3)

axes[1].plot(sol.t * 1e3, i * 1e3, 'g-', linewidth=2, label='Current')
axes[1].set_xlabel('Time (ms)')
axes[1].set_ylabel('Current (mA)')
axes[1].set_title('RLC Circuit Current')
axes[1].legend()
axes[1].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('rlc_response.png', dpi=150)
plt.show()

6.2 Analytical Solutions with SymPy

import sympy as sp

# Symbol definitions
t, s = sp.symbols('t s')
R_sym, L_sym, C_sym = sp.symbols('R L C', positive=True)
Vs_sym = sp.Symbol('Vs')

# RC circuit analytical solution (solving ODE with sympy)
v = sp.Function('v')
tau = sp.Symbol('tau', positive=True)

# dv/dt = (Vs - v) / (R*C)
ode = sp.Eq(v(t).diff(t), (Vs_sym - v(t)) / (R_sym * C_sym))

# Initial condition v(0) = 0
solution = sp.dsolve(ode, v(t), ics={v(0): 0})
print("RC circuit solution:")
print(solution)
sp.pprint(solution)

# Laplace transform calculation
f = sp.exp(-t) * sp.sin(2*t)
F = sp.laplace_transform(f, t, s, noconds=True)
print(f"\nL{{e^(-t) sin(2t)}} = {F}")

# Inverse Laplace transform
F_inv = sp.inverse_laplace_transform(
    1 / (s**2 + 2*s + 5), s, t
)
print(f"L^(-1){{1/(s^2+2s+5)}} = {F_inv}")

6.3 Comparing Different Integrators

import numpy as np
import matplotlib.pyplot as plt

# Euler method vs RK4 comparison
# Example: y' = -y, y(0) = 1, exact solution y = e^(-t)

def f(t, y):
    return -y

def euler_method(f, t0, y0, h, n):
    t = np.zeros(n+1)
    y = np.zeros(n+1)
    t[0], y[0] = t0, y0
    for i in range(n):
        y[i+1] = y[i] + h * f(t[i], y[i])
        t[i+1] = t[i] + h
    return t, y

def rk4_method(f, t0, y0, h, n):
    t = np.zeros(n+1)
    y = np.zeros(n+1)
    t[0], y[0] = t0, y0
    for i in range(n):
        k1 = h * f(t[i], y[i])
        k2 = h * f(t[i] + h/2, y[i] + k1/2)
        k3 = h * f(t[i] + h/2, y[i] + k2/2)
        k4 = h * f(t[i] + h, y[i] + k3)
        y[i+1] = y[i] + (k1 + 2*k2 + 2*k3 + k4) / 6
        t[i+1] = t[i] + h
    return t, y

h = 0.5
n = 10
t0, y0 = 0.0, 1.0

t_euler, y_euler = euler_method(f, t0, y0, h, n)
t_rk4, y_rk4 = rk4_method(f, t0, y0, h, n)
t_exact = np.linspace(0, n*h, 200)
y_exact = np.exp(-t_exact)

plt.figure(figsize=(10, 5))
plt.plot(t_exact, y_exact, 'k-', linewidth=2, label='Exact solution')
plt.plot(t_euler, y_euler, 'r--o', label=f'Euler method (h={h})')
plt.plot(t_rk4, y_rk4, 'b--s', label=f'RK4 (h={h})')
plt.xlabel('t')
plt.ylabel("y(t)")
plt.title("Euler vs RK4 Comparison: y' = -y")
plt.legend()
plt.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

# Error comparison
print(f"Euler method error at t=5: {abs(y_euler[-1] - np.exp(-5)):.6f}")
print(f"RK4 error at t=5: {abs(y_rk4[-1] - np.exp(-5)):.6f}")

6.4 Automating Laplace Transforms with Python

import sympy as sp

s, t = sp.symbols('s t', real=True)
a, omega = sp.symbols('a omega', positive=True)

# Auto-generate Laplace transform table
functions = {
    '1': sp.Integer(1),
    't': t,
    't^2': t**2,
    't^n': t**3,
    'e^(at)': sp.exp(a*t),
    'sin(wt)': sp.sin(omega*t),
    'cos(wt)': sp.cos(omega*t),
    't*e^(at)': t * sp.exp(a*t),
    'e^(at)*sin(wt)': sp.exp(a*t) * sp.sin(omega*t),
}

print("Laplace Transform Table:")
print("-" * 50)
for name, func in functions.items():
    F = sp.laplace_transform(func, t, s, noconds=True)
    print(f"L{{{name}}} = {sp.simplify(F)}")

7. Comprehensive Engineering Application Examples

7.1 Transient Analysis of RLC Circuit with Power Source

Problem: A series RLC circuit with $R = 50\,\Omega$, $L = 0.2\,\text{H}$, $C = 100\,\mu\text{F}$ is connected to a $V_s = 100\,\text{V}$ DC source at $t=0$. Initial conditions are $v_C(0) = 0$, $i(0) = 0$.

Solving with the Laplace Transform:

Kirchhoff's voltage law: $L\frac{di}{dt} + Ri + \frac{1}{C}\int i\,dt = V_s u(t)$

Substituting $i = C dv_C/dt$ and writing in terms of charge $q$: $L\ddot{q} + R\dot{q} + \frac{q}{C} = V_s$

Laplace transform: $\left(Ls^2 + Rs + \frac{1}{C}\right)Q(s) = \frac{V_s}{s}$

$Q(s) = \frac{V_s}{s\left(Ls^2 + Rs + \frac{1}{C}\right)} = \frac{V_s/L}{s\left(s^2 + \frac{R}{L}s + \frac{1}{LC}\right)}$

Substituting values: $\alpha = R/(2L) = 125\,\text{s}^{-1}$, $\omega_0 = 1/\sqrt{LC} \approx 223.6\,\text{s}^{-1}$

Since $\omega_0 > \alpha$, this is underdamped.

$\omega_d = \sqrt{\omega_0^2 - \alpha^2} \approx 185.0\,\text{s}^{-1}$

After partial fraction decomposition and inverse Laplace transform:

$v_C(t) = V_s\left[1 - e^{-\alpha t}\left(\cos\omega_d t + \frac{\alpha}{\omega_d}\sin\omega_d t\right)\right]$

$= 100\left[1 - e^{-125t}\left(\cos 185t + 0.676\sin 185t\right)\right]\,\text{V}$

7.2 Resonance Analysis of Spring-Mass-Damper System

Problem: $m = 1\,\text{kg}$, $k = 100\,\mathrm{N/m}$, $c = 4\,\mathrm{N \cdot s/m}$, $F(t) = 10\cos\omega t\,\mathrm{N}$

Natural frequency: $\omega_n = \sqrt{k/m} = 10\,\text{rad/s}$

Damping ratio: $\zeta = c/(2\sqrt{mk}) = 0.2$ (underdamped)

Steady-state response amplitude: $X = \frac{F_0/k}{\sqrt{(1-r^2)^2 + (2\zeta r)^2}}$

At $r = \omega/\omega_n = 1$ (resonance): $X_{resonance} = \frac{F_0/k}{2\zeta} = \frac{10/100}{2 \times 0.2} = 0.25\,\text{m}$

This is 2.5 times the static deflection $F_0/k = 0.1\,\text{m}$ without resonance. This is the danger of resonance.

Summary and Next Steps

Topics covered in this article:

1. Differential equation basics: Classification, order, linearity, types of solutions
2. First-order ODE methods: Separation of variables, integrating factor, exact ODE, Bernoulli equation
3. Second-order linear ODE: Characteristic equation, three cases, complete RLC circuit analysis
4. Non-homogeneous equations: Undetermined coefficients, variation of parameters, forced vibration
5. Systems of ODEs: Matrix representation, eigenvalue/eigenvector solution
6. Laplace transform: Definition, transform pairs, properties, inverse transform, ODE solving
7. Python implementation: Numerical and symbolic solutions with scipy and sympy

The next article will cover Fourier series/transform and partial differential equations (PDEs) — the core mathematical tools for signal processing and electromagnetics.

References

• Kreyszig, E. "Advanced Engineering Mathematics", 10th Edition, Wiley
• Boyce, W. & DiPrima, R. "Elementary Differential Equations", 11th Edition, Wiley
• Simmons, G. "Differential Equations with Applications and Historical Notes", 3rd Edition
• SciPy ODE Solver Documentation
• SymPy Laplace Transform Documentation

Quiz