Chapter 14: Observability Tests and Duality

Lesson 4: Observability in Canonical (Observable) Forms

This lesson studies observable canonical forms as structured realizations in which observability is built directly into the coordinates. We connect the Kalman observability matrix, PBH reasoning, output-derivative state reconstruction, and duality with controllable canonical form. The emphasis is theoretical: canonical structure is not merely a convenient notation; it exposes why the output channel carries enough independent information about every internal state.

1. Motivation: Why an Observable Canonical Form?

In previous lessons, observability was tested by the Kalman matrix and the PBH rank condition. Canonical forms provide an opposite viewpoint: instead of testing a general pair \( (A,C) \), we choose coordinates in which the output equation and state matrix reveal observability almost by inspection. For a continuous-time SISO LTI system,

\[ \dot{\mathbf{x}}(t)=A\mathbf{x}(t)+B u(t), \qquad y(t)=C\mathbf{x}(t)+D u(t), \]

the zero-input output derivatives satisfy \( y^{(k)}(0)=C A^k\mathbf{x}(0) \). Therefore the vector of the first \( n \) output derivatives is

\[ \boldsymbol{\eta}(0)= \begin{bmatrix} y(0)\\ \dot{y}(0)\\ \vdots\\ y^{(n-1)}(0) \end{bmatrix} =\underbrace{\begin{bmatrix} C\\ CA\\ \vdots\\ CA^{n-1} \end{bmatrix}}_{\mathbf{\mathcal{O}}(A,C)}\mathbf{x}(0). \]

Hence \( \mathbf{x}(0) \) is uniquely reconstructible from ideal output-derivative data exactly when \( \operatorname{rank}\mathbf{\mathcal{O}}(A,C)=n \). An observable canonical form makes this matrix nonsingular by structure, not by numerical accident.

flowchart TD
  A["Start with transfer function or observable pair (A,C)"] --> B["Choose observable coordinates"]
  B --> C["Build structured pair (Ao, Co)"]
  C --> D["Form O = [Co; Co Ao; ...; Co Ao^(n-1)]"]
  D --> E["O has full rank by triangular structure"]
  E --> F["Initial state is reconstructible from output history"]
        

2. Transfer-Function Observable Companion Form

Consider a strictly proper SISO transfer function with monic denominator

\[ G(s)=\frac{b_{n-1}s^{n-1}+b_{n-2}s^{n-2}+\cdots+b_1s+b_0} {s^n+a_{n-1}s^{n-1}+a_{n-2}s^{n-2}+\cdots+a_1s+a_0}. \]

One common observable companion realization is

\[ A_o=\begin{bmatrix} -a_{n-1} & 1 & 0 & \cdots & 0\\ -a_{n-2} & 0 & 1 & \cdots & 0\\ \vdots & \vdots & \vdots & \ddots & \vdots\\ -a_1 & 0 & 0 & \cdots & 1\\ -a_0 & 0 & 0 & \cdots & 0 \end{bmatrix},\quad B_o=\begin{bmatrix} b_{n-1}\\ b_{n-2}\\ \vdots\\ b_1\\ b_0 \end{bmatrix},\quad C_o=\begin{bmatrix}1&0&\cdots&0\end{bmatrix}. \]

The first state is measured directly, while the remaining states appear through successive dynamic propagation. The numerator coefficients enter \( B_o \), while the denominator coefficients appear in the first column of \( A_o \). This convention is the transpose-dual of a standard controllable companion convention.

To verify that this realization produces the stated transfer function, use the resolvent identity \( G(s)=C_o(sI-A_o)^{-1}B_o \). For this companion structure, the first row of \( (sI-A_o)^{-1} \) has polynomial numerator pattern

\[ C_o(sI-A_o)^{-1} =\frac{1}{p(s)}\begin{bmatrix}s^{n-1}&s^{n-2}&\cdots&s&1 \end{bmatrix}, \]

where \( p(s)=s^n+a_{n-1}s^{n-1}+\cdots+a_0 \). Multiplication by \( B_o \) gives exactly the numerator polynomial.

3. Direct Observability Proof by Triangular Structure

The observability matrix of \( (A_o,C_o) \) is

\[ \mathbf{\mathcal{O}}_o=\begin{bmatrix} C_o\\ C_oA_o\\ C_oA_o^2\\ \vdots\\ C_oA_o^{n-1} \end{bmatrix}. \]

Because \( A_o \) has ones on the superdiagonal, each multiplication by \( A_o \) shifts one new independent coordinate into the next column. More precisely, the \( k \)-th row block \( C_oA_o^k \) has a unit coefficient in column \( k+1 \) and zero coefficients in all columns to its right. Thus

\[ \mathbf{\mathcal{O}}_o= \begin{bmatrix} 1&0&0&\cdots&0\\ \star&1&0&\cdots&0\\ \star&\star&1&\cdots&0\\ \vdots&\vdots&\vdots&\ddots&\vdots\\ \star&\star&\star&\cdots&1 \end{bmatrix},\qquad \det(\mathbf{\mathcal{O}}_o)=1. \]

Therefore every observable companion realization of this form is observable for every choice of denominator coefficients. The coefficients \( a_i \) affect the eigenvalues and transient modal shape, but they do not destroy the full-rank observation chain.

Proof sketch. For \( k=0 \), \( C_o=[1,0,\ldots,0] \). Suppose \( C_oA_o^k \) has zero entries to the right of column \( k+1 \) and has unit entry at column \( k+1 \). Multiplication by \( A_o \) shifts this unit entry through the superdiagonal into column \( k+2 \), while the first column may collect terms involving the coefficients \( a_i \). Hence the next row keeps the triangular pattern. By induction, the diagonal entries of \( \mathbf{\mathcal{O}}_o \) are all one.

4. Output-Derivative Reconstruction in Observable Coordinates

In zero-input operation, the measured output and its first \( n-1 \) derivatives determine the initial state by

\[ \mathbf{x}(0)=\mathbf{\mathcal{O}}_o^{-1} \begin{bmatrix}y(0)\\\dot{y}(0)\\\vdots\\y^{(n-1)}(0) \end{bmatrix}. \]

In practice, direct differentiation of noisy measurements is usually avoided. However, the mathematical reconstruction formula is fundamental because it shows what information is present in the output signal. Later observer design replaces raw differentiation by a stable dynamical estimator, but the existence of that estimator rests on the same observability property.

For \( n=3 \), the canonical pair is

\[ A_o=\begin{bmatrix}-a_2&1&0\\-a_1&0&1\\-a_0&0&0\end{bmatrix}, \qquad C_o=\begin{bmatrix}1&0&0\end{bmatrix}. \]

The observability matrix becomes

\[ \mathbf{\mathcal{O}}_o=\begin{bmatrix} 1&0&0\\ -a_2&1&0\\ a_2^2-a_1&-a_2&1 \end{bmatrix},\qquad \det(\mathbf{\mathcal{O}}_o)=1. \]

Thus the state is recovered by solving a lower-triangular system. The formula also clarifies why canonical coordinates are useful for hand-calculation: the output reveals states sequentially rather than all at once.

5. Coordinate Transformation from a General Observable Pair

Suppose \( (A,C) \) is observable and define the observability matrix \( T=\mathbf{\mathcal{O}}(A,C) \). Since \( T \) is nonsingular, one may introduce observation coordinates

\[ \mathbf{z}=T\mathbf{x}=\begin{bmatrix} C\mathbf{x}\\ CA\mathbf{x}\\ \vdots\\ CA^{n-1}\mathbf{x} \end{bmatrix}. \]

In these coordinates, the output is simply \( y=\begin{bmatrix}1&0&\cdots&0\end{bmatrix}\mathbf{z} \) for zero direct feedthrough. The transformed state matrix is

\[ A_z=TAT^{-1},\qquad C_z=CT^{-1}=\begin{bmatrix}1&0&\cdots&0\end{bmatrix}. \]

The rows of \( T \) form a Krylov sequence generated by repeated multiplication by \( A \). Cayley-Hamilton closes this sequence after \( n \) rows. If \( p_A(s)=s^n+\alpha_{n-1}s^{n-1}+\cdots+ \alpha_0 \), then

\[ A^n+\alpha_{n-1}A^{n-1}+\cdots+\alpha_1A+ \alpha_0I=0. \]

Consequently, the transformed dynamics have companion-like closure. This demonstrates the key point: a fully observable system admits coordinates in which the output and its generated derivatives form a basis for the state space.

6. Duality with Controllable Canonical Form

Observability of \( (A,C) \) is controllability of the dual pair \( (A^T,C^T) \). Therefore, every result about controllable companion structure has a transposed observable statement. If

\[ \mathbf{\mathcal{C}}(A^T,C^T)= \begin{bmatrix}C^T&A^TC^T&\cdots&(A^T)^{n-1}C^T\end{bmatrix}, \]

then

\[ \mathbf{\mathcal{C}}(A^T,C^T) =\mathbf{\mathcal{O}}(A,C)^T. \]

Hence \( \operatorname{rank}\mathbf{\mathcal{O}}(A,C)=n \) if and only if \( \operatorname{rank}\mathbf{\mathcal{C}}(A^T,C^T)=n \). Observable companion form is therefore not a separate theory; it is the dual expression of the controllability canonical machinery already seen in the course.

flowchart TD
  A["Original pair (A, C)"] --> B["Observability matrix O(A,C)"]
  A --> C["Dual pair (A^T, C^T)"]
  C --> D["Controllability matrix Ctrb(A^T,C^T)"]
  B --> E["rank O = n"]
  D --> F["rank Ctrb = n"]
  E <--> F
  F --> G["Controllable companion logic transposes to observable form"]
        

7. PBH Interpretation in Observable Companion Form

The PBH test states that \( (A,C) \) is observable if and only if

\[ \operatorname{rank}\begin{bmatrix}\lambda I-A\\C\end{bmatrix}=n \quad \text{for every eigenvalue }\lambda\text{ of }A. \]

In canonical observable form, no eigenvector of \( A_o \) can lie in \( \ker C_o \). Indeed, if \( A_o\mathbf{v}=\lambda\mathbf{v} \) and \( C_o\mathbf{v}=0 \), then the first component \( v_1=0 \). The eigenvector equations and the superdiagonal chain force successively \( v_2=0,v_3=0,\ldots,v_n=0 \), which contradicts \( \mathbf{v}\neq\mathbf{0} \). Thus PBH also passes structurally.

This argument is conceptually important: observability fails precisely when an internal mode can oscillate, decay, or grow while remaining invisible at the output. The observable companion form prevents such a hidden eigenvector by connecting every state coordinate to the output chain.

8. Numerical Conditioning and Canonical Forms

Although canonical forms are theoretically clean, they are not always numerically ideal. The transformation \( T=\mathbf{\mathcal{O}}(A,C) \) may be poorly conditioned when output measurements distinguish some modes only weakly. The condition number

\[ \kappa(T)=\|T\|\,\|T^{-1}\| \]

measures sensitivity of coordinate conversion. If \( \kappa(T) \) is large, small errors in output data, coefficients, or floating-point operations can create large errors in reconstructed states. Therefore canonical forms should be understood as structural models and symbolic analysis tools; numerically robust implementations often use orthogonal transformations, SVD-based rank tests, or balanced coordinates.

9. Python Implementation: Constructing and Testing OCF

Chapter14_Lesson4.py


# Chapter14_Lesson4.py
# Observability in canonical (observable) forms
# Requirements: numpy, scipy (optional for eigenvalues)

import numpy as np


def observable_companion(den_coeffs, num_coeffs=None):
    """
    Build the SISO observable companion form for

        G(s) = (b_{n-1}s^{n-1} + ... + b_0) /
               (s^n + a_{n-1}s^{n-1} + ... + a_0)

    with convention:
        A_o = [[-a_{n-1}, 1, 0, ..., 0],
               [-a_{n-2}, 0, 1, ..., 0],
               ...,
               [-a_0,     0, 0, ..., 0]],
        C_o = [1, 0, ..., 0].

    den_coeffs must be [a0, a1, ..., a_{n-1}].
    num_coeffs, if supplied, must be [b0, b1, ..., b_{n-1}].
    The returned B_o is [b_{n-1}, b_{n-2}, ..., b0]^T.
    """
    den = np.asarray(den_coeffs, dtype=float).reshape(-1)
    n = den.size
    A = np.zeros((n, n), dtype=float)
    A[:, 0] = -den[::-1]
    for i in range(n - 1):
        A[i, i + 1] = 1.0
    C = np.zeros((1, n), dtype=float)
    C[0, 0] = 1.0
    if num_coeffs is None:
        B = np.zeros((n, 1), dtype=float)
    else:
        num = np.asarray(num_coeffs, dtype=float).reshape(-1)
        if num.size != n:
            raise ValueError("num_coeffs must have length n: [b0, ..., b_{n-1}]")
        B = num[::-1].reshape(n, 1)
    return A, B, C


def observability_matrix(A, C):
    """Return O = [C; C A; ...; C A^{n-1}]."""
    A = np.asarray(A, dtype=float)
    C = np.asarray(C, dtype=float)
    n = A.shape[0]
    blocks = []
    Ak = np.eye(n)
    for _ in range(n):
        blocks.append(C @ Ak)
        Ak = Ak @ A
    return np.vstack(blocks)


def matrix_rank(M, tol=1e-10):
    """Numerical rank computed from singular values."""
    s = np.linalg.svd(M, compute_uv=False)
    return int(np.sum(s > tol))


def pbh_observability_test(A, C, tol=1e-9):
    """
    PBH test for observability:
        rank([lambda I - A; C]) = n for every eigenvalue lambda of A.
    """
    A = np.asarray(A, dtype=complex)
    C = np.asarray(C, dtype=complex)
    n = A.shape[0]
    eigvals = np.linalg.eigvals(A)
    results = []
    for lam in eigvals:
        pbh = np.vstack((lam * np.eye(n) - A, C))
        r = matrix_rank(pbh, tol=tol)
        results.append((lam, r, r == n))
    return results


def reconstruct_initial_state_from_derivatives(A, C, y_derivatives_at_0):
    """
    For zero input, y^{(k)}(0) = C A^k x(0).
    Given [y(0), ydot(0), ..., y^{(n-1)}(0)], solve O x0 = eta.
    """
    O = observability_matrix(A, C)
    eta = np.asarray(y_derivatives_at_0, dtype=float).reshape(-1, 1)
    if O.shape[0] != eta.shape[0]:
        raise ValueError("Need exactly n output derivatives for an n-state SISO system.")
    return np.linalg.solve(O, eta)


if __name__ == "__main__":
    # Example: G(s) = (2 s^2 + 3 s + 4)/(s^3 + 6 s^2 + 11 s + 6)
    # denominator coefficients are [a0, a1, a2]
    den = [6.0, 11.0, 6.0]
    num = [4.0, 3.0, 2.0]

    A, B, C = observable_companion(den, num)
    O = observability_matrix(A, C)

    print("A_o =\n", A)
    print("B_o =\n", B)
    print("C_o =\n", C)
    print("Observability matrix O =\n", O)
    print("det(O) =", np.linalg.det(O))
    print("rank(O) =", matrix_rank(O), "out of", A.shape[0])

    print("\nPBH observability test:")
    for lam, r, passed in pbh_observability_test(A, C):
        print(f"lambda={lam:.6g}, rank={r}, observable_at_lambda={passed}")

    # Zero-input reconstruction example
    x0_true = np.array([[1.0], [-2.0], [0.5]])
    eta = O @ x0_true
    x0_hat = reconstruct_initial_state_from_derivatives(A, C, eta)
    print("\nTrue x0 =\n", x0_true)
    print("Reconstructed x0 =\n", x0_hat)

      

10. C++ Implementation: Observability Matrix from Scratch

Chapter14_Lesson4.cpp


// Chapter14_Lesson4.cpp
// Observability in canonical (observable) forms from scratch.
// Compile: g++ -std=c++17 Chapter14_Lesson4.cpp -O2 -o Chapter14_Lesson4

#include <cmath>
#include <iomanip>
#include <iostream>
#include <stdexcept>
#include <vector>

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

Matrix zeros(int r, int c) {
    return Matrix(r, std::vector<double>(c, 0.0));
}

Matrix identity(int n) {
    Matrix I = zeros(n, n);
    for (int i = 0; i < n; ++i) I[i][i] = 1.0;
    return I;
}

Matrix multiply(const Matrix& A, const Matrix& B) {
    int r = static_cast<int>(A.size());
    int m = static_cast<int>(A[0].size());
    int c = static_cast<int>(B[0].size());
    Matrix C = zeros(r, c);
    for (int i = 0; i < r; ++i)
        for (int k = 0; k < m; ++k)
            for (int j = 0; j < c; ++j)
                C[i][j] += A[i][k] * B[k][j];
    return C;
}

Matrix vstack(const std::vector<Matrix>& blocks) {
    int totalRows = 0;
    int cols = static_cast<int>(blocks[0][0].size());
    for (const auto& M : blocks) totalRows += static_cast<int>(M.size());
    Matrix R = zeros(totalRows, cols);
    int row = 0;
    for (const auto& M : blocks) {
        for (int i = 0; i < static_cast<int>(M.size()); ++i) {
            R[row++] = M[i];
        }
    }
    return R;
}

int rankGaussian(Matrix A, double tol = 1e-10) {
    int m = static_cast<int>(A.size());
    int n = static_cast<int>(A[0].size());
    int rank = 0;
    for (int col = 0; col < n && rank < m; ++col) {
        int pivot = rank;
        for (int i = rank + 1; i < m; ++i) {
            if (std::fabs(A[i][col]) > std::fabs(A[pivot][col])) pivot = i;
        }
        if (std::fabs(A[pivot][col]) <= tol) continue;
        std::swap(A[pivot], A[rank]);
        double div = A[rank][col];
        for (int j = col; j < n; ++j) A[rank][j] /= div;
        for (int i = 0; i < m; ++i) {
            if (i == rank) continue;
            double factor = A[i][col];
            for (int j = col; j < n; ++j) A[i][j] -= factor * A[rank][j];
        }
        ++rank;
    }
    return rank;
}

Matrix observableCompanion(const std::vector<double>& den) {
    // den = [a0, a1, ..., a_{n-1}]
    int n = static_cast<int>(den.size());
    Matrix A = zeros(n, n);
    for (int i = 0; i < n; ++i) A[i][0] = -den[n - 1 - i];
    for (int i = 0; i < n - 1; ++i) A[i][i + 1] = 1.0;
    return A;
}

Matrix observabilityMatrix(const Matrix& A, const Matrix& C) {
    int n = static_cast<int>(A.size());
    std::vector<Matrix> rows;
    Matrix Ak = identity(n);
    for (int k = 0; k < n; ++k) {
        rows.push_back(multiply(C, Ak));
        Ak = multiply(Ak, A);
    }
    return vstack(rows);
}

void printMatrix(const std::string& name, const Matrix& M) {
    std::cout << name << " =\n";
    for (const auto& row : M) {
        for (double v : row) std::cout << std::setw(12) << v << " ";
        std::cout << "\n";
    }
}

int main() {
    // Denominator: s^3 + 6s^2 + 11s + 6
    std::vector<double> den = {6.0, 11.0, 6.0};
    Matrix A = observableCompanion(den);
    Matrix C = {{1.0, 0.0, 0.0}};
    Matrix O = observabilityMatrix(A, C);

    printMatrix("A_o", A);
    printMatrix("C_o", C);
    printMatrix("O", O);
    std::cout << "rank(O) = " << rankGaussian(O) << " out of " << A.size() << "\n";

    return 0;
}

      

11. Java Implementation: Rank Test for OCF

Chapter14_Lesson4.java


// Chapter14_Lesson4.java
// Observability in canonical (observable) forms from scratch.
// Compile: javac Chapter14_Lesson4.java
// Run:     java Chapter14_Lesson4

public class Chapter14_Lesson4 {
    static double[][] zeros(int r, int c) {
        return new double[r][c];
    }

    static double[][] identity(int n) {
        double[][] I = zeros(n, n);
        for (int i = 0; i < n; i++) I[i][i] = 1.0;
        return I;
    }

    static double[][] multiply(double[][] A, double[][] B) {
        int r = A.length;
        int m = A[0].length;
        int c = B[0].length;
        double[][] C = zeros(r, c);
        for (int i = 0; i < r; i++) {
            for (int k = 0; k < m; k++) {
                for (int j = 0; j < c; j++) {
                    C[i][j] += A[i][k] * B[k][j];
                }
            }
        }
        return C;
    }

    static double[][] observableCompanion(double[] den) {
        // den = [a0, a1, ..., a_{n-1}]
        int n = den.length;
        double[][] A = zeros(n, n);
        for (int i = 0; i < n; i++) A[i][0] = -den[n - 1 - i];
        for (int i = 0; i < n - 1; i++) A[i][i + 1] = 1.0;
        return A;
    }

    static double[][] observabilityMatrix(double[][] A, double[][] C) {
        int n = A.length;
        double[][] O = zeros(n, n);
        double[][] Ak = identity(n);
        for (int block = 0; block < n; block++) {
            double[][] row = multiply(C, Ak);
            for (int j = 0; j < n; j++) O[block][j] = row[0][j];
            Ak = multiply(Ak, A);
        }
        return O;
    }

    static int rank(double[][] input, double tol) {
        int m = input.length;
        int n = input[0].length;
        double[][] A = new double[m][n];
        for (int i = 0; i < m; i++) System.arraycopy(input[i], 0, A[i], 0, n);

        int rank = 0;
        for (int col = 0; col < n && rank < m; col++) {
            int pivot = rank;
            for (int i = rank + 1; i < m; i++) {
                if (Math.abs(A[i][col]) > Math.abs(A[pivot][col])) pivot = i;
            }
            if (Math.abs(A[pivot][col]) <= tol) continue;
            double[] temp = A[pivot];
            A[pivot] = A[rank];
            A[rank] = temp;

            double div = A[rank][col];
            for (int j = col; j < n; j++) A[rank][j] /= div;
            for (int i = 0; i < m; i++) {
                if (i == rank) continue;
                double factor = A[i][col];
                for (int j = col; j < n; j++) A[i][j] -= factor * A[rank][j];
            }
            rank++;
        }
        return rank;
    }

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

    public static void main(String[] args) {
        // Denominator: s^3 + 6s^2 + 11s + 6
        double[] den = {6.0, 11.0, 6.0};
        double[][] A = observableCompanion(den);
        double[][] C = {{1.0, 0.0, 0.0}};
        double[][] O = observabilityMatrix(A, C);

        printMatrix("A_o", A);
        printMatrix("C_o", C);
        printMatrix("O", O);
        System.out.println("rank(O) = " + rank(O, 1e-10) + " out of " + A.length);
    }
}

      

12. MATLAB/Simulink Implementation

MATLAB users with Control System Toolbox can compare the custom observability matrix with obsv(A,C). In Simulink, the same \( A_o,B_o,C_o,D \) matrices can be inserted into a State-Space block. Set the block matrices to Ao, Bo, Co, and D=0, drive the input with a test signal, and log the output to compare the canonical realization with the transfer-function model.

Chapter14_Lesson4.m


% Chapter14_Lesson4.m
% Observability in canonical (observable) forms.
% Run in MATLAB or GNU Octave.

clear; clc;

% Example transfer function:
% G(s) = (2 s^2 + 3 s + 4)/(s^3 + 6 s^2 + 11 s + 6)
% den = [a0 a1 ... a_{n-1}], num = [b0 b1 ... b_{n-1}]
den = [6 11 6];
num = [4 3 2];

[Ao, Bo, Co] = observable_companion(den, num);
O = observability_matrix(Ao, Co);

fprintf('A_o =\n'); disp(Ao);
fprintf('B_o =\n'); disp(Bo);
fprintf('C_o =\n'); disp(Co);
fprintf('Observability matrix O =\n'); disp(O);
fprintf('rank(O) = %d out of %d\n', rank(O), size(Ao,1));
fprintf('det(O) = %.6g\n', det(O));

% MATLAB Control System Toolbox alternative:
% Obuiltin = obsv(Ao, Co);
% rank(Obuiltin)

% Reconstruct x(0) from output derivatives under zero input:
x0_true = [1; -2; 0.5];
eta = O*x0_true;              % eta = [y(0); ydot(0); yddot(0)]
x0_hat = O\eta;
fprintf('x0_true =\n'); disp(x0_true);
fprintf('x0_hat =\n'); disp(x0_hat);

function [A, B, C] = observable_companion(den, num)
    n = length(den);
    A = zeros(n,n);
    A(:,1) = -flipud(den(:));
    for i = 1:n-1
        A(i,i+1) = 1;
    end
    C = zeros(1,n);
    C(1) = 1;
    if nargin < 2
        B = zeros(n,1);
    else
        B = flipud(num(:));
    end
end

function O = observability_matrix(A, C)
    n = size(A,1);
    O = zeros(n,n);
    Ak = eye(n);
    for k = 1:n
        O(k,:) = C*Ak;
        Ak = Ak*A;
    end
end

      

13. Wolfram Mathematica Implementation

Chapter14_Lesson4.nb


Notebook[{
Cell["Chapter14_Lesson4.nb", "Title"],
Cell["Observability in canonical (observable) forms", "Subtitle"],
Cell[BoxData[RowBox[{"den", "=", RowBox[{"{", RowBox[{"6", ",", "11", ",", "6"}], "}"}]}]], "Input"],
Cell[BoxData[RowBox[{"num", "=", RowBox[{"{", RowBox[{"4", ",", "3", ",", "2"}], "}"}]}]], "Input"],
Cell[BoxData[RowBox[{"n", "=", RowBox[{"Length", "[", "den", "]"}]}]], "Input"],
Cell[BoxData[RowBox[{"Ao", "=", RowBox[{"ConstantArray", "[", RowBox[{"0", ",", RowBox[{"{", RowBox[{"n", ",", "n"}], "}"}]}], "]"}]}]], "Input"],
Cell[BoxData[RowBox[{RowBox[{"Ao", "[[", RowBox[{"All", ",", "1"}], "]]"}], "=", RowBox[{"-", RowBox[{"Reverse", "[", "den", "]"}]}]}]], "Input"],
Cell[BoxData[RowBox[{"Do", "[", RowBox[{RowBox[{RowBox[{"Ao", "[[", RowBox[{"i", ",", RowBox[{"i", "+", "1"}]}], "]]"}], "=", "1"}], ",", RowBox[{"{", RowBox[{"i", ",", "1", ",", RowBox[{"n", "-", "1"}]}], "}"}]}], "]"}]], "Input"],
Cell[BoxData[RowBox[{"Bo", "=", RowBox[{"Transpose", "[", RowBox[{"{", RowBox[{"Reverse", "[", "num", "]"}], "}"}], "]"}]}]], "Input"],
Cell[BoxData[RowBox[{"Co", "=", RowBox[{"{", RowBox[{"Join", "[", RowBox[{RowBox[{"{", "1", "}"}], ",", RowBox[{"ConstantArray", "[", RowBox[{"0", ",", RowBox[{"n", "-", "1"}]}], "]"}]}], "]"}], "}"}]}]], "Input"],
Cell[BoxData[RowBox[{"Omat", "=", RowBox[{"Table", "[", RowBox[{RowBox[{"First", "[", RowBox[{"Co", ".", RowBox[{"MatrixPower", "[", RowBox[{"Ao", ",", "k"}], "]"}]}], "]"}], ",", RowBox[{"{", RowBox[{"k", ",", "0", ",", RowBox[{"n", "-", "1"}]}], "}"}]}], "]"}]}]], "Input"],
Cell[BoxData[RowBox[{"MatrixForm", "/@", RowBox[{"{", RowBox[{"Ao", ",", "Bo", ",", "Co", ",", "Omat"}], "}"}]}]], "Input"],
Cell[BoxData[RowBox[{RowBox[{"{", RowBox[{RowBox[{"MatrixRank", "[", "Omat", "]"}], ",", RowBox[{"Det", "[", "Omat", "]"}], ",", RowBox[{"Eigenvalues", "[", "Ao", "]"}]}], "}"}]}]], "Input"],
Cell[BoxData[RowBox[{"x0True", "=", RowBox[{"{", RowBox[{"1", ",", RowBox[{"-", "2"}], ",", "0.5"}], "}"}]}]], "Input"],
Cell[BoxData[RowBox[{"eta", "=", RowBox[{"Omat", ".", "x0True"}]}]], "Input"],
Cell[BoxData[RowBox[{"x0Hat", "=", RowBox[{"LinearSolve", "[", RowBox[{"Omat", ",", "eta"}], "]"}]}]], "Input"]
}]

      

14. Problems and Solutions

Problem 1 (Third-Order Observable Companion Rank): Let \( A_o=\begin{bmatrix}-a_2&1&0\\-a_1&0&1\\-a_0&0&0 \end{bmatrix} \) and \( C_o=\begin{bmatrix}1&0&0\end{bmatrix} \). Compute the observability matrix and prove that the pair is observable for all real \( a_0,a_1,a_2 \).

Solution:

\[ C_oA_o=\begin{bmatrix}-a_2&1&0\end{bmatrix},\qquad C_oA_o^2=\begin{bmatrix}a_2^2-a_1&-a_2&1\end{bmatrix}. \]

\[ \mathbf{\mathcal{O}}_o=\begin{bmatrix} 1&0&0\\ -a_2&1&0\\ a_2^2-a_1&-a_2&1 \end{bmatrix}. \]

This matrix is lower triangular with diagonal entries equal to one, so \( \det(\mathbf{\mathcal{O}}_o)=1 \). Therefore \( \operatorname{rank}\mathbf{\mathcal{O}}_o=3 \) for every choice of the coefficients.

Problem 2 (Constructing OCF from a Transfer Function): Construct an observable companion realization for \( G(s)=\frac{2s^2+3s+4}{s^3+6s^2+11s+6} \).

Solution: Here \( a_2=6,a_1=11,a_0=6 \) and \( b_2=2,b_1=3,b_0=4 \). Thus

\[ A_o=\begin{bmatrix}-6&1&0\\-11&0&1\\-6&0&0\end{bmatrix}, \quad B_o=\begin{bmatrix}2\\3\\4\end{bmatrix},\quad C_o=\begin{bmatrix}1&0&0\end{bmatrix},\quad D=0. \]

Since \( C_o(sI-A_o)^{-1}=p(s)^{-1} \begin{bmatrix}s^2&s&1\end{bmatrix} \), multiplication by \( B_o \) gives the numerator \( 2s^2+3s+4 \).

Problem 3 (Initial-State Reconstruction): For the system in Problem 2, suppose zero-input measurements give \( y(0)=1 \), \( \dot y(0)=-8 \), and \( \ddot y(0)=37 \). Find \( \mathbf{x}(0) \).

Solution: The observability matrix is

\[ \mathbf{\mathcal{O}}_o=\begin{bmatrix} 1&0&0\\ -6&1&0\\ 25&-6&1 \end{bmatrix}. \]

Solve \( \mathbf{\mathcal{O}}_o\mathbf{x}(0)= \begin{bmatrix}1&-8&37\end{bmatrix}^T \). The equations are

\[ x_1(0)=1,\qquad -6x_1(0)+x_2(0)=-8, \qquad 25x_1(0)-6x_2(0)+x_3(0)=37. \]

Hence \( x_2(0)=-2 \) and \( x_3(0)=0 \). Therefore \( \mathbf{x}(0)=\begin{bmatrix}1&-2&0\end{bmatrix}^T \).

Problem 4 (Duality Argument): Prove that if \( (A,C) \) is observable, then \( (A^T,C^T) \) is controllable.

Solution: The controllability matrix of the dual pair is

\[ \mathbf{\mathcal{C}}(A^T,C^T)= \begin{bmatrix}C^T&A^TC^T&\cdots&(A^T)^{n-1}C^T\end{bmatrix}. \]

Taking the transpose of the observability matrix gives

\[ \mathbf{\mathcal{O}}(A,C)^T= \begin{bmatrix}C^T&A^TC^T&\cdots&(A^T)^{n-1}C^T\end{bmatrix}. \]

Thus \( \mathbf{\mathcal{C}}(A^T,C^T)= \mathbf{\mathcal{O}}(A,C)^T \). Since a matrix and its transpose have the same rank, observability of \( (A,C) \) is equivalent to controllability of the dual pair.

Problem 5 (PBH Hidden-Mode Interpretation): Assume there exists a nonzero vector \( \mathbf{v} \) such that \( A\mathbf{v}=\lambda\mathbf{v} \) and \( C\mathbf{v}=0 \). Explain why the system is not observable.

Solution: With initial state \( \mathbf{x}(0)=\mathbf{v} \) and zero input, the state evolves as \( \mathbf{x}(t)=e^{\lambda t} \mathbf{v} \) when \( \mathbf{v} \) is an eigenvector. The output is \( y(t)=C\mathbf{x}(t)=e^{\lambda t}C\mathbf{v}=0 \) for all \( t\geq 0 \). A nonzero internal motion is therefore invisible from the output, so the initial state cannot be uniquely reconstructed. This is precisely the PBH obstruction to observability.

15. Summary

Observable canonical forms place the measurement channel at the center of the state description. In the observable companion realization, the observability matrix is triangular with unit diagonal, so observability follows immediately. The output-derivative interpretation shows how the initial state is encoded in the output history, while duality explains why observable canonical form is the transpose counterpart of controllable canonical form. The main caution is numerical: canonical coordinates may be poorly conditioned even when the theoretical rank condition is satisfied exactly.

16. References

  1. Kalman, R.E. (1960). On the general theory of control systems. Proceedings of the First International Congress on Automatic Control, 481–492.
  2. Gilbert, E.G. (1963). Controllability and observability in multivariable control systems. SIAM Journal on Control, 1(2), 128–151.
  3. Luenberger, D.G. (1964). Observing the state of a linear system. IEEE Transactions on Military Electronics, 8(2), 74–80.
  4. Ho, B.L., & Kalman, R.E. (1966). Effective construction of linear state-variable models from input/output functions. Regelungstechnik, 14, 545–548.
  5. Popov, V.M. (1964). Hyperstability of control systems. Automation and Remote Control, 25, 1285–1305.
  6. Rosenbrock, H.H. (1965). State-space and multivariable theory. Proceedings of the Institution of Electrical Engineers, 112(8), 1527–1534.