Chapter 18: Jordan Canonical Form and General Modal Decomposition

Lesson 2: Construction of Jordan Canonical Form (Conceptual)

This lesson develops the conceptual construction of Jordan canonical form from generalized eigenspaces, nullity growth, and Jordan chains. The goal is not to promote numerical Jordan computation, but to understand what the similarity transformation means for state-space dynamics and modal decomposition.

1. Why Construction Matters in Modern Control

In earlier lessons, diagonal modal form was introduced for matrices with enough independent eigenvectors. Jordan form is the extension needed when a repeated eigenvalue does not have enough eigenvectors. For a continuous-time LTI system \( \dot{\mathbf{x} }=\mathbf{A}\mathbf{x}+\mathbf{B}\mathbf{u} \), a coordinate change \( \mathbf{x}=\mathbf{P}\mathbf{z} \) gives

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

The columns of \( \mathbf{P} \) are not arbitrary: they are ordered eigenvectors and generalized eigenvectors. If \( \mathbf{J} \) is diagonal, each state component is a pure exponential mode. If \( \mathbf{J} \) has Jordan blocks, a mode contains polynomial factors multiplying the exponential, which will be studied in Lesson 3.

\[ \mathbf{A}\mathbf{P}=\mathbf{P}\mathbf{J},\qquad \mathbf{A}\sim\mathbf{J}. \]

2. Construction Workflow

The Jordan construction is best understood as an algebraic workflow: determine the eigenvalues, measure how many independent chains exist, determine chain lengths from nullity growth, then assemble a basis of chain vectors.

flowchart TD
  A["Start with matrix A"] --> B["Find eigenvalues and algebraic multiplicities"]
  B --> C["For each lambda compute kernels of (A-lambda I)^k"]
  C --> D["Use nullity growth to infer Jordan block sizes"]
  D --> E["Choose chain-head generalized eigenvectors"]
  E --> F["Generate each chain backward to eigenvectors"]
  F --> G["Assemble P from chain vectors"]
  G --> H["Compute J = inv(P) A P"]
        

\[ \mathbf{J}= \operatorname{diag}\left( \mathbf{J}_{s_1}(\lambda_1),\ldots, \mathbf{J}_{s_r}(\lambda_r) \right). \]

3. Jordan Blocks and Generalized Eigenspaces

A Jordan block of size \( s \) associated with \( \lambda \) is

\[ \mathbf{J}_s(\lambda)= \begin{bmatrix} \lambda & 1 & 0 & \cdots & 0 \\ 0 & \lambda & 1 & \cdots & 0 \\ \vdots & \vdots & \ddots & \ddots & \vdots \\ 0 & 0 & \cdots & \lambda & 1 \\ 0 & 0 & \cdots & 0 & \lambda \end{bmatrix} =\lambda\mathbf{I}_s+\mathbf{N}_s,\qquad \mathbf{N}_s^s=\mathbf{0},\quad \mathbf{N}_s^{s-1}\ne\mathbf{0}. \]

The generalized eigenspace associated with \( \lambda \) is the stabilized kernel

\[ \mathcal{G}_\lambda =\ker\left((\mathbf{A}-\lambda\mathbf{I})^q\right),\qquad q\text{ sufficiently large}. \]

In finite dimension, taking \( q=n \) is always sufficient. The ordinary eigenspace is only the first kernel:

\[ \mathcal{E}_\lambda =\ker(\mathbf{A}-\lambda\mathbf{I}),\qquad \mathcal{E}_\lambda\subseteq\mathcal{G}_\lambda. \]

\[ 1\le g_\lambda\le a_\lambda,\qquad g_\lambda=a_\lambda\quad\Longleftrightarrow\quad \text{all Jordan blocks for }\lambda\text{ have size }1. \]

4. Jordan Chains and the Similarity Matrix

A Jordan chain of length \( s \) for \( \lambda \) is an ordered set of nonzero vectors \( \mathbf{v}_1,\ldots,\mathbf{v}_s \) satisfying

\[ (\mathbf{A}-\lambda\mathbf{I})\mathbf{v}_1=\mathbf{0}, \qquad (\mathbf{A}-\lambda\mathbf{I})\mathbf{v}_{j+1} =\mathbf{v}_j,\quad j=1,\ldots,s-1. \]

Equivalently, \( \mathbf{v}_1 \) is an eigenvector and \( \mathbf{v}_2,\ldots,\mathbf{v}_s \) are generalized eigenvectors. In this ordered basis,

\[ \mathbf{A}\mathbf{v}_1=\lambda\mathbf{v}_1,\qquad \mathbf{A}\mathbf{v}_{j+1}= \mathbf{v}_j+\lambda\mathbf{v}_{j+1}. \]

flowchart TD
  V1["v1: eigenvector"] --> V2["v2: generalized vector"]
  V2 --> V3["v3: generalized vector"]
  V3 --> VS["vs: chain head"]
  VS --> P["Put columns into P in chain order"]
  P --> J["Block Js(lambda) appears in J"]
        

\[ \mathbf{A}\mathbf{P}_\lambda =\mathbf{P}_\lambda\mathbf{J}_s(\lambda),\qquad \mathbf{P}_\lambda=[\mathbf{v}_1\ \mathbf{v}_2\ \cdots\ \mathbf{v}_s]. \]

\[ \mathbf{A}\mathbf{P}=\mathbf{P}\mathbf{J} \quad\Longrightarrow\quad \mathbf{P}^{-1}\mathbf{A}\mathbf{P}=\mathbf{J}. \]

5. Determining Block Sizes from Nullity Growth

Let \( n_k(\lambda)=\dim\ker((\mathbf{A}-\lambda\mathbf{I})^k) \) and define \( n_0(\lambda)=0 \). For a fixed eigenvalue, the number of Jordan blocks of size at least \( k \) is

\[ b_k(\lambda)=n_k(\lambda)-n_{k-1}(\lambda). \]

Therefore, the number of Jordan blocks of size exactly \( k \) is

\[ c_k(\lambda)=b_k(\lambda)-b_{k+1}(\lambda). \]

Proof. A single Jordan block of size \( s \) contributes \( \min(k,s) \) dimensions to \( \ker((\mathbf{J}_s(\lambda)-\lambda\mathbf{I})^k) \). Hence increasing \( k-1 \) to \( k \) adds one new null-vector dimension exactly for those blocks whose size is at least \( k \). This proves the formula for \( b_k \). Subtracting \( b_{k+1} \) from \( b_k \) leaves only blocks whose size is exactly \( k \).

\[ q_\lambda= \min\left\{k:\ker((\mathbf{A}-\lambda\mathbf{I})^k) =\ker((\mathbf{A}-\lambda\mathbf{I})^{k+1})\right\}. \]

\[ m_\mathbf{A}(s)=\prod_\lambda (s-\lambda)^{q_\lambda}. \]

6. Worked Example

Consider the matrix

\[ \mathbf{A}= \begin{bmatrix} 2&1&0&0&0&0\\ 0&2&1&-1&1&-1\\ 0&0&2&0&0&0\\ 0&0&0&2&-3&4\\ 0&0&0&0&-1&1\\ 0&0&0&0&0&-1 \end{bmatrix}. \]

\[ \chi_\mathbf{A}(s)=(s-2)^4(s+1)^2. \]

For \( \lambda=2 \), the nullity sequence is

\[ n_0=0,\qquad n_1=2,\qquad n_2=3,\qquad n_3=4. \]

Thus \( b_1=2 \), \( b_2=1 \), \( b_3=1 \), and \( b_4=0 \). Hence there is one block of size 3 and one block of size 1:

\[ \lambda=2:\qquad \mathbf{J}_3(2)\oplus\mathbf{J}_1(2). \]

For \( \lambda=-1 \), the nullity sequence is

\[ n_0=0,\qquad n_1=1,\qquad n_2=2. \]

\[ \lambda=-1:\qquad \mathbf{J}_2(-1). \]

\[ \mathbf{J}= \operatorname{diag}\left(\mathbf{J}_3(2), \mathbf{J}_1(2),\mathbf{J}_2(-1)\right),\qquad m_\mathbf{A}(s)=(s-2)^3(s+1)^2. \]

7. Control Interpretation

In state-space control, Jordan form is mostly an analytical tool. It clarifies how internal modes combine, how repeated poles can generate nontrivial transient terms, and why geometric multiplicity matters. If \( \mathbf{A} \) is transformed to \( \mathbf{J} \), the input and output matrices become

\[ \mathbf{B}_J=\mathbf{P}^{-1}\mathbf{B},\qquad \mathbf{C}_J=\mathbf{C}\mathbf{P}. \]

\[ \dot{\mathbf{z} }=\mathbf{J}\mathbf{z}+ \mathbf{B}_J\mathbf{u},\qquad \mathbf{y}=\mathbf{C}_J\mathbf{z}+ \mathbf{D}\mathbf{u}. \]

\[ \mathbf{G}(s)= \mathbf{C}(s\mathbf{I}-\mathbf{A})^{-1}\mathbf{B}+\mathbf{D} = \mathbf{C}_J(s\mathbf{I}-\mathbf{J})^{-1}\mathbf{B}_J+\mathbf{D}. \]

8. Numerical Caution

Jordan form is extremely sensitive to perturbations. A matrix with a repeated eigenvalue may split into distinct eigenvalues under arbitrarily small numerical perturbations. Therefore, practical floating-point control software usually prefers Schur decompositions, real Schur forms, or directly computed state-space transformations.

\[ \mathbf{A}=\mathbf{Q}\mathbf{T}\mathbf{Q}^*,\qquad \mathbf{T}\text{ upper triangular, numerically stable.} \]

9. Python Implementation — Exact Symbolic Jordan Diagnostics

Chapter18_Lesson2.py uses sympy for exact algebra and computes block sizes from nullity growth.

# Chapter18_Lesson2.py
"""
Conceptual construction of Jordan canonical form for a controlled example.

Libraries:
    sympy: exact symbolic linear algebra
    numpy/scipy/control can be used later for numerical state-space work,
    but exact Jordan construction should not be treated as a floating-point
    design algorithm.
"""

import sympy as sp


def algebraic_multiplicity(A: sp.Matrix, lam: sp.Expr) -> int:
    """Return the algebraic multiplicity of lam in det(sI - A)."""
    s = sp.symbols("s")
    factors = sp.factor_list(A.charpoly(s).as_expr())[1]
    for factor, exponent in factors:
        if sp.simplify(factor.subs(s, lam)) == 0:
            return int(exponent)
    return 0


def nullity(M: sp.Matrix) -> int:
    """Compute dim ker(M) exactly."""
    return M.cols - M.rank()


def nullity_sequence(A: sp.Matrix, lam: sp.Expr, max_power: int) -> list[int]:
    """Return [n_0, n_1, ..., n_max_power], where n_k = dim ker((A-lam I)^k)."""
    N = A - lam * sp.eye(A.rows)
    values = [0]
    for k in range(1, max_power + 1):
        values.append(nullity(N**k))
    return values


def block_sizes_from_nullities(nullities: list[int], algebraic_mult: int) -> list[int]:
    """
    If n_k = dim ker((A-lam I)^k), then
        b_k = n_k - n_{k-1}
    is the number of Jordan blocks of size at least k.
    The number of blocks of size exactly k is b_k - b_{k+1}.
    """
    trimmed = [nullities[0]]
    for value in nullities[1:]:
        trimmed.append(value)
        if value == algebraic_mult:
            break

    b_at_least = [trimmed[k] - trimmed[k - 1] for k in range(1, len(trimmed))]
    b_at_least.append(0)

    sizes = []
    for k in range(1, len(b_at_least)):
        exact_count = b_at_least[k - 1] - b_at_least[k]
        sizes.extend([k] * exact_count)

    return sorted(sizes, reverse=True)


def jordan_block(lam: sp.Expr, size: int) -> sp.Matrix:
    """Create one Jordan block J_size(lam)."""
    B = lam * sp.eye(size)
    for i in range(size - 1):
        B[i, i + 1] = 1
    return B


def build_example() -> sp.Matrix:
    """Create A = P J P^{-1}, so the expected Jordan structure is known."""
    J = sp.diag(jordan_block(2, 3), jordan_block(2, 1), jordan_block(-1, 2))
    P = sp.eye(6)
    for i in range(5):
        P[i, i + 1] = 1
    return sp.simplify(P * J * P.inv())


def main() -> None:
    A = build_example()
    print("A =")
    sp.print_latex(A)
    print(A)

    s = sp.symbols("s")
    print("\nCharacteristic polynomial:")
    print(sp.factor(A.charpoly(s).as_expr()))

    print("\nMinimal-polynomial prediction from block sizes:")
    print("(s - 2)^3 (s + 1)^2")

    for lam in [sp.Integer(2), sp.Integer(-1)]:
        alg = algebraic_multiplicity(A, lam)
        nullities = nullity_sequence(A, lam, A.rows)
        sizes = block_sizes_from_nullities(nullities, alg)

        print(f"\nlambda = {lam}")
        print("algebraic multiplicity =", alg)
        print("nullities n_k =", nullities[: alg + 1])
        print("Jordan block sizes =", sizes)

    Pj, Jj = A.jordan_form()
    print("\nSymPy Jordan form J:")
    print(Jj)
    print("\nVerification P^{-1} A P - J =")
    print(sp.simplify(Pj.inv() * A * Pj - Jj))


if __name__ == "__main__":
    main()

10. C++ Implementation — From-Scratch Nullity Diagnostics

Chapter18_Lesson2.cpp computes the nullity sequence for a known Jordan-structured example using Gaussian elimination.

// Chapter18_Lesson2.cpp
// From-scratch rank/nullity diagnostics for Jordan block sizes.
// This is educational code, not a numerically robust Jordan-form algorithm.

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

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

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

Matrix subtractLambdaI(const Matrix& A, double lambda) {
    Matrix B = A;
    for (int i = 0; i < static_cast<int>(A.size()); ++i) {
        B[i][i] -= lambda;
    }
    return B;
}

Matrix multiply(const Matrix& A, const Matrix& B) {
    int n = static_cast<int>(A.size());
    int m = static_cast<int>(B[0].size());
    int p = static_cast<int>(B.size());
    Matrix C(n, std::vector<double>(m, 0.0));

    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;
}

int rankGaussian(Matrix A, double tol = 1e-10) {
    int rows = static_cast<int>(A.size());
    int cols = static_cast<int>(A[0].size());
    int 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(A[i][c]) > std::fabs(A[pivot][c])) {
                pivot = i;
            }
        }
        if (std::fabs(A[pivot][c]) <= tol) continue;

        std::swap(A[r], A[pivot]);
        double div = A[r][c];
        for (int j = c; j < cols; ++j) A[r][j] /= div;

        for (int i = 0; i < rows; ++i) {
            if (i == r) continue;
            double factor = A[i][c];
            for (int j = c; j < cols; ++j) {
                A[i][j] -= factor * A[r][j];
            }
        }
        ++r;
    }
    return r;
}

std::vector<int> nullitySequence(const Matrix& A, double lambda, int maxPower) {
    int n = static_cast<int>(A.size());
    Matrix N = subtractLambdaI(A, lambda);
    Matrix Nk = identity(n);

    std::vector<int> nullities;
    nullities.push_back(0);

    for (int k = 1; k <= maxPower; ++k) {
        Nk = multiply(Nk, N);
        int rank = rankGaussian(Nk);
        nullities.push_back(n - rank);
    }
    return nullities;
}

std::vector<int> blockSizesFromNullities(const std::vector<int>& nullities,
                                         int algebraicMultiplicity) {
    std::vector<int> trimmed;
    trimmed.push_back(nullities[0]);

    for (size_t i = 1; i < nullities.size(); ++i) {
        trimmed.push_back(nullities[i]);
        if (nullities[i] == algebraicMultiplicity) break;
    }

    std::vector<int> blocksAtLeast;
    for (size_t k = 1; k < trimmed.size(); ++k) {
        blocksAtLeast.push_back(trimmed[k] - trimmed[k - 1]);
    }
    blocksAtLeast.push_back(0);

    std::vector<int> sizes;
    for (size_t k = 1; k < blocksAtLeast.size(); ++k) {
        int exactCount = blocksAtLeast[k - 1] - blocksAtLeast[k];
        for (int c = 0; c < exactCount; ++c) {
            sizes.push_back(static_cast<int>(k));
        }
    }

    std::sort(sizes.rbegin(), sizes.rend());
    return sizes;
}

void printVector(const std::vector<int>& v) {
    std::cout << "[";
    for (size_t i = 0; i < v.size(); ++i) {
        std::cout << v[i];
        if (i + 1 < v.size()) std::cout << ", ";
    }
    std::cout << "]";
}

int main() {
    Matrix J = {
        { 2, 1, 0, 0,  0, 0},
        { 0, 2, 1, 0,  0, 0},
        { 0, 0, 2, 0,  0, 0},
        { 0, 0, 0, 2,  0, 0},
        { 0, 0, 0, 0, -1, 1},
        { 0, 0, 0, 0,  0,-1}
    };

    std::vector<double> eigenvalues = {2.0, -1.0};
    std::vector<int> algebraicMultiplicities = {4, 2};

    for (size_t idx = 0; idx < eigenvalues.size(); ++idx) {
        double lambda = eigenvalues[idx];
        int alg = algebraicMultiplicities[idx];

        std::vector<int> nullities = nullitySequence(J, lambda, 6);
        std::vector<int> sizes = blockSizesFromNullities(nullities, alg);

        std::cout << "lambda = " << lambda << "\n";
        std::cout << "nullities n_k = ";
        printVector(nullities);
        std::cout << "\nJordan block sizes = ";
        printVector(sizes);
        std::cout << "\n\n";
    }

    std::cout << "Expected: lambda 2 has block sizes [3, 1]; lambda -1 has [2].\n";
    return 0;
}

11. Java Implementation — Rank and Nullity from Scratch

Chapter18_Lesson2.java mirrors the C++ implementation.

// Chapter18_Lesson2.java
// From-scratch rank/nullity diagnostics for Jordan block sizes.
// This is educational code, not a numerically robust Jordan-form algorithm.

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

public class Chapter18_Lesson2 {
    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[][] subtractLambdaI(double[][] A, double lambda) {
        int n = A.length;
        double[][] B = new double[n][n];
        for (int i = 0; i < n; i++) {
            System.arraycopy(A[i], 0, B[i], 0, n);
            B[i][i] -= lambda;
        }
        return B;
    }

    static double[][] multiply(double[][] A, double[][] B) {
        int n = A.length;
        int m = B[0].length;
        int p = B.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 int rankGaussian(double[][] input, double tolerance) {
        int rows = input.length;
        int cols = input[0].length;
        double[][] A = new double[rows][cols];

        for (int i = 0; i < rows; i++) {
            System.arraycopy(input[i], 0, A[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(A[i][c]) > Math.abs(A[pivot][c])) {
                    pivot = i;
                }
            }

            if (Math.abs(A[pivot][c]) <= tolerance) continue;

            double[] temp = A[r];
            A[r] = A[pivot];
            A[pivot] = temp;

            double divisor = A[r][c];
            for (int j = c; j < cols; j++) A[r][j] /= divisor;

            for (int i = 0; i < rows; i++) {
                if (i == r) continue;
                double factor = A[i][c];
                for (int j = c; j < cols; j++) {
                    A[i][j] -= factor * A[r][j];
                }
            }
            r++;
        }
        return r;
    }

    static List<Integer> nullitySequence(double[][] A, double lambda, int maxPower) {
        int n = A.length;
        double[][] N = subtractLambdaI(A, lambda);
        double[][] Nk = identity(n);

        List<Integer> nullities = new ArrayList<>();
        nullities.add(0);

        for (int k = 1; k <= maxPower; k++) {
            Nk = multiply(Nk, N);
            int rank = rankGaussian(Nk, 1e-10);
            nullities.add(n - rank);
        }
        return nullities;
    }

    static List<Integer> blockSizesFromNullities(List<Integer> nullities, int algebraicMultiplicity) {
        List<Integer> trimmed = new ArrayList<>();
        trimmed.add(nullities.get(0));

        for (int i = 1; i < nullities.size(); i++) {
            trimmed.add(nullities.get(i));
            if (nullities.get(i) == algebraicMultiplicity) break;
        }

        List<Integer> blocksAtLeast = new ArrayList<>();
        for (int k = 1; k < trimmed.size(); k++) {
            blocksAtLeast.add(trimmed.get(k) - trimmed.get(k - 1));
        }
        blocksAtLeast.add(0);

        List<Integer> sizes = new ArrayList<>();
        for (int k = 1; k < blocksAtLeast.size(); k++) {
            int exactCount = blocksAtLeast.get(k - 1) - blocksAtLeast.get(k);
            for (int c = 0; c < exactCount; c++) {
                sizes.add(k);
            }
        }

        sizes.sort(Collections.reverseOrder());
        return sizes;
    }

    public static void main(String[] args) {
        double[][] J = {
            { 2, 1, 0, 0,  0, 0},
            { 0, 2, 1, 0,  0, 0},
            { 0, 0, 2, 0,  0, 0},
            { 0, 0, 0, 2,  0, 0},
            { 0, 0, 0, 0, -1, 1},
            { 0, 0, 0, 0,  0,-1}
        };

        double[] eigenvalues = {2.0, -1.0};
        int[] algebraicMultiplicities = {4, 2};

        for (int i = 0; i < eigenvalues.length; i++) {
            List<Integer> nullities = nullitySequence(J, eigenvalues[i], 6);
            List<Integer> sizes = blockSizesFromNullities(nullities, algebraicMultiplicities[i]);

            System.out.println("lambda = " + eigenvalues[i]);
            System.out.println("nullities n_k = " + nullities);
            System.out.println("Jordan block sizes = " + sizes);
            System.out.println();
        }

        System.out.println("Expected: lambda 2 has block sizes [3, 1]; lambda -1 has [2].");
    }
}

12. MATLAB/Simulink Implementation

Chapter18_Lesson2.m uses the Symbolic Math Toolbox for jordan, the Control System Toolbox for ss, and optionally creates a Simulink State-Space block.

% Chapter18_Lesson2.m
% Exact Jordan canonical form construction plus a state-space/Simulink check.
% Required for jordan: Symbolic Math Toolbox.
% Required for ss/sim: Control System Toolbox and Simulink.

clear; clc;

A = sym([ ...
     2  1  0  0  0  0; ...
     0  2  1 -1  1 -1; ...
     0  0  2  0  0  0; ...
     0  0  0  2 -3  4; ...
     0  0  0  0 -1  1; ...
     0  0  0  0  0 -1]);

disp('A ='); disp(A);

[V,J] = jordan(A);
disp('Jordan form J ='); disp(J);
disp('Verification inv(V)*A*V - J ='); disp(simplify(inv(V)*A*V - J));

lambdaValues = [sym(2), sym(-1)];
algebraicMultiplicities = [4, 2];

for idx = 1:numel(lambdaValues)
    lambda = lambdaValues(idx);
    alg = algebraicMultiplicities(idx);
    N = A - lambda*eye(size(A));

    nullities = zeros(1, size(A,1));
    for k = 1:size(A,1)
        nullities(k) = size(A,1) - rank(N^k);
    end

    fprintf('\nlambda = %s\n', char(lambda));
    fprintf('algebraic multiplicity = %d\n', alg);
    fprintf('nullities n_k = ');
    disp(nullities);

    nWithZero = [0, nullities];
    bAtLeast = diff(nWithZero);
    bAtLeast = bAtLeast(1:find(nullities == alg, 1, 'first'));
    bAtLeast = [bAtLeast, 0];

    blockSizes = [];
    for k = 1:(numel(bAtLeast)-1)
        exactCount = bAtLeast(k) - bAtLeast(k+1);
        blockSizes = [blockSizes, k*ones(1, exactCount)]; %#ok<AGROW>
    end
    blockSizes = sort(blockSizes, 'descend');
    fprintf('Jordan block sizes = ');
    disp(blockSizes);
end

try
    sys = ss(double(A), zeros(6,1), eye(6), zeros(6,1));
    disp('Continuous-time state-space model sys = ss(A,0,I,0):');
    disp(sys);
catch ME
    warning('Control System Toolbox state-space object was not created: %s', ME.message);
end

try
    model = 'Chapter18_Lesson2_SimulinkModel';
    if bdIsLoaded(model)
        close_system(model, 0);
    end
    new_system(model);
    add_block('simulink/Continuous/State-Space', [model '/State-Space']);
    set_param([model '/State-Space'], ...
        'A', mat2str(double(A)), ...
        'B', mat2str(zeros(6,1)), ...
        'C', mat2str(eye(6)), ...
        'D', mat2str(zeros(6,1)));
    save_system(model);
    fprintf('Saved optional Simulink model: %s.slx\n', model);
catch ME
    warning('Simulink model was not created: %s', ME.message);
end

13. Wolfram Mathematica Implementation

Chapter18_Lesson2.nb uses exact symbolic capabilities through JordanDecomposition and null-space computations.

Notebook[{
Cell["Chapter18_Lesson2.nb", "Title"],
Cell["Construction of Jordan canonical form using exact symbolic linear algebra.", "Text"],
Cell[BoxData[
"ClearAll[\"Global`*\"]\n\nA = { {2, 1, 0, 0, 0, 0},\n     {0, 2, 1, -1, 1, -1},\n     {0, 0, 2, 0, 0, 0},\n     {0, 0, 0, 2, -3, 4},\n     {0, 0, 0, 0, -1, 1},\n     {0, 0, 0, 0, 0, -1} };\n\n{P, J} = JordanDecomposition[A];\nJ // MatrixForm\nSimplify[Inverse[P].A.P - J] // MatrixForm"
], "Input"],
Cell[BoxData[
"nullity[m_] := Length[NullSpace[m]];\nnullitySequence[A_, lambda_, maxPower_] :=\n  Prepend[Table[nullity[MatrixPower[A - lambda IdentityMatrix[Length[A]], k]], {k, 1, maxPower}], 0];\n\nblockSizesFromNullities[nullities_, algebraicMultiplicity_] := Module[\n  {trimmed, bAtLeast, exactCounts, sizes},\n  trimmed = Take[nullities, FirstPosition[nullities, algebraicMultiplicity][[1]]];\n  bAtLeast = Differences[trimmed];\n  bAtLeast = Append[bAtLeast, 0];\n  exactCounts = Table[bAtLeast[[k]] - bAtLeast[[k + 1]], {k, 1, Length[bAtLeast] - 1}];\n  sizes = Flatten[Table[ConstantArray[k, exactCounts[[k]]], {k, 1, Length[exactCounts]}]];\n  Reverse[Sort[sizes]]\n];\n\nDo[\n  ns = nullitySequence[A, lambda, Length[A]];\n  Print[\"lambda = \", lambda];\n  Print[\"nullities n_k = \", ns];\n  Print[\"Jordan block sizes = \", blockSizesFromNullities[ns, alg]];\n  , { {lambda, alg}, { {2, 4}, {-1, 2} } }\n]"
], "Input"]
}]

14. Problems and Solutions

Problem 1 (Block Sizes from Nullities): For an eigenvalue \( \lambda \), suppose \( n_0=0 \), \( n_1=3 \), \( n_2=5 \), \( n_3=6 \), and \( n_4=6 \). Determine the Jordan block sizes for \( \lambda \).

Solution:

\[ b_1=3,\qquad b_2=2,\qquad b_3=1,\qquad b_4=0. \]

Therefore \( c_1=1 \), \( c_2=1 \), and \( c_3=1 \). The block sizes are \( 3,2,1 \).

Problem 2 (Geometric Multiplicity): A matrix has characteristic polynomial \( \chi_\mathbf{A}(s)=(s-4)^5 \) and \( \dim\ker(\mathbf{A}-4\mathbf{I})=2 \). How many Jordan blocks are associated with \( \lambda=4 \)?

Solution:

\[ \#\text{ Jordan blocks for }\lambda =\dim\ker(\mathbf{A}-\lambda\mathbf{I}). \]

Hence there are exactly two Jordan blocks. Their sizes must add to 5, but the given information alone does not determine whether the sizes are \( 4,1 \) or \( 3,2 \). Additional nullity information is required.

Problem 3 (Constructing a Chain): Suppose \( (\mathbf{A}-\lambda\mathbf{I})\mathbf{v}_1=\mathbf{0} \) and \( (\mathbf{A}-\lambda\mathbf{I})\mathbf{v}_2=\mathbf{v}_1 \). Show that \( \mathbf{v}_1,\mathbf{v}_2 \) generate a Jordan block of size 2.

Solution:

\[ \mathbf{A}\mathbf{v}_1=\lambda\mathbf{v}_1,\qquad \mathbf{A}\mathbf{v}_2=\mathbf{v}_1+\lambda\mathbf{v}_2. \]

\[ \begin{bmatrix}\lambda&1\\0&\lambda\end{bmatrix} =\mathbf{J}_2(\lambda). \]

Problem 4 (Minimal Polynomial): A Jordan form has blocks \( \mathbf{J}_2(1)\oplus\mathbf{J}_4(1) \oplus\mathbf{J}_3(-2) \). Find the characteristic and minimal polynomials.

Solution:

\[ \chi_\mathbf{A}(s)=(s-1)^6(s+2)^3,\qquad m_\mathbf{A}(s)=(s-1)^4(s+2)^3. \]

Problem 5 (Diagonalizability Test): Show that a matrix is diagonalizable over \( \mathbb{C} \) if and only if every Jordan block has size 1.

Solution: If every Jordan block has size 1, then \( \mathbf{J} \) is diagonal, so \( \mathbf{A}=\mathbf{P}\mathbf{J}\mathbf{P}^{-1} \) is diagonalizable. Conversely, if \( \mathbf{A} \) is diagonalizable, there is a basis of eigenvectors. In such a basis, no generalized vector of height greater than 1 is needed. Hence all Jordan blocks must have size 1.

15. Summary

Jordan canonical form is constructed by identifying generalized eigenspaces, measuring nullity growth, and assembling Jordan chains into a similarity matrix. The nullity sequence determines block sizes exactly, while the chain equations determine the columns of \( \mathbf{P} \). For modern control, Jordan form explains repeated-mode structure and prepares the ground for the polynomial-exponential terms studied in the next lesson.

16. References

  1. Abo, H., Eklund, D., Kahle, T., & Peterson, C. (2016). Eigenschemes and the Jordan canonical form. Linear Algebra and Its Applications, 496, 121–151.
  2. Corless, R.M., Moreno Maza, M., & Thornton, S.E. (2017). Jordan canonical form with parameters from Frobenius form with parameters. Lecture Notes in Computer Science, 10693, 179–194.
  3. Mehl, C., Mehrmann, V., Ran, A.C.M., & Rodman, L. (2013). Jordan forms of real and complex matrices under rank one perturbations. Operators and Matrices, 7(2), 381–398.
  4. Kalogeropoulos, G. (2004). On the computation of the Jordan canonical form of regular matrix polynomials. Linear Algebra and Its Applications, 376, 185–210.
  5. Wörz, S. (2022). An efficient Jordan basis algorithm. Journal of Computational and Applied Mathematics, 401, 113782.
  6. Higham, N.J., & Lin, L. (2013). Matrix functions: A short course. Lecture Notes in Mathematics, 2046, 1–27.