Chapter 9: Stability of Linear Systems in State Space

Lesson 5: Examples of Stable, Marginal, and Unstable State Matrices

This lesson consolidates the state-space stability criteria from the previous lessons by studying explicit state matrices. We classify continuous-time autonomous systems \( \dot{\mathbf{x} }=\mathbf{A}\mathbf{x} \) into asymptotically stable, marginally stable, and unstable cases, with emphasis on eigenvalues, Jordan blocks, semisimplicity, and representative trajectories.

1. Why Examples Matter in State-Matrix Stability

For the autonomous continuous-time LTI system \( \dot{\mathbf{x} }(t)=\mathbf{A}\mathbf{x}(t) \), the solution is governed by the matrix exponential:

\[ \mathbf{x}(t)=e^{\mathbf{A}t}\mathbf{x}(0). \]

Thus, stability is not a property of a single trajectory; it is a property of the family of maps \( e^{\mathbf{A}t} \). The examples below illustrate three key phenomena:

  1. negative real eigenvalue parts force all modes to decay;
  2. imaginary-axis eigenvalues may be harmless only when they are semisimple;
  3. positive real parts or defective imaginary-axis blocks create unbounded solutions.
flowchart TD
  A["Given state matrix A"] --> B["Compute eigenvalues of A"]
  B --> C["All real parts negative?"]
  C -->|yes| S["Asymptotically stable"]
  C -->|no| D["Any real part positive?"]
  D -->|yes| U["Unstable"]
  D -->|no| E["Some eigenvalues on imaginary axis"]
  E --> F["Check Jordan blocks for boundary eigenvalues"]
  F -->|all semisimple| M["Marginally stable"]
  F -->|any defective block| U
        

2. General Classification Theorem for Continuous-Time LTI Systems

Let \( \lambda_i(\mathbf{A}) \) denote the eigenvalues of \( \mathbf{A} \), and define the spectral abscissa

\[ \alpha(\mathbf{A})=\max_i \operatorname{Re}\left(\lambda_i(\mathbf{A})\right). \]

The classification for \( \dot{\mathbf{x} }=\mathbf{A}\mathbf{x} \) is:

\[ \begin{array}{lll} \alpha(\mathbf{A})<0 & \text{asymptotically stable},\\[2mm] \alpha(\mathbf{A})\leq 0 \text{ and all eigenvalues with } \operatorname{Re}(\lambda)=0 \text{ are semisimple} & \text{marginally stable},\\[2mm] \alpha(\mathbf{A})>0 \text{ or a defective boundary eigenvalue exists} & \text{unstable}. \end{array} \]

A boundary eigenvalue is an eigenvalue with \( \operatorname{Re}(\lambda)=0 \). It includes zero eigenvalues and pure imaginary eigenvalues.

2.1 Proof Sketch via Jordan Blocks

Suppose \( \mathbf{A} \) is transformed into Jordan form:

\[ \mathbf{A}=\mathbf{T}\mathbf{J}\mathbf{T}^{-1}, \qquad e^{\mathbf{A}t}=\mathbf{T}e^{\mathbf{J}t}\mathbf{T}^{-1}. \]

For a Jordan block \( \mathbf{J}_m(\lambda)=\lambda\mathbf{I}+\mathbf{N} \), where \( \mathbf{N}^m=\mathbf{0} \), the exponential is

\[ e^{\mathbf{J}_m(\lambda)t} = e^{\lambda t} \left( \mathbf{I}+t\mathbf{N}+\frac{t^2}{2!}\mathbf{N}^2+ \cdots+ \frac{t^{m-1} }{(m-1)!}\mathbf{N}^{m-1} \right). \]

If \( \operatorname{Re}(\lambda)<0 \), the exponential decay dominates every polynomial factor. If \( \operatorname{Re}(\lambda)>0 \), the mode grows exponentially. If \( \operatorname{Re}(\lambda)=0 \), boundedness requires the nilpotent part \( \mathbf{N} \) to be zero for that block, which means the boundary eigenvalue is semisimple.

3. Examples of Asymptotically Stable State Matrices

3.1 Stable Diagonal Matrix

Consider \( \mathbf{A}_1=\begin{bmatrix}-2&0\\0&-5\end{bmatrix} \). The solution is

\[ e^{\mathbf{A}_1t} = \begin{bmatrix} e^{-2t} & 0\\ 0 & e^{-5t} \end{bmatrix}, \qquad \mathbf{x}(t)= \begin{bmatrix} e^{-2t}x_1(0)\\ e^{-5t}x_2(0) \end{bmatrix}. \]

Both modes decay to zero. Therefore, the origin is asymptotically stable.

3.2 Stable Damped Oscillator Matrix

The matrix \( \mathbf{A}_2=\begin{bmatrix}0&1\\-4&-2\end{bmatrix} \) has characteristic polynomial

\[ p(s)=\det(s\mathbf{I}-\mathbf{A}_2) = s^2+2s+4. \]

Hence

\[ \lambda_{1,2}=-1\pm j\sqrt{3}, \qquad \operatorname{Re}(\lambda_{1,2})=-1<0. \]

This is a spiral sink: trajectories oscillate while their envelope decays like \( e^{-t} \).

3.3 Stable but Nonnormal Matrix

The matrix \( \mathbf{A}_3=\begin{bmatrix}-1&50\\0&-2\end{bmatrix} \) has eigenvalues \( -1 \) and \( -2 \), so it is asymptotically stable. However, it is nonnormal:

\[ \mathbf{A}_3\mathbf{A}_3^T \neq \mathbf{A}_3^T\mathbf{A}_3. \]

Therefore, some initial conditions can exhibit large transient growth before eventual decay. This does not contradict asymptotic stability: stability concerns the behavior as \( t \) becomes large, not monotonic decay at every instant.

4. Examples of Marginally Stable State Matrices

4.1 Undamped Oscillator

Consider \( \mathbf{A}_4=\begin{bmatrix}0&1\\-1&0\end{bmatrix} \). The eigenvalues are \( \lambda_{1,2}=\pm j \). Since the eigenvalues are distinct, they are semisimple. The matrix exponential is a rotation:

\[ e^{\mathbf{A}_4t} = \begin{bmatrix} \cos t & \sin t\\ -\sin t & \cos t \end{bmatrix}. \]

Hence \( \|\mathbf{x}(t)\|_2=\|\mathbf{x}(0)\|_2 \). The origin is stable but not asymptotically stable because trajectories do not converge to zero unless \( \mathbf{x}(0)=\mathbf{0} \).

4.2 Semisimple Zero Eigenvalue

Consider \( \mathbf{A}_5=\begin{bmatrix}0&0\\0&-2\end{bmatrix} \). The eigenvalues are \( 0 \) and \( -2 \). The zero eigenvalue is semisimple because its algebraic and geometric multiplicities are both one. The solution is

\[ \mathbf{x}(t)= \begin{bmatrix} x_1(0)\\ e^{-2t}x_2(0) \end{bmatrix}. \]

The state remains bounded, but it does not generally converge to the origin because \( x_1(t)=x_1(0) \). Therefore, the system is marginally stable.

flowchart LR
  S["Stable: \nall modes decay to zero"] --> A["Example: \ndiagonal negative matrix"]
  M["Marginal: \nmodes stay bounded but may not decay"] --> B["Example: \noscillator or semisimple zero mode"]
  U["Unstable: \nat least one mode becomes unbounded"] --> C["Example: \npositive eigenvalue or defective boundary block"]
        

5. Examples of Unstable State Matrices

5.1 Positive Real Eigenvalue

Let \( \mathbf{A}_6=\begin{bmatrix}1&0\\0&-2\end{bmatrix} \). Then

\[ \mathbf{x}(t)= \begin{bmatrix} e^t x_1(0)\\ e^{-2t}x_2(0) \end{bmatrix}. \]

The first component grows exponentially whenever \( x_1(0)\neq 0 \). Therefore, the origin is unstable.

5.2 Defective Zero Eigenvalue

Consider the nilpotent matrix \( \mathbf{A}_7=\begin{bmatrix}0&1\\0&0\end{bmatrix} \). Its only eigenvalue is \( 0 \), but it has one eigenvector, not two. It is a Jordan block:

\[ \mathbf{A}_7^2=\mathbf{0}, \qquad e^{\mathbf{A}_7t} = \mathbf{I}+t\mathbf{A}_7 = \begin{bmatrix} 1&t\\ 0&1 \end{bmatrix}. \]

Thus \( x_1(t)=x_1(0)+t x_2(0) \). If \( x_2(0)\neq 0 \), the state grows linearly and the system is unstable. This is the most important warning about marginal eigenvalues: the eigenvalue \( 0 \) is not enough for marginal stability; semisimplicity is required.

5.3 Defective Imaginary-Axis Block

A real defective imaginary-axis example can be represented by

\[ \mathbf{A}_8= \begin{bmatrix} 0&1&1&0\\ -1&0&0&1\\ 0&0&0&1\\ 0&0&-1&0 \end{bmatrix} = \begin{bmatrix} \mathbf{R}&\mathbf{I}\\ \mathbf{0}&\mathbf{R} \end{bmatrix}, \qquad \mathbf{R}= \begin{bmatrix} 0&1\\ -1&0 \end{bmatrix}. \]

The eigenvalues are \( \pm j \), but the block is defective. The exponential contains terms proportional to \( t\cos t \) and \( t\sin t \), so the response is unbounded even though all eigenvalues lie on the imaginary axis.

6. Comparison of Representative Matrices

Matrix Eigenvalues Jordan Boundary Condition Classification Typical Behavior
\( \begin{bmatrix}-2&0\\0&-5\end{bmatrix} \) \( -2,-5 \) No boundary eigenvalue Asymptotically stable monotone modal decay
\( \begin{bmatrix}0&1\\-4&-2\end{bmatrix} \) \( -1\pm j\sqrt{3} \) No boundary eigenvalue Asymptotically stable decaying oscillation
\( \begin{bmatrix}0&1\\-1&0\end{bmatrix} \) \( \pm j \) semisimple Marginally stable bounded periodic motion
\( \begin{bmatrix}0&0\\0&-2\end{bmatrix} \) \( 0,-2 \) zero eigenvalue semisimple Marginally stable one constant mode, one decaying mode
\( \begin{bmatrix}1&0\\0&-2\end{bmatrix} \) \( 1,-2 \) not relevant Unstable exponential growth in one mode
\( \begin{bmatrix}0&1\\0&0\end{bmatrix} \) \( 0,0 \) defective zero block Unstable linear growth

7. Python Implementation

Chapter9_Lesson5.py uses numpy and scipy.linalg to compute eigenvalues, detect boundary eigenvalue defects, and simulate exact trajectories using expm.

"""
Chapter9_Lesson5.py

Modern Control — Chapter 9, Lesson 5
Examples of Stable, Marginal, and Unstable State Matrices

This script classifies continuous-time LTI state matrices A for x_dot = A x
using eigenvalue/Jordan-block logic, then simulates representative systems.

Libraries:
    numpy
    scipy
    matplotlib

Install:
    pip install numpy scipy matplotlib
"""

import numpy as np
from scipy.linalg import expm, eigvals
import matplotlib.pyplot as plt


def classify_matrix(A: np.ndarray, tol: float = 1e-9) -> str:
    """
    Classify x_dot = A x for continuous-time stability.

    This routine uses eigenvalue real parts for the main classification and
    adds a conservative algebraic/geometric multiplicity check for eigenvalues
    on the imaginary axis. For teaching examples, this is sufficient.

    Returns:
        "asymptotically stable", "marginally stable", or "unstable"
    """
    lam = eigvals(A)
    spectral_abscissa = np.max(np.real(lam))

    if spectral_abscissa < -tol:
        return "asymptotically stable"

    if spectral_abscissa > tol:
        return "unstable"

    # Boundary case: every eigenvalue has nonpositive real part.
    # Check whether eigenvalues on the imaginary axis are semisimple.
    n = A.shape[0]
    used = np.zeros(len(lam), dtype=bool)

    for i, value in enumerate(lam):
        if used[i]:
            continue

        close = np.abs(lam - value) < 1e-7
        used[close] = True
        algebraic_mult = int(np.sum(close))

        if abs(np.real(value)) <= 1e-7:
            rank = np.linalg.matrix_rank(A - value * np.eye(n), tol=1e-7)
            geometric_mult = n - rank
            if geometric_mult < algebraic_mult:
                return "unstable"

    return "marginally stable"


def simulate(A: np.ndarray, x0: np.ndarray, t_grid: np.ndarray) -> np.ndarray:
    """Exact simulation x(t) = exp(A t) x0 over a time grid."""
    X = np.zeros((len(t_grid), len(x0)))
    for k, t in enumerate(t_grid):
        X[k, :] = np.real(expm(A * t) @ x0)
    return X


def print_report(name: str, A: np.ndarray) -> None:
    """Print eigenvalues and stability classification."""
    lam = eigvals(A)
    print(f"\n{name}")
    print("-" * len(name))
    print("A =")
    print(A)
    print("eigenvalues =", lam)
    print("classification =", classify_matrix(A))


def main() -> None:
    examples = {
        "Stable diagonal": np.array([[-2.0, 0.0], [0.0, -5.0]]),
        "Stable damped oscillator": np.array([[0.0, 1.0], [-4.0, -2.0]]),
        "Stable but nonnormal": np.array([[-1.0, 50.0], [0.0, -2.0]]),
        "Marginal oscillator": np.array([[0.0, 1.0], [-1.0, 0.0]]),
        "Marginal with semisimple zero": np.array([[0.0, 0.0], [0.0, -2.0]]),
        "Unstable positive eigenvalue": np.array([[1.0, 0.0], [0.0, -2.0]]),
        "Unstable defective zero": np.array([[0.0, 1.0], [0.0, 0.0]]),
    }

    for name, A in examples.items():
        print_report(name, A)

    # Compare representative trajectories.
    t = np.linspace(0.0, 10.0, 500)
    x0 = np.array([1.0, 1.0])

    selected = [
        ("Stable damped oscillator", examples["Stable damped oscillator"]),
        ("Marginal oscillator", examples["Marginal oscillator"]),
        ("Unstable defective zero", examples["Unstable defective zero"]),
    ]

    plt.figure(figsize=(9, 5))
    for label, A in selected:
        X = simulate(A, x0, t)
        norm_x = np.linalg.norm(X, axis=1)
        plt.plot(t, norm_x, label=label)

    plt.xlabel("time t")
    plt.ylabel("state norm ||x(t)||")
    plt.title("Stable, Marginal, and Unstable State-Matrix Examples")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.show()


if __name__ == "__main__":
    main()

8. C++ Implementation

Chapter9_Lesson5.cpp implements closed-form eigenvalue classification for real \( 2\times 2 \) matrices without external numerical libraries. This is useful for seeing the determinant-trace connection directly.

/*
Chapter9_Lesson5.cpp

Modern Control — Chapter 9, Lesson 5
Examples of Stable, Marginal, and Unstable State Matrices

This compact C++ program classifies 2x2 real continuous-time matrices
A for x_dot = A x using the closed-form eigenvalues of a 2x2 matrix.

Compile:
    g++ -std=c++17 Chapter9_Lesson5.cpp -o Chapter9_Lesson5
Run:
    ./Chapter9_Lesson5
*/

#include <cmath>
#include <complex>
#include <iomanip>
#include <iostream>
#include <string>
#include <utility>
#include <vector>

struct Matrix2 {
    double a11, a12, a21, a22;
};

std::pair<std::complex<double>, std::complex<double>> eigenvalues2x2(const Matrix2& A) {
    double tr = A.a11 + A.a22;
    double det = A.a11 * A.a22 - A.a12 * A.a21;
    std::complex<double> disc = std::complex<double>(tr * tr - 4.0 * det, 0.0);
    std::complex<double> root = std::sqrt(disc);
    return {
        0.5 * (std::complex<double>(tr, 0.0) + root),
        0.5 * (std::complex<double>(tr, 0.0) - root)
    };
}

double determinant(const Matrix2& A) {
    return A.a11 * A.a22 - A.a12 * A.a21;
}

int rank2x2(const Matrix2& A, double tol = 1e-10) {
    if (std::fabs(A.a11) < tol && std::fabs(A.a12) < tol &&
        std::fabs(A.a21) < tol && std::fabs(A.a22) < tol) {
        return 0;
    }
    if (std::fabs(determinant(A)) > tol) {
        return 2;
    }
    return 1;
}

std::string classify2x2(const Matrix2& A, double tol = 1e-9) {
    auto [lambda1, lambda2] = eigenvalues2x2(A);
    double max_real = std::max(lambda1.real(), lambda2.real());

    if (max_real < -tol) {
        return "asymptotically stable";
    }
    if (max_real > tol) {
        return "unstable";
    }

    // Boundary case for 2x2:
    // repeated eigenvalue on the imaginary axis is dangerous if it is defective.
    if (std::abs(lambda1 - lambda2) < 1e-8 && std::fabs(lambda1.real()) < 1e-8) {
        Matrix2 shifted{
            A.a11 - lambda1.real(), A.a12,
            A.a21, A.a22 - lambda1.real()
        };
        int rank = rank2x2(shifted);
        int geometric_multiplicity = 2 - rank;
        if (geometric_multiplicity < 2) {
            return "unstable";
        }
    }

    return "marginally stable";
}

void report(const std::string& name, const Matrix2& A) {
    auto [lambda1, lambda2] = eigenvalues2x2(A);

    std::cout << "\n" << name << "\n";
    std::cout << std::string(name.size(), '-') << "\n";
    std::cout << "A = [[" << A.a11 << ", " << A.a12 << "], ["
              << A.a21 << ", " << A.a22 << "]]\n";
    std::cout << "lambda1 = " << lambda1 << "\n";
    std::cout << "lambda2 = " << lambda2 << "\n";
    std::cout << "classification = " << classify2x2(A) << "\n";
}

int main() {
    std::vector<std::pair<std::string, Matrix2>> examples = {
        {"Stable diagonal", {-2.0, 0.0, 0.0, -5.0} },
        {"Stable damped oscillator", {0.0, 1.0, -4.0, -2.0} },
        {"Stable but nonnormal", {-1.0, 50.0, 0.0, -2.0} },
        {"Marginal oscillator", {0.0, 1.0, -1.0, 0.0} },
        {"Marginal with semisimple zero", {0.0, 0.0, 0.0, -2.0} },
        {"Unstable positive eigenvalue", {1.0, 0.0, 0.0, -2.0} },
        {"Unstable defective zero", {0.0, 1.0, 0.0, 0.0} }
    };

    for (const auto& item : examples) {
        report(item.first, item.second);
    }

    return 0;
}

9. Java Implementation

Chapter9_Lesson5.java mirrors the C++ logic using a simple internal complex-number class. It classifies representative second-order systems using trace, determinant, and boundary-defect logic.

/*
Chapter9_Lesson5.java

Modern Control — Chapter 9, Lesson 5
Examples of Stable, Marginal, and Unstable State Matrices

This Java program classifies 2x2 real continuous-time matrices
A for x_dot = A x using the closed-form eigenvalues of a 2x2 matrix.

Compile:
    javac Chapter9_Lesson5.java
Run:
    java Chapter9_Lesson5
*/

public class Chapter9_Lesson5 {
    static class Matrix2 {
        final double a11, a12, a21, a22;

        Matrix2(double a11, double a12, double a21, double a22) {
            this.a11 = a11;
            this.a12 = a12;
            this.a21 = a21;
            this.a22 = a22;
        }
    }

    static class Complex {
        final double re, im;

        Complex(double re, double im) {
            this.re = re;
            this.im = im;
        }

        Complex plus(Complex other) {
            return new Complex(this.re + other.re, this.im + other.im);
        }

        Complex minus(Complex other) {
            return new Complex(this.re - other.re, this.im - other.im);
        }

        Complex scale(double s) {
            return new Complex(s * this.re, s * this.im);
        }

        double abs() {
            return Math.hypot(re, im);
        }

        static Complex sqrtReal(double x) {
            if (x >= 0.0) {
                return new Complex(Math.sqrt(x), 0.0);
            }
            return new Complex(0.0, Math.sqrt(-x));
        }

        @Override
        public String toString() {
            if (im >= 0.0) {
                return String.format("%.6f + %.6fi", re, im);
            }
            return String.format("%.6f - %.6fi", re, -im);
        }
    }

    static Complex[] eigenvalues2x2(Matrix2 A) {
        double tr = A.a11 + A.a22;
        double det = A.a11 * A.a22 - A.a12 * A.a21;
        double disc = tr * tr - 4.0 * det;
        Complex root = Complex.sqrtReal(disc);
        Complex center = new Complex(tr, 0.0).scale(0.5);
        return new Complex[] {
            center.plus(root.scale(0.5)),
            center.minus(root.scale(0.5))
        };
    }

    static double determinant(Matrix2 A) {
        return A.a11 * A.a22 - A.a12 * A.a21;
    }

    static int rank2x2(Matrix2 A, double tol) {
        if (Math.abs(A.a11) < tol && Math.abs(A.a12) < tol &&
            Math.abs(A.a21) < tol && Math.abs(A.a22) < tol) {
            return 0;
        }
        if (Math.abs(determinant(A)) > tol) {
            return 2;
        }
        return 1;
    }

    static String classify2x2(Matrix2 A) {
        double tol = 1e-9;
        Complex[] lambda = eigenvalues2x2(A);
        double maxReal = Math.max(lambda[0].re, lambda[1].re);

        if (maxReal < -tol) {
            return "asymptotically stable";
        }
        if (maxReal > tol) {
            return "unstable";
        }

        // Boundary case: repeated eigenvalue on the imaginary axis can be defective.
        if (lambda[0].minus(lambda[1]).abs() < 1e-8 && Math.abs(lambda[0].re) < 1e-8) {
            Matrix2 shifted = new Matrix2(
                A.a11 - lambda[0].re, A.a12,
                A.a21, A.a22 - lambda[0].re
            );
            int geometricMultiplicity = 2 - rank2x2(shifted, 1e-10);
            if (geometricMultiplicity < 2) {
                return "unstable";
            }
        }

        return "marginally stable";
    }

    static void report(String name, Matrix2 A) {
        Complex[] lambda = eigenvalues2x2(A);
        System.out.println("\n" + name);
        System.out.println("-".repeat(name.length()));
        System.out.printf("A = [[%.3f, %.3f], [%.3f, %.3f]]%n", A.a11, A.a12, A.a21, A.a22);
        System.out.println("lambda1 = " + lambda[0]);
        System.out.println("lambda2 = " + lambda[1]);
        System.out.println("classification = " + classify2x2(A));
    }

    public static void main(String[] args) {
        String[] names = {
            "Stable diagonal",
            "Stable damped oscillator",
            "Stable but nonnormal",
            "Marginal oscillator",
            "Marginal with semisimple zero",
            "Unstable positive eigenvalue",
            "Unstable defective zero"
        };

        Matrix2[] matrices = {
            new Matrix2(-2.0, 0.0, 0.0, -5.0),
            new Matrix2(0.0, 1.0, -4.0, -2.0),
            new Matrix2(-1.0, 50.0, 0.0, -2.0),
            new Matrix2(0.0, 1.0, -1.0, 0.0),
            new Matrix2(0.0, 0.0, 0.0, -2.0),
            new Matrix2(1.0, 0.0, 0.0, -2.0),
            new Matrix2(0.0, 1.0, 0.0, 0.0)
        };

        for (int i = 0; i < names.length; i++) {
            report(names[i], matrices[i]);
        }
    }
}

10. MATLAB/Simulink Implementation

Chapter9_Lesson5.m uses eig and expm. In Simulink, each example can be represented by a State-Space block with matrices \( \mathbf{A} \), \( \mathbf{B}=\mathbf{0} \), \( \mathbf{C}=\mathbf{I} \), and \( \mathbf{D}=\mathbf{0} \). The initial condition is set inside the State-Space block to visualize free response.

% Chapter9_Lesson5.m
%
% Modern Control — Chapter 9, Lesson 5
% Examples of Stable, Marginal, and Unstable State Matrices
%
% This MATLAB script classifies representative continuous-time LTI matrices
% for x_dot = A x, computes eigenvalues, and plots selected state norms.

clear; clc; close all;

examples = {
    'Stable diagonal',              [-2 0; 0 -5];
    'Stable damped oscillator',     [0 1; -4 -2];
    'Stable but nonnormal',         [-1 50; 0 -2];
    'Marginal oscillator',          [0 1; -1 0];
    'Marginal with semisimple zero',[0 0; 0 -2];
    'Unstable positive eigenvalue', [1 0; 0 -2];
    'Unstable defective zero',      [0 1; 0 0]
};

for k = 1:size(examples,1)
    name = examples{k,1};
    A = examples{k,2};
    fprintf('\n%s\n', name);
    fprintf('%s\n', repmat('-',1,length(name)));
    disp('A ='); disp(A);
    disp('eigenvalues ='); disp(eig(A).');
    fprintf('classification = %s\n', classifyContinuousLTI(A));
end

% Exact simulation using x(t) = expm(A t) x0.
t = linspace(0,10,500);
x0 = [1; 1];

selected = {
    'Stable damped oscillator', [0 1; -4 -2];
    'Marginal oscillator',      [0 1; -1 0];
    'Unstable defective zero',  [0 1; 0 0]
};

figure;
hold on;
for k = 1:size(selected,1)
    label = selected{k,1};
    A = selected{k,2};
    normx = zeros(size(t));
    for i = 1:numel(t)
        x = expm(A*t(i))*x0;
        normx(i) = norm(x,2);
    end
    plot(t,normx,'DisplayName',label,'LineWidth',1.5);
end
grid on;
xlabel('time t');
ylabel('state norm ||x(t)||_2');
title('Stable, Marginal, and Unstable State-Matrix Examples');
legend('Location','best');

function cls = classifyContinuousLTI(A)
    tol = 1e-9;
    lambda = eig(A);
    alpha = max(real(lambda));

    if alpha < -tol
        cls = 'asymptotically stable';
        return;
    elseif alpha > tol
        cls = 'unstable';
        return;
    end

    % Boundary case: check semisimplicity of eigenvalues on imaginary axis.
    n = size(A,1);
    used = false(numel(lambda),1);

    for i = 1:numel(lambda)
        if used(i)
            continue;
        end

        close = abs(lambda - lambda(i)) < 1e-7;
        used(close) = true;
        algebraicMultiplicity = sum(close);

        if abs(real(lambda(i))) <= 1e-7
            geometricMultiplicity = n - rank(A - lambda(i)*eye(n), 1e-7);
            if geometricMultiplicity < algebraicMultiplicity
                cls = 'unstable';
                return;
            end
        end
    end

    cls = 'marginally stable';
end

11. Wolfram Mathematica Implementation

Chapter9_Lesson5.nb uses Eigenvalues, MatrixRank, and MatrixExp to classify matrices symbolically or numerically and plot exact state norms.

(* Chapter9_Lesson5.nb *)

(* Modern Control — Chapter 9, Lesson 5
   Examples of Stable, Marginal, and Unstable State Matrices *)

ClearAll[classifyContinuousLTI];

classifyContinuousLTI[A_, tol_: 10^-9] := Module[
  {lambda, alpha, n, clusters, algebraicMultiplicity, geometricMultiplicity},
  lambda = N[Eigenvalues[A]];
  alpha = Max[Re[lambda]];

  If[alpha < -tol, Return["asymptotically stable"]];
  If[alpha > tol, Return["unstable"]];

  n = Length[A];

  clusters = Gather[lambda, Abs[#1 - #2] < 10^-7 &];

  Do[
    If[Abs[Re[cluster[[1]]]] <= 10^-7,
      algebraicMultiplicity = Length[cluster];
      geometricMultiplicity = n - MatrixRank[N[A - cluster[[1]] IdentityMatrix[n]], Tolerance -> 10^-7];
      If[geometricMultiplicity < algebraicMultiplicity, Return["unstable"]];
    ],
    {cluster, clusters}
  ];

  "marginally stable"
];

examples = {
  "Stable diagonal" -> { {-2, 0}, {0, -5} },
  "Stable damped oscillator" -> { {0, 1}, {-4, -2} },
  "Stable but nonnormal" -> { {-1, 50}, {0, -2} },
  "Marginal oscillator" -> { {0, 1}, {-1, 0} },
  "Marginal with semisimple zero" -> { {0, 0}, {0, -2} },
  "Unstable positive eigenvalue" -> { {1, 0}, {0, -2} },
  "Unstable defective zero" -> { {0, 1}, {0, 0} }
};

Table[
  {
    name,
    MatrixForm[A],
    N[Eigenvalues[A]],
    classifyContinuousLTI[A]
  },
  {item, examples},
  {name, {item[[1]]} },
  {A, {item[[2]]} }
] // Flatten[#, 1] & // TableForm[
  TableHeadings -> {None, {"Name", "A", "Eigenvalues", "Classification"} }
]

(* Exact state trajectories for selected examples *)
tmax = 10;
x0 = {1, 1};

selected = {
  "Stable damped oscillator" -> { {0, 1}, {-4, -2} },
  "Marginal oscillator" -> { {0, 1}, {-1, 0} },
  "Unstable defective zero" -> { {0, 1}, {0, 0} }
};

Plot[
  Evaluate[
    Table[
      Norm[MatrixExp[A t].x0, 2],
      {item, selected},
      {A, {item[[2]]} }
    ] // Flatten
  ],
  {t, 0, tmax},
  PlotLegends -> (selected[[All, 1]]),
  AxesLabel -> {"t", "||x(t)||2"},
  PlotLabel -> "Stable, Marginal, and Unstable State-Matrix Examples",
  GridLines -> Automatic
]

12. Problems and Solutions

Problem 1: Classify \( \mathbf{A}=\begin{bmatrix}-3&0\\0&-1\end{bmatrix} \).

Solution: The eigenvalues are \( -3 \) and \( -1 \). Hence \( \alpha(\mathbf{A})=-1<0 \), so the origin is asymptotically stable. The exact solution is

\[ \mathbf{x}(t)= \begin{bmatrix} e^{-3t}x_1(0)\\ e^{-t}x_2(0) \end{bmatrix}, \qquad \lim_{t\to\infty}\mathbf{x}(t)=\mathbf{0}. \]

Problem 2: Classify \( \mathbf{A}=\begin{bmatrix}0&2\\-2&0\end{bmatrix} \).

Solution: The characteristic polynomial is

\[ p(s)=s^2+4, \qquad \lambda_{1,2}=\pm 2j. \]

The eigenvalues are distinct and lie on the imaginary axis, hence they are semisimple. The system is marginally stable but not asymptotically stable. The trajectories are periodic rotations with angular frequency \( 2 \).

Problem 3: Classify \( \mathbf{A}=\begin{bmatrix}0&1\\0&0\end{bmatrix} \) and explain why eigenvalues alone are not sufficient on the stability boundary.

Solution: The only eigenvalue is \( \lambda=0 \) with algebraic multiplicity two. However,

\[ \operatorname{rank}(\mathbf{A})=1, \qquad \dim\ker(\mathbf{A})=1. \]

The geometric multiplicity is one, so the zero eigenvalue is defective. Since \( \mathbf{A}^2=\mathbf{0} \),

\[ e^{\mathbf{A}t}=\mathbf{I}+t\mathbf{A} = \begin{bmatrix} 1&t\\ 0&1 \end{bmatrix}. \]

The response is unbounded for initial states with \( x_2(0)\neq 0 \), so the system is unstable.

Problem 4: Let \( \mathbf{A}=\begin{bmatrix}-1&100\\0&-2\end{bmatrix} \). Is the system asymptotically stable? Can the state norm initially grow?

Solution: The eigenvalues are \( -1 \) and \( -2 \), so the system is asymptotically stable. Because the matrix is highly nonnormal, the state norm may initially grow for some initial conditions. The exponential can be computed from the triangular form:

\[ e^{\mathbf{A}t} = \begin{bmatrix} e^{-t} & 100\left(e^{-t}-e^{-2t}\right)\\ 0 & e^{-2t} \end{bmatrix}. \]

The off-diagonal term may be large for intermediate time, although all entries eventually decay to zero.

Problem 5: Prove that if \( \mathbf{A} \) is diagonalizable and all eigenvalues satisfy \( \operatorname{Re}(\lambda_i)\leq 0 \), then the system is stable provided every eigenvalue on the imaginary axis is included in the diagonalization.

Solution: If \( \mathbf{A}=\mathbf{T}\boldsymbol{\Lambda}\mathbf{T}^{-1} \), then

\[ e^{\mathbf{A}t} = \mathbf{T} \operatorname{diag}\left(e^{\lambda_1t},\dots,e^{\lambda_nt}\right) \mathbf{T}^{-1}. \]

For each eigenvalue with negative real part, the term decays. For each eigenvalue with zero real part, the term has constant magnitude. Hence there exists a constant \( M \) such that \( \|e^{\mathbf{A}t}\|\leq M \) for all \( t\geq 0 \). Therefore, the origin is Lyapunov stable. It is asymptotically stable only if no zero-real-part eigenvalue exists.

13. Summary

Stable, marginal, and unstable state matrices can be distinguished by how \( e^{\mathbf{A}t} \) behaves. Negative real eigenvalue parts imply decay. Imaginary-axis eigenvalues imply marginal behavior only if they are semisimple. Positive real parts or defective boundary eigenvalues imply instability. Nonnormal matrices can still be asymptotically stable while showing substantial transient growth.

14. References

  1. Lyapunov, A.M. (1892). The general problem of the stability of motion. Kharkov Mathematical Society.
  2. Massera, J.L. (1949). On Liapounoff's conditions of stability. Annals of Mathematics, 50(3), 705–721.
  3. LaSalle, J.P. (1960). Some extensions of Liapunov's second method. IRE Transactions on Circuit Theory, 7(4), 520–527.
  4. Perron, O. (1930). Die Stabilitätsfrage bei Differentialgleichungen. Mathematische Zeitschrift, 32, 703–728.
  5. Kreiss, H.O. (1962). Über die Stabilitätsdefinition für Differenzengleichungen die partielle Differentialgleichungen approximieren. BIT Numerical Mathematics, 2, 153–181.
  6. Moler, C., & Van Loan, C. (1978). Nineteen dubious ways to compute the exponential of a matrix. SIAM Review, 20(4), 801–836.
  7. Van Loan, C. (1977). The sensitivity of the matrix exponential. SIAM Journal on Numerical Analysis, 14(6), 971–981.
  8. Trefethen, L.N., Trefethen, A.E., Reddy, S.C., & Driscoll, T.A. (1993). Hydrodynamic stability without eigenvalues. Science, 261(5121), 578–584.