After completing this lesson, the student should be able to:

• separate nominal stabilization, adaptive compensation, and robust protection into distinct design roles;
• derive the tracking-error dynamics for a baseline controller with adaptive augmentation;
• construct a leakage- and projection-modified adaptive law;
• prove uniform ultimate boundedness in the presence of bounded matched and unmatched disturbances;
• explain why actuator limits, filters, dead zones, and supervisory logic must be included in practical implementations; and
• implement a representative combined architecture in Python, C++, Java, MATLAB/Simulink, and Wolfram Mathematica.

The fundamental engineering principle is role separation. The nominal feedback loop must stabilize the known plant model without relying on parameter adaptation. The adaptive channel should compensate only for uncertainty that is represented by the chosen regressor. The robust channel should attenuate bounded residual effects that the parametric model cannot cancel exactly. Finally, a supervisor should prevent adaptation from reacting to invalid data, actuator saturation, or operating conditions outside the model's validity region.

flowchart TD
  R["Command r"] --> RM["Reference model"]
  R --> NOM["Nominal stabilizing controller"]
  X["Measured state x"] --> NOM
  X --> REG["Regressor phi"]
  E["Tracking error e"] --> ADAPT["Adaptive law with \nleakage and projection"]
  REG --> ADAPT
  ADAPT --> AUG["Adaptive \ncancellation"]
  E --> ROB["Robust bounded \ninjection"]
  NOM --> SUM["Control sum"]
  AUG --> SUM
  ROB --> SUM
  SUP["Supervisor: limits, \nvalidity, freeze, reset"] --> ADAPT
  SUP --> SUM
  SUM --> ACT["Actuator saturation and rate limits"]
  ACT --> PLANT["Uncertain plant"]
  PLANT --> X
  RM --> E
  X --> E
        

Consider the state-space plant \( \mathbf{x}\in\mathbb{R}^{n} \) with a scalar control input:

\[ \dot{\mathbf{x} } = \mathbf{A}\mathbf{x} + \mathbf{B}\left[u + {\boldsymbol{\theta}^{*} }^{T}\boldsymbol{\phi}(\mathbf{x},t) + d_m(t)\right] + \mathbf{B}_{\perp}d_u(t). \]

The unknown constant vector \( \boldsymbol{\theta}^{*}\in\mathbb{R}^{p} \) multiplies the known regressor \( \boldsymbol{\phi}\in\mathbb{R}^{p} \). This term and \(d_m\) are matched because they enter through the control direction \(\mathbf{B}\). The term \(\mathbf{B}_{\perp}d_u\) is unmatched and cannot, in general, be canceled by the scalar input. Assume known bounds

\[ |d_m(t)|\leq \bar d_m, \qquad |d_u(t)|\leq \bar d_u, \qquad \|\boldsymbol{\theta}^{*}\|\leq \bar\theta. \]

The desired dynamics are specified by the stable reference model

\[ \dot{\mathbf{x} }_m = \mathbf{A}_m\mathbf{x}_m + \mathbf{B}_m r, \qquad \mathbf{A}_m \text{ Hurwitz}. \]

Select nominal gains \(\mathbf{K}_x\) and \(K_r\) satisfying the matching equations

\[ \mathbf{A}+\mathbf{B}\mathbf{K}_x = \mathbf{A}_m, \qquad \mathbf{B}K_r = \mathbf{B}_m. \]

These equations make the nominal closed loop reproduce the reference model when the uncertainty and disturbances are absent. The tracking error is \(\mathbf{e}=\mathbf{x}-\mathbf{x}_m\).

Decompose the input into nominal, adaptive, and robust components:

\[ u = u_{nom}+u_{ad}+u_{rob}, \]

\[ u_{nom}=\mathbf{K}_x\mathbf{x}+K_r r, \qquad u_{ad}=-\hat{\boldsymbol{\theta} }^{T}\boldsymbol{\phi}, \qquad u_{rob}=-k_s\operatorname{sat}\!\left(\frac{s}{\varepsilon}\right). \]

Let \(\mathbf{P}=\mathbf{P}^{T}>0\) solve the Lyapunov equation

\[ \mathbf{A}_m^{T}\mathbf{P}+\mathbf{P}\mathbf{A}_m=-\mathbf{Q}, \qquad \mathbf{Q}=\mathbf{Q}^{T}>0, \]

and define the scalar robust/adaptive error signal

\[ s=\mathbf{e}^{T}\mathbf{P}\mathbf{B}. \]

The continuous saturation function is

\[ \operatorname{sat}(z)= \begin{cases} -1, & z<-1,\\ z, & |z|\leq 1,\\ 1, & z>1. \end{cases} \]

The boundary-layer width \(\varepsilon>0\) avoids an ideal discontinuous sign function. The price is a nonzero residual tracking set. The robust gain \(k_s\) is selected from a disturbance bound and actuator constraints rather than increased without limit.

With \(\tilde{\boldsymbol{\theta} }=\boldsymbol{\theta}^{*}-\hat{\boldsymbol{\theta} }\), the error dynamics become

\[ \dot{\mathbf{e} } = \mathbf{A}_m\mathbf{e} + \mathbf{B}\left[ \tilde{\boldsymbol{\theta} }^{T}\boldsymbol{\phi} + d_m - k_s\operatorname{sat}\!\left(\frac{s}{\varepsilon}\right) \right] + \mathbf{B}_{\perp}d_u. \]

A practical continuous-time law is

\[ \dot{\hat{\boldsymbol{\theta} } } = \operatorname{Proj}_{\Omega}\!\left( \hat{\boldsymbol{\theta} }, \boldsymbol{\Gamma}\boldsymbol{\phi}s -\sigma\hat{\boldsymbol{\theta} } \right), \qquad \boldsymbol{\Gamma}=\boldsymbol{\Gamma}^{T}>0, \quad \sigma>0. \]

The compact convex set \(\Omega\) contains the physically plausible parameter values. Projection prevents the estimate from leaving this set. For Lyapunov analysis, the projection operator is chosen to satisfy

\[ \tilde{\boldsymbol{\theta} }^{T}\boldsymbol{\Gamma}^{-1} \left[ \operatorname{Proj}_{\Omega}(\hat{\boldsymbol{\theta} },\mathbf{y})-\mathbf{y} \right]\geq 0. \]

This inequality ensures that projection does not make the standard Lyapunov derivative worse. The leakage term \(-\sigma\hat{\boldsymbol{\theta} }\) suppresses drift when disturbances continuously excite the update law. A dead-zone implementation replaces the gradient term by zero whenever \(|s|\leq s_0\):

\[ \dot{\hat{\boldsymbol{\theta} } } = \operatorname{Proj}_{\Omega}\!\left( \hat{\boldsymbol{\theta} }, \chi_{\{|s|>s_0\} }\boldsymbol{\Gamma}\boldsymbol{\phi}s -\sigma\hat{\boldsymbol{\theta} } \right). \]

Dead zones reduce noise-driven adaptation but enlarge the residual set. Leakage improves boundedness but biases the estimate toward zero. Projection enforces physical limits but does not identify the true parameter unless the regressor is sufficiently exciting.

Use the candidate

\[ V(\mathbf{e},\tilde{\boldsymbol{\theta} }) = \mathbf{e}^{T}\mathbf{P}\mathbf{e} + \tilde{\boldsymbol{\theta} }^{T}\boldsymbol{\Gamma}^{-1} \tilde{\boldsymbol{\theta} }. \]

Ignoring only the nonpositive contribution supplied by projection, differentiation gives

\[ \begin{aligned} \dot V \leq {}& -\mathbf{e}^{T}\mathbf{Q}\mathbf{e} +2s\left[d_m-k_s\operatorname{sat}\!\left(\frac{s}{\varepsilon}\right)\right]\\ &+2\mathbf{e}^{T}\mathbf{P}\mathbf{B}_{\perp}d_u +2\sigma\tilde{\boldsymbol{\theta} }^{T} \boldsymbol{\Gamma}^{-1}\hat{\boldsymbol{\theta} }. \end{aligned} \]

The adaptive cross term \(2s\tilde{\boldsymbol{\theta} }^{T}\boldsymbol{\phi}\) has canceled exactly. Since \(\hat{\boldsymbol{\theta} }=\boldsymbol{\theta}^{*}-\tilde{\boldsymbol{\theta} }\),

\[ \begin{aligned} 2\sigma\tilde{\boldsymbol{\theta} }^{T}\boldsymbol{\Gamma}^{-1} \hat{\boldsymbol{\theta} } &=2\sigma\tilde{\boldsymbol{\theta} }^{T}\boldsymbol{\Gamma}^{-1} \boldsymbol{\theta}^{*} -2\sigma\|\tilde{\boldsymbol{\theta} }\|_{\boldsymbol{\Gamma}^{-1} }^{2}\\ &\leq -\sigma\|\tilde{\boldsymbol{\theta} }\|_{\boldsymbol{\Gamma}^{-1} }^{2} +\sigma\|\boldsymbol{\theta}^{*}\|_{\boldsymbol{\Gamma}^{-1} }^{2}. \end{aligned} \]

For \(|s|>\varepsilon\), the saturation term equals \(\operatorname{sgn}(s)\). If \(k_s\geq\bar d_m+\delta\) for some \(\delta>0\), then

\[ 2s\left[d_m-k_s\operatorname{sgn}(s)\right] \leq -2\delta |s|. \]

Inside the boundary layer, a conservative bound is

\[ 2s\left[d_m-k_s\operatorname{sat}\!\left(\frac{s}{\varepsilon}\right)\right] \leq 2(\bar d_m+k_s)\varepsilon. \]

For the unmatched disturbance, Young's inequality yields

\[ 2\mathbf{e}^{T}\mathbf{P}\mathbf{B}_{\perp}d_u \leq \frac{\lambda_{min}(\mathbf{Q})}{2}\|\mathbf{e}\|^{2} +\frac{2\|\mathbf{P}\mathbf{B}_{\perp}\|^{2} }{\lambda_{min}(\mathbf{Q})} \bar d_u^{2}. \]

Combining the bounds gives

\[ \dot V \leq -\frac{\lambda_{min}(\mathbf{Q})}{2}\|\mathbf{e}\|^{2} -\sigma\|\tilde{\boldsymbol{\theta} }\|_{\boldsymbol{\Gamma}^{-1} }^{2} +C, \]

\[ C=2(\bar d_m+k_s)\varepsilon +\frac{2\|\mathbf{P}\mathbf{B}_{\perp}\|^{2} }{\lambda_{min}(\mathbf{Q})} \bar d_u^{2} +\sigma\|\boldsymbol{\theta}^{*}\|_{\boldsymbol{\Gamma}^{-1} }^{2}. \]

Using the eigenvalue bounds on \(V\), define

\[ \alpha= \min\!\left\{ \frac{\lambda_{min}(\mathbf{Q})}{2\lambda_{max}(\mathbf{P})}, \sigma \right\}. \]

Then

\[ \dot V\leq -\alpha V+C, \qquad V(t)\leq e^{-\alpha t}V(0)+\frac{C}{\alpha}\left(1-e^{-\alpha t}\right). \]

Therefore, the tracking and parameter errors are uniformly ultimately bounded and converge to the invariant residual set

\[ \mathcal{S}_{UUB}= \left\{(\mathbf{e},\tilde{\boldsymbol{\theta} }): V\leq \frac{C}{\alpha}\right\}. \]

Exact convergence is not claimed because persistent disturbances, a finite boundary layer, leakage, unmatched uncertainty, and imperfect excitation produce nonzero residual terms. The proof clarifies the design trade-off: reducing \(\varepsilon\) and disturbance bounds shrinks the residual set, but excessive robust gain or adaptation gain may violate actuator and bandwidth limits.

6.1 Nominal stabilizing layer

The nominal controller determines the desired poles, command response, and basic stability margin for the model used in design. It must remain stabilizing when the adaptive term is temporarily frozen. Adaptation should augment a valid controller, not replace the need for one.

6.2 Adaptive layer

The adaptive term is effective for structured uncertainty represented by \({\boldsymbol{\theta}^{*} }^{T}\boldsymbol{\phi}\). It cannot generally cancel arbitrary unmatched disturbances, sensor dynamics, delays, or neglected high-frequency modes. Parameter convergence also requires excitation; bounded tracking alone does not imply \(\hat{\boldsymbol{\theta} }\to\boldsymbol{\theta}^{*}\).

6.3 Robust layer

The robust injection suppresses bounded residual matched disturbances and makes the tracking behavior less dependent on perfect parameter estimates. Continuous saturation replaces ideal switching to reduce chattering. Its gain is constrained by actuator magnitude, actuator rate, measurement noise, and unmodeled dynamics.

6.4 Supervisory layer

A supervisor detects conditions under which the adaptive model is unreliable. Typical actions include freezing the gradient term during actuator saturation, retaining leakage to relax estimates, resetting estimates after a validated mode change, limiting parameter rates, and reverting to a certified baseline controller.

Several architectures share the same conceptual separation while differing in where filtering and supervision are placed.

1. Baseline feedback plus adaptive feedforward: the simplest structure used in the proof above.
2. Baseline feedback plus robust injection plus adaptive cancellation: directly handles a known bound on matched residual disturbances.
3. Filtered adaptive augmentation: a stable low-pass filter limits the bandwidth of the adaptive control signal. This can reduce excitation of neglected high-frequency dynamics, but the filter must be included in the stability analysis.
4. Supervised adaptive augmentation: the baseline controller remains continuously available and the adaptive channel is enabled only inside a validated envelope.
5. Multiple-model or scheduled robust-adaptive control: a supervisor selects a nominal model/controller region while adaptation compensates local parametric mismatch.

A filtered adaptive channel may be written as

\[ \tau_f\dot u_{ad,f}+u_{ad,f} =-\hat{\boldsymbol{\theta} }^{T}\boldsymbol{\phi}, \qquad u=u_{nom}+u_{ad,f}+u_{rob}. \]

The filter separates rapid parameter estimation from the actuator command bandwidth. However, it also introduces additional dynamics and destroys the exact cross-term cancellation used in the preceding proof unless those dynamics are explicitly included in the Lyapunov model. Thus, filtering is an architectural choice, not a free modification.

flowchart TD
  START["Start with validated nominal controller"] --> CHECK1["Is uncertainty structured and matched?"]
  CHECK1 -->|"no"| ROBONLY["Use robust baseline; \ndo not expect adaptation \nto cancel it"]
  CHECK1 -->|"yes"| ADAPT["Add adaptive cancellation"]
  ADAPT --> CHECK2["Are disturbances, \nnoise, or unmodeled \ndynamics significant?"]
  CHECK2 -->|"yes"| MODS["Add leakage, projection, \ndead zone, normalization"]
  CHECK2 -->|"no"| BASIC["Retain bounded parameter \nset and monitoring"]
  MODS --> CHECK3["Can adaptive signal \nexcite fast dynamics \nor actuators?"]
  BASIC --> CHECK3
  CHECK3 -->|"yes"| FILTER["Add analyzed command \nfilter and rate limits"]
  CHECK3 -->|"no"| SUP["Add supervisor \nand fallback logic"]
  FILTER --> SUP
  SUP --> VERIFY["Verify bounds, saturation cases, \nsampling, and failure modes"]
        

8.1 Tune the baseline before adaptation

Verify nominal stability, command tracking, and actuator margins with \(u_{ad}=0\). Reserve control authority for robust and adaptive augmentation.

8.2 Choose parameter bounds and regressor scaling

Define \(\Omega\) from physical knowledge. Normalize regressors when their amplitudes differ substantially. Poor scaling makes one adaptation gain dominate the others.

8.3 Select adaptation gain

Larger \(\boldsymbol{\Gamma}\) accelerates parameter motion but also increases sensitivity to noise, sampling, delay, and neglected dynamics. Adaptation speed must remain compatible with the measurement and actuator bandwidths.

8.4 Select leakage and dead zone

Increase \(\sigma\) when drift is observed, but recognize the resulting parameter bias. Choose the dead-zone threshold above the effective noise floor of \(s\), not from numerical convenience alone.

8.5 Select robust gain and boundary layer

A theoretical starting point is \(k_s>\bar d_m\). Then verify that the associated input and input-rate demands remain feasible. Decreasing \(\varepsilon\) tightens tracking but raises control activity and sensitivity to measurement noise.

8.6 Test invalid operating conditions

Required tests include command steps, parameter jumps within the assumed set, disturbance pulses, sensor noise, actuator saturation, actuator rate limiting, frozen adaptation, delayed measurements, loss of excitation, and operation near projection boundaries.

The Python example uses only the standard library for simulation and CSV export, plus matplotlib for plotting. It simulates a scalar uncertain plant, robust saturation injection, leakage, interval projection, a dead zone, and a supervisor that freezes the gradient term while the actuator is saturated.

Chapter22_Lesson5.py

"""Chapter22_Lesson5.py
Combined robust + adaptive control for a scalar uncertain plant.

Plant:
    x_dot = -a_nom*x + u + theta_true*x + d(t)
Reference model:
    xm_dot = -a_m*xm + b_m*r(t)
Controller:
    u = (a_nom-a_m)*x + b_m*r - theta_hat*x
        - k_s*sat((x-xm)/epsilon)
Adaptive law:
    theta_hat_dot = gamma*e*x - sigma*theta_hat
with projection and a saturation-aware supervisory freeze.
"""

from __future__ import annotations

import csv
import math
from dataclasses import dataclass
from pathlib import Path

import matplotlib.pyplot as plt


@dataclass(frozen=True)
class Config:
    dt: float = 0.001
    final_time: float = 30.0
    a_nom: float = 1.0
    a_m: float = 1.5
    b_m: float = 1.5
    theta_true: float = 0.75
    gamma: float = 8.0
    sigma: float = 0.12
    k_s: float = 0.90
    epsilon: float = 0.08
    dead_zone: float = 0.01
    theta_min: float = -2.0
    theta_max: float = 2.0
    u_max: float = 6.0


def clamp(value: float, lower: float, upper: float) -> float:
    return max(lower, min(upper, value))


def sat(value: float) -> float:
    return clamp(value, -1.0, 1.0)


def reference(t: float) -> float:
    if t < 2.0:
        return 0.0
    return 1.0 + 0.25 * math.sin(0.45 * t)


def disturbance(t: float) -> float:
    pulse = 0.45 if 12.0 <= t <= 14.0 else 0.0
    return 0.30 * math.sin(4.0 * t) + pulse


def projected_derivative(theta_hat: float, raw_dot: float, cfg: Config) -> float:
    """Tangent-cone projection for a scalar interval."""
    if theta_hat >= cfg.theta_max and raw_dot > 0.0:
        return 0.0
    if theta_hat <= cfg.theta_min and raw_dot < 0.0:
        return 0.0
    return raw_dot


def simulate(cfg: Config) -> dict[str, list[float]]:
    x = 0.0
    xm = 0.0
    theta_hat = 0.0

    history = {
        "time": [],
        "reference": [],
        "x": [],
        "xm": [],
        "error": [],
        "theta_hat": [],
        "control": [],
        "disturbance": [],
        "adaptation_enabled": [],
    }

    steps = int(round(cfg.final_time / cfg.dt))
    for step in range(steps + 1):
        t = step * cfg.dt
        r = reference(t)
        d = disturbance(t)
        e = x - xm

        u_nom = (cfg.a_nom - cfg.a_m) * x + cfg.b_m * r
        u_ad = -theta_hat * x
        u_rob = -cfg.k_s * sat(e / cfg.epsilon)
        u_unsat = u_nom + u_ad + u_rob
        u = clamp(u_unsat, -cfg.u_max, cfg.u_max)
        actuator_saturated = abs(u - u_unsat) > 1.0e-12

        gradient = cfg.gamma * e * x if abs(e) > cfg.dead_zone else 0.0
        if actuator_saturated:
            gradient = 0.0
        raw_theta_dot = gradient - cfg.sigma * theta_hat
        theta_dot = projected_derivative(theta_hat, raw_theta_dot, cfg)

        x_dot = -cfg.a_nom * x + u + cfg.theta_true * x + d
        xm_dot = -cfg.a_m * xm + cfg.b_m * r

        history["time"].append(t)
        history["reference"].append(r)
        history["x"].append(x)
        history["xm"].append(xm)
        history["error"].append(e)
        history["theta_hat"].append(theta_hat)
        history["control"].append(u)
        history["disturbance"].append(d)
        history["adaptation_enabled"].append(0.0 if actuator_saturated else 1.0)

        x += cfg.dt * x_dot
        xm += cfg.dt * xm_dot
        theta_hat = clamp(
            theta_hat + cfg.dt * theta_dot,
            cfg.theta_min,
            cfg.theta_max,
        )

    return history


def write_csv(history: dict[str, list[float]], path: Path) -> None:
    columns = list(history.keys())
    with path.open("w", newline="", encoding="utf-8") as handle:
        writer = csv.writer(handle)
        writer.writerow(columns)
        writer.writerows(zip(*(history[column] for column in columns)))


def rms(values: list[float]) -> float:
    return math.sqrt(sum(value * value for value in values) / len(values))


def make_plot(history: dict[str, list[float]], path: Path) -> None:
    t = history["time"]

    fig, axes = plt.subplots(4, 1, figsize=(10, 11), sharex=True)
    axes[0].plot(t, history["x"], label="plant state x")
    axes[0].plot(t, history["xm"], label="reference-model state xm")
    axes[0].plot(t, history["reference"], "--", label="command r")
    axes[0].set_ylabel("state")
    axes[0].legend()
    axes[0].grid(True)

    axes[1].plot(t, history["error"])
    axes[1].axhline(0.0, linewidth=0.8)
    axes[1].set_ylabel("e = x - xm")
    axes[1].grid(True)

    axes[2].plot(t, history["theta_hat"], label="theta_hat")
    axes[2].axhline(0.75, linestyle="--", label="theta_true")
    axes[2].set_ylabel("parameter")
    axes[2].legend()
    axes[2].grid(True)

    axes[3].plot(t, history["control"], label="u")
    axes[3].plot(t, history["disturbance"], label="d(t)")
    axes[3].set_xlabel("time [s]")
    axes[3].set_ylabel("input")
    axes[3].legend()
    axes[3].grid(True)

    fig.tight_layout()
    fig.savefig(path, dpi=180)
    plt.close(fig)


def main() -> None:
    cfg = Config()
    history = simulate(cfg)

    output_dir = Path(__file__).resolve().parent
    csv_path = output_dir / "Chapter22_Lesson5_results.csv"
    plot_path = output_dir / "Chapter22_Lesson5_results.png"
    write_csv(history, csv_path)
    make_plot(history, plot_path)

    error_rms = rms(history["error"])
    max_error = max(abs(value) for value in history["error"])
    max_control = max(abs(value) for value in history["control"])
    final_theta = history["theta_hat"][-1]

    print(f"RMS tracking error: {error_rms:.6f}")
    print(f"Maximum absolute tracking error: {max_error:.6f}")
    print(f"Maximum absolute control: {max_control:.6f}")
    print(f"Final theta_hat: {final_theta:.6f}")
    print(f"Wrote: {csv_path}")
    print(f"Wrote: {plot_path}")


if __name__ == "__main__":
    main()

The C++17 implementation has no external dependencies and writes a CSV file for plotting in any numerical package.

Chapter22_Lesson5.cpp

// Chapter22_Lesson5.cpp
// Combined robust + adaptive control for a scalar uncertain plant.

#include <algorithm>
#include <cmath>
#include <fstream>
#include <iomanip>
#include <iostream>
#include <stdexcept>
#include <string>

struct Config {
    double dt = 0.001;
    double finalTime = 30.0;
    double aNom = 1.0;
    double aM = 1.5;
    double bM = 1.5;
    double thetaTrue = 0.75;
    double gamma = 8.0;
    double sigma = 0.12;
    double kS = 0.90;
    double epsilon = 0.08;
    double deadZone = 0.01;
    double thetaMin = -2.0;
    double thetaMax = 2.0;
    double uMax = 6.0;
};

double clampValue(double value, double lower, double upper) {
    return std::max(lower, std::min(upper, value));
}

double sat(double value) {
    return clampValue(value, -1.0, 1.0);
}

double reference(double t) {
    return (t < 2.0) ? 0.0 : 1.0 + 0.25 * std::sin(0.45 * t);
}

double disturbance(double t) {
    const double pulse = (t >= 12.0 && t <= 14.0) ? 0.45 : 0.0;
    return 0.30 * std::sin(4.0 * t) + pulse;
}

double projectedDerivative(double thetaHat, double rawDot, const Config& cfg) {
    if (thetaHat >= cfg.thetaMax && rawDot > 0.0) {
        return 0.0;
    }
    if (thetaHat <= cfg.thetaMin && rawDot < 0.0) {
        return 0.0;
    }
    return rawDot;
}

int main() {
    try {
        const Config cfg;
        std::ofstream csv("Chapter22_Lesson5_results.csv");
        if (!csv) {
            throw std::runtime_error("Cannot open Chapter22_Lesson5_results.csv");
        }

        csv << "time,reference,x,xm,error,theta_hat,control,disturbance,adaptation_enabled\n";
        csv << std::setprecision(12);

        double x = 0.0;
        double xm = 0.0;
        double thetaHat = 0.0;
        double sumErrorSquared = 0.0;
        double maxError = 0.0;
        double maxControl = 0.0;
        const int steps = static_cast<int>(std::llround(cfg.finalTime / cfg.dt));

        for (int step = 0; step <= steps; ++step) {
            const double t = step * cfg.dt;
            const double r = reference(t);
            const double d = disturbance(t);
            const double e = x - xm;

            const double uNom = (cfg.aNom - cfg.aM) * x + cfg.bM * r;
            const double uAd = -thetaHat * x;
            const double uRob = -cfg.kS * sat(e / cfg.epsilon);
            const double uUnsat = uNom + uAd + uRob;
            const double u = clampValue(uUnsat, -cfg.uMax, cfg.uMax);
            const bool actuatorSaturated = std::abs(u - uUnsat) > 1.0e-12;

            double gradient = (std::abs(e) > cfg.deadZone) ? cfg.gamma * e * x : 0.0;
            if (actuatorSaturated) {
                gradient = 0.0;
            }
            const double rawThetaDot = gradient - cfg.sigma * thetaHat;
            const double thetaDot = projectedDerivative(thetaHat, rawThetaDot, cfg);

            const double xDot = -cfg.aNom * x + u + cfg.thetaTrue * x + d;
            const double xmDot = -cfg.aM * xm + cfg.bM * r;

            csv << t << ',' << r << ',' << x << ',' << xm << ',' << e << ','
                << thetaHat << ',' << u << ',' << d << ','
                << (actuatorSaturated ? 0 : 1) << '\n';

            sumErrorSquared += e * e;
            maxError = std::max(maxError, std::abs(e));
            maxControl = std::max(maxControl, std::abs(u));

            x += cfg.dt * xDot;
            xm += cfg.dt * xmDot;
            thetaHat = clampValue(thetaHat + cfg.dt * thetaDot,
                                  cfg.thetaMin, cfg.thetaMax);
        }

        const double rmsError = std::sqrt(sumErrorSquared / static_cast<double>(steps + 1));
        std::cout << std::fixed << std::setprecision(6)
                  << "RMS tracking error: " << rmsError << '\n'
                  << "Maximum absolute tracking error: " << maxError << '\n'
                  << "Maximum absolute control: " << maxControl << '\n'
                  << "Final theta_hat: " << thetaHat << '\n'
                  << "Wrote: Chapter22_Lesson5_results.csv\n";
        return 0;
    } catch (const std::exception& error) {
        std::cerr << "Error: " << error.what() << '\n';
        return 1;
    }
}

The Java implementation uses the standard JDK, performs the same fixed-step simulation, and exports matching result columns.

Chapter22_Lesson5.java

// Chapter22_Lesson5.java
// Combined robust + adaptive control for a scalar uncertain plant.

import java.io.BufferedWriter;
import java.io.IOException;
import java.nio.charset.StandardCharsets;
import java.nio.file.Files;
import java.nio.file.Path;
import java.util.Locale;

public final class Chapter22_Lesson5 {
    private static final class Config {
        final double dt = 0.001;
        final double finalTime = 30.0;
        final double aNom = 1.0;
        final double aM = 1.5;
        final double bM = 1.5;
        final double thetaTrue = 0.75;
        final double gamma = 8.0;
        final double sigma = 0.12;
        final double kS = 0.90;
        final double epsilon = 0.08;
        final double deadZone = 0.01;
        final double thetaMin = -2.0;
        final double thetaMax = 2.0;
        final double uMax = 6.0;
    }

    private Chapter22_Lesson5() {
    }

    private static double clamp(double value, double lower, double upper) {
        return Math.max(lower, Math.min(upper, value));
    }

    private static double sat(double value) {
        return clamp(value, -1.0, 1.0);
    }

    private static double reference(double t) {
        return t < 2.0 ? 0.0 : 1.0 + 0.25 * Math.sin(0.45 * t);
    }

    private static double disturbance(double t) {
        double pulse = (t >= 12.0 && t <= 14.0) ? 0.45 : 0.0;
        return 0.30 * Math.sin(4.0 * t) + pulse;
    }

    private static double projectedDerivative(double thetaHat, double rawDot, Config cfg) {
        if (thetaHat >= cfg.thetaMax && rawDot > 0.0) {
            return 0.0;
        }
        if (thetaHat <= cfg.thetaMin && rawDot < 0.0) {
            return 0.0;
        }
        return rawDot;
    }

    public static void main(String[] args) {
        Locale.setDefault(Locale.US);
        Config cfg = new Config();
        Path output = Path.of("Chapter22_Lesson5_results.csv");

        double x = 0.0;
        double xm = 0.0;
        double thetaHat = 0.0;
        double sumErrorSquared = 0.0;
        double maxError = 0.0;
        double maxControl = 0.0;
        int steps = (int) Math.round(cfg.finalTime / cfg.dt);

        try (BufferedWriter writer = Files.newBufferedWriter(
                output, StandardCharsets.UTF_8)) {
            writer.write("time,reference,x,xm,error,theta_hat,control,disturbance,adaptation_enabled\n");

            for (int step = 0; step <= steps; step++) {
                double t = step * cfg.dt;
                double r = reference(t);
                double d = disturbance(t);
                double e = x - xm;

                double uNom = (cfg.aNom - cfg.aM) * x + cfg.bM * r;
                double uAd = -thetaHat * x;
                double uRob = -cfg.kS * sat(e / cfg.epsilon);
                double uUnsat = uNom + uAd + uRob;
                double u = clamp(uUnsat, -cfg.uMax, cfg.uMax);
                boolean actuatorSaturated = Math.abs(u - uUnsat) > 1.0e-12;

                double gradient = Math.abs(e) > cfg.deadZone ? cfg.gamma * e * x : 0.0;
                if (actuatorSaturated) {
                    gradient = 0.0;
                }
                double rawThetaDot = gradient - cfg.sigma * thetaHat;
                double thetaDot = projectedDerivative(thetaHat, rawThetaDot, cfg);

                double xDot = -cfg.aNom * x + u + cfg.thetaTrue * x + d;
                double xmDot = -cfg.aM * xm + cfg.bM * r;

                writer.write(String.format(Locale.US,
                        "%.12f,%.12f,%.12f,%.12f,%.12f,%.12f,%.12f,%.12f,%d%n",
                        t, r, x, xm, e, thetaHat, u, d,
                        actuatorSaturated ? 0 : 1));

                sumErrorSquared += e * e;
                maxError = Math.max(maxError, Math.abs(e));
                maxControl = Math.max(maxControl, Math.abs(u));

                x += cfg.dt * xDot;
                xm += cfg.dt * xmDot;
                thetaHat = clamp(thetaHat + cfg.dt * thetaDot,
                        cfg.thetaMin, cfg.thetaMax);
            }
        } catch (IOException error) {
            System.err.println("I/O error: " + error.getMessage());
            System.exit(1);
        }

        double rmsError = Math.sqrt(sumErrorSquared / (steps + 1.0));
        System.out.printf(Locale.US, "RMS tracking error: %.6f%n", rmsError);
        System.out.printf(Locale.US, "Maximum absolute tracking error: %.6f%n", maxError);
        System.out.printf(Locale.US, "Maximum absolute control: %.6f%n", maxControl);
        System.out.printf(Locale.US, "Final theta_hat: %.6f%n", thetaHat);
        System.out.println("Wrote: " + output.toAbsolutePath());
    }
}

The MATLAB script simulates the combined controller, writes a table, and produces diagnostic plots. In Simulink, implement the same equations with a Plant Integrator, Reference-Model Integrator, Sum blocks for \(u_{nom}\), \(u_{ad}\), and \(u_{rob}\), a Saturation block, and a MATLAB Function block for projection and supervisory adaptation. Feed the actuator saturation flag back to the adaptive-law block so the gradient term is frozen whenever the commanded input cannot be realized.

Chapter22_Lesson5.m

% Chapter22_Lesson5.m
% Combined robust + adaptive control for a scalar uncertain plant.

clear; clc; close all;

cfg.dt = 0.001;
cfg.finalTime = 30.0;
cfg.aNom = 1.0;
cfg.aM = 1.5;
cfg.bM = 1.5;
cfg.thetaTrue = 0.75;
cfg.gamma = 8.0;
cfg.sigma = 0.12;
cfg.kS = 0.90;
cfg.epsilon = 0.08;
cfg.deadZone = 0.01;
cfg.thetaMin = -2.0;
cfg.thetaMax = 2.0;
cfg.uMax = 6.0;

time = (0:cfg.dt:cfg.finalTime)';
numberOfSamples = numel(time);
referenceSignal = zeros(numberOfSamples, 1);
x = zeros(numberOfSamples, 1);
xm = zeros(numberOfSamples, 1);
errorSignal = zeros(numberOfSamples, 1);
thetaHat = zeros(numberOfSamples, 1);
control = zeros(numberOfSamples, 1);
disturbanceSignal = zeros(numberOfSamples, 1);
adaptationEnabled = ones(numberOfSamples, 1);

for k = 1:numberOfSamples - 1
    t = time(k);
    if t < 2.0
        r = 0.0;
    else
        r = 1.0 + 0.25 * sin(0.45 * t);
    end

    pulse = 0.45 * double(t >= 12.0 && t <= 14.0);
    d = 0.30 * sin(4.0 * t) + pulse;
    e = x(k) - xm(k);

    satValue = max(-1.0, min(1.0, e / cfg.epsilon));
    uNom = (cfg.aNom - cfg.aM) * x(k) + cfg.bM * r;
    uAd = -thetaHat(k) * x(k);
    uRob = -cfg.kS * satValue;
    uUnsat = uNom + uAd + uRob;
    u = max(-cfg.uMax, min(cfg.uMax, uUnsat));
    actuatorSaturated = abs(u - uUnsat) > 1.0e-12;

    if abs(e) > cfg.deadZone && ~actuatorSaturated
        gradient = cfg.gamma * e * x(k);
    else
        gradient = 0.0;
    end
    rawThetaDot = gradient - cfg.sigma * thetaHat(k);

    if thetaHat(k) >= cfg.thetaMax && rawThetaDot > 0.0
        thetaDot = 0.0;
    elseif thetaHat(k) <= cfg.thetaMin && rawThetaDot < 0.0
        thetaDot = 0.0;
    else
        thetaDot = rawThetaDot;
    end

    xDot = -cfg.aNom * x(k) + u + cfg.thetaTrue * x(k) + d;
    xmDot = -cfg.aM * xm(k) + cfg.bM * r;

    referenceSignal(k) = r;
    errorSignal(k) = e;
    control(k) = u;
    disturbanceSignal(k) = d;
    adaptationEnabled(k) = double(~actuatorSaturated);

    x(k + 1) = x(k) + cfg.dt * xDot;
    xm(k + 1) = xm(k) + cfg.dt * xmDot;
    thetaHat(k + 1) = max(cfg.thetaMin, min(cfg.thetaMax, ...
        thetaHat(k) + cfg.dt * thetaDot));
end

referenceSignal(end) = 1.0 + 0.25 * sin(0.45 * time(end));
errorSignal(end) = x(end) - xm(end);
disturbanceSignal(end) = 0.30 * sin(4.0 * time(end));
control(end) = control(end - 1);

results = table(time, referenceSignal, x, xm, errorSignal, thetaHat, ...
    control, disturbanceSignal, adaptationEnabled);
writetable(results, 'Chapter22_Lesson5_results.csv');

rmsError = sqrt(mean(errorSignal .^ 2));
fprintf('RMS tracking error: %.6f\n', rmsError);
fprintf('Maximum absolute tracking error: %.6f\n', max(abs(errorSignal)));
fprintf('Maximum absolute control: %.6f\n', max(abs(control)));
fprintf('Final theta_hat: %.6f\n', thetaHat(end));

figure('Name', 'Chapter 22 Lesson 5');
subplot(4, 1, 1);
plot(time, x, time, xm, '--', time, referenceSignal, ':', 'LineWidth', 1.1);
grid on; ylabel('state'); legend('x', 'x_m', 'r', 'Location', 'best');

subplot(4, 1, 2);
plot(time, errorSignal, 'LineWidth', 1.1);
grid on; ylabel('e');

subplot(4, 1, 3);
plot(time, thetaHat, 'LineWidth', 1.1); hold on;
yline(cfg.thetaTrue, '--');
grid on; ylabel('parameter'); legend('theta hat', 'theta true');

subplot(4, 1, 4);
plot(time, control, time, disturbanceSignal, '--', 'LineWidth', 1.1);
grid on; xlabel('time [s]'); ylabel('input'); legend('u', 'd');

saveas(gcf, 'Chapter22_Lesson5_results.png');

The Wolfram Language source uses an explicit discrete simulation so that projection, saturation logic, and supervisory switching remain transparent.

Chapter22_Lesson5.nb

(* Chapter22_Lesson5.nb
   Wolfram Language source that can be pasted into a Mathematica notebook.
   Combined robust + adaptive control for a scalar uncertain plant. *)

ClearAll["Global`*"];

cfg = <|
   "dt" -> 0.001, "finalTime" -> 30.0,
   "aNom" -> 1.0, "aM" -> 1.5, "bM" -> 1.5,
   "thetaTrue" -> 0.75, "gamma" -> 8.0, "sigma" -> 0.12,
   "kS" -> 0.90, "epsilon" -> 0.08, "deadZone" -> 0.01,
   "thetaMin" -> -2.0, "thetaMax" -> 2.0, "uMax" -> 6.0
   |>;

clip[value_, lower_, upper_] := Min[upper, Max[lower, value]];
sat[value_] := clip[value, -1.0, 1.0];
reference[t_] := If[t < 2.0, 0.0, 1.0 + 0.25 Sin[0.45 t]];
disturbance[t_] := 0.30 Sin[4.0 t] + If[12.0 <= t <= 14.0, 0.45, 0.0];

projectedDerivative[thetaHat_, rawDot_] := Which[
   thetaHat >= cfg["thetaMax"] && rawDot > 0.0, 0.0,
   thetaHat <= cfg["thetaMin"] && rawDot < 0.0, 0.0,
   True, rawDot
   ];

steps = Round[cfg["finalTime"]/cfg["dt"]];
x = 0.0; xm = 0.0; thetaHat = 0.0;
records = Reap[
     Do[
      t = step cfg["dt"];
      r = reference[t]; d = disturbance[t]; e = x - xm;
      uNom = (cfg["aNom"] - cfg["aM"]) x + cfg["bM"] r;
      uAd = -thetaHat x;
      uRob = -cfg["kS"] sat[e/cfg["epsilon"]];
      uUnsat = uNom + uAd + uRob;
      u = clip[uUnsat, -cfg["uMax"], cfg["uMax"]];
      actuatorSaturated = Abs[u - uUnsat] > 10^-12;
      gradient = If[Abs[e] > cfg["deadZone"] && ! actuatorSaturated,
        cfg["gamma"] e x, 0.0];
      rawThetaDot = gradient - cfg["sigma"] thetaHat;
      thetaDot = projectedDerivative[thetaHat, rawThetaDot];
      xDot = -cfg["aNom"] x + u + cfg["thetaTrue"] x + d;
      xmDot = -cfg["aM"] xm + cfg["bM"] r;
      Sow[{t, r, x, xm, e, thetaHat, u, d,
        If[actuatorSaturated, 0, 1]}];
      x = x + cfg["dt"] xDot;
      xm = xm + cfg["dt"] xmDot;
      thetaHat = clip[thetaHat + cfg["dt"] thetaDot,
        cfg["thetaMin"], cfg["thetaMax"]],
      {step, 0, steps}
      ]
     ][[2, 1]];

headers = {"time", "reference", "x", "xm", "error", "theta_hat",
   "control", "disturbance", "adaptation_enabled"};
Export["Chapter22_Lesson5_results.csv", Prepend[records, headers]];

rmsError = Sqrt[Mean[records[[All, 5]]^2]];
maxError = Max[Abs[records[[All, 5]]]];
maxControl = Max[Abs[records[[All, 7]]]];
Print["RMS tracking error: ", NumberForm[rmsError, {8, 6}]];
Print["Maximum absolute tracking error: ", NumberForm[maxError, {8, 6}]];
Print["Maximum absolute control: ", NumberForm[maxControl, {8, 6}]];
Print["Final theta_hat: ", NumberForm[Last[records][[6]], {8, 6}]];

statePlot = ListLinePlot[
   {records[[All, {1, 3}]], records[[All, {1, 4}]],
    records[[All, {1, 2}]]},
   PlotLegends -> {"x", "xm", "r"}, Frame -> True,
   FrameLabel -> {"time [s]", "state"}, ImageSize -> Large];
errorPlot = ListLinePlot[records[[All, {1, 5}]], Frame -> True,
   FrameLabel -> {"time [s]", "tracking error"}, ImageSize -> Large];
parameterPlot = Show[
   ListLinePlot[records[[All, {1, 6}]], Frame -> True,
    FrameLabel -> {"time [s]", "theta hat"}, ImageSize -> Large],
   Plot[cfg["thetaTrue"], {tt, 0, cfg["finalTime"]}, PlotStyle -> Dashed]
   ];
controlPlot = ListLinePlot[
   {records[[All, {1, 7}]], records[[All, {1, 8}]]},
   PlotLegends -> {"u", "d"}, Frame -> True,
   FrameLabel -> {"time [s]", "input"}, ImageSize -> Large];

GraphicsGrid[{ {statePlot}, {errorPlot}, {parameterPlot}, {controlPlot} }]

Problem 1 (Error Dynamics): For the plant and controller in Sections 2 and 3, derive the tracking-error dynamics and identify the term canceled by the adaptive law in the Lyapunov derivative.

Solution: Subtract the reference model from the plant:

\[ \begin{aligned} \dot{\mathbf{e} } &=\mathbf{A}\mathbf{x}+\mathbf{B}u +\mathbf{B}{\boldsymbol{\theta}^{*} }^{T}\boldsymbol{\phi} +\mathbf{B}d_m+\mathbf{B}_{\perp}d_u -\mathbf{A}_m\mathbf{x}_m-\mathbf{B}_m r. \end{aligned} \]

Substitute the matching equations and the three-part input. The nominal terms reduce to \(\mathbf{A}_m\mathbf{e}\), while the true and estimated parametric terms combine into \(\mathbf{B}\tilde{\boldsymbol{\theta} }^{T}\boldsymbol{\phi}\). Therefore,

\[ \dot{\mathbf{e} } =\mathbf{A}_m\mathbf{e} +\mathbf{B}\left[\tilde{\boldsymbol{\theta} }^{T}\boldsymbol{\phi} +d_m-k_s\operatorname{sat}(s/\varepsilon)\right] +\mathbf{B}_{\perp}d_u. \]

In \(\dot V\), the state-error part produces \(+2s\tilde{\boldsymbol{\theta} }^{T}\boldsymbol{\phi}\). Because \(\dot{\tilde{\boldsymbol{\theta} } }=-\dot{\hat{\boldsymbol{\theta} } }\), the gradient adaptation term produces its negative, so the cross term cancels.

Problem 2 (Ultimate Bound): Suppose \(\dot V\leq-\alpha V+C\), where \(\alpha>0\) and \(C\geq0\). Derive a bound on \(V(t)\) and state the limiting value.

Solution: Multiply by the integrating factor \(e^{\alpha t}\):

\[ \frac{d}{dt}\left(e^{\alpha t}V(t)\right) \leq Ce^{\alpha t}. \]

Integrating from zero to \(t\) gives

\[ V(t)\leq e^{-\alpha t}V(0) +\frac{C}{\alpha}\left(1-e^{-\alpha t}\right). \]

Hence \(\limsup_{t\to\infty}V(t)\leq C/\alpha\).

Problem 3 (Leakage Bias): Assume perfect tracking \(s=0\), no projection activity, and a constant true parameter. What does the leakage-modified estimate do?

Solution: The adaptive law reduces to

\[ \dot{\hat{\boldsymbol{\theta} } }=-\sigma\hat{\boldsymbol{\theta} }, \qquad \hat{\boldsymbol{\theta} }(t)=e^{-\sigma t} \hat{\boldsymbol{\theta} }(0). \]

Thus the estimate decays toward zero rather than remaining at the true parameter. Leakage is introduced for robustness and boundedness, not unbiased identification. Tracking may remain satisfactory because the feedback and robust layers compensate the residual mismatch.

Problem 4 (Projection for a Scalar Parameter): Let \(\hat\theta\in[\theta_{min},\theta_{max}]\) and let \(y\) be the unconstrained update. Define a scalar projection law.

Solution: A tangent-cone projection is

\[ \operatorname{Proj}(\hat\theta,y)= \begin{cases} 0, & \hat\theta=\theta_{max}\text{ and }y>0,\\ 0, & \hat\theta=\theta_{min}\text{ and }y<0,\\ y, & \text{otherwise}. \end{cases} \]

The law removes only outward-pointing velocity at the interval boundary; inward motion remains permitted.

Problem 5 (Matched versus Unmatched Disturbance): Why can a robust term applied through \(\mathbf{B}\) directly dominate \(d_m\) but not generally eliminate \(\mathbf{B}_{\perp}d_u\)?

Solution: The matched disturbance and control share the same input direction:

\[ \mathbf{B}(u+d_m). \]

Therefore, a control component can oppose the scalar disturbance in that channel. An unmatched disturbance acts in a state-space direction outside the range of \(\mathbf{B}\). The available input cannot instantaneously create an equal and opposite vector in that direction. Its effect must be attenuated by the stable closed-loop dynamics, and it appears as a nonzero term in the ultimate bound.

Problem 6 (Actuator-Saturation Supervisor): Explain why the parameter-gradient term should often be frozen while the actuator is saturated.

Solution: During saturation, the plant receives

\[ u_{act}=\operatorname{sat}_{u_{max} }(u_{cmd})\neq u_{cmd}. \]

The error generated by this input deficiency is not evidence of an incorrect plant parameter. If adaptation continues to interpret it as parametric mismatch, the estimate may drift toward a projection boundary and command still larger inputs. Freezing the gradient term breaks this positive feedback mechanism. Leakage may remain active so estimates relax rather than stay at an extreme value.

Combined robust-adaptive control is best understood as a layered architecture. A validated nominal controller supplies baseline stability and desired dynamics. An adaptive term compensates structured matched uncertainty. Robust injection, leakage, projection, and dead zones limit the effects of disturbances, noise, unmodeled dynamics, and parameter drift. A supervisor enforces the validity of the adaptive mechanism under actuator limits and abnormal conditions. For the representative model-reference system, the Lyapunov analysis established uniform ultimate boundedness and made the residual-set dependencies explicit. The main practical conclusion is that adaptation gain, robust gain, filtering, and actuator protection must be designed together rather than tuned independently.