A continuous-time LTI state-space model (introduced in Chapter 4) is described by \( \mathbf{x}(t) \): the state, \( \mathbf{u}(t) \): the input, and \( \mathbf{y}(t) \): the output:

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

The transfer function (or transfer function matrix in MIMO) is an input–output object: it encodes how \( \mathbf{u}(t) \): influences \( \mathbf{y}(t) \): under the standard transfer-function convention of zero initial conditions \( \mathbf{x}(0)=\mathbf{0} \):.

flowchart TD
  A0["Given: A, B, C, D"] --> A1["Write state model: xdot = A x + B u; y = C x + D u"]
  A1 --> A2["Take Laplace with x(0)=0"]
  A2 --> A3["Solve: (sI - A) X = B U"]
  A3 --> A4["Compute: Y = C X + D U"]
  A4 --> A5["Identify: G(s) = Y / U = C*(sI-A)^(-1)*B + D"]
        

Throughout this lesson, \( A \in \mathbb{R}^{n\times n} \):, \( B \in \mathbb{R}^{n\times m} \):, \( C \in \mathbb{R}^{p\times n} \):, \( D \in \mathbb{R}^{p\times m} \):. The result is the MIMO transfer function matrix \( G(s) \in \mathbb{C}^{p\times m} \):.

We use Laplace transforms as in classical linear control, but now applied to the vector–matrix dynamics. Let \( \mathcal{L}\{\cdot\} \): denote the Laplace transform, and define \( \mathbf{X}(s) \): \(=\mathcal{L}\{\mathbf{x}(t)\} \), \( \mathbf{U}(s) \): \(=\mathcal{L}\{\mathbf{u}(t)\} \), \( \mathbf{Y}(s) \): \(=\mathcal{L}\{\mathbf{y}(t)\} \).

Theorem 1 (Transfer Function from State Space):

For the LTI system \( \dot{\mathbf{x}}=A\mathbf{x}+B\mathbf{u} \):, \( \mathbf{y}=C\mathbf{x}+D\mathbf{u} \):, the Laplace-domain input–output map is

\[ \mathbf{Y}(s) = C(sI-A)^{-1}\mathbf{x}(0) \;+\; \Big(C(sI-A)^{-1}B + D\Big)\mathbf{U}(s), \]

and under zero initial conditions \( \mathbf{x}(0)=\mathbf{0} \):, the transfer function matrix is

\[ G(s) \;=\; \frac{\mathbf{Y}(s)}{\mathbf{U}(s)} \;=\; C(sI-A)^{-1}B + D. \]

Proof (step-by-step):

Start from the state equation and take Laplace transforms. Using the standard identity \( \mathcal{L}\{\dot{\mathbf{x}}(t)\}:\) \(= s\mathbf{X}(s) - \mathbf{x}(0)\), we get

\[ s\mathbf{X}(s) - \mathbf{x}(0) = A\mathbf{X}(s) + B\mathbf{U}(s). \]

Rearrange terms:

\[ (sI - A)\mathbf{X}(s) = \mathbf{x}(0) + B\mathbf{U}(s). \]

For any complex \( s \): such that \( (sI-A) \): is invertible, we solve:

\[ \mathbf{X}(s) = (sI-A)^{-1}\mathbf{x}(0) + (sI-A)^{-1}B\mathbf{U}(s). \]

Now take Laplace transforms of the output equation \( \mathbf{y}=C\mathbf{x}+D\mathbf{u} \)::

\[ \mathbf{Y}(s) = C\mathbf{X}(s) + D\mathbf{U}(s). \]

Substitute the expression for \( \mathbf{X}(s) \)::

\[ \mathbf{Y}(s) = C(sI-A)^{-1}\mathbf{x}(0) + \Big(C(sI-A)^{-1}B + D\Big)\mathbf{U}(s). \]

Under the transfer-function convention \( \mathbf{x}(0)=\mathbf{0} \):, the initial-condition term vanishes, so \( \mathbf{Y}(s)=G(s)\mathbf{U}(s) \): with \( G(s)=C(sI-A)^{-1}B + D \):. This completes the proof.

Interpretation: the matrix \( (sI-A)^{-1} \): is the resolvent of \( A \):. Multiplying by \( B \): injects the input into the state dynamics, multiplying by \( C \): reads out the state to the output, and \( D \): accounts for direct (instantaneous) feedthrough.

Since \( (sI-A) \): is a matrix whose entries are polynomials in \( s \):, its inverse (when it exists) has entries that are rational functions of \( s \):. A purely algebraic identity from linear algebra gives an explicit form.

Lemma 1 (Adjugate Identity): For any square matrix \( M \): with \( \det(M)\neq 0 \):,

\[ M^{-1} = \frac{\operatorname{adj}(M)}{\det(M)}. \]

Apply Lemma 1 to \( M = (sI-A) \)::

\[ (sI-A)^{-1} = \frac{\operatorname{adj}(sI-A)}{\det(sI-A)}. \]

Substituting into the transfer function formula yields the rational representation

\[ G(s) = C\frac{\operatorname{adj}(sI-A)}{\det(sI-A)}B + D = \frac{C\,\operatorname{adj}(sI-A)\,B + D\,\det(sI-A)}{\det(sI-A)}. \]

In particular, for strictly proper systems \( D=\mathbf{0} \): the transfer matrix is

\[ G(s) = \frac{C\,\operatorname{adj}(sI-A)\,B}{\det(sI-A)}. \]

This expression is useful theoretically: it shows every entry of \( G(s) \): is a rational function whose denominator divides \( \det(sI-A) \):, a monic polynomial of degree \( n \):. (A deeper discussion of poles and cancellations is deferred to Lesson 2.)

A second representation connects state-space matrices to high-frequency behavior. For large \( |s| \):, factor out \( s \)::

\[ (sI-A)^{-1} = \frac{1}{s}\Big(I-\frac{A}{s}\Big)^{-1}. \]

Using the matrix geometric series (from Chapter 3), when \( \|A/s\| \): is sufficiently small,

\[ \Big(I-\frac{A}{s}\Big)^{-1} = \sum_{k=0}^{\infty}\Big(\frac{A}{s}\Big)^k, \qquad (sI-A)^{-1} = \sum_{k=0}^{\infty}\frac{A^k}{s^{k+1}}. \]

Therefore,

\[ G(s) = D + \sum_{k=0}^{\infty}\frac{CA^kB}{s^{k+1}}. \]

A key immediate corollary is the classical limit:

\[ \lim_{s → \infty} G(s) = D. \]

Hence, a system is strictly proper (transfer goes to zero at high frequency) precisely when \( D=\mathbf{0} \): (formalized again in the Problems section).

In numerical work, explicitly forming \( (sI-A)^{-1} \): is often inferior to solving a linear system. For any chosen complex frequency \( s \): (e.g., \( s=j\omega \):), compute \( X_s \): as the solution of \( (sI-A)X_s = B \):, then \( G(s)=CX_s + D \):.

\[ \text{Solve } (sI-A)X_s = B \quad \Rightarrow \quad G(s)=CX_s + D. \]

flowchart TD
  S0["Pick s (e.g., s = j*w)"] --> S1["Form M = s*I - A"]
  S1 --> S2["Solve linear system: M * X = B"]
  S2 --> S3["Compute G(s) = C*X + D"]
  S3 --> S4["Repeat for grid of s values (frequency response)"]
        

This approach generalizes seamlessly to MIMO: \( X_s \): is an \( n\times m \): matrix, solved via multiple right-hand sides, and \( G(s) \): is \( p\times m \):.

Recommended libraries: python-control (state-space/transfer tools), NumPy/SciPy (linear algebra), and SymPy (symbolic derivations when desired).

import numpy as np

# Example (SISO, n=2)
A = np.array([[0.0, 1.0],
              [-2.0, -3.0]])
B = np.array([[0.0],
              [1.0]])
C = np.array([[1.0, 0.0]])
D = np.array([[0.0]])

# 1) Direct evaluation of G(s) = C (sI-A)^{-1} B + D at a chosen s
def G_of_s(s, A, B, C, D):
    n = A.shape[0]
    M = (s * np.eye(n) - A).astype(complex)
    X = np.linalg.solve(M, B.astype(complex))  # solves (sI-A)X = B
    return C.astype(complex) @ X + D.astype(complex)

s = 1j * 2.0
print("G(j2) =", G_of_s(s, A, B, C, D))

# 2) Convert ss -> tf using python-control
# pip install control
import control
sys_ss = control.ss(A, B, C, D)

# Transfer function form (for SISO, yields a scalar tf; for MIMO, a tf matrix)
sys_tf = control.tf(sys_ss)
print(sys_tf)

# 3) Frequency response check
w = np.logspace(-2, 2, 200)
mag, phase, omega = control.bode(sys_ss, w, Plot=False)  # uses stable numerical routines
print("Computed bode arrays:", mag.shape, phase.shape)
      

The function \( G(s) \): above uses the numerically preferred linear-solve route rather than explicitly computing \( (sI-A)^{-1} \):.

In C++, a common pattern is to evaluate \( G(s) \): on a frequency grid using complex linear algebra. Recommended library: Eigen (matrix factorizations and linear solves).

#include <iostream>
#include <complex>
#include <Eigen/Dense>

int main() {
  using cd = std::complex<double>;
  using Mat = Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>;
  using MatC = Eigen::Matrix<cd, Eigen::Dynamic, Eigen::Dynamic>;

  // Example: n=2, SISO
  Mat A(2,2); A << 0.0, 1.0,
                -2.0, -3.0;
  Mat B(2,1); B << 0.0,
                1.0;
  Mat C(1,2); C << 1.0, 0.0;
  Mat D(1,1); D << 0.0;

  cd s(0.0, 2.0); // s = j*2

  MatC I = MatC::Identity(2,2);
  MatC M = s*I - A.cast<cd>();
  MatC X = M.partialPivLu().solve(B.cast<cd>()); // (sI-A)X = B
  MatC G = C.cast<cd>()*X + D.cast<cd>();

  std::cout << "G(j2) = " << G(0,0) << std::endl;
  return 0;
}
      

For MIMO, make B have multiple columns; Eigen will solve all right-hand sides simultaneously, producing \( X_s \): as \( n\times m \):, then \( G(s)=CX_s + D \):.

Recommended Java library: EJML, which supports complex matrices using ZMatrixRMaj. Strategy: form \( M=sI-A \):, solve \( MX=B \):, then compute \( G=CX+D \):.

import org.ejml.data.ZMatrixRMaj;
import org.ejml.dense.row.CommonOps_ZDRM;

public class TransferFromStateSpace {
  public static void main(String[] args) {
    // Example: n=2, SISO
    // A, B, C, D are real but stored in complex form (imag part = 0)
    ZMatrixRMaj A = new ZMatrixRMaj(2,2);
    A.set(0,0, 0.0, 0.0);  A.set(0,1, 1.0, 0.0);
    A.set(1,0,-2.0, 0.0);  A.set(1,1,-3.0, 0.0);

    ZMatrixRMaj B = new ZMatrixRMaj(2,1);
    B.set(0,0, 0.0, 0.0);
    B.set(1,0, 1.0, 0.0);

    ZMatrixRMaj C = new ZMatrixRMaj(1,2);
    C.set(0,0, 1.0, 0.0);
    C.set(0,1, 0.0, 0.0);

    ZMatrixRMaj D = new ZMatrixRMaj(1,1);
    D.set(0,0, 0.0, 0.0);

    // s = j*2
    double sRe = 0.0;
    double sIm = 2.0;

    // M = sI - A
    ZMatrixRMaj M = new ZMatrixRMaj(2,2);
    CommonOps_ZDRM.setIdentity(M);
    // scale I by s
    CommonOps_ZDRM.scale(sRe, sIm, M);
    // M = M - A
    CommonOps_ZDRM.subtractEquals(M, A);

    // Solve M X = B
    ZMatrixRMaj X = new ZMatrixRMaj(2,1);
    if (!CommonOps_ZDRM.solve(M, B, X)) {
      throw new RuntimeException("Solve failed (matrix may be singular for this s).");
    }

    // G = C X + D
    ZMatrixRMaj G = new ZMatrixRMaj(1,1);
    CommonOps_ZDRM.mult(C, X, G);
    CommonOps_ZDRM.addEquals(G, D);

    System.out.println("G(j2) = " + G.getReal(0,0) + " + j" + G.getImag(0,0));
  }
}
      

This evaluation-first approach is often preferable in real-time or embedded contexts where one needs frequency responses or transfer values without symbolic polynomial construction.

MATLAB’s Control System Toolbox provides direct conversion and evaluation. For MIMO, tf(sys) returns a transfer matrix object.

% Example matrices
A = [0 1; -2 -3];
B = [0; 1];
C = [1 0];
D = 0;

% Build state-space model
sys = ss(A,B,C,D);

% Convert to transfer function
G = tf(sys)

% Evaluate G(s) at s = j*2
s = 1j*2;
Gj2 = evalfr(sys, s)

% For SISO, you can also extract numerator/denominator explicitly
[num, den] = ss2tf(A,B,C,D);
num, den
      

Simulink note (conceptual): Use the State-Space block to implement \( (A,B,C,D) \):. To obtain a transfer function from a Simulink model, linearize around an operating point and convert the resulting LTI object to tf. Typical workflows use linear analysis tools (e.g., model linearization functions) and then apply tf to the linearized system.

Mathematica supports state-space and transfer-function models and can perform symbolic manipulations.

(* Example matrices *)
A = { {0, 1}, {-2, -3} };
B = { {0}, {1} };
C = { {1, 0} };
D = { {0} };

(* State-space model *)
ss = StateSpaceModel[{A, B, C, D}];

(* Convert to transfer function model *)
tf = TransferFunctionModel[ss, s]
(* Alternatively: TransferFunctionModel[ss] may infer variables in some versions *)

(* Evaluate at s = I*2 *)
Gj2 = tf /. s -> I*2 // Simplify

(* Direct formula: G(s)=C.(sI-A)^(-1).B + D *)
Gdirect[s_] := C . Inverse[s IdentityMatrix[Length[A]] - A] . B + D
Simplify[Gdirect[s]]
      

Problem 1 (Direct Derivation with Initial Conditions): Consider \( \dot{\mathbf{x}}=A\mathbf{x}+B\mathbf{u} \):, \( \mathbf{y}=C\mathbf{x}+D\mathbf{u} \):. Starting from Laplace transforms, derive \( \mathbf{Y}(s) \): in terms of \( \mathbf{U}(s) \): and \( \mathbf{x}(0) \):.

Solution:

Taking Laplace transforms gives \( s\mathbf{X}(s)-\mathbf{x}(0)=A\mathbf{X}(s)+B\mathbf{U}(s) \):, hence \( (sI-A)\mathbf{X}(s)=\mathbf{x}(0)+B\mathbf{U}(s) \):, so

\[ \mathbf{X}(s)=(sI-A)^{-1}\mathbf{x}(0) + (sI-A)^{-1}B\mathbf{U}(s). \]

Then \( \mathbf{Y}(s)=C\mathbf{X}(s)+D\mathbf{U}(s) \):, so

\[ \mathbf{Y}(s) = C(sI-A)^{-1}\mathbf{x}(0) + \Big(C(sI-A)^{-1}B + D\Big)\mathbf{U}(s). \]

Under \( \mathbf{x}(0)=\mathbf{0} \):, the transfer function matrix is \( G(s)=C(sI-A)^{-1}B+D \):.

Problem 2 (Compute a SISO Transfer Function): Let \( A=\begin{bmatrix}0 & 1\\ -2 & -3\end{bmatrix} \):, \( B=\begin{bmatrix}0\\1\end{bmatrix} \):, \( C=\begin{bmatrix}1 & 0\end{bmatrix} \):, \( D=0 \):. Compute \( G(s) \): explicitly as a rational function.

Solution:

First compute \( sI-A \)::

\[ sI-A = \begin{bmatrix} s & -1\\ 2 & s+3 \end{bmatrix}. \]

Its determinant is \( \det(sI-A)=s(s+3)+2=s^2+3s+2 \):. The inverse is

\[ (sI-A)^{-1} = \frac{1}{s^2+3s+2} \begin{bmatrix} s+3 & 1\\ -2 & s \end{bmatrix}. \]

Now compute \( (sI-A)^{-1}B \)::

\[ (sI-A)^{-1}B = \frac{1}{s^2+3s+2} \begin{bmatrix} s+3 & 1\\ -2 & s \end{bmatrix} \begin{bmatrix} 0\\1 \end{bmatrix} = \frac{1}{s^2+3s+2} \begin{bmatrix} 1\\ s \end{bmatrix}. \]

Finally, \( G(s)=C(sI-A)^{-1}B + D \): gives

\[ G(s)= \begin{bmatrix}1 & 0\end{bmatrix} \frac{1}{s^2+3s+2} \begin{bmatrix}1\\ s\end{bmatrix} = \frac{1}{s^2+3s+2}. \]

Problem 3 (Strictly Proper iff \( D=0 \)): Prove that \( \lim_{s → \infty} G(s) = D \):, and conclude that the system is strictly proper iff \( D=\mathbf{0} \):.

Solution:

Start from \( G(s)=C(sI-A)^{-1}B + D \):. Factor \( s \)::

\[ (sI-A)^{-1} = \frac{1}{s}\Big(I-\frac{A}{s}\Big)^{-1}. \]

For large \( |s| \):, the inverse admits the convergent series \( \Big(I-\frac{A}{s}\Big)^{-1} = \sum_{k=0}^{\infty}\Big(\frac{A}{s}\Big)^k \):, hence

\[ C(sI-A)^{-1}B = C\left(\sum_{k=0}^{\infty}\frac{A^k}{s^{k+1}}\right)B = \sum_{k=0}^{\infty}\frac{CA^kB}{s^{k+1}}. \]

Every term in the sum tends to \( 0 \): as \( s → \infty \):, therefore \( \lim_{s → \infty} C(sI-A)^{-1}B = 0 \):, implying

\[ \lim_{s → \infty} G(s) = \lim_{s → \infty}\Big(C(sI-A)^{-1}B + D\Big) = D. \]

By definition, “strictly proper” means \( \lim_{s → \infty} G(s)=0 \):, which holds iff \( D=\mathbf{0} \):.

Problem 4 (Efficient Computation at a Frequency): Show that evaluating \( G(s) \): at a single \( s \): can be done by solving \( (sI-A)X=B \): and computing \( G=CX+D \):, without explicitly forming \( (sI-A)^{-1} \):. Explain why this is numerically preferable.

Solution:

Algebraically, \( X=(sI-A)^{-1}B \): is exactly the matrix solution of \( (sI-A)X=B \):. Substituting into \( G=C(sI-A)^{-1}B + D \): yields \( G=CX+D \):.

Numerically, computing an explicit inverse is typically less stable and more expensive than a factorization-and-solve: linear solvers (LU/QR) reduce rounding error amplification and reuse factorizations across multiple right-hand sides (the columns of \( B \):) and across many frequency points in a sweep.

We proved the fundamental conversion from state-space to transfer-function form: \( G(s): \(=C(sI-A)^{-1}B + D \), with a clear separation between the zero-initial-condition transfer mapping and the additional initial-state response term \( C(sI-A)^{-1}\mathbf{x}(0) \):. We also derived rational and series representations and presented practical numerical evaluation via linear solves, along with implementations in Python, C++, Java, MATLAB/Simulink, and Mathematica.