Consider the continuous-time LTI system \( \dot{\mathbf{x} }=\mathbf{A}\mathbf{x}+\mathbf{B}\mathbf{u} \), where \( \mathbf{x}\in\mathbb{R}^{n} \) and \( \mathbf{u}\in\mathbb{R}^{m} \). Controllability asks whether every state direction can be influenced by the input. System decomposition asks a more precise question: which subspace is influenced, and which coordinate directions are outside the input authority?

flowchart TD
  X["State space R^n"] --> W["Build Wc = [B, AB, ..., A^(n-1)B]"]
  W --> R["Column space of Wc = controllable subspace"]
  R --> Q["Choose basis Qc"]
  X --> U["Complete with complement Qu"]
  Q --> T["T = [Qc Qu]"]
  U --> T
  T --> D["Coordinates: controllable part plus complement"]
        

Starting from zero initial condition, the finite-time state reached at \( T > 0 \) is

\[ \mathbf{x}(T)=\int_{0}^{T}e^{\mathbf{A}(T-s)}\mathbf{B}\mathbf{u}(s)\,ds. \]

The reachable or controllable subspace is

\[ \mathcal{R}(\mathbf{A},\mathbf{B})= \operatorname{Im}(\mathcal{C}),\qquad \mathcal{C}=\begin{bmatrix} \mathbf{B} & \mathbf{AB} & \mathbf{A}^{2}\mathbf{B} & \cdots & \mathbf{A}^{n-1}\mathbf{B} \end{bmatrix}. \]

Therefore the system is controllable exactly when

\[ \operatorname{rank}(\mathcal{C})=n. \]

The controllable subspace is unique. A complementary subspace \( \mathcal{U} \) satisfying \( \mathbb{R}^{n}=\mathcal{R}\oplus\mathcal{U} \) is not unique; it is a coordinate choice. The quotient dynamics on \( \mathbb{R}^{n}/\mathcal{R} \) is the intrinsic uncontrollable part.

Theorem: For a finite-dimensional LTI system, the set of all states reachable from the origin equals \( \operatorname{Im}(\mathcal{C}) \).

Proof, first inclusion. Use the exponential series

\[ e^{\mathbf{A}(T-s)}\mathbf{B}=\sum_{k=0}^{\infty} \frac{(T-s)^{k} }{k!}\mathbf{A}^{k}\mathbf{B}. \]

By Cayley-Hamilton, every power \( \mathbf{A}^{k} \) with \( k\ge n \) is a linear combination of lower powers \( \mathbf{I},\mathbf{A},\ldots,\mathbf{A}^{n-1} \). Hence the integrand and its integral lie in \( \operatorname{Im}(\mathcal{C}) \).

Proof, reverse inclusion. The finite-time controllability Gramian is

\[ \mathbf{W}_{c}(T)=\int_{0}^{T}e^{\mathbf{A}s}\mathbf{B}\mathbf{B}^{T} e^{\mathbf{A}^{T}s}\,ds. \]

From the Gramian result established earlier, \( \operatorname{Im}(\mathbf{W}_{c}(T))=\operatorname{Im}(\mathcal{C}) \) for \( T > 0 \). The Gramian image is the reachable set; hence the reverse inclusion follows.

The controllable subspace is \( \mathbf{A} \)-invariant:

\[ \mathbf{A}\mathcal{R}(\mathbf{A},\mathbf{B})\subseteq \mathcal{R}(\mathbf{A},\mathbf{B}). \]

For any \( \mathbf{v}\in\mathcal{R} \), write

\[ \mathbf{v}=\sum_{k=0}^{n-1}\mathbf{A}^{k}\mathbf{B}\boldsymbol{\alpha}_{k}. \]

Then

\[ \mathbf{A}\mathbf{v}=\sum_{k=0}^{n-1}\mathbf{A}^{k+1}\mathbf{B}\boldsymbol{\alpha}_{k}. \]

The only apparently new term is \( \mathbf{A}^{n}\mathbf{B} \), and Cayley-Hamilton reduces it to a combination of the columns already present in \( \mathcal{C} \). Therefore \( \mathbf{A}\mathbf{v}\in\mathcal{R} \).

Let \( r=\operatorname{rank}(\mathcal{C}) \). Choose \( \mathbf{T}_{c} \) with columns spanning \( \mathcal{R}(\mathbf{A},\mathbf{B}) \), then complete it to a nonsingular matrix

\[ \mathbf{T}=\begin{bmatrix}\mathbf{T}_{c} & \mathbf{T}_{u}\end{bmatrix}, \qquad \mathbf{x}=\mathbf{T}\mathbf{z},\qquad \mathbf{z}=\begin{bmatrix}\mathbf{z}_{c}\\\mathbf{z}_{u}\end{bmatrix}. \]

The transformed system is

\[ \dot{\mathbf{z} }=\bar{\mathbf{A} }\mathbf{z}+\bar{\mathbf{B} }\mathbf{u}, \qquad \bar{\mathbf{A} }=\mathbf{T}^{-1}\mathbf{A}\mathbf{T}, \qquad \bar{\mathbf{B} }=\mathbf{T}^{-1}\mathbf{B}. \]

Because \( \mathcal{R} \) is invariant and \( \operatorname{Im}(\mathbf{B})\subseteq\mathcal{R} \),

\[ \bar{\mathbf{A} }= \begin{bmatrix}\mathbf{A}_{c} & \mathbf{A}_{12}\\ \mathbf{0} & \mathbf{A}_{u}\end{bmatrix}, \qquad \bar{\mathbf{B} }=\begin{bmatrix}\mathbf{B}_{c}\\\mathbf{0}\end{bmatrix}. \]

Thus the uncontrollable complement obeys

\[ \dot{\mathbf{z} }_{u}=\mathbf{A}_{u}\mathbf{z}_{u}. \]

No input appears in this equation, so these modes cannot be assigned by state feedback.

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

\[ \mathcal{C}=\begin{bmatrix}\mathbf{B}&\mathbf{AB}&\mathbf{A}^{2}\mathbf{B}\end{bmatrix} =\begin{bmatrix}0&1&0\\1&0&0\\0&0&0\end{bmatrix}, \qquad \operatorname{rank}(\mathcal{C})=2. \]

Hence \( x_{1} \) and \( x_{2} \) span the reachable subspace, while \( x_{3} \) is an autonomous uncontrollable coordinate with \( \dot{x}_{3}=-2x_{3} \).

A mode \( \lambda \) is uncontrollable if there is a nonzero left eigenvector satisfying

\[ \mathbf{q}^{T}\mathbf{A}=\lambda\mathbf{q}^{T}, \qquad \mathbf{q}^{T}\mathbf{B}=\mathbf{0}. \]

Equivalently,

\[ \operatorname{rank}\begin{bmatrix}\lambda\mathbf{I}-\mathbf{A} & \mathbf{B}\end{bmatrix}<n. \]

With transformed feedback \( \mathbf{u}=-\bar{\mathbf{K} }\mathbf{z} \),

\[ \bar{\mathbf{A} }-\bar{\mathbf{B} }\bar{\mathbf{K} }= \begin{bmatrix} \mathbf{A}_{c}-\mathbf{B}_{c}\mathbf{K}_{c} & \mathbf{A}_{12}-\mathbf{B}_{c}\mathbf{K}_{u}\\ \mathbf{0} & \mathbf{A}_{u} \end{bmatrix}. \]

The block \( \mathbf{A}_{u} \) is unchanged; its eigenvalues are unassignable.

Recommended tools include NumPy/SciPy and python-control in Python, Eigen or Armadillo in C++, EJML or Apache Commons Math in Java, MATLAB Control System Toolbox, Simulink State-Space blocks, and Wolfram Mathematica for exact symbolic verification.

File: Chapter19_Lesson1.py

"""
Chapter19_Lesson1.py
Modern Control - Chapter 19, Lesson 1
Controllable/Uncontrollable Subspaces for continuous-time LTI systems.

Required packages:
    pip install numpy
Optional packages for broader control workflows:
    pip install scipy control slycot
"""

import numpy as np


def controllability_matrix(A: np.ndarray, B: np.ndarray) -> np.ndarray:
    """Return Wc = [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 numerical_rank(M: np.ndarray, tol: float | None = None) -> tuple[int, np.ndarray, float]:
    """Return numerical rank, singular values, and tolerance."""
    s = np.linalg.svd(M, compute_uv=False)
    if tol is None:
        tol = max(M.shape) * np.finfo(float).eps * (s[0] if len(s) else 1.0)
    return int(np.sum(s > tol)), s, tol


def controllability_decomposition(A: np.ndarray, B: np.ndarray, tol: float | None = None):
    """
    Construct an orthonormal coordinate matrix T = [Qc Qu].

    Qc spans the reachable/controllable subspace R(A,B).
    Qu spans the orthogonal complement used here as one possible coordinate
    complement. The complement itself is not unique.
    """
    A = np.asarray(A, dtype=float)
    B = np.asarray(B, dtype=float)
    n = A.shape[0]
    Wc = controllability_matrix(A, B)
    U, s, _ = np.linalg.svd(Wc, full_matrices=True)
    if tol is None:
        tol = max(Wc.shape) * np.finfo(float).eps * (s[0] if len(s) else 1.0)
    r = int(np.sum(s > tol))
    Qc = U[:, :r]
    Qu = U[:, r:n]
    T = np.hstack((Qc, Qu))

    # T is orthogonal because it was built from left singular vectors.
    Abar = T.T @ A @ T
    Bbar = T.T @ B
    return Wc, r, s, tol, Qc, Qu, T, Abar, Bbar


def print_matrix(name: str, M: np.ndarray, digits: int = 6) -> None:
    print(f"\n{name} =")
    print(np.array2string(M, precision=digits, suppress_small=True))


def main() -> None:
    # Example: x3 is an autonomous uncontrollable state.
    A = np.array([
        [0.0, 1.0, 0.0],
        [0.0, 0.0, 0.0],
        [0.0, 0.0, -2.0],
    ])
    B = np.array([
        [0.0],
        [1.0],
        [0.0],
    ])

    Wc, r, s, tol, Qc, Qu, T, Abar, Bbar = controllability_decomposition(A, B)
    n = A.shape[0]

    print_matrix("A", A)
    print_matrix("B", B)
    print_matrix("Wc = [B AB ... A^(n-1)B]", Wc)
    print(f"\nSingular values: {s}")
    print(f"Numerical rank: {r} out of n = {n}; tolerance = {tol:.3e}")
    print("System is controllable." if r == n else "System is not controllable.")

    print_matrix("Qc: basis for controllable subspace", Qc)
    print_matrix("Qu: orthogonal complement basis", Qu)
    print_matrix("T = [Qc Qu]", T)
    print_matrix("Abar = T^T A T", Abar)
    print_matrix("Bbar = T^T B", Bbar)

    lower_left = Abar[r:, :r]
    lower_B = Bbar[r:, :]
    print_matrix("lower-left block Abar_uc", lower_left)
    print_matrix("lower block Bbar_u", lower_B)
    print("\nFor a correct reachable decomposition, both blocks should be numerically zero.")


if __name__ == "__main__":
    main()

File: Chapter19_Lesson1.cpp

/*
Chapter19_Lesson1.cpp
Modern Control - Chapter 19, Lesson 1
Controllable/Uncontrollable Subspaces using Eigen.

Compile example:
    g++ -std=c++17 Chapter19_Lesson1.cpp -I /path/to/eigen -O2 -o Chapter19_Lesson1
*/

#include <Eigen/Dense>
#include <iostream>
#include <iomanip>
#include <limits>
#include <algorithm>

using Eigen::MatrixXd;
using Eigen::VectorXd;

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

int numericalRankFromSingularValues(const VectorXd& s, double tol) {
    int r = 0;
    for (int i = 0; i < s.size(); ++i) {
        if (s(i) > tol) {
            ++r;
        }
    }
    return r;
}

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

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

    MatrixXd Wc = controllabilityMatrix(A, B);

    Eigen::JacobiSVD<MatrixXd> svd(Wc, Eigen::ComputeFullU | Eigen::ComputeFullV);
    VectorXd s = svd.singularValues();
    double sigmaMax = (s.size() > 0) ? s(0) : 1.0;
    double eps = std::numeric_limits<double>::epsilon();
    double tol = std::max(Wc.rows(), Wc.cols()) * eps * sigmaMax;
    int r = numericalRankFromSingularValues(s, tol);
    int n = static_cast<int>(A.rows());

    MatrixXd U = svd.matrixU();
    MatrixXd Qc = U.leftCols(r);
    MatrixXd Qu = U.rightCols(n - r);
    MatrixXd T(n, n);
    T << Qc, Qu;

    MatrixXd Abar = T.transpose() * A * T;
    MatrixXd Bbar = T.transpose() * B;

    std::cout << std::fixed << std::setprecision(6);
    std::cout << "A =\n" << A << "\n\n";
    std::cout << "B =\n" << B << "\n\n";
    std::cout << "Wc = [B AB ... A^(n-1)B] =\n" << Wc << "\n\n";
    std::cout << "Singular values = " << s.transpose() << "\n";
    std::cout << "Rank = " << r << " out of n = " << n << "\n";
    std::cout << ((r == n) ? "System is controllable.\n\n" : "System is not controllable.\n\n");

    std::cout << "Qc basis =\n" << Qc << "\n\n";
    std::cout << "Qu complement basis =\n" << Qu << "\n\n";
    std::cout << "T = [Qc Qu] =\n" << T << "\n\n";
    std::cout << "Abar = T^T A T =\n" << Abar << "\n\n";
    std::cout << "Bbar = T^T B =\n" << Bbar << "\n\n";

    if (r < n) {
        std::cout << "Abar lower-left block =\n" << Abar.block(r, 0, n - r, r) << "\n\n";
        std::cout << "Bbar lower block =\n" << Bbar.block(r, 0, n - r, B.cols()) << "\n";
    }

    return 0;
}

File: Chapter19_Lesson1.java

/*
Chapter19_Lesson1.java
Modern Control - Chapter 19, Lesson 1
Controllable/Uncontrollable Subspaces from scratch.

Compile and run:
    javac Chapter19_Lesson1.java
    java Chapter19_Lesson1

For larger numerical control projects, use EJML, Apache Commons Math, or jblas.
*/

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

public class Chapter19_Lesson1 {
    static final double TOL = 1e-10;

    static double[][] multiply(double[][] A, double[][] B) {
        int n = A.length;
        int p = A[0].length;
        int m = B[0].length;
        double[][] C = new double[n][m];
        for (int i = 0; i < n; i++) {
            for (int k = 0; k < p; k++) {
                for (int j = 0; j < m; j++) {
                    C[i][j] += A[i][k] * B[k][j];
                }
            }
        }
        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[][] hstack(List<double[][]> blocks) {
        int rows = blocks.get(0).length;
        int cols = 0;
        for (double[][] M : blocks) {
            cols += M[0].length;
        }
        double[][] H = new double[rows][cols];
        int c0 = 0;
        for (double[][] M : blocks) {
            for (int i = 0; i < rows; i++) {
                for (int j = 0; j < M[0].length; j++) {
                    H[i][c0 + j] = M[i][j];
                }
            }
            c0 += M[0].length;
        }
        return H;
    }

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

    static int rank(double[][] M) {
        double[][] A = copy(M);
        int rows = A.length;
        int cols = A[0].length;
        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(A[i][c]) > Math.abs(A[pivot][c])) {
                    pivot = i;
                }
            }
            if (Math.abs(A[pivot][c]) <= TOL) {
                continue;
            }
            double[] temp = A[r];
            A[r] = A[pivot];
            A[pivot] = temp;
            double pv = A[r][c];
            for (int j = c; j < cols; j++) {
                A[r][j] /= pv;
            }
            for (int i = 0; i < rows; i++) {
                if (i != r) {
                    double factor = A[i][c];
                    for (int j = c; j < cols; j++) {
                        A[i][j] -= factor * A[r][j];
                    }
                }
            }
            r++;
        }
        return r;
    }

    static double[][] independentColumnBasis(double[][] M) {
        int rows = M.length;
        List<double[]> cols = new ArrayList<>();
        double[][] current = new double[rows][0];
        int currentRank = 0;
        for (int j = 0; j < M[0].length; j++) {
            double[][] candidate = appendColumn(current, column(M, j));
            int newRank = rank(candidate);
            if (newRank > currentRank) {
                cols.add(column(M, j));
                current = candidate;
                currentRank = newRank;
            }
        }
        return columnsToMatrix(cols, rows);
    }

    static double[][] completeBasis(double[][] basis, int n) {
        double[][] T = basis;
        int currentRank = rank(T);
        for (int j = 0; j < n; j++) {
            double[] e = new double[n];
            e[j] = 1.0;
            double[][] candidate = appendColumn(T, e);
            int newRank = rank(candidate);
            if (newRank > currentRank) {
                T = candidate;
                currentRank = newRank;
            }
            if (currentRank == n) {
                break;
            }
        }
        return T;
    }

    static double[][] inverse(double[][] M) {
        int n = M.length;
        double[][] A = new double[n][2 * n];
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                A[i][j] = M[i][j];
            }
            A[i][n + i] = 1.0;
        }
        for (int c = 0; c < n; c++) {
            int pivot = c;
            for (int i = c + 1; i < n; i++) {
                if (Math.abs(A[i][c]) > Math.abs(A[pivot][c])) {
                    pivot = i;
                }
            }
            if (Math.abs(A[pivot][c]) <= TOL) {
                throw new IllegalArgumentException("Matrix is singular.");
            }
            double[] temp = A[c];
            A[c] = A[pivot];
            A[pivot] = temp;
            double pv = A[c][c];
            for (int j = 0; j < 2 * n; j++) {
                A[c][j] /= pv;
            }
            for (int i = 0; i < n; i++) {
                if (i != c) {
                    double factor = A[i][c];
                    for (int j = 0; j < 2 * n; j++) {
                        A[i][j] -= factor * A[c][j];
                    }
                }
            }
        }
        double[][] Inv = new double[n][n];
        for (int i = 0; i < n; i++) {
            for (int j = 0; j < n; j++) {
                Inv[i][j] = A[i][n + j];
            }
        }
        return Inv;
    }

    static double[][] copy(double[][] M) {
        double[][] C = new double[M.length][M[0].length];
        for (int i = 0; i < M.length; i++) {
            C[i] = Arrays.copyOf(M[i], M[i].length);
        }
        return C;
    }

    static double[] column(double[][] M, int j) {
        double[] c = new double[M.length];
        for (int i = 0; i < M.length; i++) {
            c[i] = M[i][j];
        }
        return c;
    }

    static double[][] appendColumn(double[][] M, double[] col) {
        int rows = col.length;
        int cols = M[0].length;
        double[][] N = new double[rows][cols + 1];
        for (int i = 0; i < rows; i++) {
            for (int j = 0; j < cols; j++) {
                N[i][j] = M[i][j];
            }
            N[i][cols] = col[i];
        }
        return N;
    }

    static double[][] columnsToMatrix(List<double[]> cols, int rows) {
        double[][] M = new double[rows][cols.size()];
        for (int j = 0; j < cols.size(); j++) {
            for (int i = 0; i < rows; i++) {
                M[i][j] = cols.get(j)[i];
            }
        }
        return M;
    }

    static void printMatrix(String name, double[][] M) {
        System.out.println("\n" + name + " =");
        for (double[] row : M) {
            for (double v : row) {
                System.out.printf("%10.5f ", v);
            }
            System.out.println();
        }
    }

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

        double[][] Wc = controllabilityMatrix(A, B);
        int r = rank(Wc);
        int n = A.length;
        double[][] Qc = independentColumnBasis(Wc);
        double[][] T = completeBasis(Qc, n);
        double[][] Tinv = inverse(T);
        double[][] Abar = multiply(multiply(Tinv, A), T);
        double[][] Bbar = multiply(Tinv, B);

        printMatrix("A", A);
        printMatrix("B", B);
        printMatrix("Wc = [B AB ... A^(n-1)B]", Wc);
        System.out.println("\nRank = " + r + " out of n = " + n);
        System.out.println(r == n ? "System is controllable." : "System is not controllable.");
        printMatrix("Basis for controllable subspace", Qc);
        printMatrix("T = [controllable basis, complement]", T);
        printMatrix("Abar = inv(T) A T", Abar);
        printMatrix("Bbar = inv(T) B", Bbar);
    }
}

File: Chapter19_Lesson1.m

% Chapter19_Lesson1.m
% Modern Control - Chapter 19, Lesson 1
% Controllable/Uncontrollable Subspaces.
%
% Related MATLAB tools:
%   ctrb(A,B)          - controllability matrix
%   rank(Wc)           - numerical rank
%   orth(Wc)           - orthonormal basis for column space
%   null(Qc')          - orthogonal complement
%   ss(A,B,C,D)        - state-space model object
%   canon(sys,'modal') - canonical-form utility, when applicable

clear; clc;

A = [0 1 0;
     0 0 0;
     0 0 -2];
B = [0; 1; 0];

n = size(A,1);

% Use Control System Toolbox if available; otherwise build manually.
if exist('ctrb','file') == 2
    Wc = ctrb(A,B);
else
    Wc = [];
    Ak = eye(n);
    for k = 1:n
        Wc = [Wc, Ak*B]; %#ok<AGROW>
        Ak = Ak*A;
    end
end

r = rank(Wc);
fprintf('Rank(Wc) = %d out of n = %d\n', r, n);
if r == n
    fprintf('System is controllable.\n');
else
    fprintf('System is not controllable.\n');
end

Qc = orth(Wc);       % basis for controllable subspace
Qu = null(Qc');      % one orthogonal complement
T  = [Qc Qu];        % orthogonal coordinate transform when [Qc Qu] is square

Abar = T' * A * T;
Bbar = T' * B;

disp('Wc = [B AB ... A^(n-1)B]'); disp(Wc);
disp('Qc: basis for controllable subspace'); disp(Qc);
disp('Qu: orthogonal complement basis'); disp(Qu);
disp('Abar = T'' A T'); disp(Abar);
disp('Bbar = T'' B'); disp(Bbar);

if r < n
    disp('Lower-left block of Abar (should be numerically zero):');
    disp(Abar(r+1:end,1:r));
    disp('Lower block of Bbar (should be numerically zero):');
    disp(Bbar(r+1:end,:));
end

% Optional visualization of subspace dimension for this example
figure;
bar([r, n-r]);
set(gca,'XTickLabel',{'controllable','uncontrollable complement'});
ylabel('dimension');
title('Reachable Subspace Dimension Decomposition');
grid on;

File: Chapter19_Lesson1.nb

(* Chapter19_Lesson1.nb *)
ClearAll["Global`*"]

A = { {0, 1, 0}, {0, 0, 0}, {0, 0, -2} };
B = { {0}, {1}, {0} };
n = Length[A];

Wc = ArrayFlatten[{Table[MatrixPower[A, k].B, {k, 0, n - 1}]}];
r = MatrixRank[Wc];

Qc = Transpose[Orthogonalize[Transpose[Wc]]];
Qu = Transpose[NullSpace[Transpose[Qc]]];
T = ArrayFlatten[{ {Qc, Qu} }];

Abar = Chop[Inverse[T].A.T];
Bbar = Chop[Inverse[T].B];

Print["Wc = ", MatrixForm[Wc]];
Print["Rank(Wc) = ", r, " out of n = ", n];
Print["Qc = ", MatrixForm[Qc]];
Print["Qu = ", MatrixForm[Qu]];
Print["Abar = ", MatrixForm[Abar]];
Print["Bbar = ", MatrixForm[Bbar]];

Problem 1: For \( \mathbf{A}=\begin{bmatrix}0&1\\0&0\end{bmatrix} \) and \( \mathbf{B}=\begin{bmatrix}0\\1\end{bmatrix} \), compute the reachable subspace.

Solution:

\[ \mathcal{C}=\begin{bmatrix}0&1\\1&0\end{bmatrix},\qquad \operatorname{rank}(\mathcal{C})=2. \]

Thus \( \mathcal{R}=\mathbb{R}^{2} \).

Problem 2: For \( \mathbf{A}=\operatorname{diag}(1,-3) \), \( \mathbf{B}=[1\;0]^T \), identify the uncontrollable mode.

Solution:

\[ \mathcal{C}=\begin{bmatrix}1&1\\0&0\end{bmatrix},\qquad \operatorname{rank}(\mathcal{C})=1. \]

The reachable subspace is \( \operatorname{span}\{[1\;0]^T\} \), and the uncontrollable eigenvalue is \( -3 \).

Problem 3: Prove that \( \mathcal{R}(\mathbf{A},\mathbf{B}) \) is \( \mathbf{A} \)-invariant.

Solution: If \( \mathbf{v}=\sum_{k=0}^{n-1}\mathbf{A}^{k}\mathbf{B}\boldsymbol{\alpha}_{k} \), then \( \mathbf{A}\mathbf{v} \) contains powers up to \( \mathbf{A}^{n}\mathbf{B} \). Cayley-Hamilton reduces that final term to lower powers, so \( \mathbf{A}\mathbf{v}\in\mathcal{R} \).

Problem 4: Show that transformed \( \bar{\mathbf{B} } \) has zero lower block.

Solution: Since \( \operatorname{Im}(\mathbf{B})\subseteq\mathcal{R} \), every input vector has coordinates only in the basis \( \mathbf{T}_{c} \). Therefore the complement coordinates of \( \bar{\mathbf{B} } \) are zero.

Problem 5: Explain why uncontrollable modes cannot be changed by state feedback.

Solution: The closed-loop transformed matrix remains block upper triangular with lower-right block \( \mathbf{A}_{u} \). Therefore the eigenvalues of \( \mathbf{A}_{u} \) are always part of the closed-loop spectrum.

The controllable subspace is \( \operatorname{Im}[\mathbf{B}\;\mathbf{AB}\;\cdots\;\mathbf{A}^{n-1}\mathbf{B}] \). It is unique and \( \mathbf{A} \)-invariant. If its dimension is less than \( n \), a coordinate transformation separates the controllable coordinates from an uncontrollable complement, whose modes are unassignable by state feedback.