Chapter 8: State-Space Modeling
Lesson 2: Deriving State-Space Models from Differential Equations and Physical Laws
This lesson develops rigorous, reproducible procedures for deriving state-space models from (i) higher-order differential equations and (ii) first-principles physical laws (Newton/Kirchhoff). Emphasis is placed on selecting physically meaningful state variables, proving equivalence between descriptions, and producing simulation-ready models for multi-input multi-output systems.
1. Motivation and Learning Objectives
In Lesson 1, we introduced the idea that the state is a minimal set of variables that summarizes the system’s “memory.” Here we operationalize that concept: given either a differential equation or governing physical laws, we will construct a state model of the form \( \dot{\mathbf{x}}(t) = \mathbf{f}(\mathbf{x}(t), \mathbf{u}(t), t) \) and \( \mathbf{y}(t) = \mathbf{g}(\mathbf{x}(t), \mathbf{u}(t), t) \).
For linear time-invariant (LTI) systems, the model becomes \( \dot{\mathbf{x}}(t) = \mathbf{A}\mathbf{x}(t) + \mathbf{B}\mathbf{u}(t) \), \( \mathbf{y}(t) = \mathbf{C}\mathbf{x}(t) + \mathbf{D}\mathbf{u}(t) \). A key engineering benefit is that the state directly supports simulation, interconnection, and later analysis (e.g., eigenvalue-based stability in Chapter 11).
flowchart TD
S["Start: given ODEs or physical laws"] --> I["Identify inputs u(t) and outputs y(t)"]
I --> E["Find energy-storage variables (candidates for state)"]
E --> C["Choose state vector x(t): minimal + physically meaningful"]
C --> R["Rewrite dynamics as first-order: xdot = f(x,u,t)"]
R --> O["Write output law: y = g(x,u,t)"]
O --> V["Check units + initial conditions + consistency"]
V --> L["If linear LTI: form A,B,C,D"]
L --> SIM["Simulate and validate (numerically)"]
2. Scalar n-th Order ODE → First-Order State Model (Constructive Method)
Consider a scalar input-output ODE with input \( u(t) \) and output \( y(t) \). A standard linear form is
\[ y^{(n)}(t) + a_{n-1} y^{(n-1)}(t) + \cdots + a_1 \dot{y}(t) + a_0 y(t) = b_m u^{(m)}(t) + \cdots + b_1 \dot{u}(t) + b_0 u(t), \quad m \le n. \]
In many control-oriented derivations, we avoid input derivatives in the state model by assuming the physical input enters without differentiation (or by augmenting states, which we avoid here). So we focus on the common case
\[ y^{(n)}(t) + a_{n-1} y^{(n-1)}(t) + \cdots + a_1 \dot{y}(t) + a_0 y(t) = b_0 u(t). \]
Define state variables (a constructive choice): \( x_1(t) = y(t) \), \( x_2(t) = \dot{y}(t) \), …, \( x_n(t) = y^{(n-1)}(t) \). Then by definition, \( \dot{x}_1(t) = x_2(t) \), …, \( \dot{x}_{n-1}(t) = x_n(t) \). The last derivative comes from the ODE:
\[ \dot{x}_n(t) = y^{(n)}(t) = -a_{n-1} x_n(t) - a_{n-2} x_{n-1}(t) - \cdots - a_1 x_2(t) - a_0 x_1(t) + b_0 u(t). \]
Stacking into \( \mathbf{x}(t) = [x_1(t)\; x_2(t)\; \cdots\; x_n(t)]^\top \), we obtain an LTI state model:
\[ \dot{\mathbf{x}}(t) = \underbrace{\begin{bmatrix} 0 & 1 & 0 & \cdots & 0 \\ 0 & 0 & 1 & \cdots & 0 \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & 0 & \cdots & 1 \\ -a_0 & -a_1 & -a_2 & \cdots & -a_{n-1} \end{bmatrix}}_{\mathbf{A}} \mathbf{x}(t) + \underbrace{\begin{bmatrix} 0\\0\\ \vdots\\0\\ b_0 \end{bmatrix}}_{\mathbf{B}} u(t). \]
With output \( y(t)=x_1(t) \) we have \( \mathbf{C} = [1\; 0\; \cdots\; 0] \), \( \mathbf{D}=0 \).
2.1 Proof of Equivalence (ODE solution ⇔ State trajectory)
Proposition: Let \( y(t) \) be any sufficiently differentiable function. Define \( x_1,\ldots,x_n \) by \( x_k(t)=y^{(k-1)}(t) \). Then \( y(t) \) satisfies the n-th order ODE if and only if \( \mathbf{x}(t) \) satisfies the constructed first-order state equation with the same input.
Proof:
(Only if) If \( y(t) \) satisfies the ODE, then by definition \( \dot{x}_k(t)=x_{k+1}(t) \) for \( k=1,\ldots,n-1 \). The ODE gives an explicit expression for \( y^{(n)}(t)=\dot{x}_n(t) \) as a linear combination of \( x_1,\ldots,x_n \) and \( u(t) \), which matches the last row of \( \mathbf{A}\mathbf{x}+\mathbf{B}u \). Hence \( \dot{\mathbf{x}}(t)=\mathbf{A}\mathbf{x}(t)+\mathbf{B}u(t) \).
(If) Conversely, if \( \mathbf{x}(t) \) satisfies the state model and we define \( y(t)=x_1(t) \), then repeated differentiation yields \( \dot{y}(t)=x_2(t) \), …, \( y^{(n-1)}(t)=x_n(t) \). Differentiating one more time and using the last state equation gives \( y^{(n)}(t) = -a_{n-1}y^{(n-1)}(t) - \cdots - a_0 y(t) + b_0 u(t) \), which is exactly the ODE. ∎
2.2 Initial Conditions Mapping
An n-th order ODE requires initial data \( y(0), \dot{y}(0), \ldots, y^{(n-1)}(0) \). Under the construction, these are precisely the state initial condition \( \mathbf{x}(0) = [y(0)\; \dot{y}(0)\; \cdots\; y^{(n-1)}(0)]^\top \).
3. Coupled Differential Equations → Vector State Model
Many physical systems naturally yield coupled ODEs. Suppose we have \( q \) scalar variables collected into \( \mathbf{z}(t) \in \mathbb{R}^q \) governed by a first-order form already: \( \dot{\mathbf{z}}(t) = \mathbf{h}(\mathbf{z}(t), \mathbf{u}(t), t) \). Then we can directly set \( \mathbf{x}(t) = \mathbf{z}(t) \).
More commonly, we obtain mixed orders (e.g., some second-order, some first-order). A systematic method is:
- For each second-order coordinate \( q_i(t) \), define \( x_{q_i}(t)=q_i(t) \) and \( x_{\dot{q}_i}(t)=\dot{q}_i(t) \).
- For each first-order energy variable (e.g., capacitor voltage, inductor current), include it directly.
- Use the governing equations to write \( \dot{\mathbf{x}}(t) \) solely in terms of \( \mathbf{x}(t) \), \( \mathbf{u}(t) \), and possibly \( t \).
For linear models, this yields \( \dot{\mathbf{x}}(t)=\mathbf{A}\mathbf{x}(t)+\mathbf{B}\mathbf{u}(t) \). For nonlinear models, it yields the general form \( \dot{\mathbf{x}}(t)=\mathbf{f}(\mathbf{x}(t),\mathbf{u}(t),t) \).
3.1 Well-Posedness (Why State Determines the Future)
The “state” claim is meaningful only if the model is well-posed. A standard sufficient condition is local Lipschitz continuity in \( \mathbf{x} \).
Theorem (Existence–Uniqueness, classical): If \( \mathbf{f}(\mathbf{x},\mathbf{u},t) \) is continuous in \( t \) and locally Lipschitz in \( \mathbf{x} \) (for fixed \( \mathbf{u}(t) \)), then for any initial condition \( \mathbf{x}(t_0)=\mathbf{x}_0 \) there exists a unique trajectory in some interval around \( t_0 \).
This formally justifies the statement: the pair \( (\mathbf{x}(t_0), \mathbf{u}(\cdot)) \) uniquely determines future evolution (locally in time).
4. Derivation from Physical Laws: A State-Selection Principle
When deriving from first principles (Newton/Kirchhoff), a reliable guideline is: choose as states the variables associated with energy storage. In the chapters on mechanical and electrical modeling you have already seen typical storage elements:
- Mechanical: mass (kinetic energy via velocity), spring (potential energy via displacement/extension).
- Electrical: inductor (magnetic energy via current), capacitor (electric energy via voltage/charge).
- Thermal/Fluid (later expansions): capacitances store energy/mass, resistances dissipate.
This state choice is not arbitrary: it tends to produce minimal models without algebraic loops, and yields states with clear units and measurable meaning.
flowchart TD
A["Physical laws: Newton/Kirchhoff"] --> B["Write balance equations"]
B --> C["Identify storage elements"]
C --> D["Pick storage variables as states"]
D --> E["Express derivatives using laws"]
E --> F["Obtain xdot = f(x,u,t)"]
F --> G["Add output equation y = g(x,u,t)"]
4.1 Example A: Mass–Spring–Damper (Force Input)
Consider the translational system (Chapter 3) with mass \( m \), damping \( c \), stiffness \( k \), displacement \( x(t) \), and force input \( u(t)=F(t) \). Newton’s law gives
\[ m \ddot{x}(t) + c \dot{x}(t) + k x(t) = u(t). \]
Choose states (energy storage variables): \( x_1(t)=x(t) \) (spring storage coordinate), \( x_2(t)=\dot{x}(t) \) (mass velocity). Then \( \dot{x}_1(t)=x_2(t) \) and
\[ \dot{x}_2(t) = \ddot{x}(t) = -\frac{k}{m}x_1(t) - \frac{c}{m}x_2(t) + \frac{1}{m}u(t). \]
So the state-space matrices (with output as displacement \( y(t)=x(t) \)) are:
\[ \dot{\mathbf{x}}(t) = \underbrace{\begin{bmatrix} 0 & 1\\ -\frac{k}{m} & -\frac{c}{m} \end{bmatrix}}_{\mathbf{A}} \mathbf{x}(t) + \underbrace{\begin{bmatrix} 0\\ \frac{1}{m} \end{bmatrix}}_{\mathbf{B}} u(t), \qquad y(t) = \underbrace{\begin{bmatrix} 1 & 0 \end{bmatrix}}_{\mathbf{C}} \mathbf{x}(t) + \underbrace{[0]}_{\mathbf{D}}u(t). \]
4.2 Example B: Series RLC Circuit (Voltage Input, Capacitor Voltage Output)
For a series RLC driven by input voltage \( u(t)=v_{in}(t) \), with inductor current \( i(t) \) and capacitor voltage \( v_C(t) \), Kirchhoff’s voltage law yields
\[ L\frac{di(t)}{dt} + R i(t) + v_C(t) = u(t), \qquad C\frac{dv_C(t)}{dt} = i(t). \]
Choose states as storage variables: \( x_1(t)=i(t) \) (inductor storage), \( x_2(t)=v_C(t) \) (capacitor storage). Then
\[ \dot{x}_1(t) = -\frac{R}{L}x_1(t) - \frac{1}{L}x_2(t) + \frac{1}{L}u(t), \qquad \dot{x}_2(t) = \frac{1}{C}x_1(t). \]
With output \( y(t)=v_C(t)=x_2(t) \), the matrices are:
\[ \mathbf{A}= \begin{bmatrix} -\frac{R}{L} & -\frac{1}{L}\\ \frac{1}{C} & 0 \end{bmatrix}, \quad \mathbf{B}= \begin{bmatrix} \frac{1}{L}\\ 0 \end{bmatrix}, \quad \mathbf{C}= \begin{bmatrix} 0 & 1 \end{bmatrix}, \quad \mathbf{D}=[0]. \]
5. Output Equations and Direct Feedthrough
The output equation \( \mathbf{y}(t)=\mathbf{g}(\mathbf{x}(t),\mathbf{u}(t),t) \) encodes what is measured or of interest. For LTI systems, \( \mathbf{y}(t)=\mathbf{C}\mathbf{x}(t)+\mathbf{D}\mathbf{u}(t) \).
Interpretation of \( \mathbf{D} \): If \( \mathbf{D}\ne 0 \), the input affects the output “instantaneously.” For many mechanical systems, displacement does not depend instantaneously on force (so \( \mathbf{D}=0 \)). In some electrical measurements (e.g., output chosen as resistor voltage \( v_R=Ri \)), you may still have \( \mathbf{D}=0 \) because \( i \) is a state; but in other interconnections, algebraic relations can produce direct terms.
5.1 Consistency Check: Units
A practical derivation check is unit consistency: \( [\dot{\mathbf{x}}] = [\mathbf{A}][\mathbf{x}] \) and \( [\dot{\mathbf{x}}] = [\mathbf{B}][\mathbf{u}] \). For the mass–spring–damper, \( [x_1]=\mathrm{m} \), \( [x_2]=\mathrm{m/s} \), so \( A_{12} \) must have units \( \mathrm{s}^{-1} \) (it equals 1 in our coordinate choice), and \( k/m \) has units \( \mathrm{s}^{-2} \), as expected.
6. Python Lab: Derive, Build, Simulate (NumPy/SciPy + python-control)
Recommended Python libraries for this lesson: \( \texttt{numpy} \), \( \texttt{scipy.integrate} \) for simulation; \( \texttt{control} \) (python-control) for LTI state-space objects; \( \texttt{sympy} \) for symbolic checks.
import numpy as np
from scipy.integrate import solve_ivp
# Mass-spring-damper parameters
m = 1.0
c = 0.4
k = 4.0
A = np.array([[0.0, 1.0],
[-k/m, -c/m]])
B = np.array([[0.0],
[1.0/m]])
C = np.array([[1.0, 0.0]])
D = np.array([[0.0]])
def u_of_t(t):
# Step input u(t) = 1 for t >= 0
return 1.0
def f(t, x):
x = x.reshape((-1, 1))
dx = A @ x + B * u_of_t(t)
return dx.flatten()
t_span = (0.0, 10.0)
x0 = np.array([0.0, 0.0]) # initial displacement and velocity
sol = solve_ivp(f, t_span, x0, max_step=0.01, dense_output=True)
# Output y = C x + D u
t_eval = np.linspace(t_span[0], t_span[1], 1001)
x_eval = sol.sol(t_eval).T
y_eval = (x_eval @ C.T).flatten()
print("Final state:", x_eval[-1])
print("Final output:", y_eval[-1])
# Optional: create an LTI state-space object (python-control)
try:
import control
sys = control.ss(A, B, C, D)
print(sys)
except Exception as e:
print("python-control not available in this environment:", e)
Engineering note: the numerical solver does not require explicit construction of \( \mathbf{A},\mathbf{B} \) if the model is nonlinear; it only requires a function computing \( \dot{\mathbf{x}}=\mathbf{f}(\mathbf{x},\mathbf{u},t) \).
7. C++ Lab: Simulation from State Equations (Eigen + Boost.Odeint)
Recommended C++ libraries: Eigen for matrices/vectors and Boost.Odeint for ODE integration. Below uses a compact state type and integrates the mass–spring–damper model.
#include <array>
#include <iostream>
#include <boost/numeric/odeint.hpp>
using state_type = std::array<double, 2>;
using boost::numeric::odeint::runge_kutta4;
struct MSD {
double m, c, k;
double u(double /*t*/) const {
return 1.0; // step input
}
void operator()(const state_type &x, state_type &dxdt, double t) const {
// x[0] = displacement, x[1] = velocity
dxdt[0] = x[1];
dxdt[1] = -(k/m) * x[0] - (c/m) * x[1] + (1.0/m) * u(t);
}
};
int main() {
MSD sys{1.0, 0.4, 4.0};
state_type x{{0.0, 0.0}};
runge_kutta4<state_type> stepper;
double t0 = 0.0, tf = 10.0, dt = 0.01;
double t = t0;
while (t < tf) {
stepper.do_step(sys, x, t, dt);
t += dt;
}
std::cout << "Final state: x = [" << x[0] << ", " << x[1] << "]\n";
std::cout << "Final output (displacement): y = " << x[0] << "\n";
return 0;
}
This code is a direct translation of \( \dot{x}_1=x_2 \), \( \dot{x}_2=-(k/m)x_1-(c/m)x_2+(1/m)u \).
8. Java Lab: State-Space Simulation (Apache Commons Math)
Recommended Java libraries: Apache Commons Math for ODE solvers and EJML for matrix computations (optional). The following integrates the same mass–spring–damper state model.
import org.apache.commons.math3.ode.FirstOrderDifferentialEquations;
import org.apache.commons.math3.ode.nonstiff.DormandPrince54Integrator;
public class MassSpringDamperSim {
static class MSD implements FirstOrderDifferentialEquations {
private final double m, c, k;
MSD(double m, double c, double k) {
this.m = m; this.c = c; this.k = k;
}
public int getDimension() {
return 2;
}
private double u(double t) {
return 1.0; // step input
}
public void computeDerivatives(double t, double[] x, double[] xDot) {
// x[0] = displacement, x[1] = velocity
xDot[0] = x[1];
xDot[1] = -(k/m) * x[0] - (c/m) * x[1] + (1.0/m) * u(t);
}
}
public static void main(String[] args) {
MSD ode = new MSD(1.0, 0.4, 4.0);
double t0 = 0.0;
double tf = 10.0;
double[] x0 = new double[] {0.0, 0.0};
DormandPrince54Integrator integrator =
new DormandPrince54Integrator(1.0e-6, 1.0e-2, 1.0e-10, 1.0e-10);
double[] x = x0.clone();
integrator.integrate(ode, t0, x, tf, x);
System.out.println("Final state: x = [" + x[0] + ", " + x[1] + "]");
System.out.println("Final output (displacement): y = " + x[0]);
}
}
9. MATLAB/Simulink Lab: ss() and a State-Space Block
Recommended MATLAB toolboxes:
Control System Toolbox (for ss,
lsim) and Simulink (for graphical
simulation).
9.1 MATLAB Script
% Mass-spring-damper parameters
m = 1.0; c = 0.4; k = 4.0;
A = [0 1; -k/m -c/m];
B = [0; 1/m];
C = [1 0];
D = 0;
sys = ss(A,B,C,D);
t = linspace(0,10,1001);
u = ones(size(t)); % step input
x0 = [0; 0];
[y, t_out, x_out] = lsim(sys, u, t, x0);
disp("Final state:");
disp(x_out(end,:));
disp("Final output:");
disp(y(end));
9.2 Simulink (Conceptual Build Steps)
- Open Simulink and place a Step block for \( u(t) \).
- Place a State-Space block and set A,B,C,D to the matrices above.
- Connect Step → State-Space input; connect output to a Scope.
- Set initial conditions in the block to \( \mathbf{x}(0)=[0\;0]^\top \).
This Simulink model is a direct numerical realization of the derived continuous-time state equations.
10. Wolfram Mathematica Lab: StateSpaceModel and NDSolve
Mathematica provides native support for state-space objects and responses. Below are two approaches: (i) explicit ODE integration and (ii) a state-space model object.
(* Mass-spring-damper parameters *)
m = 1.0; c = 0.4; k = 4.0;
A = { {0, 1}, {-k/m, -c/m} };
B = { {0}, {1/m} };
Cmat = { {1, 0} };
Dmat = { {0} };
u[t_] := 1; (* step input *)
(* (i) Direct ODE integration *)
eqns = {
x1'[t] == x2[t],
x2'[t] == -(k/m) x1[t] - (c/m) x2[t] + (1/m) u[t],
x1[0] == 0, x2[0] == 0
};
sol = NDSolve[eqns, {x1, x2}, {t, 0, 10}];
y[t_] := x1[t] /. sol[[1]];
{y[10], (x1[10] /. sol[[1]]), (x2[10] /. sol[[1]])}
(* (ii) State-space model object *)
Needs["ControlSystems`"];
sys = StateSpaceModel[{A, B, Cmat, Dmat}];
(* Response to unit step input *)
resp = OutputResponse[sys, 1, {t, 0, 10}];
resp[10]
11. Problems and Solutions
The following problems reinforce the derivation process and the discipline of selecting states and outputs.
Problem 1 (Third-Order ODE Realization): Given \( y^{(3)}(t) + 6 y^{(2)}(t) + 11 \dot{y}(t) + 6 y(t) = 2 u(t) \), construct a state-space model with \( y(t)=x_1(t) \).
Solution: Let \( x_1=y \), \( x_2=\dot{y} \), \( x_3=y^{(2)} \). Then \( \dot{x}_1=x_2 \), \( \dot{x}_2=x_3 \), and
\[ \dot{x}_3 = y^{(3)} = -6x_3 - 11x_2 - 6x_1 + 2u. \]
Hence
\[ \mathbf{A}= \begin{bmatrix} 0 & 1 & 0\\ 0 & 0 & 1\\ -6 & -11 & -6 \end{bmatrix}, \quad \mathbf{B}= \begin{bmatrix} 0\\0\\2 \end{bmatrix}, \quad \mathbf{C}=\begin{bmatrix} 1 & 0 & 0 \end{bmatrix}, \quad \mathbf{D}=[0]. \]
Problem 2 (Mass–Spring–Damper with Velocity Output): For \( m\ddot{x}+c\dot{x}+kx=u \), choose \( x_1=x \), \( x_2=\dot{x} \), but now let the output be \( y=\dot{x} \). Find \( \mathbf{C},\mathbf{D} \).
Solution: Since \( y=x_2 \), \( \mathbf{C}=[0\;1] \) and \( \mathbf{D}=0 \). The state matrices \( \mathbf{A},\mathbf{B} \) are unchanged from Section 4.1.
Problem 3 (Series RLC with Current Output): Using the series RLC equations \( L\dot{i}+Ri+v_C=u \), \( C\dot{v}_C=i \), choose states \( x_1=i \), \( x_2=v_C \). Let output be \( y=i \). Determine \( \mathbf{C},\mathbf{D} \).
Solution: Since \( y=x_1 \), \( \mathbf{C}=[1\;0] \), \( \mathbf{D}=0 \). The matrices \( \mathbf{A},\mathbf{B} \) are those in Section 4.2.
Problem 4 (Input–Output ODE from a Given State Model): Consider the 2-state LTI system \( \dot{x}_1=x_2 \), \( \dot{x}_2=-4x_1-2x_2+u \), \( y=x_1 \). Eliminate the state to obtain a single ODE in \( y \).
Solution: Since \( y=x_1 \) and \( \dot{y}=x_2 \), differentiate again: \( \ddot{y}=\dot{x}_2=-4x_1-2x_2+u \). Substitute \( x_1=y \), \( x_2=\dot{y} \):
\[ \ddot{y}(t) + 2\dot{y}(t) + 4y(t) = u(t). \]
Problem 5 (State Choice Justification): Explain (mathematically) why choosing \( \mathbf{x}=[x\;\dot{x}]^\top \) for the mass–spring–damper system yields a first-order system in \( \mathbf{x} \) without requiring higher derivatives.
Solution: The governing law is second order in \( x \). Defining \( x_1=x \), \( x_2=\dot{x} \) makes \( \dot{x}_1 \) equal to an existing state and expresses \( \dot{x}_2=\ddot{x} \) directly from the physical law as a function of \( x_1,x_2,u \). Thus the vector field closes: \( \dot{\mathbf{x}}=\mathbf{f}(\mathbf{x},u) \).
12. Summary
We derived state-space models using two complementary perspectives: (i) constructive rewriting of higher-order ODEs into first-order vector form, and (ii) first-principles derivations guided by energy-storage variables (mass/spring; inductor/capacitor). We proved equivalence between ODE and state trajectories and implemented the resulting models in Python, C++, Java, MATLAB/Simulink, and Mathematica for simulation and validation.
13. References
- Kalman, R.E. (1960). On the general theory of control systems. Proceedings of the First International Congress of Automatic Control, 481–492.
- Kalman, R.E., & Bertram, J.E. (1960). Control system analysis and design via the “second method” of Lyapunov: Part I: Continuous-time systems. Transactions of the ASME, Journal of Basic Engineering, 82, 371–393.
- Gilbert, E.G. (1963). Controllability and observability in multivariable control systems. SIAM Journal on Control, 1(2), 128–151.
- Ho, B.L., & Kalman, R.E. (1966). Effective construction of linear state-variable models from input/output functions. Regelungstechnik, 14, 545–548.
- Wonham, W.M. (1967). On pole assignment in multi-input controllable linear systems. IEEE Transactions on Automatic Control, 12(6), 660–665.
- Luenberger, D.G. (1964). Observing the state of a linear system. IEEE Transactions on Military Electronics, 8(2), 74–80.
- Brockett, R.W. (1970). Finite Dimensional Linear Systems (selected foundational journal developments summarized therein). Wiley.