In lumped-parameter modeling, a small number of coordinates yields ODEs. In distributed-parameter modeling, the state is a function \( u(x,t) \) evolving by a PDE. Separation of variables and modal decomposition replace the PDE by an (infinite) collection of ODEs in time, one per spatial eigenmode.

The key idea is: choose an eigenfunction basis \( \{\phi_n(x)\}_{n=1}^{\infty} \) adapted to the boundary conditions, expand \( u(x,t) = \sum_{n=1}^{\infty} q_n(t)\phi_n(x) \), and project the PDE onto each mode to obtain evolution laws for \( q_n(t) \).

flowchart TD
  A["Given PDE + boundary conditions"] --> B["Assume product form u(x,t)=X(x)T(t) (prototype)"]
  B --> C["Derive spatial eigenproblem (Sturm-Liouville)"]
  C --> D["Compute eigenpairs (lambda_n, phi_n)"]
  D --> E["Expand solution: u = sum q_n(t) phi_n(x)"]
  E --> F["Project PDE onto each phi_n (inner product)"]
  F --> G["Obtain decoupled modal ODEs for q_n(t)"]
  G --> H["Solve modal ODEs + apply initial conditions"]
  H --> I["Reconstruct u(x,t) from truncated or full series"]
        

We start from canonical 1D PDEs on a spatial domain \( x \in (0,L) \). Consider the homogeneous heat equation with constant diffusivity \( \alpha \):

\[ \frac{\partial u}{\partial t}(x,t) = \alpha \frac{\partial^2 u}{\partial x^2}(x,t). \]

Suppose Dirichlet boundary conditions: \( u(0,t)=0 \) and \( u(L,t)=0 \). Seek a separable solution \( u(x,t) = X(x)T(t) \). Substitution yields:

\[ X(x)\,T'(t) = \alpha X''(x)\,T(t). \]

Divide by \( \alpha X(x)T(t) \) (for nontrivial solutions):

\[ \frac{1}{\alpha}\frac{T'(t)}{T(t)} = \frac{X''(x)}{X(x)} = -\lambda, \]

where \( -\lambda \) is a separation constant (the sign is chosen to match typical boundary-value spectra). This gives two ODEs:

\[ \begin{aligned} X''(x) + \lambda X(x) &= 0, \\ T'(t) + \alpha \lambda T(t) &= 0. \end{aligned} \]

The boundary conditions impose \( X(0)=0 \), \( X(L)=0 \), producing a discrete spectrum \( \lambda_n \) and eigenfunctions \( X_n(x) \). For each \( \lambda_n \), the temporal factor satisfies \( T_n(t) = \exp(-\alpha \lambda_n t) \) up to a constant.

A parallel calculation applies to the wave equation:

\[ \frac{\partial^2 u}{\partial t^2}(x,t) = c^2 \frac{\partial^2 u}{\partial x^2}(x,t), \]

where separation yields the same spatial eigenproblem but a different time ODE:

\[ T_n''(t) + c^2 \lambda_n T_n(t) = 0. \]

In system dynamics, boundary conditions are not a side detail: they define the operator domain and therefore the modal basis. A broad class of PDEs leads to a Sturm–Liouville (S–L) problem:

\[ \mathcal{L}\phi := -\frac{d}{dx}\!\left(p(x)\frac{d\phi}{dx}\right) + q(x)\phi = \lambda w(x)\phi,\quad x\in(a,b), \]

with weight \( w(x) > 0 \) and coefficients \( p(x) > 0 \), \( q(x)\ge 0 \). Typical homogeneous boundary conditions are separated (Robin form):

\[ \begin{aligned} \beta_1 \phi(a) + \beta_2 p(a)\phi'(a) &= 0, \\ \gamma_1 \phi(b) + \gamma_2 p(b)\phi'(b) &= 0, \end{aligned} \]

where Dirichlet is \( \phi(a)=0 \), Neumann is \( \phi'(a)=0 \), and Robin mixes both. Under mild conditions, \( \mathcal{L} \) is self-adjoint on the weighted inner product space: \( \langle f,g\rangle_w := \int_a^b f(x)g(x)w(x)\,dx \).

For the heat equation with constant coefficients and Dirichlet BC on \( (0,L) \), we have: \( p(x)=1 \), \( q(x)=0 \), \( w(x)=1 \), giving \( -\phi'' = \lambda \phi \) with \( \phi(0)=\phi(L)=0 \).

The modal decomposition hinges on orthogonality of eigenfunctions and the ability to expand admissible initial conditions in that basis.

4.1 Orthogonality (core proof)

Let \( \phi_m \) and \( \phi_n \) be eigenfunctions with distinct eigenvalues \( \lambda_m \neq \lambda_n \):

\[ \mathcal{L}\phi_n = \lambda_n w \phi_n,\quad \mathcal{L}\phi_m = \lambda_m w \phi_m. \]

Multiply the first equation by \( \phi_m \), the second by \( \phi_n \), subtract, and integrate:

\[ \int_a^b \left(\phi_m \mathcal{L}\phi_n - \phi_n \mathcal{L}\phi_m\right)\,dx = (\lambda_n-\lambda_m)\int_a^b w(x)\phi_n(x)\phi_m(x)\,dx. \]

For a self-adjoint \( \mathcal{L} \), the left-hand side is a boundary term that vanishes under the homogeneous boundary conditions, hence:

\[ (\lambda_n-\lambda_m)\langle \phi_n,\phi_m\rangle_w = 0 \quad \Longrightarrow \quad \langle \phi_n,\phi_m\rangle_w = 0 \;\; (m\neq n). \]

Therefore eigenfunctions associated with distinct eigenvalues are orthogonal in the weighted inner product.

4.2 Expansion and coefficients

Assume \( \{\phi_n\} \) is complete in an appropriate function space (typically \( L^2_w(a,b) \)). For an initial condition \( u(x,0)=u_0(x) \), expand:

\[ u_0(x) = \sum_{n=1}^{\infty} q_n(0)\phi_n(x), \quad q_n(0) = \frac{\langle u_0,\phi_n\rangle_w}{\langle \phi_n,\phi_n\rangle_w}. \]

In the Dirichlet heat equation on \( (0,L) \), \( \phi_n(x)=\sin\!\left(\frac{n\pi x}{L}\right) \) and \( \langle \phi_n,\phi_n\rangle = \frac{L}{2} \), giving the classical sine-series coefficient: \( q_n(0) = \frac{2}{L}\int_0^L u_0(x)\sin\!\left(\frac{n\pi x}{L}\right)\,dx \).

Consider a forced diffusion/heat-type PDE with distributed actuation:

\[ \frac{\partial u}{\partial t}(x,t) = \alpha \frac{\partial^2 u}{\partial x^2}(x,t) + b(x)\,u_{\text{in}}(t), \quad u(0,t)=u(L,t)=0. \]

Expand \( u(x,t)=\sum_{n=1}^{\infty} q_n(t)\phi_n(x) \) with \( \phi_n(x)=\sin\!\left(\frac{n\pi x}{L}\right) \). Using orthogonality, project onto mode \( n \) by taking the inner product with \( \phi_n \). The spatial operator satisfies \( \phi_n'' = -\lambda_n \phi_n \) with \( \lambda_n=\left(\frac{n\pi}{L}\right)^2 \).

This yields the decoupled modal dynamics:

\[ \dot{q}_n(t) = -\alpha \lambda_n q_n(t) + B_n u_{\text{in}}(t), \quad B_n = \frac{\langle b,\phi_n\rangle}{\langle \phi_n,\phi_n\rangle}. \]

For the wave equation with distributed forcing:

\[ \frac{\partial^2 u}{\partial t^2}(x,t) = c^2 \frac{\partial^2 u}{\partial x^2}(x,t) + b(x)\,u_{\text{in}}(t), \]

the modal ODE becomes a second-order oscillator per mode:

\[ \ddot{q}_n(t) + c^2 \lambda_n q_n(t) = B_n u_{\text{in}}(t). \]

This is a major conceptual bridge: distributed dynamics become an infinite cascade of stable first-order modes (heat) or oscillatory modes (wave). Later, we will design controllers on a truncated subset of dominant modes.

For the (unforced) heat equation with Dirichlet BC, the energy-like quantity \( E(t)=\frac{1}{2}\int_0^L u(x,t)^2\,dx \) is non-increasing. Differentiate and use integration by parts:

\[ \frac{dE}{dt} = \int_0^L u\,u_t\,dx = \alpha\int_0^L u\,u_{xx}\,dx = -\alpha\int_0^L (u_x)^2\,dx \le 0, \]

where boundary terms vanish because \( u(0,t)=u(L,t)=0 \). In modal coordinates, \( u=\sum q_n\phi_n \) and orthogonality implies:

\[ E(t) = \frac{1}{2}\sum_{n=1}^{\infty} \langle \phi_n,\phi_n\rangle\, q_n(t)^2, \quad \dot{q}_n(t) = -\alpha\lambda_n q_n(t). \]

Thus each mode decays exponentially at rate \( \alpha\lambda_n \). Higher spatial frequencies (larger \( n \)) decay faster, explaining spatial smoothing in diffusion processes.

In practice, we approximate \( u(x,t) \) with the first \( N \) modes: \( u_N(x,t)=\sum_{n=1}^{N} q_n(t)\phi_n(x) \). For the forced heat equation, the truncated model is:

\[ \dot{\mathbf{q}}(t) = \mathbf{A}\mathbf{q}(t) + \mathbf{B} u_{\text{in}}(t), \quad \mathbf{A} = -\alpha\,\mathrm{diag}(\lambda_1,\dots,\lambda_N), \quad \mathbf{B} = [B_1,\dots,B_N]^{\mathsf{T}}. \]

Measurements are typically spatial samples or weighted integrals. If a sensor measures \( y(t)=\int_0^L c(x)u(x,t)\,dx \), then:

\[ y(t) \approx \mathbf{C}\mathbf{q}(t), \quad C_n = \langle c,\phi_n\rangle. \]

flowchart TD
  P0["Distributed PDE model"] --> P1["Choose eigenbasis from boundary conditions"]
  P1 --> P2["Project onto first N modes"]
  P2 --> P3["Finite ODEs: qdot = A q + B u_in"]
  P3 --> P4["Outputs: y = C q (sensors/averages)"]
  P4 --> P5["Simulate / identify / (later) design control on truncated system"]
  P5 --> P6["Validate truncation: increase N until error acceptable"]
        

Truncation note. For diffusion, higher modes decay rapidly; low-order truncations can be accurate for moderate/large times. For wave dynamics, high modes do not necessarily decay (without damping), so truncation must be validated more carefully.

We implement the Dirichlet 1D heat equation and its forced variant using a sine eigenfunction basis: \( \phi_n(x)=\sin(n\pi x/L) \), \( \lambda_n=(n\pi/L)^2 \). The goal is to compute coefficients \( q_n(0) \) from \( u_0(x) \), integrate modal ODEs, and reconstruct \( u(x,t) \).

8.1 Python (SciPy) — modal ODE simulation + reconstruction

File: Chapter19_Lesson2.py

# Chapter19_Lesson2.py
# Separation of Variables & Modal Decomposition (1D Heat Equation, Dirichlet BC)
#
# PDE: u_t = alpha * u_xx + b(x) * u_in(t),   x in (0, L)
# BC:  u(0,t)=0, u(L,t)=0
# IC:  u(x,0)=u0(x)
#
# Modal expansion: u(x,t) = sum_{n=1..N} q_n(t) * sin(n*pi*x/L)
# ODEs: q_n'(t) = -alpha*(n*pi/L)^2 * q_n(t) + B_n * u_in(t)
#
# Requires: numpy, scipy, matplotlib

import numpy as np
from math import pi
from dataclasses import dataclass

try:
    from scipy.integrate import quad
    from scipy.integrate import solve_ivp
except Exception as e:
    raise RuntimeError("This script requires SciPy (scipy.integrate). Install with: pip install scipy") from e

@dataclass
class Heat1DModal:
    L: float = 1.0
    alpha: float = 0.02
    N: int = 30  # number of modes

    def phi(self, n, x):
        return np.sin(n * pi * x / self.L)

    def lam(self, n):
        return (n * pi / self.L) ** 2

    def project(self, func, n):
        # a_n = (2/L) * ∫_0^L func(x) * sin(n*pi*x/L) dx
        L = self.L
        integrand = lambda x: func(x) * np.sin(n * pi * x / L)
        val, _ = quad(integrand, 0.0, L, limit=200)
        return (2.0 / L) * val

    def build_matrices(self, b_func):
        # diagonal Lambda, and input projection vector B
        n = np.arange(1, self.N + 1)
        Lambda = np.diag((n * pi / self.L) ** 2)
        B = np.array([self.project(b_func, int(k)) for k in n])
        return Lambda, B

    def simulate(self, u0_func, b_func, u_in, t_span=(0.0, 5.0), t_eval=None):
        # initial modal coordinates q(0)
        q0 = np.array([self.project(u0_func, n) for n in range(1, self.N + 1)])

        Lambda, B = self.build_matrices(b_func)

        def ode(t, q):
            return -self.alpha * (Lambda @ q) + B * u_in(t)

        if t_eval is None:
            t_eval = np.linspace(t_span[0], t_span[1], 300)

        sol = solve_ivp(ode, t_span, q0, t_eval=t_eval, method="RK45")
        if not sol.success:
            raise RuntimeError("ODE solve failed: " + str(sol.message))
        return sol.t, sol.y, q0, Lambda, B

    def reconstruct(self, x_grid, q):
        # u(x) from modal coefficients q at a given time
        x = np.asarray(x_grid)
        u = np.zeros_like(x, dtype=float)
        for n in range(1, self.N + 1):
            u += q[n - 1] * np.sin(n * pi * x / self.L)
        return u


def main():
    model = Heat1DModal(L=1.0, alpha=0.02, N=40)

    # Initial condition and distributed input shape
    u0 = lambda x: x * (1.0 - x)          # smooth, satisfies Dirichlet at x=0,1
    b  = lambda x: 4.0 * x * (1.0 - x)    # actuation profile (also Dirichlet-compatible)

    # Input signal
    omega = 2.0
    u_in = lambda t: 0.8 * np.sin(omega * t)

    t, q_hist, q0, Lambda, B = model.simulate(u0, b, u_in, t_span=(0.0, 6.0))

    # Reconstruct solution at a few time instants
    x = np.linspace(0.0, model.L, 200)
    sample_times = [0.0, 1.0, 2.5, 5.0, 6.0]
    idx = [int(np.argmin(np.abs(t - ts))) for ts in sample_times]

    # Plot
    import matplotlib.pyplot as plt

    for k, ts in zip(idx, sample_times):
        u_xt = model.reconstruct(x, q_hist[:, k])
        plt.plot(x, u_xt, label=f"t={t[k]:.2f}")

    plt.xlabel("x")
    plt.ylabel("u(x,t)")
    plt.title("1D Heat equation via modal decomposition (Dirichlet BC)")
    plt.legend()
    plt.grid(True)
    plt.show()

    # Print first few modal time constants
    n = np.arange(1, 6)
    tau = 1.0 / (model.alpha * (n * pi / model.L) ** 2)
    print("First 5 modal decay time constants:", tau)

    # Example: compare truncation levels at final time
    t_final_idx = -1
    u_true = model.reconstruct(x, q_hist[:, t_final_idx])

    for Ntest in [3, 6, 12, 24, model.N]:
        modelN = Heat1DModal(L=model.L, alpha=model.alpha, N=Ntest)
        # reuse the first Ntest coefficients from the simulated solution
        qN = q_hist[:Ntest, t_final_idx]
        uN = modelN.reconstruct(x, qN)
        err = np.linalg.norm(u_true - uN) / np.linalg.norm(u_true)
        print(f"Relative L2 error with N={Ntest:2d}: {err:.3e}")

if __name__ == "__main__":
    main()
      

8.2 C++ — coefficient computation + series evaluation

File: Chapter19_Lesson2.cpp

// Chapter19_Lesson2.cpp
// Separation of Variables & Modal Decomposition (1D Heat Equation, Dirichlet BC)
//
// PDE: u_t = alpha * u_xx,  x in (0,L),  u(0,t)=u(L,t)=0
// Modal series: u(x,t) = sum_{n=1..N} a_n * exp(-alpha*(n*pi/L)^2 * t) * sin(n*pi*x/L)
//
// This demo computes a_n numerically (trapezoidal rule) for u0(x)=x(L-x),
// then evaluates u(x,t) on a grid for a user-specified time.

#include <iostream>
#include <vector>
#include <cmath>
#include <iomanip>

static double u0(double x, double L) {
    return x * (L - x);
}

static double phi(int n, double x, double L) {
    return std::sin(n * M_PI * x / L);
}

static double trapz_coeff(int n, double L, int M) {
    // a_n = (2/L) * ∫_0^L u0(x)*sin(n*pi*x/L) dx
    double h = L / (M - 1);
    double sum = 0.0;
    for (int i = 0; i < M; ++i) {
        double x = i * h;
        double w = (i == 0 || i == M - 1) ? 0.5 : 1.0;
        sum += w * u0(x, L) * phi(n, x, L);
    }
    double integral = h * sum;
    return (2.0 / L) * integral;
}

int main() {
    const double L = 1.0;
    const double alpha = 0.02;
    const int N = 40;      // modes
    const int M = 5001;    // integration points
    const int Nx = 200;    // spatial grid points

    double t;
    std::cout << "Enter time t (e.g., 2.5): ";
    std::cin >> t;

    // Compute coefficients
    std::vector<double> a(N+1, 0.0);
    for (int n = 1; n <= N; ++n) {
        a[n] = trapz_coeff(n, L, M);
    }

    // Evaluate u(x,t)
    std::cout << "# x u(x,t)\\n";
    for (int i = 0; i < Nx; ++i) {
        double x = (L * i) / (Nx - 1);
        double u = 0.0;
        for (int n = 1; n <= N; ++n) {
            double lam = (n * M_PI / L) * (n * M_PI / L);
            u += a[n] * std::exp(-alpha * lam * t) * phi(n, x, L);
        }
        std::cout << std::fixed << std::setprecision(8) << x << " " << u << "\\n";
    }

    std::cerr << "Computed with N=" << N << " modes.\\n";
    return 0;
}
      

8.3 Java — same computation (CSV output)

File: Chapter19_Lesson2.java

// Chapter19_Lesson2.java
// Separation of Variables & Modal Decomposition (1D Heat Equation, Dirichlet BC)
//
// PDE: u_t = alpha * u_xx,  u(0,t)=u(L,t)=0
// Solution: u(x,t) = sum_{n=1..N} a_n exp(-alpha (n*pi/L)^2 t) sin(n*pi x/L)
//
// This program:
// 1) computes a_n numerically using the trapezoidal rule for u0(x)=x(L-x)
// 2) prints u(x,t) on a grid as CSV (x,u)

import java.io.BufferedReader;
import java.io.InputStreamReader;
import java.util.Locale;

public class Chapter19_Lesson2 {
    static double u0(double x, double L) { return x * (L - x); }
    static double phi(int n, double x, double L) { return Math.sin(n * Math.PI * x / L); }

    static double trapzCoeff(int n, double L, int M) {
        double h = L / (M - 1);
        double sum = 0.0;
        for (int i = 0; i < M; i++) {
            double x = i * h;
            double w = (i == 0 || i == M - 1) ? 0.5 : 1.0;
            sum += w * u0(x, L) * phi(n, x, L);
        }
        double integral = h * sum;
        return (2.0 / L) * integral;
    }

    public static void main(String[] args) throws Exception {
        Locale.setDefault(Locale.US);

        final double L = 1.0;
        final double alpha = 0.02;
        final int N = 40;     // modes
        final int M = 5001;   // integration points
        final int Nx = 200;   // spatial grid points

        BufferedReader br = new BufferedReader(new InputStreamReader(System.in));
        System.out.print("Enter time t (e.g., 2.5): ");
        double t = Double.parseDouble(br.readLine().trim());

        double[] a = new double[N + 1];
        for (int n = 1; n <= N; n++) {
            a[n] = trapzCoeff(n, L, M);
        }

        System.out.println("x,u");
        for (int i = 0; i < Nx; i++) {
            double x = (L * i) / (Nx - 1);
            double u = 0.0;
            for (int n = 1; n <= N; n++) {
                double lam = (n * Math.PI / L) * (n * Math.PI / L);
                u += a[n] * Math.exp(-alpha * lam * t) * phi(n, x, L);
            }
            System.out.printf("%.8f,%.10f%n", x, u);
        }

        System.err.println("Computed with N=" + N + " modes.");
    }
}
      

8.4 MATLAB — modal coefficients + plotting (and Simulink hint)

File: Chapter19_Lesson2.m

% Chapter19_Lesson2.m
% Separation of Variables & Modal Decomposition (1D Heat Equation, Dirichlet BC)
%
% PDE: u_t = alpha * u_xx, x in (0,L), u(0,t)=u(L,t)=0
% Series: u(x,t) = sum_{n=1..N} a_n exp(-alpha*(n*pi/L)^2 t) sin(n*pi x/L)
%
% This script computes coefficients a_n numerically for u0(x)=x(L-x),
% then plots u(x,t) at several times.

clear; clc;

L = 1.0;
alpha = 0.02;
N = 40;             % modes
M = 5001;           % integration points for coefficients
x_int = linspace(0, L, M);

u0 = @(x) x .* (L - x);

% Compute coefficients: a_n = (2/L) * int_0^L u0(x)*sin(n*pi*x/L) dx
a = zeros(N,1);
for n = 1:N
    phi = sin(n*pi*x_int/L);
    a(n) = (2/L) * trapz(x_int, u0(x_int) .* phi);
end

% Evaluate on spatial grid
x = linspace(0, L, 200);

times = [0.0, 1.0, 2.5, 5.0, 6.0];
figure; hold on; grid on;
for k = 1:numel(times)
    t = times(k);
    u = zeros(size(x));
    for n = 1:N
        lam = (n*pi/L)^2;
        u = u + a(n) * exp(-alpha*lam*t) * sin(n*pi*x/L);
    end
    plot(x, u, 'DisplayName', sprintf('t=%.2f', t));
end
xlabel('x'); ylabel('u(x,t)');
title('1D Heat equation via modal decomposition (Dirichlet BC)');
legend('show');

% ---- Simulink hint (modal truncation as finite-dimensional LTI) ----
% A = -alpha * diag((1:N).^2*pi^2/L^2)
% qdot = A*q + B*u_in(t),  y = C*q  (choose C to reconstruct u at sensors)
% Use "State-Space" block with matrices A,B,C,D.
      

8.5 Wolfram Mathematica — symbolic coefficients + truncated series

File: Chapter19_Lesson2.nb

(* Chapter19_Lesson2.nb (Wolfram Mathematica Notebook content)
   Separation of Variables & Modal Decomposition (1D Heat Equation, Dirichlet BC)

   PDE: u_t = alpha u_xx, u(0,t)=u(L,t)=0
   u(x,t) = Sum_{n=1..Infinity} a_n Exp[-alpha (n Pi/L)^2 t] Sin[n Pi x/L]
*)

L = 1;
alpha = 0.02;

u0[x_] := x (L - x);

a[n_] := (2/L) Integrate[u0[x] Sin[n Pi x/L], {x, 0, L},
                        Assumptions -> n ∈ Integers && n >= 1];

uSeries[x_, t_, N_] := Sum[a[n] Exp[-alpha (n Pi/L)^2 t] Sin[n Pi x/L], {n, 1, N}];

Plot[Evaluate@Table[uSeries[x, t, 30], {t, {0, 1, 2.5, 5, 6}}], {x, 0, L},
     PlotLegends -> {"t=0", "t=1", "t=2.5", "t=5", "t=6"},
     PlotRange -> All, GridLines -> Automatic,
     AxesLabel -> {"x", "u(x,t)"}]

(* Truncation error indicator (heuristic) *)
With[{t0 = 3},
  Table[{N, NIntegrate[(uSeries[x, t0, 120] - uSeries[x, t0, N])^2, {x, 0, L}]},
        {N, {3, 6, 12, 24, 48}}]
]
      

Problem 1 (Dirichlet spectrum): Solve the eigenvalue problem \( X''(x) + \lambda X(x)=0 \) on \( (0,L) \) with \( X(0)=X(L)=0 \). Find \( \lambda_n \) and \( X_n(x) \).

Solution: For \( \lambda = 0 \), solutions are linear and cannot satisfy both Dirichlet conditions unless trivial. For \( \lambda < 0 \), solutions are hyperbolic and again cannot satisfy both endpoints nontrivially. For \( \lambda > 0 \), write \( \lambda = k^2 \), solutions \( X(x)=A\sin(kx)+B\cos(kx) \). From \( X(0)=0 \) get \( B=0 \). From \( X(L)=0 \) get \( \sin(kL)=0 \), so \( kL=n\pi \).

\[ \lambda_n = \left(\frac{n\pi}{L}\right)^2,\quad X_n(x)=\sin\!\left(\frac{n\pi x}{L}\right),\quad n=1,2,\dots \]

Problem 2 (Orthogonality check): Show that for \( m\neq n \), \( \int_0^L \sin\!\left(\frac{m\pi x}{L}\right)\sin\!\left(\frac{n\pi x}{L}\right)\,dx = 0 \).

Solution: Use the product-to-sum identity: \( \sin(a)\sin(b) = \frac{1}{2}[\cos(a-b)-\cos(a+b)] \). Then both cosine terms integrate to zero over \( (0,L) \) because their arguments are integer multiples of \( \pi \) at endpoints.

Problem 3 (Modal ODE derivation with distributed input): For \( u_t=\alpha u_{xx}+b(x)u_{\text{in}}(t) \) with Dirichlet BC, derive \( \dot{q}_n = -\alpha\lambda_n q_n + B_n u_{\text{in}} \) and give \( B_n \).

Solution: Expand \( u=\sum q_n\phi_n \) in sine basis. Multiply PDE by \( \phi_n \), integrate over \( (0,L) \), and use orthogonality plus \( \phi_n''=-\lambda_n\phi_n \). The forcing projects as:

\[ B_n = \frac{\int_0^L b(x)\phi_n(x)\,dx}{\int_0^L \phi_n(x)^2\,dx} = \frac{2}{L}\int_0^L b(x)\sin\!\left(\frac{n\pi x}{L}\right)\,dx. \]

Problem 4 (Energy decay): For the unforced heat equation with Dirichlet BC, prove \( \frac{d}{dt}\frac{1}{2}\int_0^L u^2\,dx = -\alpha\int_0^L (u_x)^2\,dx \).

Solution: Differentiate under the integral sign: \( \frac{d}{dt}\frac{1}{2}\int u^2 dx = \int u u_t dx \). Substitute \( u_t=\alpha u_{xx} \) and integrate by parts: \( \int u u_{xx} dx = [u u_x]_0^L - \int (u_x)^2 dx \). The boundary term vanishes because \( u(0,t)=u(L,t)=0 \).

Problem 5 (Truncation sanity check): Let the exact heat solution be \( u(x,t)=\sum_{n=1}^{\infty} a_n e^{-\alpha\lambda_n t}\phi_n(x) \). Show that the tail energy satisfies:

\[ \int_0^L \left(u(x,t)-u_N(x,t)\right)^2 dx = \sum_{n=N+1}^{\infty} \langle \phi_n,\phi_n\rangle\, a_n^2 e^{-2\alpha\lambda_n t}. \]

Solution: The error is the omitted tail \( e_N=\sum_{n=N+1}^{\infty} a_n e^{-\alpha\lambda_n t}\phi_n \). Square and integrate. Cross-terms vanish by orthogonality, leaving the stated series. This identity explains why diffusion truncations improve as \( t \) grows: the exponential factors suppress high \( n \).

We derived separation-of-variables for heat and wave PDEs, showed how boundary conditions induce Sturm–Liouville eigenproblems, proved orthogonality, and constructed modal expansions with explicit coefficient formulas. By projecting forced PDEs onto eigenmodes, we obtained decoupled modal ODEs and introduced modal truncation as a principled path to finite-dimensional state-space models for simulation and (later) control.