In Lesson 1, the state transition matrix was introduced as the matrix that maps an initial state to a later state for the homogeneous LTI system

\[ \dot{\mathbf{x} }(t)=\mathbf{A}\mathbf{x}(t), \qquad \mathbf{x}(0)=\mathbf{x}_0. \]

The solution is \( \mathbf{x}(t)=\Phi(t)\mathbf{x}_0 \), where \( \Phi(t)=e^{\mathbf{A}t} \). In this lesson, we study the internal algebra of \( \Phi(t) \). By the end of the lesson, students should be able to:

• prove the semigroup property \( \Phi(t+s)=\Phi(t)\Phi(s) \);
• prove that \( \Phi(t) \) is nonsingular for every finite \( t \);
• derive \( \Phi(t)^{-1}=\Phi(-t) \);
• interpret forward and backward propagation geometrically in state space;
• verify the identities computationally in Python, C++, Java, MATLAB, Simulink, and Wolfram Mathematica.

For a constant matrix \( \mathbf{A}\in\mathbb{R}^{n\times n} \), the state transition matrix is

\[ \Phi(t)=e^{\mathbf{A}t} =\sum_{k=0}^{\infty}\frac{(\mathbf{A}t)^k}{k!}. \]

It satisfies the matrix initial value problem

\[ \frac{d}{dt}\Phi(t)=\mathbf{A}\Phi(t), \qquad \Phi(0)=\mathbf{I}. \]

Therefore each column of \( \Phi(t) \) is a solution of the homogeneous system initialized from one of the basis vectors \( \mathbf{e}_1,\dots,\mathbf{e}_n \). If \( \mathbf{x}_0=\sum_{i=1}^{n}\alpha_i\mathbf{e}_i \), linearity gives

\[ \mathbf{x}(t)=\Phi(t)\mathbf{x}_0 =\sum_{i=1}^{n}\alpha_i\Phi(t)\mathbf{e}_i. \]

Thus \( \Phi(t) \) is not merely a computational shortcut. It is the complete linear map that carries every initial condition at time zero to its state at time \( t \).

The semigroup property states that propagation over a total interval can be decomposed into consecutive subintervals:

\[ \boxed{\Phi(t+s)=\Phi(t)\Phi(s)}. \]

In words, first propagate the state by time \( s \), then propagate the result by time \( t \). The final state is the same as propagating once by time \( t+s \).

flowchart LR
  X0["Initial state: x0"]
  XS["Intermediate state: after time s"]
  XTS["Final state: after time t plus s"]

  X0 -- "Phi_s" --> XS
  XS -- "Phi_t" --> XTS
  X0 -- "Phi_t_plus_s" --> XTS
  XTS -- "Phi_minus_t_plus_s" --> X0
  

For an initial state \( \mathbf{x}_0 \), the one-step and two-step descriptions are

\[ \mathbf{x}(t+s)=\Phi(t+s)\mathbf{x}_0, \]

\[ \mathbf{x}(t+s)=\Phi(t)\mathbf{x}(s) =\Phi(t)\Phi(s)\mathbf{x}_0. \]

Since this equality holds for every \( \mathbf{x}_0\in\mathbb{R}^{n} \), the matrices must be equal:

\[ \Phi(t+s)=\Phi(t)\Phi(s). \]

Because \( \mathbf{A}t \) and \( \mathbf{A}s \) are scalar multiples of the same matrix, they commute:

\[ (\mathbf{A}t)(\mathbf{A}s)=ts\mathbf{A}^2 =(\mathbf{A}s)(\mathbf{A}t). \]

For commuting matrices, the exponential identity \( e^{\mathbf{M}+\mathbf{N} }=e^{\mathbf{M} }e^{\mathbf{N} } \) holds. Setting \( \mathbf{M}=\mathbf{A}t \) and \( \mathbf{N}=\mathbf{A}s \) gives

\[ e^{\mathbf{A}(t+s)} =e^{\mathbf{A}t+\mathbf{A}s} =e^{\mathbf{A}t}e^{\mathbf{A}s}. \]

Hence

\[ \Phi(t+s)=\Phi(t)\Phi(s). \]

A direct series derivation is also possible. Since \( \Phi(t)=\sum_{i=0}^{\infty}(\mathbf{A}t)^i/i! \) and \( \Phi(s)=\sum_{j=0}^{\infty}(\mathbf{A}s)^j/j! \),

\[ \Phi(t)\Phi(s) =\sum_{i=0}^{\infty}\sum_{j=0}^{\infty} \frac{\mathbf{A}^{i+j}t^i s^j}{i!j!}. \]

Grouping all terms for which \( i+j=k \), we obtain

\[ \Phi(t)\Phi(s) =\sum_{k=0}^{\infty}\mathbf{A}^{k} \sum_{i=0}^{k}\frac{t^i s^{k-i} }{i!(k-i)!}. \]

By the binomial theorem,

\[ \sum_{i=0}^{k}\frac{t^i s^{k-i} }{i!(k-i)!} =\frac{(t+s)^k}{k!}. \]

Therefore,

\[ \Phi(t)\Phi(s) =\sum_{k=0}^{\infty}\frac{\mathbf{A}^{k}(t+s)^k}{k!} =e^{\mathbf{A}(t+s)}=\Phi(t+s). \]

The inverse property follows from the semigroup property by choosing \( s=-t \):

\[ \Phi(t)\Phi(-t)=\Phi(t-t)=\Phi(0)=\mathbf{I}. \]

Similarly,

\[ \Phi(-t)\Phi(t)=\Phi(-t+t)=\Phi(0)=\mathbf{I}. \]

Hence

\[ \boxed{\Phi(t)^{-1}=\Phi(-t)}. \]

This means that the homogeneous LTI state map is reversible as a mathematical mapping: if \( \mathbf{x}(t)=\Phi(t)\mathbf{x}_0 \), then the initial state can be recovered by

\[ \mathbf{x}_0=\Phi(-t)\mathbf{x}(t). \]

This reversibility does not mean that a physical system is necessarily easy to reverse in practice. If some modes decay rapidly, then \( \Phi(-t) \) amplifies those directions. Therefore the inverse can be numerically sensitive even though it exists analytically.

A second proof of invertibility uses the determinant identity for the matrix exponential:

\[ \det(e^{\mathbf{M} })=e^{\operatorname{tr}(\mathbf{M})}. \]

Taking \( \mathbf{M}=\mathbf{A}t \) yields

\[ \det(\Phi(t))=\det(e^{\mathbf{A}t}) =e^{\operatorname{tr}(\mathbf{A})t}. \]

Since \( e^{\operatorname{tr}(\mathbf{A})t}\neq 0 \) for every finite \( t \), \( \Phi(t) \) is nonsingular:

\[ \det(\Phi(t))\neq 0 \quad \Longrightarrow \quad \Phi(t)^{-1}\ \text{exists}. \]

This result has an important geometric meaning. For a small volume element in state space, the homogeneous flow scales volume by \( |\det(\Phi(t))| \). Thus

\[ |\det(\Phi(t))|=e^{\operatorname{tr}(\mathbf{A})t}. \]

If \( \operatorname{tr}(\mathbf{A})<0 \), volumes contract. If \( \operatorname{tr}(\mathbf{A})>0 \), volumes expand. If \( \operatorname{tr}(\mathbf{A})=0 \), volumes are preserved.

Although this chapter focuses on LTI systems, it is useful to introduce the notation \( \Phi(t,t_0) \), meaning the transition matrix from initial time \( t_0 \) to final time \( t \). For an LTI system,

\[ \Phi(t,t_0)=e^{\mathbf{A}(t-t_0)}. \]

The composition law becomes

\[ \boxed{\Phi(t_2,t_0)=\Phi(t_2,t_1)\Phi(t_1,t_0)}. \]

The inverse becomes

\[ \boxed{\Phi(t,t_0)^{-1}=\Phi(t_0,t)}. \]

For the LTI case, this is simply

\[ \Phi(t,t_0)^{-1} =\left(e^{\mathbf{A}(t-t_0)}\right)^{-1} =e^{-\mathbf{A}(t-t_0)} =e^{\mathbf{A}(t_0-t)} =\Phi(t_0,t). \]

This two-time notation is especially useful when later lessons discuss finite-time controllability, observability, and time-varying systems.

The identities \( \Phi(t+s)=\Phi(t)\Phi(s) \) and \( \Phi(t)^{-1}=\Phi(-t) \) are exact mathematical identities. In floating-point computation, however, the residuals

\[ \|\Phi(t+s)-\Phi(t)\Phi(s)\| \]

and

\[ \|\Phi(t)^{-1}-\Phi(-t)\| \]

are usually small but not exactly zero. Numerical error depends on the algorithm used to compute the matrix exponential, the norm of \( \mathbf{A}t \), and the conditioning of the state basis.

flowchart TD
  A["Choose state matrix A"] --> B["Compute Phi(t), Phi(s), Phi(t+s)"]
  B --> C["Check semigroup residual"]
  C --> D["Compute Phi(-t) and inverse of Phi(t)"]
  D --> E["Check inverse residual"]
  E --> F["Propagate state directly and in two steps"]
  F --> G["Compare numerical errors"]
        

A useful diagnostic is the condition number \( \kappa(\Phi(t)) \). If it is large, then reconstructing the initial state from \( \mathbf{x}(t) \) may be highly sensitive to measurement noise:

\[ \frac{\|\delta \mathbf{x}_0\|}{\|\mathbf{x}_0\|} \lesssim \kappa(\Phi(t)) \frac{\|\delta \mathbf{x}(t)\|}{\|\mathbf{x}(t)\|}. \]

Therefore the inverse exists mathematically, but inverse propagation can be numerically ill-conditioned for systems with rapidly decaying or rapidly growing modes.

The Python implementation uses scipy.linalg.expm to compute the state transition matrix and verifies semigroup and inverse residuals.

"""
Chapter8_Lesson2.py
Modern Control — Chapter 8, Lesson 2
Semigroup Property and Inverse of the State Transition Matrix Phi(t)

Requires: numpy, scipy
Install if needed: pip install numpy scipy
"""

import numpy as np
from scipy.linalg import expm, norm


def phi(A: np.ndarray, t: float) -> np.ndarray:
    """State transition matrix Phi(t) = exp(A t)."""
    return expm(A * t)


def main() -> None:
    # A stable second-order continuous-time state matrix.
    A = np.array([[0.0, 1.0],
                  [-4.0, -0.8]])

    t = 0.7
    s = 1.3
    x0 = np.array([1.0, -0.25])

    Phi_t = phi(A, t)
    Phi_s = phi(A, s)
    Phi_t_plus_s = phi(A, t + s)
    Phi_minus_t = phi(A, -t)

    # Semigroup property: Phi(t+s) = Phi(t) Phi(s).
    semigroup_error = norm(Phi_t_plus_s - Phi_t @ Phi_s, ord="fro")

    # Inverse property: Phi(t)^(-1) = Phi(-t).
    inverse_error = norm(np.linalg.inv(Phi_t) - Phi_minus_t, ord="fro")
    identity_error = norm(Phi_t @ Phi_minus_t - np.eye(A.shape[0]), ord="fro")

    # State propagation equivalence.
    x_direct = Phi_t_plus_s @ x0
    x_two_step = Phi_t @ (Phi_s @ x0)
    state_error = norm(x_direct - x_two_step)

    print("A =")
    print(A)
    print("\nPhi(t) =")
    print(Phi_t)
    print("\nPhi(s) =")
    print(Phi_s)
    print("\nPhi(t+s) =")
    print(Phi_t_plus_s)
    print("\nSemigroup error ||Phi(t+s)-Phi(t)Phi(s)||_F =", semigroup_error)
    print("Inverse error ||inv(Phi(t))-Phi(-t)||_F =", inverse_error)
    print("Identity error ||Phi(t)Phi(-t)-I||_F =", identity_error)
    print("State two-step propagation error =", state_error)


if __name__ == "__main__":
    main()
      

The C++ implementation uses a compact from-scratch scaling-and-squaring Taylor approximation for a two-state example. In production control software, one should normally use robust numerical linear algebra libraries such as Eigen, Armadillo, or vendor BLAS/LAPACK wrappers.

/*
Chapter8_Lesson2.cpp
Modern Control — Chapter 8, Lesson 2
Semigroup Property and Inverse of the State Transition Matrix Phi(t)

Compile:
  g++ -std=c++17 Chapter8_Lesson2.cpp -O2 -o Chapter8_Lesson2
*/

#include <array>
#include <cmath>
#include <iomanip>
#include <iostream>

using Matrix2 = std::array<std::array<double, 2>, 2>;
using Vector2 = std::array<double, 2>;

Matrix2 zero2() {
    return { { {0.0, 0.0}, {0.0, 0.0} } };
}

Matrix2 identity2() {
    return { { {1.0, 0.0}, {0.0, 1.0} } };
}

Matrix2 add(const Matrix2& A, const Matrix2& B) {
    Matrix2 C = zero2();
    for (int i = 0; i < 2; ++i)
        for (int j = 0; j < 2; ++j)
            C[i][j] = A[i][j] + B[i][j];
    return C;
}

Matrix2 scale(const Matrix2& A, double c) {
    Matrix2 B = zero2();
    for (int i = 0; i < 2; ++i)
        for (int j = 0; j < 2; ++j)
            B[i][j] = c * A[i][j];
    return B;
}

Matrix2 multiply(const Matrix2& A, const Matrix2& B) {
    Matrix2 C = zero2();
    for (int i = 0; i < 2; ++i)
        for (int j = 0; j < 2; ++j)
            for (int k = 0; k < 2; ++k)
                C[i][j] += A[i][k] * B[k][j];
    return C;
}

Vector2 multiply(const Matrix2& A, const Vector2& x) {
    return {A[0][0] * x[0] + A[0][1] * x[1],
            A[1][0] * x[0] + A[1][1] * x[1]};
}

double frobeniusNorm(const Matrix2& A) {
    double sum = 0.0;
    for (int i = 0; i < 2; ++i)
        for (int j = 0; j < 2; ++j)
            sum += A[i][j] * A[i][j];
    return std::sqrt(sum);
}

double vectorNorm(const Vector2& x) {
    return std::sqrt(x[0] * x[0] + x[1] * x[1]);
}

Matrix2 subtract(const Matrix2& A, const Matrix2& B) {
    return add(A, scale(B, -1.0));
}

Vector2 subtract(const Vector2& a, const Vector2& b) {
    return {a[0] - b[0], a[1] - b[1]};
}

Matrix2 inverse2(const Matrix2& A) {
    double det = A[0][0] * A[1][1] - A[0][1] * A[1][0];
    return { { { A[1][1] / det, -A[0][1] / det},
             {-A[1][0] / det,  A[0][0] / det} } };
}

Matrix2 matrixExpTaylor2(const Matrix2& A, double t) {
    // Scaling-and-squaring with Taylor series.
    // For teaching purposes: robust enough for small 2x2 examples.
    Matrix2 M = scale(A, t);
    double normM = frobeniusNorm(M);
    int scalingPower = std::max(0, static_cast<int>(std::ceil(std::log2(normM + 1e-12))));
    double divisor = std::pow(2.0, scalingPower);
    Matrix2 B = scale(M, 1.0 / divisor);

    Matrix2 result = identity2();
    Matrix2 term = identity2();
    const int terms = 40;
    for (int k = 1; k <= terms; ++k) {
        term = multiply(term, scale(B, 1.0 / static_cast<double>(k)));
        result = add(result, term);
    }

    for (int i = 0; i < scalingPower; ++i) {
        result = multiply(result, result);
    }
    return result;
}

void printMatrix(const std::string& name, const Matrix2& A) {
    std::cout << name << " =\n";
    for (const auto& row : A) {
        std::cout << "  [" << std::setw(12) << row[0] << ", "
                  << std::setw(12) << row[1] << "]\n";
    }
}

int main() {
    Matrix2 A = { { {0.0, 1.0}, {-4.0, -0.8} } };
    double t = 0.7;
    double s = 1.3;
    Vector2 x0 = {1.0, -0.25};

    Matrix2 Phi_t = matrixExpTaylor2(A, t);
    Matrix2 Phi_s = matrixExpTaylor2(A, s);
    Matrix2 Phi_t_plus_s = matrixExpTaylor2(A, t + s);
    Matrix2 Phi_minus_t = matrixExpTaylor2(A, -t);

    Matrix2 semigroupResidual = subtract(Phi_t_plus_s, multiply(Phi_t, Phi_s));
    Matrix2 inverseResidual = subtract(inverse2(Phi_t), Phi_minus_t);
    Matrix2 identityResidual = subtract(multiply(Phi_t, Phi_minus_t), identity2());

    Vector2 x_direct = multiply(Phi_t_plus_s, x0);
    Vector2 x_two_step = multiply(Phi_t, multiply(Phi_s, x0));

    std::cout << std::setprecision(10);
    printMatrix("Phi(t)", Phi_t);
    printMatrix("Phi(s)", Phi_s);
    printMatrix("Phi(t+s)", Phi_t_plus_s);
    std::cout << "Semigroup error = " << frobeniusNorm(semigroupResidual) << "\n";
    std::cout << "Inverse error = " << frobeniusNorm(inverseResidual) << "\n";
    std::cout << "Identity error = " << frobeniusNorm(identityResidual) << "\n";
    std::cout << "State two-step propagation error = " << vectorNorm(subtract(x_direct, x_two_step)) << "\n";
    return 0;
}
      

The Java implementation mirrors the C++ version and computes the matrix exponential from scratch for a two-state system.

/*
Chapter8_Lesson2.java
Modern Control — Chapter 8, Lesson 2
Semigroup Property and Inverse of the State Transition Matrix Phi(t)

Compile and run:
  javac Chapter8_Lesson2.java
  java Chapter8_Lesson2
*/

public class Chapter8_Lesson2 {
    static double[][] identity2() {
        return new double[][]{ {1.0, 0.0}, {0.0, 1.0} };
    }

    static double[][] add(double[][] A, double[][] B) {
        double[][] C = new double[2][2];
        for (int i = 0; i < 2; i++)
            for (int j = 0; j < 2; j++)
                C[i][j] = A[i][j] + B[i][j];
        return C;
    }

    static double[][] subtract(double[][] A, double[][] B) {
        return add(A, scale(B, -1.0));
    }

    static double[][] scale(double[][] A, double c) {
        double[][] B = new double[2][2];
        for (int i = 0; i < 2; i++)
            for (int j = 0; j < 2; j++)
                B[i][j] = c * A[i][j];
        return B;
    }

    static double[][] multiply(double[][] A, double[][] B) {
        double[][] C = new double[2][2];
        for (int i = 0; i < 2; i++)
            for (int j = 0; j < 2; j++)
                for (int k = 0; k < 2; k++)
                    C[i][j] += A[i][k] * B[k][j];
        return C;
    }

    static double[] multiply(double[][] A, double[] x) {
        return new double[]{A[0][0] * x[0] + A[0][1] * x[1],
                            A[1][0] * x[0] + A[1][1] * x[1]};
    }

    static double[] subtract(double[] a, double[] b) {
        return new double[]{a[0] - b[0], a[1] - b[1]};
    }

    static double frobeniusNorm(double[][] A) {
        double sum = 0.0;
        for (int i = 0; i < 2; i++)
            for (int j = 0; j < 2; j++)
                sum += A[i][j] * A[i][j];
        return Math.sqrt(sum);
    }

    static double vectorNorm(double[] x) {
        return Math.sqrt(x[0] * x[0] + x[1] * x[1]);
    }

    static double[][] inverse2(double[][] A) {
        double det = A[0][0] * A[1][1] - A[0][1] * A[1][0];
        return new double[][]{ { A[1][1] / det, -A[0][1] / det},
                              {-A[1][0] / det,  A[0][0] / det} };
    }

    static double[][] matrixExpTaylor2(double[][] A, double t) {
        // Scaling-and-squaring with Taylor series for a small 2x2 teaching example.
        double[][] M = scale(A, t);
        double normM = frobeniusNorm(M);
        int scalingPower = Math.max(0, (int)Math.ceil(Math.log(normM + 1e-12) / Math.log(2.0)));
        double divisor = Math.pow(2.0, scalingPower);
        double[][] B = scale(M, 1.0 / divisor);

        double[][] result = identity2();
        double[][] term = identity2();
        int terms = 40;
        for (int k = 1; k <= terms; k++) {
            term = multiply(term, scale(B, 1.0 / k));
            result = add(result, term);
        }
        for (int i = 0; i < scalingPower; i++) {
            result = multiply(result, result);
        }
        return result;
    }

    static void printMatrix(String name, double[][] A) {
        System.out.println(name + " =");
        for (int i = 0; i < 2; i++) {
            System.out.printf("  [%12.8f, %12.8f]%n", A[i][0], A[i][1]);
        }
    }

    public static void main(String[] args) {
        double[][] A = { {0.0, 1.0}, {-4.0, -0.8} };
        double t = 0.7;
        double s = 1.3;
        double[] x0 = {1.0, -0.25};

        double[][] Phi_t = matrixExpTaylor2(A, t);
        double[][] Phi_s = matrixExpTaylor2(A, s);
        double[][] Phi_t_plus_s = matrixExpTaylor2(A, t + s);
        double[][] Phi_minus_t = matrixExpTaylor2(A, -t);

        double semigroupError = frobeniusNorm(subtract(Phi_t_plus_s, multiply(Phi_t, Phi_s)));
        double inverseError = frobeniusNorm(subtract(inverse2(Phi_t), Phi_minus_t));
        double identityError = frobeniusNorm(subtract(multiply(Phi_t, Phi_minus_t), identity2()));

        double[] xDirect = multiply(Phi_t_plus_s, x0);
        double[] xTwoStep = multiply(Phi_t, multiply(Phi_s, x0));

        printMatrix("Phi(t)", Phi_t);
        printMatrix("Phi(s)", Phi_s);
        printMatrix("Phi(t+s)", Phi_t_plus_s);
        System.out.println("Semigroup error = " + semigroupError);
        System.out.println("Inverse error = " + inverseError);
        System.out.println("Identity error = " + identityError);
        System.out.println("State two-step propagation error = " + vectorNorm(subtract(xDirect, xTwoStep)));
    }
}
      

MATLAB provides expm for the matrix exponential. The script also attempts to build a simple Simulink model for the homogeneous state-space system.

% Chapter8_Lesson2.m
% Modern Control — Chapter 8, Lesson 2
% Semigroup Property and Inverse of the State Transition Matrix Phi(t)
%
% This script verifies Phi(t+s)=Phi(t)Phi(s) and inv(Phi(t))=Phi(-t).
% It also contains an optional Simulink construction section.

clear; clc;

A = [0 1; -4 -0.8];
t = 0.7;
s = 1.3;
x0 = [1; -0.25];

Phi_t = expm(A*t);
Phi_s = expm(A*s);
Phi_t_plus_s = expm(A*(t+s));
Phi_minus_t = expm(-A*t);

semigroup_error = norm(Phi_t_plus_s - Phi_t*Phi_s, 'fro');
inverse_error = norm(inv(Phi_t) - Phi_minus_t, 'fro');
identity_error = norm(Phi_t*Phi_minus_t - eye(size(A)), 'fro');

x_direct = Phi_t_plus_s*x0;
x_two_step = Phi_t*(Phi_s*x0);
state_error = norm(x_direct - x_two_step);

disp('Phi(t) ='); disp(Phi_t);
disp('Phi(s) ='); disp(Phi_s);
disp('Phi(t+s) ='); disp(Phi_t_plus_s);
fprintf('Semigroup error = %.3e\n', semigroup_error);
fprintf('Inverse error = %.3e\n', inverse_error);
fprintf('Identity error = %.3e\n', identity_error);
fprintf('State two-step propagation error = %.3e\n', state_error);

% Optional Simulink model generation:
% The model contains a continuous-time State-Space block with zero input.
% It is a vehicle for visualizing the homogeneous evolution xdot = A x.
try
    modelName = 'Chapter8_Lesson2_Simulink';
    if bdIsLoaded(modelName)
        close_system(modelName, 0);
    end
    new_system(modelName);
    open_system(modelName);

    add_block('simulink/Sources/Constant', [modelName '/ZeroInput'], ...
        'Value', '0', 'Position', [80 80 120 110]);
    add_block('simulink/Continuous/State-Space', [modelName '/StateSpace'], ...
        'A', mat2str(A), 'B', '[0;0]', 'C', 'eye(2)', 'D', '[0;0]', ...
        'X0', mat2str(x0), 'Position', [180 65 300 125]);
    add_block('simulink/Sinks/Scope', [modelName '/Scope'], ...
        'Position', [360 70 420 120]);
    add_line(modelName, 'ZeroInput/1', 'StateSpace/1');
    add_line(modelName, 'StateSpace/1', 'Scope/1');
    set_param(modelName, 'StopTime', num2str(t+s));

    fprintf('Created Simulink model: %s\n', modelName);
    fprintf('Run sim(modelName) in MATLAB to visualize x(t).\n');
catch ME
    fprintf('Simulink model generation skipped: %s\n', ME.message);
end
      

The Mathematica implementation uses MatrixExp to verify the same algebraic properties.

A = { {0, 1}, {-4, -0.8} };
t = 0.7; s = 1.3; x0 = {1, -0.25};

Phi[tau_] := MatrixExp[A tau];

PhiT = Phi[t];
PhiS = Phi[s];
PhiTS = Phi[t + s];
PhiMinusT = Phi[-t];

semigroupError = Norm[PhiTS - PhiT.PhiS];
inverseError = Norm[Inverse[PhiT] - PhiMinusT];
identityError = Norm[PhiT.PhiMinusT - IdentityMatrix[2]];
stateError = Norm[PhiTS.x0 - PhiT.(PhiS.x0)];

{semigroupError, inverseError, identityError, stateError}
      

Problem 1: Let \( \Phi(t)=e^{\mathbf{A}t} \). Prove that \( \Phi(t)\Phi(s)=\Phi(s)\Phi(t) \) for all real \( t,s \).

Solution: Since \( \Phi(t)\Phi(s)=e^{\mathbf{A}t}e^{\mathbf{A}s} \), and \( \mathbf{A}t \) commutes with \( \mathbf{A}s \), we have

\[ e^{\mathbf{A}t}e^{\mathbf{A}s} =e^{\mathbf{A}(t+s)} =e^{\mathbf{A}(s+t)} =e^{\mathbf{A}s}e^{\mathbf{A}t}. \]

Therefore \( \Phi(t)\Phi(s)=\Phi(s)\Phi(t) \).

Problem 2: Show that \( \Phi(0)=\mathbf{I} \) and use this fact to prove \( \Phi(t)^{-1}=\Phi(-t) \).

Solution: From the matrix exponential series,

\[ \Phi(0)=e^{\mathbf{A}0} =\sum_{k=0}^{\infty}\frac{(\mathbf{A}0)^k}{k!} =\mathbf{I}. \]

Using the semigroup property,

\[ \Phi(t)\Phi(-t)=\Phi(t-t)=\Phi(0)=\mathbf{I}. \]

Similarly, \( \Phi(-t)\Phi(t)=\mathbf{I} \). Therefore \( \Phi(-t) \) is the inverse of \( \Phi(t) \).

Problem 3: Let \( \mathbf{A}=\begin{bmatrix}0 & 1\\ -2 & -3\end{bmatrix} \). Without explicitly computing \( e^{\mathbf{A}t} \), find \( \det(\Phi(t)) \).

Solution: The trace is \( \operatorname{tr}(\mathbf{A})=0+(-3)=-3 \). Hence

\[ \det(\Phi(t))=e^{\operatorname{tr}(\mathbf{A})t}=e^{-3t}. \]

Since \( e^{-3t}\neq 0 \), the transition matrix is invertible for every finite \( t \).

Problem 4: Suppose \( \mathbf{x}(2)=\Phi(2)\mathbf{x}_0 \). Express \( \mathbf{x}_0 \) in terms of \( \mathbf{x}(2) \).

Solution: Since \( \Phi(2)^{-1}=\Phi(-2) \),

\[ \mathbf{x}_0=\Phi(2)^{-1}\mathbf{x}(2)=\Phi(-2)\mathbf{x}(2). \]

Problem 5: For an LTI system, show that \( \Phi(t_2,t_0)=\Phi(t_2,t_1)\Phi(t_1,t_0) \).

Solution: For an LTI system,

\[ \Phi(t_2,t_0)=e^{\mathbf{A}(t_2-t_0)}. \]

Also,

\[ \Phi(t_2,t_1)\Phi(t_1,t_0) =e^{\mathbf{A}(t_2-t_1)}e^{\mathbf{A}(t_1-t_0)}. \]

The exponents commute because both are scalar multiples of \( \mathbf{A} \). Hence

\[ e^{\mathbf{A}(t_2-t_1)}e^{\mathbf{A}(t_1-t_0)} =e^{\mathbf{A}(t_2-t_1+t_1-t_0)} =e^{\mathbf{A}(t_2-t_0)} =\Phi(t_2,t_0). \]

Problem 6: Explain why the existence of \( \Phi(t)^{-1} \) does not guarantee that backward reconstruction of \( \mathbf{x}_0 \) is numerically reliable.

Solution: If the system has rapidly decaying modes, then forward propagation compresses state-space directions associated with those modes. Backward propagation by \( \Phi(-t) \) expands the same directions. Therefore small measurement errors in \( \mathbf{x}(t) \) may be greatly amplified when reconstructing \( \mathbf{x}_0 \). A large condition number \( \kappa(\Phi(t)) \) indicates this sensitivity.

The state transition matrix of a continuous-time LTI system satisfies two central identities:

\[ \Phi(t+s)=\Phi(t)\Phi(s), \qquad \Phi(t)^{-1}=\Phi(-t). \]

The first identity formalizes the idea that state propagation can be split into consecutive intervals. The second identity shows that the homogeneous state map is invertible for every finite time. The determinant formula \( \det(\Phi(t))=e^{\operatorname{tr}(\mathbf{A})t} \) confirms nonsingularity and gives a geometric interpretation as state-space volume scaling. These properties are essential for later topics such as controllability, observability, Gramian analysis, and canonical representations.