After completing this lesson, students should be able to:

• derive the common error model used to compare plain and robust MRAC laws;
• explain why the nominal Lyapunov proof does not establish robustness to persistent disturbances;
• derive an ultimate-bound estimate for σ-modified MRAC;
• distinguish the mechanisms of σ-modification, e-modification, and dead-zone adaptation;
• interpret simulation metrics without treating one modification as universally superior;
• select a modification according to noise level, unmodeled dynamics, required adaptation speed, and acceptable steady-state bias.

Four questions organize the case studies. Does the tracking error remain small? Do the adaptive parameters remain bounded? Does the control input excite neglected high-frequency dynamics? How much nominal tracking performance is sacrificed to obtain robustness?

Consider the reduced-order plant model \( \dot{x}=-a x+b u+d(t) \), where \(a>0\), \(b>0\), and \(d(t)\) represents a bounded disturbance or the effect of neglected dynamics. The desired response is generated by

\[ \dot{x}_m=-a_m x_m+b_m r, \qquad a_m>0. \]

The direct MRAC control law is linearly parameterized:

\[ u=\theta_x x+\theta_r r=\boldsymbol{\theta}^{T}\boldsymbol{\phi}, \qquad \boldsymbol{\phi}=\begin{bmatrix}x & r\end{bmatrix}^{T}. \]

Under exact matching, the ideal parameters satisfy

\[ \theta_x^{*}=\frac{a-a_m}{b}, \qquad \theta_r^{*}=\frac{b_m}{b}. \]

Define the tracking error \(e=x-x_m\) and parameter error \(\widetilde{\boldsymbol{\theta} }= \boldsymbol{\theta}-\boldsymbol{\theta}^{*}\). Then

\[ \dot e=-a_m e+b\widetilde{\boldsymbol{\theta} }^{T} \boldsymbol{\phi}+w(t), \]

where \(w(t)\) collects disturbances, measurement-error effects, and reduced-order modeling error. The benchmark assumes \(|w(t)|\leq \bar w\), but does not assume that \(w(t)\) tends to zero.

The plain gradient law is

\[ \dot{\boldsymbol{\theta} }=-\mathbf{\Gamma} \boldsymbol{\phi}e, \qquad \mathbf{\Gamma}=\mathbf{\Gamma}^{T}>0. \]

Use the Lyapunov candidate

\[ V=\frac{1}{2}e^2+\frac{b}{2} \widetilde{\boldsymbol{\theta} }^{T}\mathbf{\Gamma}^{-1} \widetilde{\boldsymbol{\theta} }. \]

For the ideal model, \(w(t)=0\), differentiation gives

\[ \begin{aligned} \dot V &=e\left(-a_m e+b\widetilde{\boldsymbol{\theta} }^{T} \boldsymbol{\phi}\right) +b\widetilde{\boldsymbol{\theta} }^{T}\mathbf{\Gamma}^{-1} \dot{\boldsymbol{\theta} } \\ &=-a_m e^2. \end{aligned} \]

This proves boundedness and, with the regularity assumptions developed in earlier lessons, convergence of the tracking error. In the perturbed case,

\[ \dot V=-a_m e^2+e w(t). \]

Young's inequality yields

\[ e w\leq \frac{a_m}{2}e^2+\frac{\bar w^2}{2a_m}, \qquad \dot V\leq -\frac{a_m}{2}e^2+\frac{\bar w^2}{2a_m}. \]

The inequality supplies no negative term in \(\widetilde{\boldsymbol{\theta} }\). It therefore does not establish a compact bound on the parameter error. A persistent small error can be integrated by the update law, producing parameter drift even when the measured tracking error appears acceptable. This is the robustness gap addressed by leakage and gating modifications.

4.1 Plain MRAC

\[ \dot{\boldsymbol{\theta} }_{\text{plain} } =-\mathbf{\Gamma}\boldsymbol{\phi}e. \]

4.2 σ-Modification

\[ \dot{\boldsymbol{\theta} }_{\sigma} =-\mathbf{\Gamma}\left(\boldsymbol{\phi}e+ \sigma\boldsymbol{\theta}\right), \qquad \sigma>0. \]

The constant leakage continuously pulls the estimate toward the origin. It supplies direct damping in parameter space, but shifts the equilibrium away from the ideal parameter whenever \(\boldsymbol{\theta}^{*}\neq \mathbf{0}\).

4.3 e-Modification

\[ \dot{\boldsymbol{\theta} }_{e} =-\mathbf{\Gamma}\left(\boldsymbol{\phi}e+ \sigma_e |e|\boldsymbol{\theta}\right), \qquad \sigma_e>0. \]

Leakage becomes strong when the tracking error is large and weak when tracking is good. This generally reduces quiet-operation bias compared with constant leakage, but it also provides less parameter damping near \(e=0\).

4.4 Dead-Zone Adaptation

\[ e_d=\begin{cases} 0, & |e|\leq e_0,\\ e-e_0\operatorname{sgn}(e), & |e|>e_0, \end{cases} \qquad \dot{\boldsymbol{\theta} }_{\text{dz} } =-\mathbf{\Gamma}\boldsymbol{\phi}e_d. \]

The continuous dead-zone used here avoids a discontinuous jump in the update at \(|e|=e_0\). It prevents adaptation to errors that are plausibly explained by noise, but it does not directly damp a large parameter estimate once the error leaves the dead zone.

flowchart LR
  A["Measured tracking error \nand regressor"] --> B["Plain: integrate gradient \ncontinuously"]
  A --> C["Sigma modification: \ngradient plus constant leakage"]
  A --> D["Error modification: \ngradient plus \nerror-scaled leakage"]
  A --> E["Dead zone: stop or reduce \nadaptation near zero error"]
  B --> F["Best nominal accuracy; \nweakest drift protection"]
  C --> G["Strong parameter bound; \nconstant bias"]
  D --> H["Conditional damping; \nsmaller quiet-operation bias"]
  E --> I["Noise rejection; residual \nerror and no direct \nparameter damping"]
        

Substitute the σ-modified law into the derivative of the same Lyapunov function. Since \(\boldsymbol{\theta}= \widetilde{\boldsymbol{\theta} }+\boldsymbol{\theta}^{*}\),

\[ \dot V=-a_m e^2+e w -b\sigma\widetilde{\boldsymbol{\theta} }^{T} \boldsymbol{\theta}. \]

\[ -b\sigma\widetilde{\boldsymbol{\theta} }^{T} \boldsymbol{\theta} =-b\sigma\|\widetilde{\boldsymbol{\theta} }\|^2 -b\sigma\widetilde{\boldsymbol{\theta} }^{T} \boldsymbol{\theta}^{*}. \]

Apply Young's inequality to the disturbance and cross term:

\[ e w\leq \frac{a_m}{2}e^2+\frac{\bar w^2}{2a_m}, \]

\[ -\widetilde{\boldsymbol{\theta} }^{T} \boldsymbol{\theta}^{*} \leq \frac{1}{2}\|\widetilde{\boldsymbol{\theta} }\|^2 +\frac{1}{2}\|\boldsymbol{\theta}^{*}\|^2. \]

Therefore,

\[ \dot V\leq -\frac{a_m}{2}e^2 -\frac{b\sigma}{2}\|\widetilde{\boldsymbol{\theta} }\|^2 +\frac{\bar w^2}{2a_m} +\frac{b\sigma}{2}\|\boldsymbol{\theta}^{*}\|^2. \]

Let \(q=\min\{a_m,b\sigma\}\) and \(c=\bar w^2/(2a_m)+b\sigma \|\boldsymbol{\theta}^{*}\|^2/2\). Outside a compact set in \((e,\widetilde{\boldsymbol{\theta} })\), the negative quadratic terms dominate \(c\). Thus the tracking and parameter errors are uniformly ultimately bounded. Increasing \(\sigma\) strengthens parameter damping but also enlarges the bias term containing \(\boldsymbol{\theta}^{*}\).

6.1 e-Modification

For e-modification, the additional Lyapunov term is

\[ -b\sigma_e |e| \widetilde{\boldsymbol{\theta} }^{T}\boldsymbol{\theta} \leq -\frac{b\sigma_e|e|}{2} \|\widetilde{\boldsymbol{\theta} }\|^2 +\frac{b\sigma_e|e|}{2} \|\boldsymbol{\theta}^{*}\|^2. \]

Parameter damping is state dependent. During a large transient it can be substantial; close to perfect tracking it vanishes. Consequently, the simplest quadratic proof is less uniform than for constant σ-modification. The practical advantage is less leakage-induced bias after the transient.

6.2 Dead Zone

If the known error contribution of noise and disturbances is bounded by \(e_0\), setting the update to zero inside that band prevents the parameter integrator from reacting to indistinguishable small errors. The price is that adaptation cannot force the measured error below the selected threshold. Moreover, dead-zone logic is not a substitute for leakage when large errors are caused by neglected actuator dynamics: once \(|e|>e_0\), the gradient law may still drive the parameters to large values.

The simulations use a physical plant whose input is filtered by a first-order actuator:

\[ \dot x=-a x+b v+d(t), \qquad \tau\dot v=-v+u. \]

The adaptive controller is designed as though \(v=u\). Therefore \(\tau\) is an unmodeled high-frequency time constant. The reference model is

\[ \dot x_m=-2x_m+2r. \]

The deterministic reference and perturbations are

\[ r(t)=r_{\mathrm{sq} }(t)+0.2\sin(0.6t), \]

\[ d(t)=d_0\left(0.6\sin(2t)+0.4\sin(7t)\right), \]

\[ n(t)=n_0\left(0.7\sin(70t)+0.3\sin(113t)\right), \qquad y=x+n(t). \]

The controller uses \(y\) rather than the uncorrupted state. All implementations use forward Euler integration so the numerical algorithm is transparent and identical across languages. Production work should use a validated fixed-step or variable-step integrator and explicit sample-time analysis.

The nominal case uses a small actuator time constant \(\tau=0.005\), no disturbance, no measurement noise, adaptation gain \(\gamma=15\), and a control limit of \(20\). Representative deterministic results are:

Plain MRAC gives the smallest nominal RMS error because it contains no deliberate leakage or inactive band. Constant leakage has the largest nominal bias. The e-modification approaches plain MRAC after the transient because its leakage decreases with \(|e|\). The dead zone preserves a nonzero error band and therefore cannot be expected to minimize asymptotic tracking error.

With measurement noise and a bounded disturbance, the plain update law behaves like an integrator driven by the product \(\boldsymbol{\phi}e\). Even when the mean error is small, correlation between the regressor and disturbance can create a nonzero mean update. Over long intervals this causes drift.

The σ-modification introduces a restoring term for every parameter and is therefore the most direct defense against drift. The e-modification applies restoring action primarily when the disturbance produces a visible tracking error. The dead zone is effective when its width exceeds the noise-induced error, but an excessively large threshold suppresses useful adaptation to genuine parameter changes.

\[ e_0 \gtrsim e_{\mathrm{noise,max} }+e_{\mathrm{disturbance,max} } \]

This inequality is a design heuristic, not an exact theorem. The relevant quantities are closed-loop error contributions, not raw sensor-noise and disturbance amplitudes.

The stress case uses \(\tau=0.10\), \(\gamma=80\), disturbance amplitude \(d_0=0.03\), noise amplitude \(n_0=0.01\), and input saturation at \(|u|\leq20\). The controller still assumes \(v=u\). Representative results are:

The principal robustness metric is not only RMS error. Plain MRAC and the dead-zone law drive the parameter norm above \(74\) and hit actuator saturation. The dead zone does not help much because the parasitic dynamics generate errors well outside the inactive band. Constant and error-scaled leakage keep the parameter norm near single digits. In this tuning, σ-modification also gives the best stress-case RMS error, while e-modification is more conservative during large-error intervals.

These numbers are properties of this benchmark and tuning, not universal rankings. A different leakage coefficient, dead-zone width, adaptation gain, or actuator bandwidth can change the order of the RMS errors. The qualitative distinction—direct parameter damping versus no direct damping— is more transferable than the numerical ranking.

Suppose the true plant coefficient changes from \(a_1\) to \(a_2\). The ideal feedback parameter changes by

\[ \Delta\theta_x^{*}=\frac{a_2-a_1}{b}. \]

Plain MRAC generally responds fastest because the full gradient is always active. Constant leakage may slow convergence toward a nonzero ideal value and creates steady-state bias. e-modification applies stronger damping during the large transient, which may reduce overshoot but can slow the parameter movement. Dead-zone adaptation begins only after the tracking error exceeds the threshold; a small physical change may therefore remain undetected.

Robust modification should not be confused with fault detection. If abrupt plant changes are safety critical, a supervisory monitor should compare tracking error, parameter rate, control saturation, and model residuals, and should be able to freeze or reset adaptation.

flowchart TD
  A["Start from a verified nominal MRAC design"] --> B["Estimate sensor noise and disturbance-induced error"]
  B --> C["Check neglected actuator and plant bandwidth"]
  C --> D["Does parameter drift occur in long tests?"]
  D -->|"yes"| E["Need direct damping?"]
  E -->|"yes"| F["Use sigma modification \nas baseline"]
  E -->|"only during large errors"| G["Use error modification"]
  D -->|"no; small error is mostly noise"| H["Use a justified dead zone"]
  F --> I["Sweep leakage and adaptation gain"]
  G --> I
  H --> I
  I --> J["Test saturation, delays, noise, parameter jumps, and long duration"]
  J --> K["Accept only if tracking, parameter, and control bounds are satisfied"]
        

A practical tuning sequence is to reduce the adaptation gain until the neglected dynamics are not strongly excited, introduce the smallest robust modification that prevents drift, and then retune the baseline feedback or reference model if the resulting tracking error is unacceptable. Robust modification cannot repair a fundamentally invalid matching model or an actuator whose bandwidth is too low for the commanded reference dynamics.

The implementation uses NumPy for array operations and Matplotlib for plots. The python-control package can be added for transfer- function and frequency-response analysis, while SciPy provides higher-order integrators. The code deliberately uses explicit Euler integration to keep the cross-language comparison identical.

Chapter8_Lesson5.py

"""Chapter 8, Lesson 5: Plain MRAC versus robust modifications.

Dependencies:
    pip install numpy matplotlib

The script runs two deterministic benchmarks:
1. nominal reduced-order plant;
2. stress case with actuator lag, bounded disturbance, measurement noise,
   fast adaptation, and input saturation.
"""

from __future__ import annotations

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

import matplotlib.pyplot as plt
import numpy as np


class Method(str, Enum):
    PLAIN = "plain"
    SIGMA = "sigma"
    EMOD = "e-mod"
    DEADZONE = "dead-zone"


@dataclass(frozen=True)
class Scenario:
    name: str
    duration: float
    dt: float
    gamma: float
    actuator_tau: float
    disturbance_amplitude: float
    noise_amplitude: float
    input_limit: float


@dataclass
class Result:
    t: np.ndarray
    r: np.ndarray
    x: np.ndarray
    xm: np.ndarray
    e: np.ndarray
    u: np.ndarray
    theta_x: np.ndarray
    theta_r: np.ndarray

    @property
    def theta_norm(self) -> np.ndarray:
        return np.sqrt(self.theta_x**2 + self.theta_r**2)

    def metrics(self) -> dict[str, float]:
        return {
            "rms_error": float(np.sqrt(np.mean(self.e**2))),
            "peak_error": float(np.max(np.abs(self.e))),
            "peak_theta_norm": float(np.max(self.theta_norm)),
            "peak_control": float(np.max(np.abs(self.u))),
        }


def reference(t: float) -> float:
    square = 1.0 if (t % 6.0) < 3.0 else -1.0
    return square + 0.2 * math.sin(0.6 * t)


def simulate(method: Method, scenario: Scenario) -> Result:
    # Reduced-order plant: x_dot = -a*x + b*v + d.
    # Unmodeled actuator: tau*v_dot = -v + u.
    # Reference model: xm_dot = -am*xm + bm*r.
    a, b = 1.0, 1.0
    am, bm = 2.0, 2.0

    # Robust-modification coefficients.
    sigma = 0.08
    sigma_e = 0.20
    dead_zone = 0.04

    n = int(round(scenario.duration / scenario.dt)) + 1
    t_vec = np.linspace(0.0, scenario.duration, n)
    r_vec = np.zeros(n)
    x_vec = np.zeros(n)
    xm_vec = np.zeros(n)
    e_vec = np.zeros(n)
    u_vec = np.zeros(n)
    tx_vec = np.zeros(n)
    tr_vec = np.zeros(n)

    x = 0.0
    v = 0.0
    xm = 0.0
    theta = np.zeros(2)
    gamma = np.array([scenario.gamma, scenario.gamma])

    for k, t in enumerate(t_vec):
        r = reference(t)
        disturbance = scenario.disturbance_amplitude * (
            0.6 * math.sin(2.0 * t) + 0.4 * math.sin(7.0 * t)
        )
        noise = scenario.noise_amplitude * (
            0.7 * math.sin(70.0 * t) + 0.3 * math.sin(113.0 * t)
        )

        y = x + noise
        e = y - xm
        phi = np.array([y, r])
        u = float(theta @ phi)
        u = float(np.clip(u, -scenario.input_limit, scenario.input_limit))

        if method is Method.PLAIN:
            theta_dot = -gamma * phi * e
        elif method is Method.SIGMA:
            theta_dot = -gamma * (phi * e + sigma * theta)
        elif method is Method.EMOD:
            theta_dot = -gamma * (phi * e + sigma_e * abs(e) * theta)
        else:
            if abs(e) <= dead_zone:
                e_dead = 0.0
            else:
                e_dead = e - math.copysign(dead_zone, e)
            theta_dot = -gamma * phi * e_dead

        x_dot = -a * x + b * v + disturbance
        v_dot = (-v + u) / scenario.actuator_tau
        xm_dot = -am * xm + bm * r

        r_vec[k] = r
        x_vec[k] = x
        xm_vec[k] = xm
        e_vec[k] = e
        u_vec[k] = u
        tx_vec[k] = theta[0]
        tr_vec[k] = theta[1]

        # Forward Euler is intentionally used so every language implementation
        # follows exactly the same transparent numerical algorithm.
        x += scenario.dt * x_dot
        v += scenario.dt * v_dot
        xm += scenario.dt * xm_dot
        theta += scenario.dt * theta_dot

    return Result(t_vec, r_vec, x_vec, xm_vec, e_vec, u_vec, tx_vec, tr_vec)


def write_csv(path: Path, result: Result) -> None:
    with path.open("w", newline="", encoding="utf-8") as handle:
        writer = csv.writer(handle)
        writer.writerow(
            ["t", "r", "x", "xm", "e", "u", "theta_x", "theta_r", "theta_norm"]
        )
        for row in zip(
            result.t,
            result.r,
            result.x,
            result.xm,
            result.e,
            result.u,
            result.theta_x,
            result.theta_r,
            result.theta_norm,
        ):
            writer.writerow(row)


def plot_results(scenario: Scenario, results: dict[Method, Result], output: Path) -> None:
    fig, axes = plt.subplots(3, 1, figsize=(10, 10), sharex=True)
    for method, result in results.items():
        axes[0].plot(result.t, result.e, label=method.value)
        axes[1].plot(result.t, result.theta_norm, label=method.value)
        axes[2].plot(result.t, result.u, label=method.value)

    axes[0].set_ylabel("tracking error")
    axes[1].set_ylabel("parameter norm")
    axes[2].set_ylabel("control input")
    axes[2].set_xlabel("time [s]")
    axes[0].set_title(f"MRAC comparison: {scenario.name}")
    for axis in axes:
        axis.grid(True)
        axis.legend()
    fig.tight_layout()
    fig.savefig(output, dpi=180)
    plt.close(fig)


def main() -> None:
    scenarios = [
        Scenario(
            name="nominal",
            duration=30.0,
            dt=0.001,
            gamma=15.0,
            actuator_tau=0.005,
            disturbance_amplitude=0.0,
            noise_amplitude=0.0,
            input_limit=20.0,
        ),
        Scenario(
            name="stress",
            duration=30.0,
            dt=0.001,
            gamma=80.0,
            actuator_tau=0.10,
            disturbance_amplitude=0.03,
            noise_amplitude=0.01,
            input_limit=20.0,
        ),
    ]

    output_dir = Path("Chapter8_Lesson5_output")
    output_dir.mkdir(exist_ok=True)

    for scenario in scenarios:
        results: dict[Method, Result] = {}
        print(f"\nScenario: {scenario.name}")
        print("method       RMS(e)    peak|e|   peak||theta||   peak|u|")
        print("--------------------------------------------------------")
        for method in Method:
            result = simulate(method, scenario)
            results[method] = result
            metrics = result.metrics()
            print(
                f"{method.value:10s} "
                f"{metrics['rms_error']:9.4f} "
                f"{metrics['peak_error']:10.4f} "
                f"{metrics['peak_theta_norm']:15.4f} "
                f"{metrics['peak_control']:10.4f}"
            )
            write_csv(output_dir / f"{scenario.name}_{method.value}.csv", result)

        plot_results(
            scenario,
            results,
            output_dir / f"Chapter8_Lesson5_{scenario.name}.png",
        )


if __name__ == "__main__":
    main()

This version uses only the C++17 standard library. Eigen is useful when the scalar benchmark is generalized to matrix-valued MIMO MRAC, and libraries such as Boost.Odeint can replace the transparent Euler step with validated ODE solvers.

Chapter8_Lesson5.cpp

// Chapter 8, Lesson 5: Plain MRAC versus robust modifications.
// Build: g++ -std=c++17 -O2 Chapter8_Lesson5.cpp -o Chapter8_Lesson5

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

constexpr double PI = 3.14159265358979323846;

enum class Method { Plain, Sigma, EMod, DeadZone };

struct Scenario {
    std::string name;
    double duration;
    double dt;
    double gamma;
    double actuator_tau;
    double disturbance_amplitude;
    double noise_amplitude;
    double input_limit;
};

struct Sample {
    double t{};
    double r{};
    double x{};
    double xm{};
    double e{};
    double u{};
    double theta_x{};
    double theta_r{};
};

struct Metrics {
    double rms_error{};
    double peak_error{};
    double peak_theta_norm{};
    double peak_control{};
};

std::string method_name(Method method) {
    switch (method) {
        case Method::Plain: return "plain";
        case Method::Sigma: return "sigma";
        case Method::EMod: return "e-mod";
        case Method::DeadZone: return "dead-zone";
    }
    throw std::logic_error("Unknown method");
}

double reference(double t) {
    const double square = std::fmod(t, 6.0) < 3.0 ? 1.0 : -1.0;
    return square + 0.2 * std::sin(0.6 * t);
}

std::vector<Sample> simulate(Method method, const Scenario& scenario) {
    const double a = 1.0;
    const double b = 1.0;
    const double am = 2.0;
    const double bm = 2.0;
    const double sigma = 0.08;
    const double sigma_e = 0.20;
    const double dead_zone = 0.04;

    const std::size_t n = static_cast<std::size_t>(
        std::llround(scenario.duration / scenario.dt)) + 1U;
    std::vector<Sample> history;
    history.reserve(n);

    double x = 0.0;
    double v = 0.0;
    double xm = 0.0;
    std::array<double, 2> theta{0.0, 0.0};

    for (std::size_t k = 0; k < n; ++k) {
        const double t = static_cast<double>(k) * scenario.dt;
        const double r = reference(t);
        const double disturbance = scenario.disturbance_amplitude *
            (0.6 * std::sin(2.0 * t) + 0.4 * std::sin(7.0 * t));
        const double noise = scenario.noise_amplitude *
            (0.7 * std::sin(70.0 * t) + 0.3 * std::sin(113.0 * t));

        const double y = x + noise;
        const double e = y - xm;
        const std::array<double, 2> phi{y, r};
        double u = theta[0] * phi[0] + theta[1] * phi[1];
        u = std::clamp(u, -scenario.input_limit, scenario.input_limit);

        std::array<double, 2> theta_dot{};
        if (method == Method::Plain) {
            theta_dot[0] = -scenario.gamma * phi[0] * e;
            theta_dot[1] = -scenario.gamma * phi[1] * e;
        } else if (method == Method::Sigma) {
            theta_dot[0] = -scenario.gamma * (phi[0] * e + sigma * theta[0]);
            theta_dot[1] = -scenario.gamma * (phi[1] * e + sigma * theta[1]);
        } else if (method == Method::EMod) {
            theta_dot[0] = -scenario.gamma *
                (phi[0] * e + sigma_e * std::abs(e) * theta[0]);
            theta_dot[1] = -scenario.gamma *
                (phi[1] * e + sigma_e * std::abs(e) * theta[1]);
        } else {
            double e_dead = 0.0;
            if (std::abs(e) > dead_zone) {
                e_dead = e - std::copysign(dead_zone, e);
            }
            theta_dot[0] = -scenario.gamma * phi[0] * e_dead;
            theta_dot[1] = -scenario.gamma * phi[1] * e_dead;
        }

        const double x_dot = -a * x + b * v + disturbance;
        const double v_dot = (-v + u) / scenario.actuator_tau;
        const double xm_dot = -am * xm + bm * r;

        history.push_back({t, r, x, xm, e, u, theta[0], theta[1]});

        x += scenario.dt * x_dot;
        v += scenario.dt * v_dot;
        xm += scenario.dt * xm_dot;
        theta[0] += scenario.dt * theta_dot[0];
        theta[1] += scenario.dt * theta_dot[1];
    }

    return history;
}

Metrics calculate_metrics(const std::vector<Sample>& history) {
    double squared_error_sum = 0.0;
    double peak_error = 0.0;
    double peak_theta_norm = 0.0;
    double peak_control = 0.0;

    for (const Sample& sample : history) {
        squared_error_sum += sample.e * sample.e;
        peak_error = std::max(peak_error, std::abs(sample.e));
        peak_theta_norm = std::max(
            peak_theta_norm,
            std::hypot(sample.theta_x, sample.theta_r));
        peak_control = std::max(peak_control, std::abs(sample.u));
    }

    return {
        std::sqrt(squared_error_sum / static_cast<double>(history.size())),
        peak_error,
        peak_theta_norm,
        peak_control
    };
}

void write_csv(
    const std::string& filename,
    const std::vector<Sample>& history) {
    std::ofstream file(filename);
    if (!file) {
        throw std::runtime_error("Cannot open output file: " + filename);
    }
    file << "t,r,x,xm,e,u,theta_x,theta_r,theta_norm\n";
    file << std::setprecision(12);
    for (const Sample& sample : history) {
        file << sample.t << ',' << sample.r << ',' << sample.x << ','
             << sample.xm << ',' << sample.e << ',' << sample.u << ','
             << sample.theta_x << ',' << sample.theta_r << ','
             << std::hypot(sample.theta_x, sample.theta_r) << '\n';
    }
}

int main() {
    const std::array<Scenario, 2> scenarios{ {
        {"nominal", 30.0, 0.001, 15.0, 0.005, 0.0, 0.0, 20.0},
        {"stress", 30.0, 0.001, 80.0, 0.10, 0.03, 0.01, 20.0}
    } };
    const std::array<Method, 4> methods{ {
        Method::Plain, Method::Sigma, Method::EMod, Method::DeadZone
    } };

    for (const Scenario& scenario : scenarios) {
        std::cout << "\nScenario: " << scenario.name << '\n';
        std::cout << "method       RMS(e)    peak|e|   peak||theta||   peak|u|\n";
        std::cout << "--------------------------------------------------------\n";

        for (Method method : methods) {
            const auto history = simulate(method, scenario);
            const Metrics metrics = calculate_metrics(history);
            const std::string name = method_name(method);

            std::cout << std::left << std::setw(11) << name
                      << std::right << std::fixed << std::setprecision(4)
                      << std::setw(9) << metrics.rms_error
                      << std::setw(11) << metrics.peak_error
                      << std::setw(16) << metrics.peak_theta_norm
                      << std::setw(11) << metrics.peak_control << '\n';

            write_csv("Chapter8_Lesson5_" + scenario.name + "_" + name + ".csv", history);
        }
    }

    return 0;
}

The Java implementation uses records, enums, and standard file I/O. EJML or Apache Commons Math can support vector and matrix operations for higher- order adaptive systems. JavaFX can be used to display the generated time histories in a desktop laboratory application.

Chapter8_Lesson5.java

// Chapter 8, Lesson 5: Plain MRAC versus robust modifications.
// Build and run:
//   javac Chapter8_Lesson5.java
//   java Chapter8_Lesson5

import java.io.BufferedWriter;
import java.io.IOException;
import java.nio.file.Files;
import java.nio.file.Path;
import java.util.ArrayList;
import java.util.List;
import java.util.Locale;

public final class Chapter8_Lesson5 {
    private enum Method {
        PLAIN("plain"),
        SIGMA("sigma"),
        EMOD("e-mod"),
        DEADZONE("dead-zone");

        private final String label;

        Method(String label) {
            this.label = label;
        }

        public String label() {
            return label;
        }
    }

    private record Scenario(
        String name,
        double duration,
        double dt,
        double gamma,
        double actuatorTau,
        double disturbanceAmplitude,
        double noiseAmplitude,
        double inputLimit
    ) {}

    private record Sample(
        double t,
        double r,
        double x,
        double xm,
        double e,
        double u,
        double thetaX,
        double thetaR
    ) {}

    private record Metrics(
        double rmsError,
        double peakError,
        double peakThetaNorm,
        double peakControl
    ) {}

    private static double reference(double t) {
        double square = (t % 6.0) < 3.0 ? 1.0 : -1.0;
        return square + 0.2 * Math.sin(0.6 * t);
    }

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

    private static List<Sample> simulate(Method method, Scenario scenario) {
        final double a = 1.0;
        final double b = 1.0;
        final double am = 2.0;
        final double bm = 2.0;
        final double sigma = 0.08;
        final double sigmaE = 0.20;
        final double deadZone = 0.04;

        int n = (int) Math.round(scenario.duration() / scenario.dt()) + 1;
        List<Sample> history = new ArrayList<>(n);

        double x = 0.0;
        double v = 0.0;
        double xm = 0.0;
        double thetaX = 0.0;
        double thetaR = 0.0;

        for (int k = 0; k < n; k++) {
            double t = k * scenario.dt();
            double r = reference(t);
            double disturbance = scenario.disturbanceAmplitude()
                * (0.6 * Math.sin(2.0 * t) + 0.4 * Math.sin(7.0 * t));
            double noise = scenario.noiseAmplitude()
                * (0.7 * Math.sin(70.0 * t) + 0.3 * Math.sin(113.0 * t));

            double y = x + noise;
            double e = y - xm;
            double phiX = y;
            double phiR = r;
            double u = thetaX * phiX + thetaR * phiR;
            u = clamp(u, -scenario.inputLimit(), scenario.inputLimit());

            double thetaXDot;
            double thetaRDot;
            switch (method) {
                case PLAIN -> {
                    thetaXDot = -scenario.gamma() * phiX * e;
                    thetaRDot = -scenario.gamma() * phiR * e;
                }
                case SIGMA -> {
                    thetaXDot = -scenario.gamma() * (phiX * e + sigma * thetaX);
                    thetaRDot = -scenario.gamma() * (phiR * e + sigma * thetaR);
                }
                case EMOD -> {
                    thetaXDot = -scenario.gamma()
                        * (phiX * e + sigmaE * Math.abs(e) * thetaX);
                    thetaRDot = -scenario.gamma()
                        * (phiR * e + sigmaE * Math.abs(e) * thetaR);
                }
                case DEADZONE -> {
                    double eDead = 0.0;
                    if (Math.abs(e) > deadZone) {
                        eDead = e - Math.copySign(deadZone, e);
                    }
                    thetaXDot = -scenario.gamma() * phiX * eDead;
                    thetaRDot = -scenario.gamma() * phiR * eDead;
                }
                default -> throw new IllegalStateException("Unknown method");
            }

            double xDot = -a * x + b * v + disturbance;
            double vDot = (-v + u) / scenario.actuatorTau();
            double xmDot = -am * xm + bm * r;

            history.add(new Sample(t, r, x, xm, e, u, thetaX, thetaR));

            x += scenario.dt() * xDot;
            v += scenario.dt() * vDot;
            xm += scenario.dt() * xmDot;
            thetaX += scenario.dt() * thetaXDot;
            thetaR += scenario.dt() * thetaRDot;
        }

        return history;
    }

    private static Metrics calculateMetrics(List<Sample> history) {
        double squaredErrorSum = 0.0;
        double peakError = 0.0;
        double peakThetaNorm = 0.0;
        double peakControl = 0.0;

        for (Sample sample : history) {
            squaredErrorSum += sample.e() * sample.e();
            peakError = Math.max(peakError, Math.abs(sample.e()));
            peakThetaNorm = Math.max(
                peakThetaNorm,
                Math.hypot(sample.thetaX(), sample.thetaR())
            );
            peakControl = Math.max(peakControl, Math.abs(sample.u()));
        }

        return new Metrics(
            Math.sqrt(squaredErrorSum / history.size()),
            peakError,
            peakThetaNorm,
            peakControl
        );
    }

    private static void writeCsv(Path path, List<Sample> history) throws IOException {
        try (BufferedWriter writer = Files.newBufferedWriter(path)) {
            writer.write("t,r,x,xm,e,u,theta_x,theta_r,theta_norm\n");
            for (Sample sample : history) {
                writer.write(String.format(
                    Locale.US,
                    "%.12f,%.12f,%.12f,%.12f,%.12f,%.12f,%.12f,%.12f,%.12f%n",
                    sample.t(),
                    sample.r(),
                    sample.x(),
                    sample.xm(),
                    sample.e(),
                    sample.u(),
                    sample.thetaX(),
                    sample.thetaR(),
                    Math.hypot(sample.thetaX(), sample.thetaR())
                ));
            }
        }
    }

    public static void main(String[] args) throws IOException {
        Locale.setDefault(Locale.US);

        Scenario[] scenarios = {
            new Scenario("nominal", 30.0, 0.001, 15.0, 0.005, 0.0, 0.0, 20.0),
            new Scenario("stress", 30.0, 0.001, 80.0, 0.10, 0.03, 0.01, 20.0)
        };

        for (Scenario scenario : scenarios) {
            System.out.println("\nScenario: " + scenario.name());
            System.out.println("method       RMS(e)    peak|e|   peak||theta||   peak|u|");
            System.out.println("--------------------------------------------------------");

            for (Method method : Method.values()) {
                List<Sample> history = simulate(method, scenario);
                Metrics metrics = calculateMetrics(history);

                System.out.printf(
                    Locale.US,
                    "%-11s%9.4f%11.4f%16.4f%11.4f%n",
                    method.label(),
                    metrics.rmsError(),
                    metrics.peakError(),
                    metrics.peakThetaNorm(),
                    metrics.peakControl()
                );

                writeCsv(
                    Path.of("Chapter8_Lesson5_" + scenario.name() + "_" + method.label() + ".csv"),
                    history
                );
            }
        }
    }
}

The script requires only base MATLAB. For Simulink, implement the plant, actuator, and reference model with Integrator blocks; compute the regressor and tracking error with Product and Sum blocks; place the selected update law in a MATLAB Function block; and integrate each adaptive parameter. Saturation should be placed on the physical control input, while parameter monitoring should occur before and after any projection or reset logic.

Chapter8_Lesson5.m

%% Chapter 8, Lesson 5: Plain MRAC versus robust modifications
% Run this script in MATLAB. It requires only base MATLAB.
% A Simulink realization can use the same equations in MATLAB Function blocks.

clear; clc; close all;

methods = ["plain", "sigma", "e-mod", "dead-zone"];

nominal = struct( ...
    'name', "nominal", ...
    'duration', 30.0, ...
    'dt', 0.001, ...
    'gamma', 15.0, ...
    'actuatorTau', 0.005, ...
    'disturbanceAmplitude', 0.0, ...
    'noiseAmplitude', 0.0, ...
    'inputLimit', 20.0);

stress = struct( ...
    'name', "stress", ...
    'duration', 30.0, ...
    'dt', 0.001, ...
    'gamma', 80.0, ...
    'actuatorTau', 0.10, ...
    'disturbanceAmplitude', 0.03, ...
    'noiseAmplitude', 0.01, ...
    'inputLimit', 20.0);

scenarios = {nominal, stress};

for s = 1:numel(scenarios)
    scenario = scenarios{s};
    fprintf('\nScenario: %s\n', scenario.name);
    fprintf('method       RMS(e)    peak|e|   peak||theta||   peak|u|\n');
    fprintf('--------------------------------------------------------\n');

    results = cell(1, numel(methods));
    for k = 1:numel(methods)
        result = simulateMRAC(methods(k), scenario);
        results{k} = result;
        metrics = calculateMetrics(result);

        fprintf('%-11s%9.4f%11.4f%16.4f%11.4f\n', ...
            methods(k), metrics.rmsError, metrics.peakError, ...
            metrics.peakThetaNorm, metrics.peakControl);

        outputTable = table( ...
            result.t, result.r, result.x, result.xm, result.e, result.u, ...
            result.thetaX, result.thetaR, result.thetaNorm, ...
            'VariableNames', ...
            {'t','r','x','xm','e','u','theta_x','theta_r','theta_norm'});
        writetable(outputTable, ...
            "Chapter8_Lesson5_" + scenario.name + "_" + methods(k) + ".csv");
    end

    figure('Name', "MRAC comparison: " + scenario.name);
    tiledlayout(3,1);

    nexttile;
    hold on;
    for k = 1:numel(methods)
        plot(results{k}.t, results{k}.e, 'DisplayName', methods(k));
    end
    grid on; ylabel('tracking error'); legend('Location','best');
    title("MRAC comparison: " + scenario.name);

    nexttile;
    hold on;
    for k = 1:numel(methods)
        plot(results{k}.t, results{k}.thetaNorm, 'DisplayName', methods(k));
    end
    grid on; ylabel('parameter norm'); legend('Location','best');

    nexttile;
    hold on;
    for k = 1:numel(methods)
        plot(results{k}.t, results{k}.u, 'DisplayName', methods(k));
    end
    grid on; ylabel('control input'); xlabel('time [s]');
    legend('Location','best');
end

function result = simulateMRAC(method, scenario)
    a = 1.0;
    b = 1.0;
    am = 2.0;
    bm = 2.0;

    sigma = 0.08;
    sigmaE = 0.20;
    deadZone = 0.04;

    n = round(scenario.duration / scenario.dt) + 1;
    t = (0:n-1)' * scenario.dt;
    r = zeros(n,1);
    xHistory = zeros(n,1);
    xmHistory = zeros(n,1);
    eHistory = zeros(n,1);
    uHistory = zeros(n,1);
    thetaXHistory = zeros(n,1);
    thetaRHistory = zeros(n,1);

    x = 0.0;
    v = 0.0;
    xm = 0.0;
    theta = [0.0; 0.0];
    Gamma = diag([scenario.gamma, scenario.gamma]);

    for k = 1:n
        tk = t(k);
        squareReference = 1.0;
        if mod(tk, 6.0) >= 3.0
            squareReference = -1.0;
        end
        rk = squareReference + 0.2 * sin(0.6 * tk);

        disturbance = scenario.disturbanceAmplitude * ...
            (0.6 * sin(2.0 * tk) + 0.4 * sin(7.0 * tk));
        noise = scenario.noiseAmplitude * ...
            (0.7 * sin(70.0 * tk) + 0.3 * sin(113.0 * tk));

        y = x + noise;
        e = y - xm;
        phi = [y; rk];
        u = theta' * phi;
        u = min(max(u, -scenario.inputLimit), scenario.inputLimit);

        switch method
            case "plain"
                thetaDot = -Gamma * phi * e;
            case "sigma"
                thetaDot = -Gamma * (phi * e + sigma * theta);
            case "e-mod"
                thetaDot = -Gamma * (phi * e + sigmaE * abs(e) * theta);
            case "dead-zone"
                if abs(e) <= deadZone
                    eDead = 0.0;
                else
                    eDead = e - sign(e) * deadZone;
                end
                thetaDot = -Gamma * phi * eDead;
            otherwise
                error('Unknown method: %s', method);
        end

        xDot = -a * x + b * v + disturbance;
        vDot = (-v + u) / scenario.actuatorTau;
        xmDot = -am * xm + bm * rk;

        r(k) = rk;
        xHistory(k) = x;
        xmHistory(k) = xm;
        eHistory(k) = e;
        uHistory(k) = u;
        thetaXHistory(k) = theta(1);
        thetaRHistory(k) = theta(2);

        x = x + scenario.dt * xDot;
        v = v + scenario.dt * vDot;
        xm = xm + scenario.dt * xmDot;
        theta = theta + scenario.dt * thetaDot;
    end

    result = struct( ...
        't', t, ...
        'r', r, ...
        'x', xHistory, ...
        'xm', xmHistory, ...
        'e', eHistory, ...
        'u', uHistory, ...
        'thetaX', thetaXHistory, ...
        'thetaR', thetaRHistory, ...
        'thetaNorm', hypot(thetaXHistory, thetaRHistory));
end

function metrics = calculateMetrics(result)
    metrics = struct( ...
        'rmsError', sqrt(mean(result.e.^2)), ...
        'peakError', max(abs(result.e)), ...
        'peakThetaNorm', max(result.thetaNorm), ...
        'peakControl', max(abs(result.u)));
end

The notebook wrapper loads a companion Wolfram Language file containing the complete simulation. For symbolic work, the same environment can be used to verify Lyapunov derivatives and manipulate ultimate-bound inequalities. NDSolve is appropriate when higher-order numerical integration is desired.

Chapter8_Lesson5.nb

Notebook[{
  Cell["Chapter 8, Lesson 5: Plain MRAC versus Robust Modifications", "Title"],
  Cell["This notebook loads and executes the companion Wolfram Language implementation.", "Text"],
  Cell[BoxData[
    RowBox[{"Get", "[", 
      RowBox[{"FileNameJoin", "[", 
        RowBox[{"{", 
          RowBox[{"NotebookDirectory", "[", "]", ",", "\"Chapter8_Lesson5.wl\""}], 
        "}"}], "]"}], "]"}]], "Input"]
}, WindowTitle -> "Chapter8_Lesson5", Saveable -> True]

Chapter8_Lesson5.wl

(* Chapter 8, Lesson 5: Plain MRAC versus robust modifications. *)

ClearAll[reference, clip, simulateMRAC, metrics, exportResult];

reference[t_] := If[Mod[t, 6.0] < 3.0, 1.0, -1.0] + 0.2 Sin[0.6 t];
clip[value_, limit_] := Min[Max[value, -limit], limit];

simulateMRAC[method_String, scenario_Association] := Module[
  {
    a = 1.0, b = 1.0, am = 2.0, bm = 2.0,
    sigma = 0.08, sigmaE = 0.20, deadZone = 0.04,
    n, dt, x = 0.0, v = 0.0, xm = 0.0, theta = {0.0, 0.0},
    gamma, rows, t, r, disturbance, noise, y, e, phi, u,
    thetaDot, eDead, xDot, vDot, xmDot
  },

  dt = scenario["dt"];
  gamma = scenario["gamma"];
  n = Round[scenario["duration"]/dt] + 1;
  rows = ConstantArray[0.0, {n, 9}];

  Do[
    t = (k - 1) dt;
    r = reference[t];
    disturbance = scenario["disturbanceAmplitude"]
      (0.6 Sin[2.0 t] + 0.4 Sin[7.0 t]);
    noise = scenario["noiseAmplitude"]
      (0.7 Sin[70.0 t] + 0.3 Sin[113.0 t]);

    y = x + noise;
    e = y - xm;
    phi = {y, r};
    u = clip[theta . phi, scenario["inputLimit"]];

    thetaDot = Switch[method,
      "plain",
        -gamma phi e,
      "sigma",
        -gamma (phi e + sigma theta),
      "e-mod",
        -gamma (phi e + sigmaE Abs[e] theta),
      "dead-zone",
        eDead = If[Abs[e] <= deadZone, 0.0, e - Sign[e] deadZone];
        -gamma phi eDead,
      _,
        Print["Unknown method: ", method]; Abort[]
    ];

    xDot = -a x + b v + disturbance;
    vDot = (-v + u)/scenario["actuatorTau"];
    xmDot = -am xm + bm r;

    rows[[k]] = {
      t, r, x, xm, e, u, theta[[1]], theta[[2]], Norm[theta]
    };

    x = x + dt xDot;
    v = v + dt vDot;
    xm = xm + dt xmDot;
    theta = theta + dt thetaDot,
    {k, 1, n}
  ];

  Association[
    "data" -> rows,
    "columns" -> {"t", "r", "x", "xm", "e", "u", "theta_x", "theta_r", "theta_norm"}
  ]
];

metrics[result_Association] := Module[
  {data = result["data"], e, u, thetaNorm},
  e = data[[All, 5]];
  u = data[[All, 6]];
  thetaNorm = data[[All, 9]];
  Association[
    "rmsError" -> Sqrt[Mean[e^2]],
    "peakError" -> Max[Abs[e]],
    "peakThetaNorm" -> Max[thetaNorm],
    "peakControl" -> Max[Abs[u]]
  ]
];

exportResult[filename_String, result_Association] := Export[
  filename,
  Prepend[result["data"], result["columns"]],
  "CSV"
];

methods = {"plain", "sigma", "e-mod", "dead-zone"};
scenarios = {
  Association[
    "name" -> "nominal",
    "duration" -> 30.0,
    "dt" -> 0.001,
    "gamma" -> 15.0,
    "actuatorTau" -> 0.005,
    "disturbanceAmplitude" -> 0.0,
    "noiseAmplitude" -> 0.0,
    "inputLimit" -> 20.0
  ],
  Association[
    "name" -> "stress",
    "duration" -> 30.0,
    "dt" -> 0.001,
    "gamma" -> 80.0,
    "actuatorTau" -> 0.10,
    "disturbanceAmplitude" -> 0.03,
    "noiseAmplitude" -> 0.01,
    "inputLimit" -> 20.0
  ]
};

allResults = Association[];

Do[
  Print["\nScenario: ", scenario["name"]];
  Print["method       RMS(e)    peak|e|   peak||theta||   peak|u|"];
  Print["--------------------------------------------------------"];

  Do[
    result = simulateMRAC[method, scenario];
    metric = metrics[result];
    allResults[scenario["name"] <> "-" <> method] = result;

    Print @ Row[{
      PaddedForm[method, {11, 0}], " ",
      NumberForm[metric["rmsError"], {7, 4}], " ",
      NumberForm[metric["peakError"], {8, 4}], " ",
      NumberForm[metric["peakThetaNorm"], {12, 4}], " ",
      NumberForm[metric["peakControl"], {8, 4}]
    }];

    exportResult[
      "Chapter8_Lesson5_" <> scenario["name"] <> "_" <> method <> ".csv",
      result
    ],
    {method, methods}
  ],
  {scenario, scenarios}
];

stressPlots = Table[
  With[{data = allResults["stress-" <> method]["data"]},
    {
      ListLinePlot[data[[All, {1, 5}]], PlotRange -> All, PlotLegends -> {method}],
      ListLinePlot[data[[All, {1, 9}]], PlotRange -> All, PlotLegends -> {method}],
      ListLinePlot[data[[All, {1, 6}]], PlotRange -> All, PlotLegends -> {method}]
    }
  ],
  {method, methods}
];

GraphicsGrid[
  Transpose[stressPlots],
  Frame -> All,
  ImageSize -> Large
]

Problem 1 (Failure of the Nominal Proof): For plain MRAC with bounded \(|w(t)|\leq\bar w\), show why the standard Lyapunov derivative does not prove bounded parameter error.

Solution:

\[ \dot V=-a_m e^2+e w \leq -\frac{a_m}{2}e^2+\frac{\bar w^2}{2a_m}. \]

The upper bound is negative only when the tracking error is sufficiently large. It contains no negative term involving \(\widetilde{\boldsymbol{\theta} }\). Hence a large parameter error with small tracking error is not forced back toward a compact set by this inequality. The result is insufficient to establish boundedness of \(V\), because its parameter component can increase while the derivative estimate remains nonnegative inside the error residual set.

Problem 2 (Ultimate Set for σ-Modification): Starting from the inequality in Section 5, identify a sufficient compact set outside which \(\dot V<0\).

Solution: Define

\[ c=\frac{\bar w^2}{2a_m} +\frac{b\sigma}{2}\|\boldsymbol{\theta}^{*}\|^2. \]

A sufficient condition for a negative derivative is

\[ \frac{a_m}{2}e^2+\frac{b\sigma}{2} \|\widetilde{\boldsymbol{\theta} }\|^2>c. \]

Therefore trajectories are driven toward the ellipsoidal set defined by the reverse inequality. The estimate is conservative because Young's inequality discards sign and correlation information. Increasing \(\sigma\) strengthens the parameter term but also changes \(c\) through the ideal-parameter bias contribution.

Problem 3 (Dead-Zone Design): Suppose closed-loop tests show that sensor noise alone produces a tracking-error amplitude of \(0.018\), while a small but important plant change produces a steady error of \(0.030\). Select a plausible dead-zone width and explain the trade-off.

Solution:

A value between the two amplitudes, such as \(e_0=0.022\) to \(0.025\), rejects most noise-only updates while allowing the plant-change error to activate adaptation. Choosing \(e_0\leq0.018\) leaves frequent noise-driven updates. Choosing \(e_0\geq0.030\) can hide the meaningful plant change. A margin is required because the measured bounds are uncertain and may vary with operating point.

Problem 4 (Why Dead Zone Did Not Solve the Stress Case): Explain why the dead-zone law produced a parameter norm near that of plain MRAC in the stress benchmark.

Solution:

The actuator lag creates tracking errors substantially larger than the threshold. Outside the dead zone, the law is essentially the plain gradient update with only a small subtraction from the error magnitude. It contains no term proportional to the current parameter estimate. Consequently, once the parasitic dynamics keep the error outside the inactive band, the parameters can continue growing until input saturation or another nonlinear limitation dominates. The result demonstrates that noise gating and parameter damping solve different robustness problems.

Problem 5 (Leakage Bias in the Disturbance-Free Equilibrium): Assume a constant sufficiently rich regressor and perfect tracking would require \(\boldsymbol{\theta}=\boldsymbol{\theta}^{*}\neq0\). Why can this point not be an equilibrium of the σ-modified update law?

Solution:

At perfect tracking, \(e=0\). The update becomes

\[ \dot{\boldsymbol{\theta} }=-\mathbf{\Gamma}\sigma \boldsymbol{\theta}. \]

At \(\boldsymbol{\theta}=\boldsymbol{\theta}^{*}\neq0\), this derivative is nonzero, so the ideal parameter is not an equilibrium. As leakage moves the estimate toward the origin, a tracking error develops; the gradient term then balances the leakage. The resulting equilibrium has nonzero tracking or prediction error, which is the structural source of leakage bias.

Plain MRAC is optimized for the matched nominal model and can provide the best nominal tracking, but its Lyapunov proof does not automatically extend to persistent disturbances or unmodeled dynamics. σ-modification adds uniform parameter damping and supports a direct ultimate-boundedness proof, at the cost of constant leakage bias. e-modification concentrates damping during large-error intervals and can reduce quiet-operation bias. Dead-zone adaptation is effective against known small-error noise but does not directly restrain parameters when neglected dynamics generate large errors. Robust design therefore requires simultaneous evaluation of tracking error, parameter magnitude, control effort, saturation, and long-duration behavior.