Consider the continuous-time LTI system \( \dot{\mathbf{x} } = \mathbf{A}\mathbf{x} + \mathbf{B}\mathbf{u} \), \( \mathbf{y} = \mathbf{C}\mathbf{x} + \mathbf{D}\mathbf{u} \). Kalman decomposition is not merely a change of notation. It separates the state vector into coordinates that have different physical roles with respect to actuators and sensors.

\[ \bar{\mathbf{x} } = \begin{bmatrix} \mathbf{x}_{co} \\ \mathbf{x}_{c\bar{o} } \\ \mathbf{x}_{\bar{c}o} \\ \mathbf{x}_{\bar{c}\bar{o} } \end{bmatrix}, \quad \mathbf{x} = \mathbf{T}\bar{\mathbf{x} }, \quad \bar{\mathbf{A} } = \mathbf{T}^{-1}\mathbf{A}\mathbf{T}, \quad \bar{\mathbf{B} } = \mathbf{T}^{-1}\mathbf{B}, \quad \bar{\mathbf{C} } = \mathbf{C}\mathbf{T}. \]

The subscript \( co \) means controllable and observable; \( c\bar{o} \) means controllable but unobservable; \( \bar{c}o \) means uncontrollable but observable; and \( \bar{c}\bar{o} \) means neither controllable nor observable. These are not four different physical devices; they are four questions about the same internal dynamics: input authority, output visibility, both, or neither.

flowchart LR
  U["Input u"] --> CO["co block: \nactuator moves it, \nsensor sees it"]
  U --> CU["c_no block: \nactuator moves it, \nsensor cannot see it"]
  NO["no_o block: \nactuator cannot move it, \nsensor sees initial motion"] --> Y["Output y"]
  CO --> Y
  H["no_no block: \nhidden internal state"] --> H
  CU --> HI["Hidden driven motion"]
  H --> LAT["Latent unmeasured \nevolution"]
        

In a Kalman-decomposed coordinate system, one useful block form is

\[ \bar{\mathbf{A} } = \begin{bmatrix} \mathbf{A}_{11} & \mathbf{0} & \mathbf{A}_{13} & \mathbf{0} \\ \mathbf{A}_{21} & \mathbf{A}_{22} & \mathbf{A}_{23} & \mathbf{A}_{24} \\ \mathbf{0} & \mathbf{0} & \mathbf{A}_{33} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \mathbf{A}_{43} & \mathbf{A}_{44} \end{bmatrix}, \quad \bar{\mathbf{B} } = \begin{bmatrix} \mathbf{B}_1 \\ \mathbf{B}_2 \\ \mathbf{0} \\ \mathbf{0} \end{bmatrix}, \quad \bar{\mathbf{C} } = \begin{bmatrix} \mathbf{C}_1 & \mathbf{0} & \mathbf{C}_3 & \mathbf{0} \end{bmatrix}. \]

The corresponding state equations are

\[ \begin{aligned} \dot{\mathbf{x} }_{co} &= \mathbf{A}_{11}\mathbf{x}_{co} + \mathbf{A}_{13}\mathbf{x}_{\bar{c}o} + \mathbf{B}_1\mathbf{u}, \\ \dot{\mathbf{x} }_{c\bar{o} } &= \mathbf{A}_{21}\mathbf{x}_{co} + \mathbf{A}_{22}\mathbf{x}_{c\bar{o} } + \mathbf{A}_{23}\mathbf{x}_{\bar{c}o} + \mathbf{A}_{24}\mathbf{x}_{\bar{c}\bar{o} } + \mathbf{B}_2\mathbf{u}, \\ \dot{\mathbf{x} }_{\bar{c}o} &= \mathbf{A}_{33}\mathbf{x}_{\bar{c}o}, \\ \dot{\mathbf{x} }_{\bar{c}\bar{o} } &= \mathbf{A}_{43}\mathbf{x}_{\bar{c}o} + \mathbf{A}_{44}\mathbf{x}_{\bar{c}\bar{o} }, \\ \mathbf{y} &= \mathbf{C}_1\mathbf{x}_{co} + \mathbf{C}_3\mathbf{x}_{\bar{c}o} + \mathbf{D}\mathbf{u}. \end{aligned} \]

This block form encodes two invariance facts. The reachable subspace contains \( \mathbf{x}_{co} \) and \( \mathbf{x}_{c\bar{o} } \), so the input cannot directly create \( \mathbf{x}_{\bar{c}o} \) or \( \mathbf{x}_{\bar{c}\bar{o} } \). The unobservable subspace contains \( \mathbf{x}_{c\bar{o} } \) and \( \mathbf{x}_{\bar{c}\bar{o} } \), so these coordinates cannot appear in \( \mathbf{y} \).

The controllable-observable block is the externally active part of the system. For zero initial condition, it is the only block that contributes to the transfer matrix. Physically, it contains modes for which the actuator has authority and the sensor has visibility.

\[ \mathbf{G}(s) = \bar{\mathbf{C} }(s\mathbf{I}-\bar{\mathbf{A} })^{-1} \bar{\mathbf{B} } + \mathbf{D} = \mathbf{C}_1(s\mathbf{I}-\mathbf{A}_{11})^{-1}\mathbf{B}_1 + \mathbf{D}. \]

The controllable-unobservable block is input-driven but sensor-blind. It can absorb control energy or contain internal oscillations without appearing in output data. If such a mode is unstable, the measured output may look acceptable while an internal variable grows.

\[ \mathbf{u} \neq \mathbf{0} \quad \Longrightarrow \quad \mathbf{x}_{c\bar{o} }(t) \text{ may change}, \qquad \bar{\mathbf{C} }\begin{bmatrix}\mathbf{0}\\\mathbf{x}_{c\bar{o} }\\\mathbf{0}\\\mathbf{0}\end{bmatrix}=\mathbf{0}. \]

The uncontrollable-observable block is sensor-visible but not input reachable. It may be excited by initial condition or disturbance, and it may appear as a natural transient in measured output, but it cannot be generated from zero initial condition by the input channel.

\[ \mathbf{x}_{\bar{c}o}(0) \neq \mathbf{0} \quad \Longrightarrow \quad \mathbf{y}_{ic}(t) = \mathbf{C}_3 e^{\mathbf{A}_{33}t} \mathbf{x}_{\bar{c}o}(0). \]

The uncontrollable-unobservable block is latent. It is absent from the zero-initial input-output map and absent from measured output. Its existence usually indicates a redundant state, an internal storage variable not connected to ports, or a coordinate introduced by a nonminimal model.

flowchart LR
  Q1["Can input reach the state?"] -->|"yes"| Q2["Can output see the state?"]
  Q1 -->|"no"| Q3["Can output see the state?"]
  Q2 -->|"yes"| CO["co: transfer-function \ndynamics"]
  Q2 -->|"no"| CU["c_no: hidden \ndriven dynamics"]
  Q3 -->|"yes"| UO["no_o: visible \ninitial-condition dynamics"]
  Q3 -->|"no"| UU["no_no: latent \nredundant dynamics"]
        

Assume \( \bar{\mathbf{x} }(0)=\mathbf{0} \). Since the last two rows of \( \bar{\mathbf{B} } \) are zero and the equations for \( \mathbf{x}_{\bar{c}o} \) and \( \mathbf{x}_{\bar{c}\bar{o} } \) have no input forcing, uniqueness of solutions gives

\[ \mathbf{x}_{\bar{c}o}(0)=\mathbf{0},\quad \mathbf{x}_{\bar{c}\bar{o} }(0)=\mathbf{0} \quad \Longrightarrow \quad \mathbf{x}_{\bar{c}o}(t)=\mathbf{0},\quad \mathbf{x}_{\bar{c}\bar{o} }(t)=\mathbf{0}. \]

Therefore the zero-initial forced equations reduce to

\[ \dot{\mathbf{x} }_{co} = \mathbf{A}_{11}\mathbf{x}_{co} + \mathbf{B}_1\mathbf{u}, \qquad \dot{\mathbf{x} }_{c\bar{o} } = \mathbf{A}_{21}\mathbf{x}_{co}+ \mathbf{A}_{22}\mathbf{x}_{c\bar{o} } + \mathbf{B}_2\mathbf{u}, \qquad \mathbf{y} = \mathbf{C}_1\mathbf{x}_{co} + \mathbf{D}\mathbf{u}. \]

The equation for \( \mathbf{x}_{co} \) does not depend on \( \mathbf{x}_{c\bar{o} } \). Thus \( \mathbf{x}_{c\bar{o} } \) may evolve internally, but it cannot affect \( \mathbf{y} \). Taking the Laplace transform with zero initial condition gives

\[ s\mathbf{X}_{co}(s) = \mathbf{A}_{11}\mathbf{X}_{co}(s) + \mathbf{B}_1\mathbf{U}(s), \quad \mathbf{Y}(s)=\mathbf{C}_1\mathbf{X}_{co}(s)+\mathbf{D}\mathbf{U}(s), \]

\[ \mathbf{Y}(s)=\left[\mathbf{C}_1(s\mathbf{I}-\mathbf{A}_{11})^{-1} \mathbf{B}_1+\mathbf{D}\right]\mathbf{U}(s). \]

Hence the poles of \( \mathbf{G}(s) \) are precisely the externally visible poles remaining after cancellation of uncontrollable or unobservable modes. In Lesson 4 this block was the minimal realization; here it receives its physical interpretation.

A full response can be split into forced response and initial-condition response. The zero-initial forced response uses only the \( co \) block, but the initial-condition response can contain the observable uncontrollable block:

\[ \mathbf{y}(t) = \underbrace{\mathbf{C}_1\int_0^t e^{\mathbf{A}_{11}(t-\tau)}\mathbf{B}_1\mathbf{u}(\tau)d\tau + \mathbf{D}\mathbf{u}(t)}_{\text{forced input-output response} } + \underbrace{\mathbf{C}_1e^{\mathbf{A}_{11}t}\mathbf{x}_{co}(0) + \mathbf{C}_3e^{\mathbf{A}_{33}t}\mathbf{x}_{\bar{c}o}(0)}_{\text{visible initial-condition response} }. \]

The hidden initial-condition contribution from \( \mathbf{x}_{c\bar{o} }(0) \) and \( \mathbf{x}_{\bar{c}\bar{o} }(0) \) is absent from output. This leads to three diagnostic rules:

Rule 1. If a mode is controllable but unobservable, it can be excited by the actuator but cannot be diagnosed from the selected sensor output. Rule 2. If a mode is uncontrollable but observable, the sensor may report transients that no input can remove by direct state steering. Rule 3. If a mode is both uncontrollable and unobservable, it is irrelevant to the external map but may still matter for internal safety if it represents a real physical energy-storage element.

\[ \mathcal{R} = \operatorname{span}\left\{\mathbf{B},\mathbf{A}\mathbf{B}, \dots,\mathbf{A}^{n-1}\mathbf{B}\right\}, \qquad \mathcal{N} = \bigcap_{k=0}^{n-1}\ker(\mathbf{C}\mathbf{A}^k). \]

The four physical sectors are generated by intersecting the reachable directions and unobservable directions with their complements in the decomposition basis. In numerical software, these sectors are usually extracted from rank-revealing factorizations of controllability and observability matrices rather than by visual inspection.

Consider the decomposed system with state order \( [x_{co,1},x_{co,2},x_{c\bar{o} },x_{\bar{c}o},x_{\bar{c}\bar{o} }] \):

\[ \bar{\mathbf{A} } = \begin{bmatrix} 0 & 1 & 0 & 0.2 & 0 \\ -2 & -3 & 0 & 0 & 0 \\ 0 & 0.3 & -4 & 0.1 & 0 \\ 0 & 0 & 0 & -0.5 & 0 \\ 0 & 0 & 0 & 0 & 0.2 \end{bmatrix}, \quad \bar{\mathbf{B} }=\begin{bmatrix}0\\1\\1\\0\\0\end{bmatrix}, \quad \bar{\mathbf{C} }=\begin{bmatrix}1&0&0&1&0\end{bmatrix}, \quad D=0. \]

The controllability and observability matrices are

\[ \mathcal{C} = \begin{bmatrix}\bar{\mathbf{B} } & \bar{\mathbf{A} }\bar{\mathbf{B} } & \cdots & \bar{\mathbf{A} }^{n-1}\bar{\mathbf{B} }\end{bmatrix}, \qquad \mathcal{O} = \begin{bmatrix}\bar{\mathbf{C} } \\ \bar{\mathbf{C} }\bar{\mathbf{A} } \\ \vdots \\ \bar{\mathbf{C} }\bar{\mathbf{A} }^{n-1}\end{bmatrix}. \]

For this example, \( \operatorname{rank}(\mathcal{C})=3 \) and \( \operatorname{rank}(\mathcal{O})=3 \). This agrees with the physical partition: two \( co \) states plus one \( c\bar{o} \) state are reachable, while two \( co \) states plus one \( \bar{c}o \) state are observable. The minimal transfer-function order is therefore \( 2 \).

File: Chapter19_Lesson5.py

The Python version uses NumPy for matrix operations and a from-scratch Euler simulator. For larger control projects, students can also use scipy.signal.StateSpace or the python-control package for state-space and transfer-function analysis.

# Chapter19_Lesson5.py
# Physical Interpretation of Decomposed Subsystems in Kalman Decomposition
# Libraries: NumPy is used for matrix operations. SciPy / python-control can be
# added for advanced state-space and transfer-function workflows.

import numpy as np


def controllability_matrix(A, B):
    """Return [B, AB, ..., A^(n-1)B]."""
    n = A.shape[0]
    blocks = []
    Ak = np.eye(n)
    for _ in range(n):
        blocks.append(Ak @ B)
        Ak = Ak @ A
    return np.hstack(blocks)


def observability_matrix(A, C):
    """Return stacked [C; CA; ...; CA^(n-1)]."""
    n = A.shape[0]
    blocks = []
    Ak = np.eye(n)
    for _ in range(n):
        blocks.append(C @ Ak)
        Ak = Ak @ A
    return np.vstack(blocks)


def euler_simulate(A, B, C, D, u_func, x0, t_final=8.0, dt=0.002):
    """Simple from-scratch forward Euler simulation."""
    t = np.arange(0.0, t_final + dt, dt)
    x = np.zeros((len(t), A.shape[0]))
    y = np.zeros((len(t), C.shape[0]))
    x[0] = x0.reshape(-1)
    for k in range(len(t) - 1):
        u = np.asarray(u_func(t[k])).reshape(B.shape[1], 1)
        xk = x[k].reshape(-1, 1)
        y[k] = (C @ xk + D @ u).reshape(-1)
        dx = A @ xk + B @ u
        x[k + 1] = (xk + dt * dx).reshape(-1)
    u = np.asarray(u_func(t[-1])).reshape(B.shape[1], 1)
    y[-1] = (C @ x[-1].reshape(-1, 1) + D @ u).reshape(-1)
    return t, x, y


def build_kalman_decomposed_example():
    """Create a system already arranged as [co, c_no, no_o, no_no]."""
    A = np.array([
        [0.0,  1.0,  0.0,  0.2,  0.0],
        [-2.0, -3.0, 0.0,  0.0,  0.0],
        [0.0,  0.3, -4.0,  0.1,  0.0],
        [0.0,  0.0,  0.0, -0.5,  0.0],
        [0.0,  0.0,  0.0,  0.0,  0.2]
    ])
    B = np.array([[0.0], [1.0], [1.0], [0.0], [0.0]])
    C = np.array([[1.0, 0.0, 0.0, 1.0, 0.0]])
    D = np.array([[0.0]])
    return A, B, C, D


if __name__ == "__main__":
    A, B, C, D = build_kalman_decomposed_example()
    n = A.shape[0]
    Wc = controllability_matrix(A, B)
    Wo = observability_matrix(A, C)

    print("A =\n", A)
    print("B =\n", B)
    print("C =\n", C)
    print("rank(Wc) =", np.linalg.matrix_rank(Wc), "out of", n)
    print("rank(Wo) =", np.linalg.matrix_rank(Wo), "out of", n)
    print("Expected reachable dimension = 3 = dim(co) + dim(c_no)")
    print("Expected observable dimension = 3 = dim(co) + dim(no_o)")
    print("Minimal input-output order = dim(co) = 2")

    # Zero-initial step response: c_no is internally excited but hidden from y.
    step = lambda t: [1.0]
    t, x, y = euler_simulate(A, B, C, D, step, np.zeros(n), t_final=6.0)
    print("Final zero-initial output under unit step:", y[-1, 0])
    print("Final hidden c_no state under unit step:", x[-1, 2])

    # Initial-condition response: no_o is visible; no_no is invisible.
    x0 = np.array([0.0, 0.0, 0.0, 2.0, 1.0])
    zero = lambda t: [0.0]
    t2, x2, y2 = euler_simulate(A, B, C, D, zero, x0, t_final=6.0)
    print("Initial-condition output y(0) from no_o state:", y2[0, 0])
    print("Final no_no state (can evolve while hidden):", x2[-1, 4])

File: Chapter19_Lesson5.cpp

The C++ version implements matrix multiplication and Gaussian-elimination rank computation from scratch, which is useful when external linear algebra libraries are not allowed. In production, libraries such as Eigen, Armadillo, or Blaze should be preferred.

// Chapter19_Lesson5.cpp
// Physical Interpretation of Decomposed Subsystems in Kalman Decomposition
// Self-contained C++17 implementation: matrix products, powers, and rank.

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

using Matrix = std::vector<std::vector<double>>;

Matrix zeros(int r, int c) { return Matrix(r, std::vector<double>(c, 0.0)); }
Matrix eye(int n) { Matrix I = zeros(n, n); for (int i = 0; i < n; ++i) I[i][i] = 1.0; return I; }

Matrix multiply(const Matrix& A, const Matrix& B) {
    int r = (int)A.size(), m = (int)A[0].size(), c = (int)B[0].size();
    Matrix C = zeros(r, c);
    for (int i = 0; i < r; ++i)
        for (int k = 0; k < m; ++k)
            for (int j = 0; j < c; ++j)
                C[i][j] += A[i][k] * B[k][j];
    return C;
}

Matrix hstack(const std::vector<Matrix>& blocks) {
    int r = (int)blocks[0].size(), ctotal = 0;
    for (const auto& M : blocks) ctotal += (int)M[0].size();
    Matrix H = zeros(r, ctotal);
    int offset = 0;
    for (const auto& M : blocks) {
        for (int i = 0; i < r; ++i)
            for (int j = 0; j < (int)M[0].size(); ++j)
                H[i][offset + j] = M[i][j];
        offset += (int)M[0].size();
    }
    return H;
}

Matrix vstack(const std::vector<Matrix>& blocks) {
    int c = (int)blocks[0][0].size(), rtotal = 0;
    for (const auto& M : blocks) rtotal += (int)M.size();
    Matrix V = zeros(rtotal, c);
    int offset = 0;
    for (const auto& M : blocks) {
        for (int i = 0; i < (int)M.size(); ++i)
            for (int j = 0; j < c; ++j)
                V[offset + i][j] = M[i][j];
        offset += (int)M.size();
    }
    return V;
}

int rank(Matrix M, double tol = 1e-9) {
    int rows = (int)M.size(), cols = (int)M[0].size(), r = 0;
    for (int c = 0; c < cols && r < rows; ++c) {
        int pivot = r;
        for (int i = r + 1; i < rows; ++i)
            if (std::fabs(M[i][c]) > std::fabs(M[pivot][c])) pivot = i;
        if (std::fabs(M[pivot][c]) < tol) continue;
        std::swap(M[r], M[pivot]);
        double div = M[r][c];
        for (int j = c; j < cols; ++j) M[r][j] /= div;
        for (int i = 0; i < rows; ++i) {
            if (i == r) continue;
            double factor = M[i][c];
            for (int j = c; j < cols; ++j) M[i][j] -= factor * M[r][j];
        }
        ++r;
    }
    return r;
}

Matrix controllabilityMatrix(const Matrix& A, const Matrix& B) {
    int n = (int)A.size();
    std::vector<Matrix> blocks;
    Matrix Ak = eye(n);
    for (int k = 0; k < n; ++k) { blocks.push_back(multiply(Ak, B)); Ak = multiply(Ak, A); }
    return hstack(blocks);
}

Matrix observabilityMatrix(const Matrix& A, const Matrix& C) {
    int n = (int)A.size();
    std::vector<Matrix> blocks;
    Matrix Ak = eye(n);
    for (int k = 0; k < n; ++k) { blocks.push_back(multiply(C, Ak)); Ak = multiply(Ak, A); }
    return vstack(blocks);
}

int main() {
    Matrix A = {
        {0.0, 1.0, 0.0, 0.2, 0.0},
        {-2.0, -3.0, 0.0, 0.0, 0.0},
        {0.0, 0.3, -4.0, 0.1, 0.0},
        {0.0, 0.0, 0.0, -0.5, 0.0},
        {0.0, 0.0, 0.0, 0.0, 0.2}
    };
    Matrix B = { {0.0}, {1.0}, {1.0}, {0.0}, {0.0} };
    Matrix C = { {1.0, 0.0, 0.0, 1.0, 0.0} };

    Matrix Wc = controllabilityMatrix(A, B);
    Matrix Wo = observabilityMatrix(A, C);
    std::cout << "rank(Wc) = " << rank(Wc) << " out of 5\n";
    std::cout << "rank(Wo) = " << rank(Wo) << " out of 5\n";
    std::cout << "Reachable states: co plus c_no. Observable states: co plus no_o.\n";
    std::cout << "Only the co block contributes to the zero-initial transfer function.\n";
    return 0;
}

File: Chapter19_Lesson5.java

The Java version is self-contained. For larger numerical projects, students may use EJML, Apache Commons Math, or ojAlgo for matrix computations and rank-revealing decompositions.

// Chapter19_Lesson5.java
// Physical Interpretation of Decomposed Subsystems in Kalman Decomposition
// Self-contained Java implementation of controllability and observability ranks.

import java.util.ArrayList;
import java.util.List;

public class Chapter19_Lesson5 {
    static double[][] zeros(int r, int c) { return new double[r][c]; }
    static double[][] eye(int n) { double[][] I = zeros(n, n); for (int i = 0; i < n; i++) I[i][i] = 1.0; return I; }

    static double[][] multiply(double[][] A, double[][] B) {
        int r = A.length, m = A[0].length, c = B[0].length;
        double[][] C = zeros(r, c);
        for (int i = 0; i < r; i++)
            for (int k = 0; k < m; k++)
                for (int j = 0; j < c; j++)
                    C[i][j] += A[i][k] * B[k][j];
        return C;
    }

    static double[][] hstack(List<double[][]> blocks) {
        int rows = blocks.get(0).length, cols = 0;
        for (double[][] M : blocks) cols += M[0].length;
        double[][] H = zeros(rows, cols);
        int offset = 0;
        for (double[][] M : blocks) {
            for (int i = 0; i < rows; i++)
                for (int j = 0; j < M[0].length; j++)
                    H[i][offset + j] = M[i][j];
            offset += M[0].length;
        }
        return H;
    }

    static double[][] vstack(List<double[][]> blocks) {
        int cols = blocks.get(0)[0].length, rows = 0;
        for (double[][] M : blocks) rows += M.length;
        double[][] V = zeros(rows, cols);
        int offset = 0;
        for (double[][] M : blocks) {
            for (int i = 0; i < M.length; i++)
                for (int j = 0; j < cols; j++)
                    V[offset + i][j] = M[i][j];
            offset += M.length;
        }
        return V;
    }

    static int rank(double[][] input, double tol) {
        int rows = input.length, cols = input[0].length;
        double[][] M = zeros(rows, cols);
        for (int i = 0; i < rows; i++) System.arraycopy(input[i], 0, M[i], 0, cols);
        int r = 0;
        for (int c = 0; c < cols && r < rows; c++) {
            int pivot = r;
            for (int i = r + 1; i < rows; i++)
                if (Math.abs(M[i][c]) > Math.abs(M[pivot][c])) pivot = i;
            if (Math.abs(M[pivot][c]) < tol) continue;
            double[] tmp = M[r]; M[r] = M[pivot]; M[pivot] = tmp;
            double div = M[r][c];
            for (int j = c; j < cols; j++) M[r][j] /= div;
            for (int i = 0; i < rows; i++) {
                if (i == r) continue;
                double factor = M[i][c];
                for (int j = c; j < cols; j++) M[i][j] -= factor * M[r][j];
            }
            r++;
        }
        return r;
    }

    static double[][] controllabilityMatrix(double[][] A, double[][] B) {
        int n = A.length;
        List<double[][]> blocks = new ArrayList<>();
        double[][] Ak = eye(n);
        for (int k = 0; k < n; k++) { blocks.add(multiply(Ak, B)); Ak = multiply(Ak, A); }
        return hstack(blocks);
    }

    static double[][] observabilityMatrix(double[][] A, double[][] C) {
        int n = A.length;
        List<double[][]> blocks = new ArrayList<>();
        double[][] Ak = eye(n);
        for (int k = 0; k < n; k++) { blocks.add(multiply(C, Ak)); Ak = multiply(Ak, A); }
        return vstack(blocks);
    }

    public static void main(String[] args) {
        double[][] A = {
            {0.0, 1.0, 0.0, 0.2, 0.0},
            {-2.0, -3.0, 0.0, 0.0, 0.0},
            {0.0, 0.3, -4.0, 0.1, 0.0},
            {0.0, 0.0, 0.0, -0.5, 0.0},
            {0.0, 0.0, 0.0, 0.0, 0.2}
        };
        double[][] B = { {0.0}, {1.0}, {1.0}, {0.0}, {0.0} };
        double[][] C = { {1.0, 0.0, 0.0, 1.0, 0.0} };

        System.out.println("rank(Wc) = " + rank(controllabilityMatrix(A, B), 1e-9) + " out of 5");
        System.out.println("rank(Wo) = " + rank(observabilityMatrix(A, C), 1e-9) + " out of 5");
        System.out.println("Reachable states: co plus c_no. Observable states: co plus no_o.");
        System.out.println("Only the co block contributes to the zero-initial transfer function.");
    }
}

File: Chapter19_Lesson5.m

MATLAB's Control System Toolbox provides ss, tf, ctrb, and obsv. The script below includes local controllability and observability builders so that the main computations remain transparent. In Simulink, the same matrices can be inserted into a State-Space block.

% Chapter19_Lesson5.m
% Physical Interpretation of Decomposed Subsystems in Kalman Decomposition
% MATLAB / Simulink-oriented script. Uses Control System Toolbox if available,
% but includes local controllability and observability matrix builders.

clear; clc;

% State order: [x_co1, x_co2, x_c_no, x_no_o, x_no_no]
A = [ 0.0  1.0  0.0  0.2  0.0;
     -2.0 -3.0  0.0  0.0  0.0;
      0.0  0.3 -4.0  0.1  0.0;
      0.0  0.0  0.0 -0.5  0.0;
      0.0  0.0  0.0  0.0  0.2 ];
B = [0; 1; 1; 0; 0];
C = [1 0 0 1 0];
D = 0;

Wc = local_ctrb(A, B);
Wo = local_obsv(A, C);

fprintf('rank(Wc) = %d out of %d\n', rank(Wc), size(A,1));
fprintf('rank(Wo) = %d out of %d\n', rank(Wo), size(A,1));
fprintf('Reachable dimension should be 3: co plus c_no.\n');
fprintf('Observable dimension should be 3: co plus no_o.\n');

% If Control System Toolbox is installed, build the state-space model.
if exist('ss','file') == 2
    sys = ss(A, B, C, D);
    disp('State-space model:');
    disp(sys);
    disp('Zero-initial transfer function:');
    disp(tf(sys));
end

% Zero-initial step response: x_c_no is internally excited, but hidden from y.
dt = 0.002; tFinal = 6; t = 0:dt:tFinal;
x = zeros(5, numel(t)); y = zeros(1, numel(t));
for k = 1:numel(t)-1
    u = 1;
    y(k) = C*x(:,k) + D*u;
    x(:,k+1) = x(:,k) + dt*(A*x(:,k) + B*u);
end
y(end) = C*x(:,end) + D;
fprintf('Final step-response output y = %.6f\n', y(end));
fprintf('Final hidden x_c_no = %.6f\n', x(3,end));

% Simulink note:
% Use a State-Space block with A, B, C, D defined above. Feed a Step block into
% the input and route y to a Scope. To inspect hidden states, enable state
% logging or replace the block with an integrator-level realization.

function Wc = local_ctrb(A, B)
    n = size(A,1);
    Wc = [];
    Ak = eye(n);
    for k = 1:n
        Wc = [Wc Ak*B]; %#ok<AGROW>
        Ak = Ak*A;
    end
end

function Wo = local_obsv(A, C)
    n = size(A,1);
    Wo = [];
    Ak = eye(n);
    for k = 1:n
        Wo = [Wo; C*Ak]; %#ok<AGROW>
        Ak = Ak*A;
    end
end

File: Chapter19_Lesson5.nb

Mathematica can compute ranks symbolically or numerically. The notebook expression below builds controllability and observability matrices and converts the state-space model into a transfer-function model.

Notebook[{
Cell["Chapter19_Lesson5.nb", "Title"],
Cell["Physical Interpretation of Decomposed Subsystems in Kalman Decomposition", "Subtitle"],
Cell[BoxData[ToBoxes[
A = { {0, 1, 0, 0.2, 0}, {-2, -3, 0, 0, 0}, {0, 0.3, -4, 0.1, 0}, {0, 0, 0, -0.5, 0}, {0, 0, 0, 0, 0.2} };
B = { {0}, {1}, {1}, {0}, {0} };
Cmat = { {1, 0, 0, 1, 0} };
n = Length[A];
controllabilityMatrix = Transpose[Flatten[Table[MatrixPower[A, k].B, {k, 0, n - 1}], 1]];
observabilityMatrix = Flatten[Table[Cmat.MatrixPower[A, k], {k, 0, n - 1}], 1];
{MatrixRank[controllabilityMatrix], MatrixRank[observabilityMatrix]}
]], "Input"],
Cell[BoxData[ToBoxes[
sys = StateSpaceModel[{A, B, Cmat, { {0} } }];
TransferFunctionModel[sys]
]], "Input"],
Cell["The ranks identify reachable dimension 3 and observable dimension 3. The zero-initial transfer function is determined only by the controllable-observable block.", "Text"]
}]

Problem 1 (Transfer Function of the Decomposed System): For the block form in Section 2, prove that the zero-initial transfer matrix is \( \mathbf{G}(s)=\mathbf{C}_1(s\mathbf{I}-\mathbf{A}_{11})^{-1}\mathbf{B}_1+\mathbf{D} \).

Solution: With zero initial condition, the uncontrollable coordinates satisfy homogeneous equations and remain zero. The output is therefore \( \mathbf{y}=\mathbf{C}_1\mathbf{x}_{co}+\mathbf{D}\mathbf{u} \). The \( co \) state satisfies \( \dot{\mathbf{x} }_{co}=\mathbf{A}_{11}\mathbf{x}_{co}+\mathbf{B}_1\mathbf{u} \). Taking the Laplace transform gives the claimed expression.

Problem 2 (Hidden Driven State): Suppose \( \mathbf{x}_{c\bar{o} } \) is scalar and \( \dot{x}_{c\bar{o} }=-x_{c\bar{o} }+u \), while the output is independent of \( x_{c\bar{o} } \). For unit-step input and zero initial condition, find \( x_{c\bar{o} }(t) \) and explain why it is hidden.

Solution: The solution is

\[ x_{c\bar{o} }(t)=\int_0^t e^{-(t-\tau)}d\tau = 1-e^{-t}. \]

It is hidden because the output matrix has zero column in the \( c\bar{o} \) coordinate, so \( y(t) \) contains no direct or indirect dependence on this state in the decomposed form.

Problem 3 (Visible but Unreachable Initial State): Let \( \dot{x}_{\bar{c}o}=-2x_{\bar{c}o} \) and \( y=x_{\bar{c}o} \). If \( x_{\bar{c}o}(0)=3 \) and \( u(t)=0 \), find \( y(t) \). Does this mode appear in the zero-initial transfer function?

Solution:

\[ y(t)=3e^{-2t}. \]

The mode is observable because it appears in the output initial-condition response. It is not controllable, so it cannot be generated by the input from zero initial condition. Therefore it does not appear in the zero-initial transfer function.

Problem 4 (Rank Interpretation): A fifth-order realization has \( \operatorname{rank}(\mathcal{C})=3 \), \( \operatorname{rank}(\mathcal{O})=4 \), and a controllable-observable block of dimension \( 2 \). Determine the dimensions of the other three sectors.

Solution: Since reachable dimension is \( 3 \),

\[ \dim(c\bar{o})=\operatorname{rank}(\mathcal{C})-\dim(co)=3-2=1. \]

Since observable dimension is \( 4 \),

\[ \dim(\bar{c}o)=\operatorname{rank}(\mathcal{O})-\dim(co)=4-2=2. \]

The total dimension is \( 5 \), so

\[ \dim(\bar{c}\bar{o})=5-2-1-2=0. \]

Problem 5 (Minimal Realization Decision): Explain why removing the \( c\bar{o} \), \( \bar{c}o \), and \( \bar{c}\bar{o} \) blocks preserves the zero-initial transfer function but may not preserve all internal initial-condition behavior.

Solution: The zero-initial transfer function depends only on the forced input-output map. By the proof in Section 4, this map is fully determined by the \( co \) block. However, \( \bar{c}o \) initial states can create measured output transients, and hidden states may represent internal variables with their own physical evolution. Thus minimal realization preserves external equivalence under zero initial condition, not necessarily every internal trajectory of the original nonminimal model.

Kalman decomposition assigns physical meaning to state coordinates by asking whether each coordinate is reachable from inputs and visible in outputs. The \( co \) block is the minimal externally active subsystem. The \( c\bar{o} \) block is driven but hidden, the \( \bar{c}o \) block is visible only through initial conditions or disturbances, and the \( \bar{c}\bar{o} \) block is latent from the input-output viewpoint. This interpretation prepares students for minimal realizations, zero dynamics, and later feedback-design limits.