A reference-tracking and disturbance-rejection problem is not completely specified by the pair \( (A,B) \). We must also state how the reference enters the controller, how the actuator enters the plant, how disturbances enter the state equation, and how the output is measured. The continuous-time LTI model used in this lesson is

\[ \begin{aligned} \dot{x}(t) &= A x(t) + B u(t) + E d(t),\\ y(t) &= C x(t) + D u(t),\\ e(t) &= r(t)-y(t). \end{aligned} \]

Here \( x\in\mathbb{R}^n \) is the state, \( u\in\mathbb{R}^m \) is the plant input, \( d\in\mathbb{R}^{q} \) is the disturbance, and \( y\in\mathbb{R}^{p} \) is the controlled output. The reference \( r \) is not generally injected into the plant directly; it enters through a command prefilter or an error integrator.

flowchart TD
  R["reference r"] --> REF["prefilter or error integrator"]
  REF --> U["control input u"]
  U --> PLANT["plant: xdot = A x + B u + E d"]
  D["disturbance d"] --> PLANT
  PLANT --> Y["output y = C x + D u"]
  Y --> ERR["tracking error e = r - y"]
  ERR --> REF
        

The central design question is whether the chosen controller can force \( e(t)\to 0 \) for the desired reference and disturbance classes while preserving internal stability.

Consider the normalized mechanical plant \( m\ddot{q}+b\dot{q}+kq=u+d \), where \( u \) is the commanded actuator force and \( d \) is an additive load disturbance. Define \( x_1=q \) and \( x_2=\dot{q} \). Then

\[ x = \begin{bmatrix} q\\ \dot{q} \end{bmatrix}, \qquad \dot{x} = \begin{bmatrix} 0 & 1\\ -k/m & -b/m \end{bmatrix}x + \begin{bmatrix} 0\\ 1/m \end{bmatrix}u + \begin{bmatrix} 0\\ 1/m \end{bmatrix}d, \qquad y = \begin{bmatrix} 1 & 0 \end{bmatrix}x. \]

The disturbance is matched because it enters through the same direction as the actuator input: \( E=B \). A matched disturbance can be cancelled directly if it is measured, but if it is unknown and constant, the controller needs an internal mechanism to estimate or reject its steady effect.

In this lesson we use \( m=1 \), \( b=0.6 \), and \( k=2 \). Thus

\[ A= \begin{bmatrix} 0 & 1\\ -2 & -0.6 \end{bmatrix}, \qquad B=E= \begin{bmatrix} 0\\ 1 \end{bmatrix}, \qquad C= \begin{bmatrix} 1 & 0 \end{bmatrix}. \]

Suppose the controller is \( u=-Kx+N_r r \), where \( K\in\mathbb{R}^{1\times n} \) places the closed-loop poles and \( N_r \) is a scalar reference prefilter. The closed-loop state equation for \( d=0 \) is

\[ \dot{x} = (A-BK)x + B N_r r, \qquad y=Cx. \]

For a constant reference \( r(t)=r_0 \), steady state satisfies

\[ 0=(A-BK)x_{ss}+B N_r r_0, \qquad y_{ss}=C x_{ss}. \]

If \( A-BK \) is nonsingular, then

\[ x_{ss}=-(A-BK)^{-1}B N_r r_0, \qquad y_{ss}=-C(A-BK)^{-1}B N_r r_0. \]

Therefore exact unit DC tracking requires

\[ \boxed{ N_r = -\left[C(A-BK)^{-1}B\right]^{-1} }. \]

This formula is meaningful only when the scalar \( C(A-BK)^{-1}B \) is nonzero. In MIMO systems, the corresponding matrix must have the appropriate rank for the desired tracking channels.

Let \( K=[K_1\;K_2] \). For the mechanical example,

\[ A-BK = \begin{bmatrix} 0 & 1\\ -(k+K_1)/m & -(b+K_2)/m \end{bmatrix}. \]

The characteristic polynomial of \( A-BK \) is

\[ \det(sI-(A-BK)) = s^2+\frac{b+K_2}{m}s+\frac{k+K_1}{m}. \]

If the desired closed-loop polynomial is \( s^2+a_1s+a_0 \), coefficient matching gives

\[ \boxed{ K_1=ma_0-k, \qquad K_2=ma_1-b. } \]

Choosing desired poles \( -2 \) and \( -3 \) gives

\[ (s+2)(s+3)=s^2+5s+6, \qquad K_1=4, \qquad K_2=4.4. \]

For this plant, the static prefilter simplifies to \( N_r=k+K_1 \), so \( N_r=6 \).

Now keep the same controller \( u=-Kx+N_r r \), but let \( d(t)=d_0 \) be a constant unknown disturbance. The closed-loop state equation becomes

\[ \dot{x}=(A-BK)x+B N_r r_0+E d_0. \]

At steady state,

\[ x_{ss}=-(A-BK)^{-1}(B N_r r_0+E d_0), \]

and therefore

\[ y_{ss} = -C(A-BK)^{-1}B N_r r_0 - C(A-BK)^{-1}E d_0. \]

Since \( N_r \) was chosen so that the first term is \( r_0 \), the steady-state tracking error is

\[ e_{ss} = r_0-y_{ss} = C(A-BK)^{-1}E d_0. \]

Thus static feedforward tracking is not, by itself, a disturbance rejection mechanism. It calibrates the reference channel but does not create an internal model for the disturbance channel.

In the present matched case \( E=B \), the error becomes

\[ e_{ss} = C(A-BK)^{-1}B d_0 = -\frac{d_0}{N_r}. \]

Hence a positive load disturbance produces a nonzero output bias unless it is measured and cancelled or rejected by integral action.

To reject constant disturbances and track constant references, augment the plant with an integral of the tracking error:

\[ z(t)=\int_0^t (r(\sigma)-y(\sigma))\,d\sigma, \qquad \dot{z}=r-Cx-Du. \]

For the common strictly proper case \( D=0 \), use

\[ \dot{z}=r-Cx, \qquad u=-K_xx+K_i z. \]

The augmented closed-loop dynamics are

\[ \begin{bmatrix} \dot{x}\\ \dot{z} \end{bmatrix} = \begin{bmatrix} A-BK_x & BK_i\\ -C & 0 \end{bmatrix} \begin{bmatrix} x\\ z \end{bmatrix} + \begin{bmatrix} 0\\ 1 \end{bmatrix}r + \begin{bmatrix} E\\ 0 \end{bmatrix}d. \]

For the mechanical plant, the closed-loop characteristic polynomial is

\[ s^3+\frac{b+K_2}{m}s^2+\frac{k+K_1}{m}s+\frac{K_i}{m}. \]

If the desired integral-servo polynomial is \( s^3+\alpha_2s^2+\alpha_1s+\alpha_0 \), then

\[ \boxed{ K_1=m\alpha_1-k, \qquad K_2=m\alpha_2-b, \qquad K_i=m\alpha_0. } \]

Selecting poles \( -2,-3,-5 \) gives

\[ (s+2)(s+3)(s+5)=s^3+10s^2+31s+30, \]

\[ K_1=29, \qquad K_2=9.4, \qquad K_i=30. \]

Assume the augmented closed-loop matrix is Hurwitz, the reference is constant, and the disturbance is constant. Since the system is stable, every state converges to a finite equilibrium. At equilibrium, \( \dot{z}_{ss}=0 \). But

\[ \dot{z}_{ss}=r_0-Cx_{ss}=r_0-y_{ss}. \]

Therefore

\[ 0=r_0-y_{ss} \qquad\Longrightarrow\qquad \boxed{e_{ss}=0.} \]

This proof is deliberately independent of the exact value of \( d_0 \). The constant disturbance changes the equilibrium value of the integrator state \( z_{ss} \) and the steady control input \( u_{ss} \), but it does not change the final tracking equality \( y_{ss}=r_0 \) as long as the augmented closed loop is internally stable.

This is the state-space form of the internal model principle for constant signals: a stable closed loop must contain an integrator in the error path if it is to reject unknown constant bias disturbances without measuring them directly.

flowchart TD
  S["Start with xdot = A x + B u + E d, y = C x + D u"] --> R["Need exact constant-reference tracking?"]
  R -->|no| STAB["Stabilize with u = -K x"]
  R -->|yes| D0["Is disturbance zero or measured?"]
  D0 -->|yes| FF["Use u = -K x + Nr r"]
  D0 -->|no| CLASS["Disturbance class known? \nconstant / ramp / sinusoid"]
  CLASS --> CONST["For constant disturbance: \nadd error integrator"]
  CONST --> AUG["Place poles of augmented pair"]
  AUG --> CHECK["Check internal stability and actuator effort"]
  FF --> CHECK
  STAB --> CHECK
        

The key practical point is that the reference channel and disturbance channel impose different algebraic requirements. A prefilter adjusts the reference-to-output DC gain. An internal model changes the dynamics so that persistent exogenous signals are rejected through feedback.

For a MIMO plant with \( y\in\mathbb{R}^{p} \) and \( r\in\mathbb{R}^{p} \), a static prefilter \( N_r\in\mathbb{R}^{m\times p} \) is chosen from

\[ y_{ss} = -C(A-BK)^{-1}B N_r r_0. \]

Exact constant-reference tracking requires

\[ -C(A-BK)^{-1}B N_r = I_p. \]

A solution exists if the matrix \( G_0=-C(A-BK)^{-1}B \) has full row rank \( p \). If \( m=p \) and \( G_0 \) is nonsingular, then

\[ N_r=G_0^{-1}. \]

For integral tracking with \( p \) controlled outputs, define \( z\in\mathbb{R}^{p} \) by \( \dot{z}=r-y \). The augmented open-loop pair is

\[ A_a= \begin{bmatrix} A & 0\\ -C & 0 \end{bmatrix}, \qquad B_a= \begin{bmatrix} B\\ -D \end{bmatrix}. \]

The augmented pair must be stabilizable. If it is not, at least one required tracking channel is incompatible with the plant input-output structure, often because of an invariant zero at the origin.

The Python implementation uses scipy.signal.place_poles for pole placement and solve_ivp for numerical integration.

Chapter27_Lesson5.py

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp
from scipy.signal import place_poles


def feedforward_gain(A, B, C, K):
    Acl = A - B @ K
    dc_gain = C @ np.linalg.solve(Acl, B)
    return float(-1.0 / dc_gain)


m, b, k = 1.0, 0.6, 2.0
A = np.array([[0.0, 1.0], [-k / m, -b / m]])
B = np.array([[0.0], [1.0 / m]])
E = B.copy()
C = np.array([[1.0, 0.0]])

K = place_poles(A, B, [-2.0, -3.0]).gain_matrix
Nbar = feedforward_gain(A, B, C, K)

A_aug = np.block([[A, np.zeros((2, 1))], [-C, np.zeros((1, 1))]])
B_aug = np.vstack([B, [[0.0]]])
K_aug = place_poles(A_aug, B_aug, [-2.0, -3.0, -5.0]).gain_matrix
Kx = K_aug[:, :2]
Ki = -K_aug[0, 2]

r_value = 1.0
d_value = 0.4


def simulate_feedforward(d):
    Acl = A - B @ K

    def rhs(t, x):
        return Acl @ x + B.flatten() * Nbar * r_value + E.flatten() * d

    sol = solve_ivp(rhs, [0, 8], np.zeros(2), max_step=0.01)
    y = (C @ sol.y).flatten()
    return sol.t, y


def simulate_integral(d):
    def rhs(t, xa):
        x = xa[:2]
        z = xa[2]
        y = float(C @ x)
        u = float(-Kx @ x + Ki * z)
        return np.r_[A @ x + B.flatten() * u + E.flatten() * d, r_value - y]

    sol = solve_ivp(rhs, [0, 8], np.zeros(3), max_step=0.01)
    y = (C @ sol.y[:2, :]).flatten()
    return sol.t, y


t0, y0 = simulate_feedforward(0.0)
t1, y1 = simulate_feedforward(d_value)
t2, y2 = simulate_integral(d_value)

print("K =", K)
print("Nbar =", Nbar)
print("Kx =", Kx)
print("Ki =", Ki)
print("final y, feedforward no disturbance =", y0[-1])
print("final y, feedforward with disturbance =", y1[-1])
print("final y, integral with disturbance =", y2[-1])

plt.figure()
plt.plot(t0, y0, label="FF, d=0")
plt.plot(t1, y1, label="FF, d=0.4")
plt.plot(t2, y2, label="Integral, d=0.4")
plt.axhline(r_value, linestyle="--", label="reference")
plt.xlabel("time [s]")
plt.ylabel("output y")
plt.title("Reference Tracking and Disturbance Rejection")
plt.legend()
plt.grid(True)
plt.show()
      

The C++ version uses a self-contained RK4 integrator. For this second-order plant, the feedback gains are obtained by coefficient matching.

Chapter27_Lesson5.cpp

#include <fstream>
#include <iostream>

struct Params {
    double m = 1.0;
    double b = 0.6;
    double k = 2.0;
    double r = 1.0;
    double d = 0.4;
};

struct State2 {
    double q;
    double v;
};

struct State3 {
    double q;
    double v;
    double z;
};

State2 ffeed(const State2& x, const Params& p, double K1, double K2, double N, double d) {
    double u = -K1 * x.q - K2 * x.v + N * p.r;
    return {x.v, (-p.k * x.q - p.b * x.v + u + d) / p.m};
}

State3 fint(const State3& x, const Params& p, double K1, double K2, double Ki, double d) {
    double u = -K1 * x.q - K2 * x.v + Ki * x.z;
    return {x.v, (-p.k * x.q - p.b * x.v + u + d) / p.m, p.r - x.q};
}

int main() {
    Params p;

    double K1 = p.m * 6.0 - p.k;
    double K2 = p.m * 5.0 - p.b;
    double Nbar = p.k + K1;

    double K1i = p.m * 31.0 - p.k;
    double K2i = p.m * 10.0 - p.b;
    double Ki = p.m * 30.0;

    std::cout << "Feedforward K1,K2,Nbar = "
              << K1 << ", " << K2 << ", " << Nbar << "\\n";
    std::cout << "Integral K1,K2,Ki = "
              << K1i << ", " << K2i << ", " << Ki << "\\n";

    return 0;
}
      

The downloadable file contains the full RK4 simulation and writes a CSV file with the three output responses.

The Java implementation mirrors the C++ version and is useful when the simulation is embedded into a larger desktop or server-side application.

Chapter27_Lesson5.java

public class Chapter27_Lesson5 {
    static class Params {
        double m = 1.0;
        double b = 0.6;
        double k = 2.0;
        double r = 1.0;
        double d = 0.4;
    }

    static double[] derivIntegral(double[] x, Params p, double K1, double K2, double Ki) {
        double u = -K1 * x[0] - K2 * x[1] + Ki * x[2];
        return new double[] {
            x[1],
            (-p.k * x[0] - p.b * x[1] + u + p.d) / p.m,
            p.r - x[0]
        };
    }

    public static void main(String[] args) {
        Params p = new Params();

        double K1 = p.m * 6.0 - p.k;
        double K2 = p.m * 5.0 - p.b;
        double Nbar = p.k + K1;

        double K1i = p.m * 31.0 - p.k;
        double K2i = p.m * 10.0 - p.b;
        double Ki = p.m * 30.0;

        System.out.println("Feedforward K1,K2,Nbar = " + K1 + ", " + K2 + ", " + Nbar);
        System.out.println("Integral K1,K2,Ki = " + K1i + ", " + K2i + ", " + Ki);
    }
}
      

The downloadable Java file includes the complete RK4 simulation and CSV export.

MATLAB directly supports state-space modeling, pole placement, and simulation. The related Simulink implementation uses a State-Space block for the plant, a Gain block for \( -K_x \), an Integrator block for \( \dot{z}=r-y \), and a Sum block for \( u=-K_xx+K_i z \).

Chapter27_Lesson5.m

clear; clc; close all;

m = 1.0; b = 0.6; k = 2.0;
A = [0 1; -k/m -b/m];
B = [0; 1/m];
E = B;
C = [1 0];

r = 1.0;
d = 0.4;

K = place(A, B, [-2 -3]);
Nbar = -1 / (C * ((A - B*K) \ B));

Aaug = [A zeros(2,1); -C 0];
Baug = [B; 0];
Kaug = place(Aaug, Baug, [-2 -3 -5]);
Kx = Kaug(1:2);
Ki = -Kaug(3);

ff_rhs_d = @(t,x) (A - B*K)*x + B*Nbar*r + E*d;
[t1, x1] = ode45(ff_rhs_d, [0 8], [0;0]);
y1 = (C*x1')';

int_rhs = @(t,xa) [
    A*xa(1:2) + B*(-Kx*xa(1:2) + Ki*xa(3)) + E*d;
    r - C*xa(1:2)
];

[t2, xa2] = ode45(int_rhs, [0 8], [0;0;0]);
y2 = (C*xa2(:,1:2)')';

figure;
plot(t1, y1, 'LineWidth', 1.5); hold on;
plot(t2, y2, 'LineWidth', 1.5);
yline(r, '--');
grid on;
xlabel('time [s]');
ylabel('output y');
legend('feedforward with disturbance', 'integral with disturbance', 'reference');
      

Mathematica can symbolically derive the characteristic polynomial and numerically simulate the nonlinear-looking closed-loop equations using NDSolveValue.

Chapter27_Lesson5.nb

ClearAll["Global`*"];

m = 1.0; b = 0.6; k = 2.0;
r = 1.0; d = 0.4;

K1 = m*6.0 - k;
K2 = m*5.0 - b;
Nbar = k + K1;

K1i = m*31.0 - k;
K2i = m*10.0 - b;
Ki = m*30.0;

integralSolution = NDSolveValue[
  {
    q'[t] == v[t],
    v'[t] == (-k*q[t] - b*v[t] - K1i*q[t] - K2i*v[t] + Ki*z[t] + d)/m,
    z'[t] == r - q[t],
    q[0] == 0,
    v[0] == 0,
    z[0] == 0
  },
  {q, v, z},
  {t, 0, 8}
];

{qint, vint, zint} = integralSolution;

Plot[
  {qint[t], r},
  {t, 0, 8},
  PlotLegends -> {"integral output", "reference"},
  AxesLabel -> {"time [s]", "output y"},
  GridLines -> Automatic
]
      

Problem 1: For the plant \( \dot{x}=Ax+Bu,\; y=Cx \) with controller \( u=-Kx+N_r r \), derive the static prefilter \( N_r \) for exact constant-reference tracking.

Solution: For constant reference, the steady state is

\[ 0=(A-BK)x_{ss}+B N_r r_0. \]

Thus

\[ x_{ss}=-(A-BK)^{-1}B N_r r_0, \qquad y_{ss}=-C(A-BK)^{-1}B N_r r_0. \]

Enforcing \( y_{ss}=r_0 \) gives

\[ \boxed{ N_r=-\left[C(A-BK)^{-1}B\right]^{-1}. } \]

Problem 2: For the mass-spring-damper plant with \( m=1,\;b=0.6,\;k=2 \), compute \( K=[K_1\;K_2] \) so that the closed-loop poles are \( -2 \) and \( -3 \).

Solution: The desired polynomial is

\[ (s+2)(s+3)=s^2+5s+6. \]

Since \( \det(sI-(A-BK))=s^2+(b+K_2)s+(k+K_1) \) for \( m=1 \), coefficient matching gives

\[ K_1=6-2=4, \qquad K_2=5-0.6=4.4. \]

Problem 3: Show that a constant unknown disturbance generally produces steady-state error under \( u=-Kx+N_r r \).

Solution: With constant disturbance,

\[ x_{ss}=-(A-BK)^{-1}(B N_r r_0+E d_0). \]

Therefore

\[ e_{ss} = r_0-y_{ss} = C(A-BK)^{-1}E d_0. \]

Unless this term is zero for structural reasons, the feedforward-only design does not reject the disturbance.

Problem 4: For the same plant, design an integral servo with desired poles \( -2,-3,-5 \).

Solution: The desired polynomial is

\[ (s+2)(s+3)(s+5)=s^3+10s^2+31s+30. \]

The augmented characteristic polynomial is \( s^3+(b+K_2)s^2+(k+K_1)s+K_i \) for \( m=1 \). Hence

\[ K_1=31-2=29, \qquad K_2=10-0.6=9.4, \qquad K_i=30. \]

Problem 5: Prove that, if the integral-servo closed loop is internally stable, then constant-reference tracking is exact despite an unknown constant disturbance.

Solution: Internal stability implies convergence to an equilibrium. At equilibrium,

\[ 0=\dot{z}_{ss}=r_0-y_{ss}. \]

Therefore \( y_{ss}=r_0 \) and \( e_{ss}=0 \). The disturbance changes the required equilibrium value of the control input but not the final tracking equality.

This lesson compared three important cases. First, static state feedback with a reference prefilter gives exact constant-reference tracking when no persistent disturbance is present and the DC gain is invertible. Second, the same prefilter does not generally reject unknown constant disturbances because it only calibrates the reference channel. Third, integral augmentation adds an internal model of constant signals and can remove steady-state error for both constant references and constant disturbances, provided the augmented closed loop is internally stable.