Consider the finite-dimensional continuous-time LTI system \( \dot{\mathbf{x} }(t)=A\mathbf{x}(t)+B\mathbf{u}(t) \), where \( A\in\mathbb{R}^{n\times n} \), \( B\in\mathbb{R}^{n\times m} \), \( \mathbf{x}(t)\in\mathbb{R}^{n} \), and \( \mathbf{u}(t)\in\mathbb{R}^{m} \). In Lesson 1, the Kalman test declared controllability through the rank of \( \mathcal{C}=[B\;AB\;\cdots\;A^{n-1}B] \). The PBH test reaches the same conclusion by checking each eigenvalue of \( A \) separately.

The central idea is modal: a mode associated with \( \lambda\in\sigma(A) \) is controllable only if no nonzero left eigenvector of \( A \) at that eigenvalue is orthogonal to every input direction in \( B \). Thus PBH asks: does the actuator matrix see every left eigenmode?

\[ \boxed{\begin{aligned} (A,B)\text{ controllable} &\iff \operatorname{rank}\left[\lambda I-A\;\;B\right]=n \quad \text{for every }\lambda\in\sigma(A) \\ &\iff \nexists\,\mathbf{q}\neq\mathbf{0}:\; \mathbf{q}^{\ast}A=\lambda\mathbf{q}^{\ast},\;\mathbf{q}^{\ast}B=\mathbf{0}. \end{aligned} } \]

The matrix \( [\lambda I-A\;\;B] \) has \( n \) rows and \( n+m \) columns. PBH requires full row rank at every eigenvalue of \( A \). If the rank drops, the lost row direction identifies an input-invisible left eigenmode.

Theorem (PBH controllability test). The pair \( (A,B) \) is controllable if and only if

\[ \operatorname{rank}\left[\lambda I-A\;\;B\right]=n \quad \forall\lambda\in\sigma(A). \]

Equivalently, using the left-eigenspace \( \mathcal{L}_{\lambda}=\ker(\lambda I-A)^{\ast} \), the PBH condition can be written as

\[ \mathcal{L}_{\lambda}\cap\ker(B^{\ast})=\{\mathbf{0}\} \quad \forall\lambda\in\sigma(A). \]

In words, no left eigenvector is allowed to annihilate all columns of \( B \). For real systems with complex conjugate eigenvalues, the rank test is evaluated over the complex field, because \( \lambda I-A \) may be complex even when \( A \) and \( B \) are real.

flowchart TD
  A0["Start with state matrices A and B"] --> A1["Compute distinct eigenvalues of A"]
  A1 --> A2["For each eigenvalue lambda"]
  A2 --> A3["Build matrix [lambda I - A, B]"]
  A3 --> A4["Compute numerical rank"]
  A4 --> A5{"Rank equals n?"}
  A5 -->|"yes"| A6["Mode at lambda is input-visible"]
  A5 -->|"no"| A7["Uncontrollable mode detected"]
  A6 --> A8{"All eigenvalues checked?"}
  A8 -->|"no"| A2
  A8 -->|"yes"| A9["System is controllable"]
  A7 --> A10["System is not controllable"]
        

We prove the theorem using the Kalman matrix \( \mathcal{C}=[B\;AB\;\cdots\;A^{n-1}B] \). The proof exposes why left eigenvectors are the correct modal objects for controllability.

Step 1: PBH failure implies Kalman rank failure.

Suppose PBH fails at \( \lambda\in\sigma(A) \). Then \( [\lambda I-A\;\;B] \) does not have full row rank, so there exists a nonzero vector \( \mathbf{q}\in\mathbb{C}^{n} \) such that

\[ \mathbf{q}^{\ast}\left[\lambda I-A\;\;B\right]=\mathbf{0}. \]

This is equivalent to the two equations

\[ \mathbf{q}^{\ast}(\lambda I-A)=\mathbf{0}, \qquad \mathbf{q}^{\ast}B=\mathbf{0}. \]

Hence \( \mathbf{q}^{\ast}A=\lambda\mathbf{q}^{\ast} \), and for every integer \( k\ge 0 \),

\[ \mathbf{q}^{\ast}A^{k}B=\lambda^{k}\mathbf{q}^{\ast}B=\mathbf{0}. \]

Therefore

\[ \mathbf{q}^{\ast}\mathcal{C}=\left[\mathbf{q}^{\ast}B\;\; \mathbf{q}^{\ast}AB\;\;\cdots\;\;\mathbf{q}^{\ast}A^{n-1}B\right] =\mathbf{0}. \]

Since a nonzero row vector annihilates \( \mathcal{C} \), \( \operatorname{rank}(\mathcal{C})<n \), so the system is not controllable by the Kalman criterion.

Step 2: Kalman rank failure implies PBH failure.

Suppose \( \operatorname{rank}(\mathcal{C})<n \). Define the reachable subspace generated by \( A \) and \( B \) as

\[ \mathcal{R}=\operatorname{range}(\mathcal{C}) =\operatorname{span}\{B,AB,\dots,A^{n-1}B\}. \]

This subspace is proper and is invariant under \( A \). Invariance follows from Cayley-Hamilton: every \( A^{n}B \) can be expressed as a linear combination of \( B,AB,\dots,A^{n-1}B \). Hence \( A\mathcal{R}\subseteq\mathcal{R} \).

Because \( \mathcal{R} \) is proper, its orthogonal complement \( \mathcal{R}^{\perp} \) is nontrivial. If \( \mathbf{v}\in\mathcal{R}^{\perp} \) and \( \mathbf{r}\in\mathcal{R} \), then

\[ \langle \mathbf{r},A^{\ast}\mathbf{v}\rangle =\langle A\mathbf{r},\mathbf{v}\rangle=0, \]

so \( \mathcal{R}^{\perp} \) is invariant under \( A^{\ast} \). A nonzero finite-dimensional invariant subspace of \( A^{\ast} \) contains an eigenvector; thus there exists \( \mathbf{q}\neq\mathbf{0} \) and \( \lambda\in\sigma(A) \) such that

\[ A^{\ast}\mathbf{q}=\overline{\lambda}\mathbf{q}, \qquad \mathbf{q}^{\ast}A=\lambda\mathbf{q}^{\ast}. \]

Since \( \operatorname{range}(B)\subseteq\mathcal{R} \) and \( \mathbf{q}\in\mathcal{R}^{\perp} \), we also have \( \mathbf{q}^{\ast}B=\mathbf{0} \). Therefore \( \mathbf{q}^{\ast}[\lambda I-A\;\;B]=\mathbf{0} \), and PBH fails. This completes the equivalence.

When \( A \) is diagonalizable with left eigenvectors \( \mathbf{q}_{i}^{\ast}A=\lambda_i\mathbf{q}_{i}^{\ast} \), PBH says that every modal row direction must satisfy \( \mathbf{q}_{i}^{\ast}B\neq\mathbf{0} \). If \( \mathbf{q}_{i}^{\ast}B=\mathbf{0} \), no input can inject motion into that modal coordinate.

\[ \dot{z}_i=\lambda_i z_i+\mathbf{q}_{i}^{\ast}B\mathbf{u}, \qquad z_i=\mathbf{q}_{i}^{\ast}\mathbf{x}. \]

For repeated eigenvalues, it is not enough to excite only one direction inside the eigenspace. PBH requires the input matrix to cover the whole left eigenspace associated with that repeated eigenvalue:

\[ \operatorname{rank}\left( Q_{\lambda}^{\ast}B\right) =\dim\ker(\lambda I-A)^{\ast}, \]

where the columns of \( Q_{\lambda} \) form a basis for the left eigenspace. This immediately explains why repeated eigenvalues may require multiple independent input directions.

flowchart LR
  U["Input channels u"] --> B0["Input matrix B"]
  B0 --> M1["Mode lambda 1: \nvisible"]
  B0 --> M2["Mode lambda 2: \nvisible"]
  B0 -. no projection .-> M3["Mode lambda 3: \ninvisible"]
  M1 --> C0["Controllable contribution"]
  M2 --> C0
  M3 --> C1["Uncontrollable subspace"]
        

Example 1: one unactuated eigenmode. Let

\[ A=\begin{bmatrix}0&0&0\\0&-1&0\\0&0&-2\end{bmatrix}, \qquad B=\begin{bmatrix}0\\1\\1\end{bmatrix}. \]

For \( \lambda=0 \),

\[ [\lambda I-A\;\;B]= \begin{bmatrix}0&0&0&0\\0&1&0&1\\0&0&2&1\end{bmatrix}, \qquad \operatorname{rank}=2<3. \]

Therefore the pair is not controllable. The left eigenvector \( \mathbf{q}^{\top}=[1\;0\;0] \) satisfies \( \mathbf{q}^{\top}A=0\mathbf{q}^{\top} \) and \( \mathbf{q}^{\top}B=0 \), meaning the first state-mode is invisible to the actuator.

Example 2: all diagonal modes actuated. With the same \( A \), choose

\[ B=\begin{bmatrix}1\\1\\1\end{bmatrix}. \]

At \( \lambda=0,-1,-2 \), the PBH matrices all have row rank \( 3 \). Hence the system is controllable. For a diagonal matrix with distinct eigenvalues and a single input, this is equivalent to requiring each component of \( B \) in the eigenbasis to be nonzero.

Example 3: repeated eigenvalue. Let

\[ A=\operatorname{diag}(1,1,2),\qquad B=\begin{bmatrix}1\\1\\1\end{bmatrix}. \]

For \( \lambda=1 \),

\[ [I-A\;\;B]= \begin{bmatrix}0&0&0&1\\0&0&0&1\\0&0&-1&1\end{bmatrix}, \qquad \operatorname{rank}=2<3. \]

Although the input touches the repeated eigenspace, it does not span enough independent directions inside that two-dimensional left eigenspace. A two-input choice such as

\[ B=\begin{bmatrix}1&0\\0&1\\1&1\end{bmatrix} \]

restores full PBH rank at \( \lambda=1 \) and \( \lambda=2 \).

Numerically, PBH is often preferable when one wants to diagnose a specific uncontrollable eigenvalue. However, rank decisions near a singularity depend on the tolerance. For a matrix \( M_{\lambda}=[\lambda I-A\;\;B] \) with singular values \( s_1\ge\cdots\ge s_n\ge 0 \), a practical rank estimate is

\[ \widehat{r}_{\lambda}=\#\{i:s_i>\epsilon\}, \qquad \epsilon\approx \max(n,n+m)\,s_1\,\varepsilon_{\mathrm{mach} }. \]

Useful libraries include NumPy/SciPy/python-control in Python, Eigen in C++, Apache Commons Math or EJML in Java, the Control System Toolbox in MATLAB, and built-in symbolic/numeric linear algebra in Wolfram Mathematica.

Chapter11_Lesson2.py

# Chapter11_Lesson2.py
"""
PBH (Popov-Belevitch-Hautus) controllability test for continuous-time LTI systems.
The pair (A, B) is controllable iff rank([lambda I - A, B]) = n
for every eigenvalue lambda of A.
"""

import numpy as np


def unique_complex_values(values, tol=1e-8):
    """Cluster numerically repeated eigenvalues."""
    unique = []
    for value in values:
        if not any(abs(value - old) <= tol for old in unique):
            unique.append(value)
    return unique


def numerical_rank(M, tol=1e-9):
    """SVD-based numerical rank for real or complex matrices."""
    s = np.linalg.svd(M, compute_uv=False)
    return int(np.sum(s > tol)), s


def controllability_matrix(A, B):
    """Kalman controllability matrix [B, AB, ..., A^(n-1)B]."""
    A = np.asarray(A, dtype=float)
    B = np.asarray(B, dtype=float)
    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 pbh_rank_test(A, B, tol=1e-9, eig_cluster_tol=1e-8):
    """Return PBH controllability verdict and per-eigenvalue diagnostics."""
    A = np.asarray(A, dtype=complex)
    B = np.asarray(B, dtype=complex)
    n = A.shape[0]
    if A.shape != (n, n):
        raise ValueError("A must be square")
    if B.shape[0] != n:
        raise ValueError("B must have the same number of rows as A")

    eigenvalues = unique_complex_values(np.linalg.eigvals(A), eig_cluster_tol)
    details = []
    ok = True
    for lam in eigenvalues:
        M = np.hstack((lam * np.eye(n, dtype=complex) - A, B))
        rank, singular_values = numerical_rank(M, tol)
        passed = (rank == n)
        ok = ok and passed
        details.append({
            "lambda": lam,
            "rank": rank,
            "required_rank": n,
            "passed": passed,
            "singular_values": singular_values,
        })
    return ok, details


def print_pbh_report(A, B, name):
    print(f"\n{name}")
    print("A =\n", np.array(A, dtype=float))
    print("B =\n", np.array(B, dtype=float))

    Ctrb = controllability_matrix(A, B)
    kalman_rank, kalman_s = numerical_rank(Ctrb)
    print("Kalman matrix rank:", kalman_rank)
    print("Kalman singular values:", kalman_s)

    ok, details = pbh_rank_test(A, B)
    print("PBH controllable:", ok)
    for d in details:
        print(
            f"  lambda={d['lambda']:.6g}, "
            f"rank={d['rank']}/{d['required_rank']}, "
            f"passed={d['passed']}, "
            f"singular_values={np.round(d['singular_values'], 6)}"
        )


if __name__ == "__main__":
    # Example 1: diagonal system with an unactuated mode at lambda = 0.
    A1 = np.diag([0.0, -1.0, -2.0])
    B1 = np.array([[0.0], [1.0], [1.0]])
    print_pbh_report(A1, B1, "Example 1: uncontrollable because the first mode is not actuated")

    # Example 2: same A, but every distinct eigenmode is directly reached.
    B2 = np.array([[1.0], [1.0], [1.0]])
    print_pbh_report(A1, B2, "Example 2: controllable diagonal system")

    # Example 3: repeated eigenvalue with one input is not enough to cover the full eigenspace.
    A3 = np.diag([1.0, 1.0, 2.0])
    B3 = np.array([[1.0], [1.0], [1.0]])
    print_pbh_report(A3, B3, "Example 3: uncontrollable repeated-eigenvalue eigenspace")

    # Example 4: two inputs cover the repeated eigenspace at lambda = 1.
    B4 = np.array([[1.0, 0.0], [0.0, 1.0], [1.0, 1.0]])
    print_pbh_report(A3, B4, "Example 4: controllable repeated-eigenvalue case with enough input directions")
      

Chapter11_Lesson2.cpp

// Chapter11_Lesson2.cpp
// PBH controllability test using the Eigen C++ library.
// Compile example:
//   g++ -std=c++17 Chapter11_Lesson2.cpp -I /path/to/eigen -O2 -o pbh_test

#include <Eigen/Dense>
#include <Eigen/Eigenvalues>
#include <Eigen/SVD>
#include <complex>
#include <iostream>
#include <vector>

using Eigen::ComplexEigenSolver;
using Eigen::JacobiSVD;
using Eigen::MatrixXd;
using Eigen::MatrixXcd;
using Eigen::VectorXcd;
using std::complex;
using std::cout;
using std::endl;
using std::vector;

int numericalRank(const MatrixXcd& M, double tol) {
    JacobiSVD<MatrixXcd> svd(M);
    auto s = svd.singularValues();
    int rank = 0;
    for (int i = 0; i < s.size(); ++i) {
        if (s(i) > tol) {
            ++rank;
        }
    }
    return rank;
}

vector<complex<double>> uniqueEigenvalues(const VectorXcd& eig, double tol) {
    vector<complex<double>> values;
    for (int i = 0; i < eig.size(); ++i) {
        bool found = false;
        for (const auto& old : values) {
            if (std::abs(eig(i) - old) <= tol) {
                found = true;
                break;
            }
        }
        if (!found) {
            values.push_back(eig(i));
        }
    }
    return values;
}

MatrixXd controllabilityMatrix(const MatrixXd& A, const MatrixXd& B) {
    const int n = A.rows();
    const int m = B.cols();
    MatrixXd C(n, n * m);
    MatrixXd Ak = MatrixXd::Identity(n, n);
    for (int k = 0; k < n; ++k) {
        C.block(0, k * m, n, m) = Ak * B;
        Ak = Ak * A;
    }
    return C;
}

bool pbhRankTest(const MatrixXd& A, const MatrixXd& B, double tol = 1e-9) {
    const int n = A.rows();
    const int m = B.cols();
    ComplexEigenSolver<MatrixXd> eig(A);
    vector<complex<double>> lambdas = uniqueEigenvalues(eig.eigenvalues(), 1e-8);

    bool controllable = true;
    for (const auto& lambda : lambdas) {
        MatrixXcd M(n, n + m);
        M.block(0, 0, n, n) = lambda * MatrixXcd::Identity(n, n) - A.cast<complex<double>>();
        M.block(0, n, n, m) = B.cast<complex<double>>();
        int r = numericalRank(M, tol);
        cout << "lambda = " << lambda << ", rank([lambda I - A, B]) = " << r << "/" << n << endl;
        if (r != n) {
            controllable = false;
        }
    }
    return controllable;
}

int main() {
    MatrixXd A(3, 3);
    A << 0.0, 0.0, 0.0,
         0.0, -1.0, 0.0,
         0.0, 0.0, -2.0;

    MatrixXd B_bad(3, 1);
    B_bad << 0.0, 1.0, 1.0;

    MatrixXd B_good(3, 1);
    B_good << 1.0, 1.0, 1.0;

    cout << "Example 1: unactuated first mode" << endl;
    cout << "Kalman rank = " << numericalRank(controllabilityMatrix(A, B_bad).cast<complex<double>>(), 1e-9) << endl;
    cout << "PBH controllable? " << (pbhRankTest(A, B_bad) ? "yes" : "no") << "\n" << endl;

    cout << "Example 2: all modes actuated" << endl;
    cout << "Kalman rank = " << numericalRank(controllabilityMatrix(A, B_good).cast<complex<double>>(), 1e-9) << endl;
    cout << "PBH controllable? " << (pbhRankTest(A, B_good) ? "yes" : "no") << endl;

    return 0;
}
      

Chapter11_Lesson2.java

// Chapter11_Lesson2.java
// PBH controllability test for real-eigenvalue examples using Apache Commons Math.
// Compile example:
//   javac -cp commons-math3-3.6.1.jar Chapter11_Lesson2.java
// Run example:
//   java -cp .:commons-math3-3.6.1.jar Chapter11_Lesson2

import java.util.ArrayList;
import java.util.List;
import org.apache.commons.math3.linear.Array2DRowRealMatrix;
import org.apache.commons.math3.linear.EigenDecomposition;
import org.apache.commons.math3.linear.RealMatrix;
import org.apache.commons.math3.linear.SingularValueDecomposition;

public class Chapter11_Lesson2 {
    static double[][] hstackPbh(double[][] A, double[][] B, double lambda) {
        int n = A.length;
        int m = B[0].length;
        double[][] M = new double[n][n + m];
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                M[i][j] = (i == j ? lambda : 0.0) - A[i][j];
            }
            for (int j = 0; j < m; j++) {
                M[i][n + j] = B[i][j];
            }
        }
        return M;
    }

    static int matrixRank(double[][] M) {
        RealMatrix R = new Array2DRowRealMatrix(M);
        return new SingularValueDecomposition(R).getRank();
    }

    static List<Double> uniqueRealEigenvalues(double[][] A, double tol) {
        RealMatrix R = new Array2DRowRealMatrix(A);
        EigenDecomposition ed = new EigenDecomposition(R);
        List<Double> values = new ArrayList<Double>();
        for (int i = 0; i < A.length; i++) {
            double real = ed.getRealEigenvalue(i);
            double imag = ed.getImagEigenvalue(i);
            if (Math.abs(imag) > tol) {
                throw new IllegalArgumentException(
                    "This compact Java demo handles real eigenvalues only. " +
                    "Use a complex linear algebra package for lambda = " + real + " + " + imag + "i."
                );
            }
            boolean found = false;
            for (double old : values) {
                if (Math.abs(real - old) <= tol) {
                    found = true;
                    break;
                }
            }
            if (!found) {
                values.add(real);
            }
        }
        return values;
    }

    static boolean pbhRankTest(double[][] A, double[][] B) {
        int n = A.length;
        boolean ok = true;
        for (double lambda : uniqueRealEigenvalues(A, 1e-8)) {
            double[][] M = hstackPbh(A, B, lambda);
            int rank = matrixRank(M);
            System.out.printf("lambda = %.6f, rank([lambda I - A, B]) = %d/%d%n", lambda, rank, n);
            if (rank != n) {
                ok = false;
            }
        }
        return ok;
    }

    static double[][] controllabilityMatrix(double[][] A, double[][] B) {
        int n = A.length;
        int m = B[0].length;
        double[][] C = new double[n][n * m];
        double[][] Ak = identity(n);
        for (int k = 0; k < n; k++) {
            double[][] block = multiply(Ak, B);
            for (int i = 0; i < n; i++) {
                for (int j = 0; j < m; j++) {
                    C[i][k * m + j] = block[i][j];
                }
            }
            Ak = multiply(Ak, A);
        }
        return C;
    }

    static double[][] identity(int n) {
        double[][] I = new double[n][n];
        for (int i = 0; i < n; i++) {
            I[i][i] = 1.0;
        }
        return I;
    }

    static double[][] multiply(double[][] X, double[][] Y) {
        int r = X.length;
        int c = Y[0].length;
        int inner = Y.length;
        double[][] Z = new double[r][c];
        for (int i = 0; i < r; i++) {
            for (int j = 0; j < c; j++) {
                double sum = 0.0;
                for (int k = 0; k < inner; k++) {
                    sum += X[i][k] * Y[k][j];
                }
                Z[i][j] = sum;
            }
        }
        return Z;
    }

    public static void main(String[] args) {
        double[][] A = {
            {0.0, 0.0, 0.0},
            {0.0, -1.0, 0.0},
            {0.0, 0.0, -2.0}
        };
        double[][] Bbad = { {0.0}, {1.0}, {1.0} };
        double[][] Bgood = { {1.0}, {1.0}, {1.0} };

        System.out.println("Example 1: unactuated first mode");
        System.out.println("Kalman rank = " + matrixRank(controllabilityMatrix(A, Bbad)));
        System.out.println("PBH controllable? " + pbhRankTest(A, Bbad));

        System.out.println("\nExample 2: all modes actuated");
        System.out.println("Kalman rank = " + matrixRank(controllabilityMatrix(A, Bgood)));
        System.out.println("PBH controllable? " + pbhRankTest(A, Bgood));
    }
}
      

In MATLAB, the PBH calculation can be implemented directly and compared with the Kalman matrix generated by ctrb(A,B). In Simulink, the same matrices are assigned in the MATLAB workspace and used inside a State-Space block; the controllability test itself is executed before simulation to verify actuator authority.

Chapter11_Lesson2.m

% Chapter11_Lesson2.m
% PBH controllability test for continuous-time LTI systems.
% Run: Chapter11_Lesson2

function Chapter11_Lesson2()
    A = diag([0, -1, -2]);
    B_bad = [0; 1; 1];
    B_good = [1; 1; 1];

    fprintf('\nExample 1: unactuated first mode\n');
    report_pbh(A, B_bad);

    fprintf('\nExample 2: all modes actuated\n');
    report_pbh(A, B_good);

    A_rep = diag([1, 1, 2]);
    B_rep_bad = [1; 1; 1];
    B_rep_good = [1 0; 0 1; 1 1];

    fprintf('\nExample 3: repeated eigenvalue, one input is insufficient\n');
    report_pbh(A_rep, B_rep_bad);

    fprintf('\nExample 4: repeated eigenvalue, two inputs cover the eigenspace\n');
    report_pbh(A_rep, B_rep_good);
end

function report_pbh(A, B)
    n = size(A, 1);
    C = local_ctrb(A, B);
    fprintf('Kalman controllability rank = %d/%d\n', rank(C), n);

    % If Control System Toolbox is available, this should match local_ctrb.
    if exist('ctrb', 'file') == 2
        fprintf('Control System Toolbox ctrb rank = %d/%d\n', rank(ctrb(A, B)), n);
    end

    [ok, table_out] = pbh_rank_test(A, B, 1e-9);
    disp(table_out);
    fprintf('PBH controllable? %d\n', ok);
end

function [ok, table_out] = pbh_rank_test(A, B, tol)
    n = size(A, 1);
    eigenvalues = eig(A);
    lambdas = unique(round(eigenvalues, 10), 'stable');
    ranks = zeros(length(lambdas), 1);
    passed = false(length(lambdas), 1);

    ok = true;
    for k = 1:length(lambdas)
        lambda = lambdas(k);
        M = [lambda * eye(n) - A, B];
        ranks(k) = rank(M, tol);
        passed(k) = (ranks(k) == n);
        ok = ok && passed(k);
    end

    table_out = table(lambdas, ranks, passed, ...
        'VariableNames', {'lambda', 'PBH_rank', 'passed'});
end

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

Chapter11_Lesson2.nb

(* Chapter11_Lesson2.nb *)
(* PBH controllability test in Wolfram Mathematica / Wolfram Language. *)

ClearAll[pbhRankTest, controllabilityMatrix, reportPBH];

controllabilityMatrix[A_, B_] := Module[{n = Length[A]},
  ArrayFlatten[{Table[MatrixPower[A, k].B, {k, 0, n - 1}]}]
];

pbhRankTest[A_, B_] := Module[{n, lambdas, rows, M, r},
  n = Length[A];
  lambdas = DeleteDuplicates[Eigenvalues[A]];
  rows = Table[
    M = ArrayFlatten[{ {lambda IdentityMatrix[n] - A, B} }];
    r = MatrixRank[M];
    <|"lambda" -> lambda, "PBHRank" -> r, "RequiredRank" -> n, "Passed" -> (r == n)|>,
    {lambda, lambdas}
  ];
  <|"Controllable" -> And @@ (rows[[All, "Passed"]]), "Details" -> rows|>
];

reportPBH[A_, B_, name_] := Module[{C, result},
  Print["\n", name];
  Print["A = ", MatrixForm[A]];
  Print["B = ", MatrixForm[B]];
  C = controllabilityMatrix[A, B];
  Print["Kalman controllability rank = ", MatrixRank[C], "/", Length[A]];
  result = pbhRankTest[A, B];
  Print[Dataset[result["Details"]]];
  Print["PBH controllable? ", result["Controllable"]];
];

A1 = DiagonalMatrix[{0, -1, -2}];
Bbad = { {0}, {1}, {1} };
Bgood = { {1}, {1}, {1} };

reportPBH[A1, Bbad, "Example 1: unactuated first mode"];
reportPBH[A1, Bgood, "Example 2: all modes actuated"];

Arep = DiagonalMatrix[{1, 1, 2}];
BrepBad = { {1}, {1}, {1} };
BrepGood = { {1, 0}, {0, 1}, {1, 1} };

reportPBH[Arep, BrepBad, "Example 3: repeated eigenvalue with insufficient input rank"];
reportPBH[Arep, BrepGood, "Example 4: repeated eigenvalue with enough input rank"];
      

Problem 1 (PBH failure from a left eigenvector): Let \( A\in\mathbb{R}^{n\times n} \) and \( B\in\mathbb{R}^{n\times m} \). Suppose there exist \( \lambda\in\sigma(A) \) and \( \mathbf{q}\neq\mathbf{0} \) such that \( \mathbf{q}^{\top}A=\lambda\mathbf{q}^{\top} \) and \( \mathbf{q}^{\top}B=\mathbf{0} \). Show that the Kalman controllability matrix cannot have full row rank.

Solution: For every \( k\ge 0 \),

\[ \mathbf{q}^{\top}A^kB=\lambda^k\mathbf{q}^{\top}B=\mathbf{0}. \]

Thus \( \mathbf{q}^{\top}[B\;AB\;\cdots\;A^{n-1}B]=\mathbf{0} \). Because \( \mathbf{q}\neq\mathbf{0} \), the rows of \( \mathcal{C} \) are linearly dependent and \( \operatorname{rank}(\mathcal{C})<n \).

Problem 2 (Diagonal system criterion): Let \( A=\operatorname{diag}(\lambda_1,\dots,\lambda_n) \) with distinct eigenvalues and \( B=[b_1\;b_2\;\cdots\;b_n]^{\top} \) be single-input. Prove that \( (A,B) \) is controllable if and only if \( b_i\neq 0 \) for all \( i \).

Solution: At \( \lambda_i \), the matrix \( \lambda_i I-A \) has a zero in row \( i \) and nonzero diagonal entries in all other rows. The only way for \( [\lambda_i I-A\;\;B] \) to recover full row rank is for the appended input column to be nonzero in row \( i \). Therefore PBH requires \( b_i\neq 0 \) for every eigenvalue index. Since the eigenvalues are distinct, these conditions are also sufficient.

Problem 3 (Repeated eigenvalue and number of inputs): Consider \( A=\operatorname{diag}(1,1,2) \) and \( B=[1\;1\;1]^{\top} \). Use PBH to determine whether the pair is controllable.

Solution: At \( \lambda=1 \),

\[ [I-A\;\;B]= \begin{bmatrix}0&0&0&1\\0&0&0&1\\0&0&-1&1\end{bmatrix}. \]

Rows 1 and 2 are identical, so the rank is \( 2 \), not \( 3 \). Therefore the pair is not controllable. The repeated eigenvalue has a two-dimensional left eigenspace, while the single input supplies only one independent direction into that eigenspace.

Problem 4 (PBH and actuator placement): For \( A=\operatorname{diag}(0,-1,-2) \), compare \( B_1=[0\;1\;1]^{\top} \) and \( B_2=[1\;1\;1]^{\top} \).

Solution: For \( B_1 \), the PBH test fails at \( \lambda=0 \), because the first row of \( [0I-A\;\;B_1] \) is zero. For \( B_2 \), each eigenvalue has a nonzero input component in its corresponding diagonal mode. Hence \( (A,B_1) \) is uncontrollable and \( (A,B_2) \) is controllable.

Problem 5 (PBH under similarity transformations): Let \( \tilde{A}=TAT^{-1} \) and \( \tilde{B}=TB \) for a nonsingular matrix \( T \). Prove that PBH controllability is invariant under this state-coordinate transformation.

Solution: Since \( \tilde{A} \) is similar to \( A \), the eigenvalues are the same. For any eigenvalue \( \lambda \),

\[ \begin{aligned} [\lambda I-\tilde{A}\;\;\tilde{B}] &= [\lambda I-TAT^{-1}\;\;TB] \\ &= T[\lambda I-A\;\;B]\begin{bmatrix}T^{-1}&0\\0&I_m\end{bmatrix}. \end{aligned} \]

Left multiplication by \( T \) and right multiplication by the block-diagonal nonsingular matrix preserve row rank. Therefore \( \operatorname{rank}[\lambda I-\tilde{A}\;\;\tilde{B}]=n \) if and only if \( \operatorname{rank}[\lambda I-A\;\;B]=n \).

The PBH test is a modal controllability criterion. It is equivalent to the Kalman rank condition, but it gives more diagnostic information: it identifies which eigenvalues and left eigenvectors are not reached by the input matrix. For simple eigenvalues, each modal input projection must be nonzero. For repeated eigenvalues, the input matrix must cover the full left eigenspace. This makes PBH central to actuator placement, controllability diagnosis, and later pole-placement feasibility.