In earlier lessons, controllability of the LTI model \( \dot{\mathbf{x} }=\mathbf{A}\mathbf{x}+\mathbf{B}\mathbf{u} \) was tested by the Kalman matrix \( \mathcal{C}=[\mathbf{B},\mathbf{A}\mathbf{B},\dots,\mathbf{A}^{n-1}\mathbf{B}] \). For a linear time-varying system, the model is

\[ \dot{\mathbf{x} }(t)=\mathbf{A}(t)\mathbf{x}(t)+\mathbf{B}(t)\mathbf{u}(t), \quad \mathbf{x}(t)\in\mathbb{R}^n, \quad \mathbf{u}(t)\in\mathbb{R}^m. \]

The main complication is that the actuator direction \( \mathbf{B}(t) \) and the natural propagation \( \mathbf{A}(t) \) can change during the interval. Therefore, controllability is an interval property, not merely a property of a single pair \( (\mathbf{A},\mathbf{B}) \). A system may be poorly actuated instantaneously but controllable over a finite interval because the dynamics rotate, stretch, or mix the input direction over time.

flowchart TD
  A["LTI controllability"] --> B["Fixed A and fixed B"]
  B --> C["Use columns: B, AB, ..., A^(n-1)B"]
  D["LTV controllability"] --> E["A(t) and B(t) vary over interval"]
  E --> F["Use Phi(tf,t) B(t) across all times"]
  F --> G["Reachable subspace depends on [t0, tf]"]
  G --> H["Check whether all state directions are generated"]
        

Conceptually, the LTV test asks: after each input impulse is injected at time \( t \), how does the remaining dynamics carry that input contribution to the final time \( t_f \)? The answer is encoded by the state transition matrix.

Let \( \boldsymbol{\Phi}(t,t_0) \) denote the state transition matrix of the homogeneous system \( \dot{\mathbf{x} }(t)=\mathbf{A}(t)\mathbf{x}(t) \). It is defined by

\[ \frac{\partial}{\partial t}\boldsymbol{\Phi}(t,t_0) =\mathbf{A}(t)\boldsymbol{\Phi}(t,t_0), \quad \boldsymbol{\Phi}(t_0,t_0)=\mathbf{I}. \]

For \( t_0 < t_f \), the solution of the nonhomogeneous LTV system is

\[ \mathbf{x}(t_f)=\boldsymbol{\Phi}(t_f,t_0)\mathbf{x}(t_0) +\int_{t_0}^{t_f}\boldsymbol{\Phi}(t_f,\tau)\mathbf{B}(\tau) \mathbf{u}(\tau)\,d\tau. \]

The first term is the free motion from the initial state. The second term is the reachable contribution due to the input. Hence, for fixed \( t_0 \) and \( t_f \), controllability is determined by the range of the linear operator

\[ \mathcal{L}_{t_0,t_f}[\mathbf{u}] =\int_{t_0}^{t_f}\boldsymbol{\Phi}(t_f,\tau)\mathbf{B}(\tau) \mathbf{u}(\tau)\,d\tau. \]

Thus, the reachable set from \( \mathbf{x}_0 \) at \( t_f \) is an affine set:

\[ \mathcal{R}(t_f;t_0,\mathbf{x}_0) =\boldsymbol{\Phi}(t_f,t_0)\mathbf{x}_0 +\operatorname{Range}(\mathcal{L}_{t_0,t_f}). \]

The LTV system is completely controllable on \( [t_0,t_f] \) if for every pair of states \( \mathbf{x}_0,\mathbf{x}_f\in\mathbb{R}^n \), there exists an admissible input \( \mathbf{u}(t) \) that transfers the system from \( \mathbf{x}(t_0)=\mathbf{x}_0 \) to \( \mathbf{x}(t_f)=\mathbf{x}_f \).

Define the finite-interval controllability matrix-like integral \( \mathbf{W}_c(t_0,t_f) \) by

\[ \mathbf{W}_c(t_0,t_f) =\int_{t_0}^{t_f}\boldsymbol{\Phi}(t_f,\tau) \mathbf{B}(\tau)\mathbf{B}^T(\tau) \boldsymbol{\Phi}^T(t_f,\tau)\,d\tau. \]

The central finite-interval criterion is

\[ \boxed{ \begin{gathered} \text{Complete controllability on }[t_0,t_f] \iff \operatorname{rank}\mathbf{W}_c(t_0,t_f)=n \\ \iff \mathbf{W}_c(t_0,t_f)\text{ is positive definite.} \end{gathered} } \]

This criterion is the LTV analogue of the Kalman rank test, but it uses all input directions after they have been propagated to the final time. The matrix \( \mathbf{W}_c \) is symmetric positive semidefinite because for any vector \( \mathbf{v} \),

\[ \mathbf{v}^T\mathbf{W}_c\mathbf{v} =\int_{t_0}^{t_f}\left\|\mathbf{B}^T(\tau) \boldsymbol{\Phi}^T(t_f,\tau)\mathbf{v}\right\|_2^2\,d\tau \geq 0. \]

flowchart TD
  S["Choose interval [t0, tf]"] --> P["Compute Phi(tf,t)"]
  P --> Q["Form propagated input directions Phi(tf,t)B(t)"]
  Q --> R["Accumulate Wc integral"]
  R --> E["Check eigenvalues or rank of Wc"]
  E -->|"full rank"| C["Controllable on this interval"]
  E -->|"rank deficient"| U["Some final-state direction is unreachable"]
        

The proof is useful because it explains why the integral test captures reachability. First, let

\[ \mathbf{r}=\mathbf{x}_f-\boldsymbol{\Phi}(t_f,t_0)\mathbf{x}_0. \]

The steering problem is equivalent to solving \( \mathcal{L}_{t_0,t_f}[\mathbf{u}]=\mathbf{r} \). If \( \mathbf{W}_c \) is nonsingular, choose

\[ \mathbf{u}(\tau)=\mathbf{B}^T(\tau) \boldsymbol{\Phi}^T(t_f,\tau) \mathbf{W}_c^{-1}(t_0,t_f)\mathbf{r}. \]

Substituting this input into the forced-response term gives

\[ \begin{aligned} \mathcal{L}_{t_0,t_f}[\mathbf{u}] &=\int_{t_0}^{t_f}\boldsymbol{\Phi}(t_f,\tau)\mathbf{B}(\tau) \mathbf{B}^T(\tau)\boldsymbol{\Phi}^T(t_f,\tau) \,d\tau\;\mathbf{W}_c^{-1}\mathbf{r}\\ &=\mathbf{W}_c\mathbf{W}_c^{-1}\mathbf{r}=\mathbf{r}. \end{aligned} \]

Therefore any target state can be reached. Conversely, suppose \( \mathbf{W}_c \) is singular. Then there exists a nonzero vector \( \mathbf{v}\neq\mathbf{0} \) such that \( \mathbf{W}_c\mathbf{v}=\mathbf{0} \). Since \( \mathbf{W}_c \) is positive semidefinite,

\[ 0=\mathbf{v}^T\mathbf{W}_c\mathbf{v} =\int_{t_0}^{t_f}\left\|\mathbf{B}^T(\tau) \boldsymbol{\Phi}^T(t_f,\tau)\mathbf{v}\right\|_2^2\,d\tau. \]

Hence \( \mathbf{B}^T(\tau)\boldsymbol{\Phi}^T(t_f,\tau)\mathbf{v}=\mathbf{0} \) almost everywhere on the interval, which implies \( \mathbf{v}^T\mathcal{L}_{t_0,t_f}[\mathbf{u}]=0 \) for every input. Therefore no input can generate a terminal displacement with a nonzero component along \( \mathbf{v} \), so the system is not completely controllable.

If \( \mathbf{A}(t)=\mathbf{A} \) and \( \mathbf{B}(t)=\mathbf{B} \) are constant, then

\[ \boldsymbol{\Phi}(t_f,\tau)=e^{\mathbf{A}(t_f-\tau)}. \]

The LTV integral criterion becomes the usual finite-horizon LTI controllability Gramian:

\[ \mathbf{W}_c(t_0,t_f)=\int_{t_0}^{t_f} e^{\mathbf{A}(t_f-\tau)}\mathbf{B}\mathbf{B}^T e^{\mathbf{A}^T(t_f-\tau)}\,d\tau. \]

For finite-dimensional LTI systems, this matrix is nonsingular on any nonzero interval exactly when \( \operatorname{rank}[\mathbf{B},\mathbf{A}\mathbf{B},\dots,\mathbf{A}^{n-1}\mathbf{B}]=n \). Thus, the LTV result is not a contradiction of the Kalman test; it is the interval-based generalization from which the LTI result follows as a special case.

A useful interpretation is that the columns of \( \boldsymbol{\Phi}(t_f,t)\mathbf{B}(t) \) are the time-propagated input directions. In the LTI case, repeated multiplication by \( \mathbf{A} \) captures the same propagation algebraically through the columns \( \mathbf{B},\mathbf{A}\mathbf{B},\dots \).

Instantaneous rank is not enough. If \( m<n \), then \( \mathbf{B}(t) \) has rank at most \( m \) at each instant, but the system may still be controllable because the propagated directions over time can span all of \( \mathbb{R}^n \).

Longer intervals usually help but do not automatically solve all problems. Extending \( [t_0,t_f] \) gives more opportunities for input directions to be propagated into independent directions. However, if the propagated input directions remain confined to a strict subspace, the rank remains deficient.

Time variation can create controllability. Consider an input vector that slowly rotates in the plane. At each time there is only one scalar actuator, but over an interval the actuator points in multiple directions. The reachable subspace can become two-dimensional.

Time variation can also destroy convenient LTI intuition. There is no single characteristic polynomial, no fixed modal basis, and no single Kalman matrix built only from constant powers of \( \mathbf{A} \). The correct object is the interval input-to-state map.

Direct quadrature of \( \boldsymbol{\Phi}(t_f,t)\mathbf{B}(t)\mathbf{B}^T(t)\boldsymbol{\Phi}^T(t_f,t) \) is possible, but a convenient forward-time method integrates the differential Lyapunov equation

\[ \dot{\mathbf{W} }(t)=\mathbf{A}(t)\mathbf{W}(t) +\mathbf{W}(t)\mathbf{A}^T(t)+\mathbf{B}(t)\mathbf{B}^T(t), \quad \mathbf{W}(t_0)=\mathbf{0}. \]

Then \( \mathbf{W}(t_f)=\mathbf{W}_c(t_0,t_f) \). In computation, controllability is usually judged by the smallest eigenvalue or singular value of \( \mathbf{W}_c \), not by exact symbolic rank:

\[ \lambda_{\min}(\mathbf{W}_c)>\varepsilon \quad \text{for a chosen numerical tolerance }\varepsilon>0. \]

A very small positive eigenvalue means the system may be mathematically controllable but numerically hard to steer in that direction. This is a bridge to the energy viewpoint developed in the next chapter.

Chapter11_Lesson5.py integrates the state transition matrix and controllability Gramian for a two-state LTV system using NumPy and SciPy. The code checks positive definiteness by the eigenvalues of \( \mathbf{W}_c \).

"""
Chapter11_Lesson5.py
Time-varying controllability for a 2-state LTV system.

This script integrates the state transition matrix Phi(t,t0) and the
finite-horizon controllability Gramian Wc(t0,t) using the differential
Lyapunov equation

    dW/dt = A(t) W + W A(t)^T + B(t) B(t)^T,  W(t0)=0.

A positive definite Wc on [t0,tf] indicates complete controllability on
that interval for this continuous-time LTV system.
"""

import numpy as np
from scipy.integrate import solve_ivp


def A(t: float) -> np.ndarray:
    """Time-varying dynamics matrix."""
    return np.array([
        [0.0, 1.0],
        [-(2.0 + 0.60 * np.sin(t)), -0.25],
    ])


def B(t: float) -> np.ndarray:
    """Time-varying input matrix."""
    return np.array([
        [0.0],
        [1.0 + 0.45 * np.cos(t)],
    ])


def augmented_ode(t: float, z: np.ndarray, n: int) -> np.ndarray:
    """ODE for Phi(t,t0) and Wc(t0,t)."""
    Phi = z[: n * n].reshape(n, n)
    W = z[n * n :].reshape(n, n)
    At = A(t)
    Bt = B(t)

    dPhi = At @ Phi
    dW = At @ W + W @ At.T + Bt @ Bt.T
    return np.concatenate([dPhi.reshape(-1), dW.reshape(-1)])


def compute_gramian(t0: float, tf: float, n: int = 2):
    Phi0 = np.eye(n)
    W0 = np.zeros((n, n))
    z0 = np.concatenate([Phi0.reshape(-1), W0.reshape(-1)])

    sol = solve_ivp(
        fun=lambda t, z: augmented_ode(t, z, n),
        t_span=(t0, tf),
        y0=z0,
        rtol=1e-9,
        atol=1e-11,
        dense_output=False,
    )
    if not sol.success:
        raise RuntimeError(sol.message)

    zf = sol.y[:, -1]
    Phi_tf_t0 = zf[: n * n].reshape(n, n)
    Wc = zf[n * n :].reshape(n, n)
    Wc = 0.5 * (Wc + Wc.T)  # remove round-off asymmetry
    return Phi_tf_t0, Wc


def is_controllable(Wc: np.ndarray, tol: float = 1e-8) -> bool:
    eigvals = np.linalg.eigvalsh(Wc)
    return bool(np.min(eigvals) > tol)


if __name__ == "__main__":
    for tf in [0.25, 0.50, 1.00, 2.00, 4.00]:
        Phi, Wc = compute_gramian(0.0, tf)
        eigvals = np.linalg.eigvalsh(Wc)
        print(f"Interval [0, {tf:.2f}]")
        print("Phi(tf,t0)=\n", Phi)
        print("Wc(t0,tf)=\n", Wc)
        print("Eigenvalues of Wc:", eigvals)
        print("Controllable on this interval?", is_controllable(Wc))
        print("-" * 60)

Chapter11_Lesson5.cpp implements the same computation from scratch with a fourth-order Runge--Kutta method for the eight scalar ODEs corresponding to \( \boldsymbol{\Phi} \) and \( \mathbf{W}_c \).

/*
Chapter11_Lesson5.cpp
Small fixed-size C++ implementation for an LTV controllability Gramian.
Compile for example:
    g++ -std=c++17 -O2 Chapter11_Lesson5.cpp -o Chapter11_Lesson5
*/

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

using Mat2 = std::array<std::array<double, 2>, 2>;
using Vec8 = std::array<double, 8>;

Mat2 A(double t) {
    return { { {0.0, 1.0}, {-(2.0 + 0.60 * std::sin(t)), -0.25} } };
}

std::array<double, 2> B(double t) {
    return { {0.0, 1.0 + 0.45 * std::cos(t)} };
}

Mat2 add(const Mat2& X, const Mat2& Y) {
    Mat2 Z{};
    for (int i = 0; i < 2; ++i)
        for (int j = 0; j < 2; ++j)
            Z[i][j] = X[i][j] + Y[i][j];
    return Z;
}

Mat2 mul(const Mat2& X, const Mat2& Y) {
    Mat2 Z{};
    for (int i = 0; i < 2; ++i)
        for (int j = 0; j < 2; ++j)
            for (int k = 0; k < 2; ++k)
                Z[i][j] += X[i][k] * Y[k][j];
    return Z;
}

Mat2 transpose(const Mat2& X) {
    return { { {X[0][0], X[1][0]}, {X[0][1], X[1][1]} } };
}

Mat2 outerBB(double t) {
    auto b = B(t);
    return { { {b[0] * b[0], b[0] * b[1]}, {b[1] * b[0], b[1] * b[1]} } };
}

Vec8 derivative(double t, const Vec8& z) {
    Mat2 Phi{ { {z[0], z[1]}, {z[2], z[3]} } };
    Mat2 W{ { {z[4], z[5]}, {z[6], z[7]} } };
    Mat2 At = A(t);
    Mat2 dPhi = mul(At, Phi);
    Mat2 dW = add(add(mul(At, W), mul(W, transpose(At))), outerBB(t));
    return { {dPhi[0][0], dPhi[0][1], dPhi[1][0], dPhi[1][1],
             dW[0][0], dW[0][1], dW[1][0], dW[1][1]} };
}

Vec8 rk4_step(double t, const Vec8& z, double h) {
    Vec8 k1 = derivative(t, z), z2{}, z3{}, z4{};
    for (int i = 0; i < 8; ++i) z2[i] = z[i] + 0.5 * h * k1[i];
    Vec8 k2 = derivative(t + 0.5 * h, z2);
    for (int i = 0; i < 8; ++i) z3[i] = z[i] + 0.5 * h * k2[i];
    Vec8 k3 = derivative(t + 0.5 * h, z3);
    for (int i = 0; i < 8; ++i) z4[i] = z[i] + h * k3[i];
    Vec8 k4 = derivative(t + h, z4);

    Vec8 zn{};
    for (int i = 0; i < 8; ++i)
        zn[i] = z[i] + h * (k1[i] + 2.0 * k2[i] + 2.0 * k3[i] + k4[i]) / 6.0;
    return zn;
}

int main() {
    const double t0 = 0.0, tf = 4.0;
    const int steps = 4000;
    const double h = (tf - t0) / steps;

    Vec8 z = { {1.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0} };
    double t = t0;
    for (int k = 0; k < steps; ++k) {
        z = rk4_step(t, z, h);
        t += h;
    }

    Mat2 W{ { {z[4], 0.5 * (z[5] + z[6])}, {0.5 * (z[5] + z[6]), z[7]} } };
    double trace = W[0][0] + W[1][1];
    double det = W[0][0] * W[1][1] - W[0][1] * W[1][0];
    double disc = std::sqrt(std::max(0.0, trace * trace - 4.0 * det));
    double lambda_min = 0.5 * (trace - disc);

    std::cout << std::setprecision(10);
    std::cout << "Wc(0,4) =\n";
    std::cout << W[0][0] << "  " << W[0][1] << "\n";
    std::cout << W[1][0] << "  " << W[1][1] << "\n";
    std::cout << "det(Wc) = " << det << "\n";
    std::cout << "lambda_min(Wc) = " << lambda_min << "\n";
    std::cout << "Controllable? " << (lambda_min > 1e-8 ? "yes" : "no") << "\n";
    return 0;
}

Chapter11_Lesson5.java mirrors the C++ implementation and is useful when students need a dependency-free implementation of the LTV Gramian test.

/*
Chapter11_Lesson5.java
Java RK4 computation of the finite-horizon LTV controllability Gramian.
Compile and run:
    javac Chapter11_Lesson5.java
    java Chapter11_Lesson5
*/

public class Chapter11_Lesson5 {
    static double[][] A(double t) {
        return new double[][]{
            {0.0, 1.0},
            {-(2.0 + 0.60 * Math.sin(t)), -0.25}
        };
    }

    static double[] B(double t) {
        return new double[]{0.0, 1.0 + 0.45 * Math.cos(t)};
    }

    static double[][] mul(double[][] X, double[][] Y) {
        double[][] Z = 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++)
                    Z[i][j] += X[i][k] * Y[k][j];
        return Z;
    }

    static double[][] transpose(double[][] X) {
        return new double[][]{ {X[0][0], X[1][0]}, {X[0][1], X[1][1]} };
    }

    static double[] derivative(double t, double[] z) {
        double[][] Phi = { {z[0], z[1]}, {z[2], z[3]} };
        double[][] W = { {z[4], z[5]}, {z[6], z[7]} };
        double[][] At = A(t);
        double[] Bt = B(t);
        double[][] dPhi = mul(At, Phi);
        double[][] AW = mul(At, W);
        double[][] WAT = mul(W, transpose(At));

        double[][] dW = new double[2][2];
        for (int i = 0; i < 2; i++) {
            for (int j = 0; j < 2; j++) {
                dW[i][j] = AW[i][j] + WAT[i][j] + Bt[i] * Bt[j];
            }
        }
        return new double[]{dPhi[0][0], dPhi[0][1], dPhi[1][0], dPhi[1][1],
                            dW[0][0], dW[0][1], dW[1][0], dW[1][1]};
    }

    static double[] rk4Step(double t, double[] z, double h) {
        double[] k1 = derivative(t, z);
        double[] z2 = new double[8], z3 = new double[8], z4 = new double[8];
        for (int i = 0; i < 8; i++) z2[i] = z[i] + 0.5 * h * k1[i];
        double[] k2 = derivative(t + 0.5 * h, z2);
        for (int i = 0; i < 8; i++) z3[i] = z[i] + 0.5 * h * k2[i];
        double[] k3 = derivative(t + 0.5 * h, z3);
        for (int i = 0; i < 8; i++) z4[i] = z[i] + h * k3[i];
        double[] k4 = derivative(t + h, z4);

        double[] zn = new double[8];
        for (int i = 0; i < 8; i++)
            zn[i] = z[i] + h * (k1[i] + 2.0 * k2[i] + 2.0 * k3[i] + k4[i]) / 6.0;
        return zn;
    }

    public static void main(String[] args) {
        double t0 = 0.0, tf = 4.0;
        int steps = 4000;
        double h = (tf - t0) / steps;
        double[] z = {1.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 0.0};
        double t = t0;
        for (int k = 0; k < steps; k++) {
            z = rk4Step(t, z, h);
            t += h;
        }

        double w00 = z[4];
        double w01 = 0.5 * (z[5] + z[6]);
        double w11 = z[7];
        double trace = w00 + w11;
        double det = w00 * w11 - w01 * w01;
        double disc = Math.sqrt(Math.max(0.0, trace * trace - 4.0 * det));
        double lambdaMin = 0.5 * (trace - disc);

        System.out.printf("Wc(0,4) =%n%.10f  %.10f%n%.10f  %.10f%n", w00, w01, w01, w11);
        System.out.printf("det(Wc) = %.10f%n", det);
        System.out.printf("lambda_min(Wc) = %.10f%n", lambdaMin);
        System.out.println("Controllable? " + (lambdaMin > 1e-8 ? "yes" : "no"));
    }
}

Chapter11_Lesson5.m uses ode45 to integrate the same differential Lyapunov equation. In Simulink, the right-hand side can be implemented by a MATLAB Function block that outputs \( \mathbf{A}(t) \), \( \mathbf{B}(t) \), and \( \dot{\mathbf{W} }(t) \), followed by Integrator blocks for the vectorized entries of \( \mathbf{W} \).

% Chapter11_Lesson5.m
% LTV controllability Gramian using MATLAB ode45.
% The same ODE can be implemented in Simulink by using a MATLAB Function
% block for A(t), B(t) and an Integrator block for the vectorized W dynamics.

clear; clc;
t0 = 0.0;
tf = 4.0;
n = 2;

Phi0 = eye(n);
W0 = zeros(n);
z0 = [Phi0(:); W0(:)];
opts = odeset('RelTol',1e-9,'AbsTol',1e-11);
[t,z] = ode45(@(t,z) augmented_ode(t,z,n), [t0 tf], z0, opts);

zf = z(end,:).';
Phi = reshape(zf(1:n*n), n, n);
Wc = reshape(zf(n*n+1:end), n, n);
Wc = 0.5*(Wc + Wc.');

disp('Phi(tf,t0) =');
disp(Phi);
disp('Wc(t0,tf) =');
disp(Wc);
disp('eig(Wc) =');
disp(eig(Wc));

if min(eig(Wc)) > 1e-8
    disp('The LTV system is controllable on this interval.');
else
    disp('The test is singular or ill-conditioned on this interval.');
end

function dz = augmented_ode(t,z,n)
    Phi = reshape(z(1:n*n), n, n);
    W = reshape(z(n*n+1:end), n, n);
    At = A_ltv(t);
    Bt = B_ltv(t);
    dPhi = At*Phi;
    dW = At*W + W*At.' + Bt*Bt.';
    dz = [dPhi(:); dW(:)];
end

function At = A_ltv(t)
    At = [0, 1; -(2.0 + 0.60*sin(t)), -0.25];
end

function Bt = B_ltv(t)
    Bt = [0; 1.0 + 0.45*cos(t)];
end

Chapter11_Lesson5.nb solves the Gramian ODE symbolically in matrix form and evaluates the finite-horizon controllability test numerically.

(* Chapter11_Lesson5.nb
   Wolfram Mathematica code for the finite-horizon LTV controllability Gramian. *)

ClearAll[A, B, t, w11, w12, w21, w22, W, t0, tf];

t0 = 0; tf = 4;
A[t_] := { {0, 1}, {-(2 + 0.60 Sin[t]), -0.25} };
B[t_] := { {0}, {1 + 0.45 Cos[t]} };
W[t_] := { {w11[t], w12[t]}, {w21[t], w22[t]} };

eqs = Thread[D[W[t], t] == A[t].W[t] + W[t].Transpose[A[t]] + B[t].Transpose[B[t]]];
ics = {w11[t0] == 0, w12[t0] == 0, w21[t0] == 0, w22[t0] == 0};

sol = NDSolveValue[Join[Flatten[eqs], ics], {w11, w12, w21, w22}, {t, t0, tf},
   WorkingPrecision -> 30, AccuracyGoal -> 15, PrecisionGoal -> 15];

Wc = N[{ {sol[[1]][tf], sol[[2]][tf]}, {sol[[3]][tf], sol[[4]][tf]} }];
WcSym = (Wc + Transpose[Wc])/2;
Eigenvalues[WcSym]
Det[WcSym]
PositiveDefiniteMatrixQ[WcSym]

Problem 1 (Reachable Set for an LTV System): For the LTV system \( \dot{\mathbf{x} }=\mathbf{A}(t)\mathbf{x}+\mathbf{B}(t)\mathbf{u} \), derive the reachable set at \( t_f \) from the initial condition \( \mathbf{x}(t_0)=\mathbf{x}_0 \).

Solution: Variation of constants gives

\[ \mathbf{x}(t_f)=\boldsymbol{\Phi}(t_f,t_0)\mathbf{x}_0 +\int_{t_0}^{t_f}\boldsymbol{\Phi}(t_f,\tau)\mathbf{B}(\tau) \mathbf{u}(\tau)\,d\tau. \]

Therefore the reachable set is the affine set

\[ \mathcal{R}(t_f;t_0,\mathbf{x}_0) =\boldsymbol{\Phi}(t_f,t_0)\mathbf{x}_0 +\operatorname{Range}(\mathcal{L}_{t_0,t_f}). \]

Problem 2 (Positive Semidefiniteness): Prove that \( \mathbf{W}_c(t_0,t_f) \) is symmetric positive semidefinite.

Solution: Symmetry follows from the integrand

\[ \boldsymbol{\Phi}\mathbf{B}\mathbf{B}^T\boldsymbol{\Phi}^T =\left(\boldsymbol{\Phi}\mathbf{B}\right) \left(\boldsymbol{\Phi}\mathbf{B}\right)^T. \]

For any \( \mathbf{v}\in\mathbb{R}^n \),

\[ \mathbf{v}^T\mathbf{W}_c\mathbf{v} =\int_{t_0}^{t_f}\left\|\mathbf{B}^T(\tau) \boldsymbol{\Phi}^T(t_f,\tau)\mathbf{v}\right\|_2^2\,d\tau\geq 0. \]

Problem 3 (Constructing a Steering Input): Assume \( \mathbf{W}_c(t_0,t_f) \) is nonsingular. Find an input that steers \( \mathbf{x}_0 \) to \( \mathbf{x}_f \).

Solution: Let

\[ \mathbf{r}=\mathbf{x}_f-\boldsymbol{\Phi}(t_f,t_0)\mathbf{x}_0. \]

One valid input is

\[ \mathbf{u}(\tau)=\mathbf{B}^T(\tau) \boldsymbol{\Phi}^T(t_f,\tau)\mathbf{W}_c^{-1}(t_0,t_f)\mathbf{r}. \]

Substitution into the forced-response integral yields \( \mathcal{L}_{t_0,t_f}[\mathbf{u}]=\mathbf{r} \), so \( \mathbf{x}(t_f)=\mathbf{x}_f \).

Problem 4 (A Scalar LTV System): Consider \( \dot{x}=a(t)x+b(t)u \) on \( [t_0,t_f] \). State the controllability condition.

Solution: Since \( n=1 \), the Gramian is

\[ W_c(t_0,t_f)=\int_{t_0}^{t_f}\Phi^2(t_f,\tau)b^2(\tau)\,d\tau. \]

Since \( \Phi(t_f,\tau)\neq0 \), the system is controllable exactly when \( b(t) \) is not zero almost everywhere on the interval.

Problem 5 (LTI Limit): Show why the finite-interval LTV test reduces to the familiar LTI Gramian when \( \mathbf{A}(t)=\mathbf{A} \) and \( \mathbf{B}(t)=\mathbf{B} \).

Solution: In the constant-coefficient case,

\[ \boldsymbol{\Phi}(t_f,\tau)=e^{\mathbf{A}(t_f-\tau)}. \]

Substitution into the LTV Gramian gives

\[ \mathbf{W}_c(t_0,t_f)=\int_{t_0}^{t_f} e^{\mathbf{A}(t_f-\tau)}\mathbf{B}\mathbf{B}^T e^{\mathbf{A}^T(t_f-\tau)}\,d\tau, \]

which is the LTI finite-horizon controllability Gramian. Its nonsingularity is equivalent to the Kalman rank condition.

Time-varying controllability is a finite-interval property. Instead of checking only the algebraic span of \( \mathbf{B},\mathbf{A}\mathbf{B},\dots \), we examine the span of all propagated input directions \( \boldsymbol{\Phi}(t_f,t)\mathbf{B}(t) \) over \( [t_0,t_f] \). The finite-interval Gramian \( \mathbf{W}_c(t_0,t_f) \) is full rank exactly when every terminal displacement can be generated. This lesson prepares the transition to Chapter 12, where the same matrix will be interpreted as a control-energy operator.