Chapter 28: Advanced Directions and Links to Other Courses

Lesson 5: Summary of Adaptive Control Foundations and Research Directions

This capstone lesson unifies the mathematical foundations developed across the course: matching conditions, model-reference error dynamics, Lyapunov adaptation, recursive estimation, robustness modifications, excitation, constraints, output feedback, nonlinear and MIMO extensions, and practical implementation. It then formulates major research directions as precise control-theoretic questions rather than as a list of fashionable topics.

1. Learning Objectives and the Adaptive-Control Knowledge Map

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

  • reduce a broad adaptive-control design to a common error-model and regressor form;
  • construct and differentiate a state-parameter Lyapunov function;
  • distinguish tracking convergence, boundedness, parameter convergence, robustness, and safety;
  • select direct, indirect, composite, robustified, output-feedback, or constrained architectures;
  • identify assumptions that are structural, assumptions that are excitation-related, and assumptions that are implementation-related;
  • formulate a research problem with a theorem, an excitation condition, a robustness claim, and reproducible numerical evidence.
flowchart TD
  A["Uncertain plant and control objective"] --> B["Choose parameterization and matching assumptions"]
  B --> C["Derive tracking or prediction error model"]
  C --> D["Select estimator or adaptive law"]
  D --> E["Construct Lyapunov or input-output certificate"]
  E --> F["Prove boundedness and tracking"]
  F --> G["Add excitation argument for parameter convergence"]
  G --> H["Add robustness, projection, normalization, and constraints"]
  H --> I["Validate sampling, noise, delay, saturation, and computation"]
  I --> J["Extend toward learning, safety, MPC, MIMO, or distributed control"]
        

The central discipline is to keep four claims separate. A controller may guarantee \( e(t)\to 0 \) without guaranteeing \( \widetilde{\theta}(t)\to 0 \); it may guarantee nominal asymptotic stability without robustness to unmodeled dynamics; and it may guarantee stability without enforcing state or input safety constraints.

2. A Unified Model for the Course

A large portion of model-reference adaptive control can be represented by the uncertain system

\[ \dot{\mathbf{x}} = \mathbf{A}\mathbf{x} + \mathbf{B} \left( u+\boldsymbol{\theta}^{*\top}\boldsymbol{\phi}(\mathbf{x},t) \right) + \mathbf{B}_d\mathbf{d}(t), \qquad \dot{\mathbf{x}}_m = \mathbf{A}_m\mathbf{x}_m+\mathbf{B}_m r. \]

The desired model matrix \( \mathbf{A}_m \) is Hurwitz. In the matched case there exist ideal gains \( \mathbf{K}_x^* \) and \( K_r^* \) satisfying

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

A direct adaptive controller is

\[ u = \widehat{\mathbf{K}}_x^{\top}\mathbf{x} + \widehat{K}_r r - \widehat{\boldsymbol{\theta}}^{\top}\boldsymbol{\phi}(\mathbf{x},t). \]

Define the tracking error \( \mathbf{e}=\mathbf{x}-\mathbf{x}_m \), collect all adjustable parameters into \( \widehat{\boldsymbol{\vartheta}} \), and define the composite regressor

\[ \boldsymbol{\omega} = \begin{bmatrix} \mathbf{x}^{\top} & r & -\boldsymbol{\phi}^{\top} \end{bmatrix}^{\top}, \qquad \widetilde{\boldsymbol{\vartheta}} = \widehat{\boldsymbol{\vartheta}}-\boldsymbol{\vartheta}^*. \]

The error dynamics become

\[ \dot{\mathbf{e}} = \mathbf{A}_m\mathbf{e} + \mathbf{B} \widetilde{\boldsymbol{\vartheta}}^{\top}\boldsymbol{\omega} + \mathbf{B}_d\mathbf{d}(t). \]

This equation is the common analytical core behind many direct MRAC, adaptive backstepping, robot adaptation, and composite-learning designs. Indirect STR methods reach a related form after a plant-parameter estimator is connected to a controller-design map.

3. The Lyapunov Skeleton: Stability and Tracking Proof

Choose symmetric matrices \( \mathbf{Q}=\mathbf{Q}^{\top}\succ 0 \) and \( \mathbf{P}=\mathbf{P}^{\top}\succ 0 \) satisfying

\[ \mathbf{A}_m^{\top}\mathbf{P} + \mathbf{P}\mathbf{A}_m = -\mathbf{Q}. \]

For a positive-definite adaptation-gain matrix \( \boldsymbol{\Gamma} \), consider

\[ V = \mathbf{e}^{\top}\mathbf{P}\mathbf{e} + \widetilde{\boldsymbol{\vartheta}}^{\top} \boldsymbol{\Gamma}^{-1} \widetilde{\boldsymbol{\vartheta}}. \]

With constant ideal parameters and no disturbance, differentiation gives

\[ \dot V = -\mathbf{e}^{\top}\mathbf{Q}\mathbf{e} + 2\mathbf{e}^{\top}\mathbf{P}\mathbf{B} \widetilde{\boldsymbol{\vartheta}}^{\top}\boldsymbol{\omega} + 2\widetilde{\boldsymbol{\vartheta}}^{\top} \boldsymbol{\Gamma}^{-1} \dot{\widehat{\boldsymbol{\vartheta}}}. \]

Select the gradient update law

\[ \dot{\widehat{\boldsymbol{\vartheta}}} = -\boldsymbol{\Gamma} \boldsymbol{\omega} \mathbf{e}^{\top}\mathbf{P}\mathbf{B}. \]

The cross terms cancel exactly:

\[ \dot V=-\mathbf{e}^{\top}\mathbf{Q}\mathbf{e} \leq -\lambda_{\min}(\mathbf{Q})\|\mathbf{e}\|^2. \]

Therefore \( \mathbf{e},\widetilde{\boldsymbol{\vartheta}}\in \mathcal{L}_{\infty} \) and \( \mathbf{e}\in\mathcal{L}_2 \). If the reference, regressor, and closed-loop derivatives remain bounded, then \( \dot{\mathbf{e}}\in\mathcal{L}_{\infty} \). Barbalat's lemma then yields

\[ \lim_{t\to\infty}\mathbf{e}(t)=\mathbf{0}. \]

This proof establishes nominal tracking and bounded adaptation. It does not by itself prove convergence of the estimated parameters to their true or ideal values.

4. Tracking Convergence Is Not Parameter Identification

The regressor is persistently exciting when there exist \( T>0 \) and \( \alpha>0 \) such that, for every \( t\geq 0 \),

\[ \int_t^{t+T} \boldsymbol{\omega}(\tau)\boldsymbol{\omega}^{\top}(\tau) \,d\tau \succeq \alpha\mathbf{I}. \]

Persistent excitation supplies information in every parameter direction. Under appropriate regularity and identifiability conditions, it upgrades bounded parameter error to asymptotic or exponential parameter convergence. Without excitation, the plant may track perfectly while the estimates settle anywhere in an indistinguishable parameter set:

\[ \mathcal{E} = \left\{ \widehat{\boldsymbol{\vartheta}} \;\middle|\; \widetilde{\boldsymbol{\vartheta}}^{\top} \boldsymbol{\omega}(t) = 0 \right\}. \]

Here, the equality must hold along the realized system trajectory.

Claim Typical sufficient ingredients What it does not automatically imply
Bounded closed-loop signals Positive Lyapunov function and nonpositive derivative Tracking convergence or correct parameters
Asymptotic tracking Boundedness, \(e\in\mathcal{L}_2\), bounded \(\dot e\) Parameter convergence
Parameter convergence PE, interval excitation plus memory, or another identifiability condition Robustness to unmodeled dynamics
Robustness Leakage, projection, filtering, normalization, margins, or ISS analysis Hard state and input safety
Safety Invariant-set, barrier, reachability, or constraint-feasibility certificate Optimality or exact learning

5. Robust Modifications and Their Mathematical Meaning

5.1 Projection

Projection constrains estimates to a known compact set \( \Omega \) and is designed to satisfy the key Lyapunov inequality

\[ \widetilde{\boldsymbol{\vartheta}}^{\top} \boldsymbol{\Gamma}^{-1} \left[ \operatorname{Proj} \left( \widehat{\boldsymbol{\vartheta}},\mathbf{y} \right)-\mathbf{y} \right] \leq 0. \]

Consequently, replacing the nominal update direction by its projection does not destroy the basic Lyapunov estimate when the ideal parameter is inside the admissible set.

5.2 Normalization

A normalized law limits update magnitude under large regressors:

\[ \dot{\widehat{\boldsymbol{\vartheta}}} = -\boldsymbol{\Gamma} \frac{ \boldsymbol{\omega}\mathbf{e}^{\top}\mathbf{P}\mathbf{B} }{ 1+\boldsymbol{\omega}^{\top}\boldsymbol{\omega} }. \]

5.3 Leakage or \(\sigma\)-Modification

Introduce \( -\sigma\widehat{\boldsymbol{\vartheta}} \) into the update law:

\[ \dot{\widehat{\boldsymbol{\vartheta}}} = -\boldsymbol{\Gamma} \boldsymbol{\omega}\mathbf{e}^{\top}\mathbf{P}\mathbf{B} - \sigma\widehat{\boldsymbol{\vartheta}}, \qquad \sigma>0. \]

Since \( \widehat{\boldsymbol{\vartheta}} =\widetilde{\boldsymbol{\vartheta}}+\boldsymbol{\vartheta}^* \), Young's inequality yields a derivative bound of the form

\[ \dot V \leq -c_e\|\mathbf{e}\|^2 - c_{\vartheta}\|\widetilde{\boldsymbol{\vartheta}}\|^2 + c_0\sigma\|\boldsymbol{\vartheta}^*\|^2. \]

Thus leakage typically exchanges exact nominal convergence for uniform ultimate boundedness and improved resistance to parameter drift.

5.4 Dead Zones and Saturation-Aware Adaptation

A dead zone stops adaptation for sufficiently small measured error, while saturation-aware laws prevent the estimator from treating actuator clipping as plant uncertainty. These modifications must be included in the proof; they are not merely software guards.

6. Architecture Selection: A Synthesis of the Course

Architecture Adjusted quantity Main analytical advantage Main limitation
Direct MRAC Controller parameters Clean tracking-error Lyapunov design Ideal gains may lack physical meaning; matching is restrictive
Indirect STR Plant estimate, then controller Interpretable model and reuse of controller synthesis Estimator-controller coupling and loss of stabilizing estimates
Direct STR Controller polynomial or gain vector Avoids explicit plant reconstruction Controller parameterization must remain valid
Adaptive backstepping Virtual and actual controls plus parameters Constructive treatment of strict-feedback nonlinear systems Complexity growth, derivatives, and robustness issues
Composite or concurrent learning Tracking and prediction-error parameters Faster learning under finite-data richness Memory management, derivative estimation, and data quality
Gain scheduling / LPV Precomputed gains indexed by scheduling variables Transparent implementation across operating regimes Not true adaptation unless unknown dynamics are learned online
Adaptive MPC Prediction model and receding-horizon policy Explicit constraints and performance optimization Recursive feasibility, computation, and dual-control effects

No architecture is universally superior. The appropriate choice depends on measurable signals, model structure, matching, control direction, excitation, constraints, computational budget, and the level of certification required.

7. Research Direction I: Learning from Finite and Stored Data

Classical gradient adaptation uses only the current tracking error. Concurrent and composite learning add prediction information from stored or filtered data. For regression samples \( y_j=\boldsymbol{\omega}_j^{\top}\boldsymbol{\vartheta}^* \), define the history Gramian

\[ \mathbf{G}_H = \sum_{j=1}^{N} \boldsymbol{\omega}_j\boldsymbol{\omega}_j^{\top}. \]

A finite-excitation condition is \( \mathbf{G}_H\succeq\alpha_H\mathbf{I} \). A representative composite update is

\[ \dot{\widehat{\boldsymbol{\vartheta}}} = -\boldsymbol{\Gamma} \boldsymbol{\omega}\mathbf{e}^{\top}\mathbf{P}\mathbf{B} - k_H\boldsymbol{\Gamma} \sum_{j=1}^{N} \boldsymbol{\omega}_j \left( \boldsymbol{\omega}_j^{\top} \widehat{\boldsymbol{\vartheta}}-y_j \right). \]

The second term contributes \( -2k_H\widetilde{\boldsymbol{\vartheta}}^{\top} \mathbf{G}_H\widetilde{\boldsymbol{\vartheta}} \) to \( \dot V \), creating a direct route to parameter-error decay without requiring the current trajectory to remain persistently exciting forever.

Open questions include data replacement under noise, derivative-free prediction errors, finite-memory robustness, time-varying parameters, and computationally efficient excitation tests.

8. Research Direction II: Dynamic Regressor Extension and Mixing

Begin with the vector regression \( y=\boldsymbol{\varphi}^{\top}\boldsymbol{\theta}^* \). Stable filtering and regressor extension produce

\[ \mathbf{Y}(t)=\boldsymbol{\Phi}(t)\boldsymbol{\theta}^*. \]

Multiplication by the adjugate matrix yields scalar regressions:

\[ \operatorname{adj}(\boldsymbol{\Phi})\mathbf{Y} = \det(\boldsymbol{\Phi})\boldsymbol{\theta}^*, \qquad \mathcal{Y}_i=\Delta\theta_i^*, \quad \Delta=\det(\boldsymbol{\Phi}). \]

With the scalar estimator

\[ \dot{\widehat{\theta}}_i = -\gamma_i\Delta \left( \Delta\widehat{\theta}_i-\mathcal{Y}_i \right). \]

the parameter error obeys

\[ \dot{\widetilde{\theta}}_i = -\gamma_i \Delta^2 \widetilde{\theta}_i, \qquad \widetilde{\theta}_i(t) = \widetilde{\theta}_i(0) \exp\left( -\gamma_i \int_0^t \Delta^2(\tau)\,d\tau \right). \]

Hence each scalar error magnitude is nonincreasing, and convergence holds when \( \Delta\notin\mathcal{L}_2 \). Research continues on filter selection, noise bias, nonlinear regressions, time-varying parameters, finite/fixed-time estimation, and integration with adaptive observers and MIMO control.

9. Research Direction III: Learning-Based Adaptive Control

Learning-based adaptive control replaces or augments a fixed regressor with a learned approximation:

\[ f(\mathbf{x}) = \mathbf{W}^{*\top}\boldsymbol{\sigma}(\mathbf{x}) + \varepsilon(\mathbf{x}), \qquad \|\varepsilon(\mathbf{x})\|\leq\bar{\varepsilon} \text{ on a certified domain}. \]

The adaptive law can estimate output weights online while basis functions, features, kernels, or neural representations are fixed or updated on a slower time scale. The Lyapunov derivative then contains an approximation residual:

\[ \dot V \leq -\lambda_Q\|\mathbf{e}\|^2 + 2\|\mathbf{e}^{\top}\mathbf{P}\mathbf{B}\| \bar{\varepsilon}. \]

This normally gives an ultimate bound rather than exact asymptotic tracking unless the approximation error vanishes or is compensated. The central research issue is not whether a neural network can fit data; it is whether the complete closed loop has verifiable stability, robustness, excitation, generalization, and computational properties.

Important interfaces include meta-learning of initial adaptive parameters, Gaussian-process uncertainty bounds, reinforcement-learning policy improvement around a stabilizing baseline, and online residual learning with a Lyapunov-certified safety filter.

10. Research Direction IV: Safety-Critical Adaptive Control

Let the safe set be \( \mathcal{C}=\{\mathbf{x}:h(\mathbf{x})\geq 0\} \). A control barrier condition for a control-affine system is

\[ L_fh(\mathbf{x}) + L_gh(\mathbf{x})u + \alpha(h(\mathbf{x})) \geq 0. \]

With parametric uncertainty, the estimated term must be supplemented by a valid uncertainty margin. A typical online safety filter solves

\[ \begin{aligned} \min_{u,\delta}\quad & \frac{1}{2}\|u-u_{\mathrm{ad}}\|^2+\rho\delta^2 \\ \text{subject to}\quad & \dot V_{\mathrm{clf}}(u,\widehat{\boldsymbol{\theta}}) \leq -cV_{\mathrm{clf}}+\delta, \\ & L_fh+L_ghu +L_{\theta}h\,\widehat{\boldsymbol{\theta}} -\mathcal{M}_{\theta} +\alpha(h) \geq 0, \\ & u_{\min}\leq u\leq u_{\max}. \end{aligned} \]

The margin \( \mathcal{M}_{\theta} \) must be justified by a parameter-error bound, set-membership estimate, confidence set, or reachability computation. The principal research challenge is preserving forward invariance while the uncertainty set changes and while the optimization remains feasible.

11. Research Direction V: Adaptive MPC, Output Feedback, and MIMO Systems

11.1 Adaptive MPC

An adaptive MPC controller uses an estimated model \( \widehat{\mathbf{x}}_{k+1} =\widehat{\mathbf{A}}_k\mathbf{x}_k+ \widehat{\mathbf{B}}_k\mathbf{u}_k \) inside a constrained finite-horizon problem. A rigorous design must connect:

  • estimator consistency or bounded error;
  • constraint tightening or robust tubes;
  • recursive feasibility;
  • closed-loop stability;
  • the dual effect, in which control actions influence future information.

11.2 Output-Feedback Adaptation

When the full state is unavailable, estimation error and parameter error interact nonlinearly. Unlike fixed-parameter linear control, a general separation principle is unavailable. Research therefore studies adaptive observers, high-gain observers, DREM-based observers, dynamic scaling, and output-feedback normal forms.

11.3 MIMO and Distributed Adaptation

MIMO systems introduce matrix high-frequency gains, interaction, noncommutativity, and richer matching conditions. Distributed adaptive control additionally requires graph connectivity, local excitation, communication delays, and resilience to inconsistent or adversarial information. A typical network Lyapunov function couples local tracking and disagreement energies:

\[ V_N = \sum_{i=1}^{N} \left( \mathbf{e}_i^{\top}\mathbf{P}_i\mathbf{e}_i + \widetilde{\boldsymbol{\theta}}_i^{\top} \boldsymbol{\Gamma}_i^{-1} \widetilde{\boldsymbol{\theta}}_i \right) + c_z\mathbf{z}^{\top} (\mathbf{L}\otimes\mathbf{I}) \mathbf{z}. \]

12. A Rigorous Research Workflow

flowchart TD
  P["Define plant class, measurements, uncertainty, and constraints"] --> A["State assumptions explicitly"]
  A --> C["Choose comparison baseline and adaptive architecture"]
  C --> M["Derive exact closed-loop error model"]
  M --> T["State theorem: stability, convergence, robustness, or safety"]
  T --> R["Prove with Lyapunov, ISS, passivity, small gain, or invariance"]
  R --> E["Specify excitation and identifiability conditions"]
  E --> N["Test noise, delay, sampling, saturation, and unmodeled dynamics"]
  N --> B["Benchmark against fixed, robust, and alternative adaptive controllers"]
  B --> X["Release parameters, code, seeds, metrics, and failure cases"]
        

A publishable contribution should identify which assumption is weakened, which guarantee is strengthened, or which implementation barrier is removed. Merely replacing a regressor by a neural network or adding an optimizer does not constitute a complete adaptive-control result.

13. Reproducibility and Evaluation Metrics

Report at least the following quantities:

\[ J_e = \sqrt{ \frac{1}{T_f-T_0} \int_{T_0}^{T_f}\|\mathbf{e}(t)\|^2dt }, \qquad J_u = \int_0^{T_f}\|u(t)\|^2dt. \]

\[ J_{\theta} = \|\widetilde{\boldsymbol{\theta}}(T_f)\|, \qquad J_{\Delta u} = \int_0^{T_f}\|\dot u(t)\|^2dt, \qquad J_{\mathrm{viol}} = \int_0^{T_f} \max(0,-h(\mathbf{x}(t)))\,dt. \]

Also report peak error, control saturation duration, minimum safety margin, adaptation-rate peaks, computation time per step, sensitivity to sampling time, and results over noise seeds. Compare with a fixed controller and with a robust nonadaptive controller; otherwise the value of adaptation is not isolated.

14. Unified Computational Lab

All five implementations simulate the scalar plant \( \dot x=-ax+bu \) and reference model \( \dot x_m=-a_mx_m+b_mr \). The controller is

\[ u=\widehat{k}_x x+\widehat{k}_r r, \qquad k_x^*=\frac{a-a_m}{b}, \qquad k_r^*=\frac{b_m}{b}. \]

The update laws use projection and leakage:

\[ \dot{\widehat{k}}_x = \operatorname{Proj} \left( \widehat{k}_x,-\gamma_xxe-\sigma\widehat{k}_x \right), \qquad \dot{\widehat{k}}_r = \operatorname{Proj} \left( \widehat{k}_r,-\gamma_rre-\sigma\widehat{k}_r \right). \]

The multi-sine command supplies useful excitation. The leakage term means the final estimates need not equal the ideal values exactly, even when tracking is accurate.

15. Python Implementation

Chapter28_Lesson5.py

Libraries: NumPy for vector integration and Matplotlib for plots. The script exports CSV data and two figures.

# Chapter28_Lesson5.py
# Projected, leakage-modified scalar MRAC used as a capstone experiment.
# Dependencies: numpy, matplotlib

from __future__ import annotations

import csv
from dataclasses import dataclass
from math import sin, sqrt
from typing import Callable

import matplotlib.pyplot as plt
import numpy as np


@dataclass(frozen=True)
class Config:
    a: float = 1.2          # unknown plant pole magnitude
    b: float = 1.0          # known positive input gain
    am: float = 2.0         # reference-model pole magnitude
    bm: float = 2.0         # reference-model command gain
    gamma_x: float = 8.0    # adaptation gain for feedback parameter
    gamma_r: float = 8.0    # adaptation gain for feedforward parameter
    sigma: float = 0.01     # leakage (sigma modification)
    lower: float = -10.0    # projection interval
    upper: float = 10.0
    dt: float = 0.001
    t_end: float = 30.0


def reference(t: float) -> float:
    """Bounded multi-sine command that supplies useful excitation."""
    return sin(0.5 * t) + 0.5 * sin(1.5 * t)


def projected_rate(value: float, raw_rate: float, lower: float, upper: float) -> float:
    """Tangent-cone projection for a scalar interval."""
    if value <= lower and raw_rate < 0.0:
        return 0.0
    if value >= upper and raw_rate > 0.0:
        return 0.0
    return raw_rate


def dynamics(state: np.ndarray, t: float, cfg: Config) -> np.ndarray:
    x, xm, kx_hat, kr_hat = state
    r = reference(t)
    e = x - xm
    u = kx_hat * x + kr_hat * r

    dx = -cfg.a * x + cfg.b * u
    dxm = -cfg.am * xm + cfg.bm * r

    raw_kx = -cfg.gamma_x * x * e - cfg.sigma * kx_hat
    raw_kr = -cfg.gamma_r * r * e - cfg.sigma * kr_hat
    dkx = projected_rate(kx_hat, raw_kx, cfg.lower, cfg.upper)
    dkr = projected_rate(kr_hat, raw_kr, cfg.lower, cfg.upper)
    return np.array([dx, dxm, dkx, dkr], dtype=float)


def rk4_step(
    f: Callable[[np.ndarray, float, Config], np.ndarray],
    state: np.ndarray,
    t: float,
    dt: float,
    cfg: Config,
) -> np.ndarray:
    k1 = f(state, t, cfg)
    k2 = f(state + 0.5 * dt * k1, t + 0.5 * dt, cfg)
    k3 = f(state + 0.5 * dt * k2, t + 0.5 * dt, cfg)
    k4 = f(state + dt * k3, t + dt, cfg)
    next_state = state + (dt / 6.0) * (k1 + 2.0 * k2 + 2.0 * k3 + k4)
    next_state[2:] = np.clip(next_state[2:], cfg.lower, cfg.upper)
    return next_state


def main() -> None:
    cfg = Config()
    n = int(round(cfg.t_end / cfg.dt)) + 1
    time = np.linspace(0.0, cfg.t_end, n)
    states = np.zeros((n, 4), dtype=float)
    controls = np.zeros(n, dtype=float)
    commands = np.zeros(n, dtype=float)

    # Deliberately non-ideal initial controller parameters.
    states[0] = np.array([0.8, 0.0, 0.0, 0.0], dtype=float)

    for k in range(n - 1):
        t = time[k]
        x, _, kx_hat, kr_hat = states[k]
        commands[k] = reference(t)
        controls[k] = kx_hat * x + kr_hat * commands[k]
        states[k + 1] = rk4_step(dynamics, states[k], t, cfg.dt, cfg)

    commands[-1] = reference(time[-1])
    controls[-1] = states[-1, 2] * states[-1, 0] + states[-1, 3] * commands[-1]

    x = states[:, 0]
    xm = states[:, 1]
    error = x - xm
    kx_hat = states[:, 2]
    kr_hat = states[:, 3]

    ideal_kx = (cfg.a - cfg.am) / cfg.b
    ideal_kr = cfg.bm / cfg.b
    second_half = time >= 0.5 * cfg.t_end
    rms_error = sqrt(float(np.mean(error[second_half] ** 2)))

    print(f"Ideal parameters: kx*={ideal_kx:.6f}, kr*={ideal_kr:.6f}")
    print(f"Final estimates:  kx={kx_hat[-1]:.6f}, kr={kr_hat[-1]:.6f}")
    print(f"RMS tracking error over second half: {rms_error:.6e}")

    with open("Chapter28_Lesson5_results.csv", "w", newline="", encoding="utf-8") as stream:
        writer = csv.writer(stream)
        writer.writerow(["t", "r", "x", "xm", "error", "u", "kx_hat", "kr_hat"])
        writer.writerows(
            zip(time, commands, x, xm, error, controls, kx_hat, kr_hat)
        )

    plt.figure()
    plt.plot(time, x, label="plant x")
    plt.plot(time, xm, "--", label="reference xm")
    plt.xlabel("time [s]")
    plt.ylabel("state")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.savefig("Chapter28_Lesson5_tracking.png", dpi=180)

    plt.figure()
    plt.plot(time, kx_hat, label="kx_hat")
    plt.plot(time, kr_hat, label="kr_hat")
    plt.axhline(ideal_kx, linestyle="--", label="kx*")
    plt.axhline(ideal_kr, linestyle="--", label="kr*")
    plt.xlabel("time [s]")
    plt.ylabel("parameter")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.savefig("Chapter28_Lesson5_parameters.png", dpi=180)
    plt.show()


if __name__ == "__main__":
    main()

16. C++ Implementation

Chapter28_Lesson5.cpp

This implementation uses only the C++17 standard library and writes a CSV file suitable for plotting in Python, MATLAB, or a spreadsheet.

// Chapter28_Lesson5.cpp
// Projected, leakage-modified scalar MRAC capstone experiment.
// Build: g++ -std=c++17 -O2 Chapter28_Lesson5.cpp -o Chapter28_Lesson5

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

struct Config {
    double a = 1.2;
    double b = 1.0;
    double am = 2.0;
    double bm = 2.0;
    double gammaX = 8.0;
    double gammaR = 8.0;
    double sigma = 0.01;
    double lower = -10.0;
    double upper = 10.0;
    double dt = 0.001;
    double tEnd = 30.0;
};

using State = std::array<double, 4>;  // x, xm, kx_hat, kr_hat

double reference(double t) {
    return std::sin(0.5 * t) + 0.5 * std::sin(1.5 * t);
}

double projectedRate(double value, double rawRate, double lower, double upper) {
    if (value <= lower && rawRate < 0.0) return 0.0;
    if (value >= upper && rawRate > 0.0) return 0.0;
    return rawRate;
}

State addScaled(const State& a, const State& b, double scale) {
    State result{};
    for (std::size_t i = 0; i < result.size(); ++i) {
        result[i] = a[i] + scale * b[i];
    }
    return result;
}

State dynamics(const State& s, double t, const Config& cfg) {
    const double x = s[0];
    const double xm = s[1];
    const double kxHat = s[2];
    const double krHat = s[3];
    const double r = reference(t);
    const double e = x - xm;
    const double u = kxHat * x + krHat * r;

    const double rawKx = -cfg.gammaX * x * e - cfg.sigma * kxHat;
    const double rawKr = -cfg.gammaR * r * e - cfg.sigma * krHat;

    return State{
        -cfg.a * x + cfg.b * u,
        -cfg.am * xm + cfg.bm * r,
        projectedRate(kxHat, rawKx, cfg.lower, cfg.upper),
        projectedRate(krHat, rawKr, cfg.lower, cfg.upper)
    };
}

State rk4Step(const State& s, double t, const Config& cfg) {
    const State k1 = dynamics(s, t, cfg);
    const State k2 = dynamics(addScaled(s, k1, 0.5 * cfg.dt), t + 0.5 * cfg.dt, cfg);
    const State k3 = dynamics(addScaled(s, k2, 0.5 * cfg.dt), t + 0.5 * cfg.dt, cfg);
    const State k4 = dynamics(addScaled(s, k3, cfg.dt), t + cfg.dt, cfg);

    State next{};
    for (std::size_t i = 0; i < next.size(); ++i) {
        next[i] = s[i] + (cfg.dt / 6.0) *
            (k1[i] + 2.0 * k2[i] + 2.0 * k3[i] + k4[i]);
    }
    next[2] = std::clamp(next[2], cfg.lower, cfg.upper);
    next[3] = std::clamp(next[3], cfg.lower, cfg.upper);
    return next;
}

int main() {
    try {
        const Config cfg;
        const long steps = static_cast<long>(std::llround(cfg.tEnd / cfg.dt));
        State state{0.8, 0.0, 0.0, 0.0};

        std::ofstream csv("Chapter28_Lesson5_results_cpp.csv");
        if (!csv) throw std::runtime_error("Cannot open output CSV file.");
        csv << "t,r,x,xm,error,u,kx_hat,kr_hat\n";
        csv << std::setprecision(12);

        double errorSquareSum = 0.0;
        long errorCount = 0;

        for (long k = 0; k <= steps; ++k) {
            const double t = k * cfg.dt;
            const double r = reference(t);
            const double error = state[0] - state[1];
            const double u = state[2] * state[0] + state[3] * r;

            csv << t << ',' << r << ',' << state[0] << ',' << state[1] << ','
                << error << ',' << u << ',' << state[2] << ',' << state[3] << '\n';

            if (t >= 0.5 * cfg.tEnd) {
                errorSquareSum += error * error;
                ++errorCount;
            }
            if (k < steps) state = rk4Step(state, t, cfg);
        }

        const double idealKx = (cfg.a - cfg.am) / cfg.b;
        const double idealKr = cfg.bm / cfg.b;
        const double rms = std::sqrt(errorSquareSum / static_cast<double>(errorCount));

        std::cout << std::fixed << std::setprecision(6);
        std::cout << "Ideal parameters: kx*=" << idealKx << ", kr*=" << idealKr << '\n';
        std::cout << "Final estimates:  kx=" << state[2] << ", kr=" << state[3] << '\n';
        std::cout << "RMS tracking error over second half: " << rms << '\n';
        return 0;
    } catch (const std::exception& ex) {
        std::cerr << "Error: " << ex.what() << '\n';
        return 1;
    }
}

17. Java Implementation

Chapter28_Lesson5.java

The Java version uses standard arrays, NIO file output, and the same RK4 and projection logic as the Python and C++ versions.

// Chapter28_Lesson5.java
// Projected, leakage-modified scalar MRAC capstone experiment.
// Build and run: javac Chapter28_Lesson5.java && java Chapter28_Lesson5

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 Chapter28_Lesson5 {
    private static final class Config {
        final double a = 1.2;
        final double b = 1.0;
        final double am = 2.0;
        final double bm = 2.0;
        final double gammaX = 8.0;
        final double gammaR = 8.0;
        final double sigma = 0.01;
        final double lower = -10.0;
        final double upper = 10.0;
        final double dt = 0.001;
        final double tEnd = 30.0;
    }

    private static double reference(double t) {
        return Math.sin(0.5 * t) + 0.5 * Math.sin(1.5 * t);
    }

    private static double projectedRate(
            double value, double rawRate, double lower, double upper) {
        if (value <= lower && rawRate < 0.0) {
            return 0.0;
        }
        if (value >= upper && rawRate > 0.0) {
            return 0.0;
        }
        return rawRate;
    }

    private static double[] dynamics(double[] state, double t, Config cfg) {
        double x = state[0];
        double xm = state[1];
        double kxHat = state[2];
        double krHat = state[3];
        double r = reference(t);
        double error = x - xm;
        double u = kxHat * x + krHat * r;

        double rawKx = -cfg.gammaX * x * error - cfg.sigma * kxHat;
        double rawKr = -cfg.gammaR * r * error - cfg.sigma * krHat;

        return new double[] {
            -cfg.a * x + cfg.b * u,
            -cfg.am * xm + cfg.bm * r,
            projectedRate(kxHat, rawKx, cfg.lower, cfg.upper),
            projectedRate(krHat, rawKr, cfg.lower, cfg.upper)
        };
    }

    private static double[] addScaled(double[] a, double[] b, double scale) {
        double[] result = new double[a.length];
        for (int i = 0; i < a.length; i++) {
            result[i] = a[i] + scale * b[i];
        }
        return result;
    }

    private static double[] rk4Step(double[] state, double t, Config cfg) {
        double[] k1 = dynamics(state, t, cfg);
        double[] k2 = dynamics(addScaled(state, k1, 0.5 * cfg.dt),
                t + 0.5 * cfg.dt, cfg);
        double[] k3 = dynamics(addScaled(state, k2, 0.5 * cfg.dt),
                t + 0.5 * cfg.dt, cfg);
        double[] k4 = dynamics(addScaled(state, k3, cfg.dt),
                t + cfg.dt, cfg);

        double[] next = new double[state.length];
        for (int i = 0; i < state.length; i++) {
            next[i] = state[i] + (cfg.dt / 6.0)
                    * (k1[i] + 2.0 * k2[i] + 2.0 * k3[i] + k4[i]);
        }
        next[2] = Math.max(cfg.lower, Math.min(cfg.upper, next[2]));
        next[3] = Math.max(cfg.lower, Math.min(cfg.upper, next[3]));
        return next;
    }

    public static void main(String[] args) {
        Locale.setDefault(Locale.US);
        Config cfg = new Config();
        long steps = Math.round(cfg.tEnd / cfg.dt);
        double[] state = {0.8, 0.0, 0.0, 0.0};

        double errorSquareSum = 0.0;
        long errorCount = 0;

        Path output = Path.of("Chapter28_Lesson5_results_java.csv");
        try (BufferedWriter writer = Files.newBufferedWriter(
                output, StandardCharsets.UTF_8)) {
            writer.write("t,r,x,xm,error,u,kx_hat,kr_hat\n");

            for (long k = 0; k <= steps; k++) {
                double t = k * cfg.dt;
                double r = reference(t);
                double error = state[0] - state[1];
                double u = state[2] * state[0] + state[3] * r;

                writer.write(String.format(Locale.US,
                        "%.12f,%.12f,%.12f,%.12f,%.12f,%.12f,%.12f,%.12f%n",
                        t, r, state[0], state[1], error, u, state[2], state[3]));

                if (t >= 0.5 * cfg.tEnd) {
                    errorSquareSum += error * error;
                    errorCount++;
                }
                if (k < steps) {
                    state = rk4Step(state, t, cfg);
                }
            }
        } catch (IOException ex) {
            System.err.println("I/O error: " + ex.getMessage());
            System.exit(1);
        }

        double idealKx = (cfg.a - cfg.am) / cfg.b;
        double idealKr = cfg.bm / cfg.b;
        double rms = Math.sqrt(errorSquareSum / errorCount);

        System.out.printf("Ideal parameters: kx*=%.6f, kr*=%.6f%n", idealKx, idealKr);
        System.out.printf("Final estimates:  kx=%.6f, kr=%.6f%n", state[2], state[3]);
        System.out.printf("RMS tracking error over second half: %.6e%n", rms);
    }
}

18. MATLAB and Simulink Implementations

Chapter28_Lesson5.m

The base-MATLAB script simulates, plots, and exports the common experiment.

% Chapter28_Lesson5.m
% Projected, leakage-modified scalar MRAC capstone experiment.
% Requires base MATLAB only.

clear; clc; close all;

cfg.a = 1.2;
cfg.b = 1.0;
cfg.am = 2.0;
cfg.bm = 2.0;
cfg.gammaX = 8.0;
cfg.gammaR = 8.0;
cfg.sigma = 0.01;
cfg.lower = -10.0;
cfg.upper = 10.0;
cfg.dt = 0.001;
cfg.tEnd = 30.0;

time = (0:cfg.dt:cfg.tEnd).';
n = numel(time);
state = zeros(n, 4);       % [x, xm, kxHat, krHat]
control = zeros(n, 1);
command = zeros(n, 1);
state(1, :) = [0.8, 0.0, 0.0, 0.0];

for k = 1:n-1
    t = time(k);
    command(k) = referenceSignal(t);
    control(k) = state(k, 3) * state(k, 1) + state(k, 4) * command(k);
    state(k + 1, :) = rk4Step(state(k, :).', t, cfg).';
end

command(end) = referenceSignal(time(end));
control(end) = state(end, 3) * state(end, 1) + state(end, 4) * command(end);

x = state(:, 1);
xm = state(:, 2);
error = x - xm;
kxHat = state(:, 3);
krHat = state(:, 4);

idealKx = (cfg.a - cfg.am) / cfg.b;
idealKr = cfg.bm / cfg.b;
secondHalf = time >= 0.5 * cfg.tEnd;
rmsError = sqrt(mean(error(secondHalf).^2));

fprintf('Ideal parameters: kx*=%.6f, kr*=%.6f\n', idealKx, idealKr);
fprintf('Final estimates:  kx=%.6f, kr=%.6f\n', kxHat(end), krHat(end));
fprintf('RMS tracking error over second half: %.6e\n', rmsError);

results = table(time, command, x, xm, error, control, kxHat, krHat, ...
    'VariableNames', {'t','r','x','xm','error','u','kx_hat','kr_hat'});
writetable(results, 'Chapter28_Lesson5_results_matlab.csv');

figure;
plot(time, x, 'LineWidth', 1.2); hold on;
plot(time, xm, '--', 'LineWidth', 1.2);
grid on; xlabel('time [s]'); ylabel('state');
legend('plant x', 'reference x_m', 'Location', 'best');
title('Projected and leakage-modified MRAC tracking');

figure;
plot(time, kxHat, 'LineWidth', 1.2); hold on;
plot(time, krHat, 'LineWidth', 1.2);
yline(idealKx, '--'); yline(idealKr, '--');
grid on; xlabel('time [s]'); ylabel('parameter');
legend('k_x estimate', 'k_r estimate', 'k_x^*', 'k_r^*', ...
    'Location', 'best');
title('Adaptive parameter trajectories');

function r = referenceSignal(t)
    r = sin(0.5 * t) + 0.5 * sin(1.5 * t);
end

function next = rk4Step(state, t, cfg)
    k1 = dynamics(state, t, cfg);
    k2 = dynamics(state + 0.5 * cfg.dt * k1, t + 0.5 * cfg.dt, cfg);
    k3 = dynamics(state + 0.5 * cfg.dt * k2, t + 0.5 * cfg.dt, cfg);
    k4 = dynamics(state + cfg.dt * k3, t + cfg.dt, cfg);
    next = state + (cfg.dt / 6.0) * (k1 + 2*k2 + 2*k3 + k4);
    next(3:4) = min(cfg.upper, max(cfg.lower, next(3:4)));
end

function ds = dynamics(state, t, cfg)
    x = state(1);
    xm = state(2);
    kxHat = state(3);
    krHat = state(4);
    r = referenceSignal(t);
    error = x - xm;
    u = kxHat * x + krHat * r;

    rawKx = -cfg.gammaX * x * error - cfg.sigma * kxHat;
    rawKr = -cfg.gammaR * r * error - cfg.sigma * krHat;

    ds = [
        -cfg.a * x + cfg.b * u;
        -cfg.am * xm + cfg.bm * r;
        projectedRate(kxHat, rawKx, cfg.lower, cfg.upper);
        projectedRate(krHat, rawKr, cfg.lower, cfg.upper)
    ];
end

function rate = projectedRate(value, rawRate, lower, upper)
    if (value <= lower && rawRate < 0) || ...
       (value >= upper && rawRate > 0)
        rate = 0;
    else
        rate = rawRate;
    end
end

Chapter28_Lesson5_Simulink.m

This script programmatically creates Chapter28_Lesson5_MRAC.slx. It requires Simulink.

% Chapter28_Lesson5_Simulink.m
% Programmatically creates a compact Simulink model for the scalar MRAC.
% Requires Simulink. The generated model uses a MATLAB Function block for
% the plant, reference model, projection, and leakage-modified update law.

modelName = 'Chapter28_Lesson5_MRAC';
if bdIsLoaded(modelName)
    close_system(modelName, 0);
end
new_system(modelName);
open_system(modelName);

add_block('simulink/Sources/Clock', [modelName '/Clock'], ...
    'Position', [40 65 70 95]);
add_block('simulink/User-Defined Functions/MATLAB Function', ...
    [modelName '/MRAC Dynamics'], 'Position', [140 35 360 125]);
add_block('simulink/Sinks/To Workspace', [modelName '/Log'], ...
    'VariableName', 'mracLog', 'SaveFormat', 'StructureWithTime', ...
    'Position', [430 55 530 105]);

functionCode = sprintf([ ...
'function y = fcn(t)\n' ...
'%%#codegen\n' ...
'persistent x xm kx kr lastT\n' ...
'if isempty(lastT)\n' ...
'    x = 0.8; xm = 0; kx = 0; kr = 0; lastT = t;\n' ...
'end\n' ...
'dt = max(t-lastT, 0); lastT = t;\n' ...
'a=1.2; b=1; am=2; bm=2; gamma=8; sigma=0.01;\n' ...
'r = sin(0.5*t) + 0.5*sin(1.5*t);\n' ...
'e = x-xm; u = kx*x + kr*r;\n' ...
'dx = -a*x + b*u; dxm = -am*xm + bm*r;\n' ...
'dkx = -gamma*x*e - sigma*kx;\n' ...
'dkr = -gamma*r*e - sigma*kr;\n' ...
'if (kx<=-10 && dkx<0) || (kx>=10 && dkx>0), dkx=0; end\n' ...
'if (kr<=-10 && dkr<0) || (kr>=10 && dkr>0), dkr=0; end\n' ...
'x=x+dt*dx; xm=xm+dt*dxm;\n' ...
'kx=min(10,max(-10,kx+dt*dkx));\n' ...
'kr=min(10,max(-10,kr+dt*dkr));\n' ...
'y=[r;x;xm;e;u;kx;kr];\n' ...
'end\n']);

block = [modelName '/MRAC Dynamics'];
chart = sfroot.find('-isa', 'Stateflow.EMChart', 'Path', block);
if isempty(chart)
    error('Could not locate the MATLAB Function chart.');
end
chart.Script = functionCode;
set_param(block, 'Position', [140 35 360 125]);

add_line(modelName, 'Clock/1', 'MRAC Dynamics/1');
add_line(modelName, 'MRAC Dynamics/1', 'Log/1');

set_param(modelName, 'StopTime', '30', 'Solver', 'ode4', ...
    'FixedStep', '0.001');
save_system(modelName, [modelName '.slx']);
sim(modelName);

disp('Generated Chapter28_Lesson5_MRAC.slx and workspace variable mracLog.');

19. Wolfram Mathematica Implementation

Chapter28_Lesson5.nb

The file is represented as a textual Wolfram Notebook expression and can be opened directly in Mathematica.

Notebook[{
  Cell["Chapter28_Lesson5.nb", "Title"],
  Cell["Projected, leakage-modified scalar MRAC capstone experiment.", "Text"],
  Cell[BoxData[
    RowBox[{
      RowBox[{"ClearAll", "[", "\"Global`*\"", "]"}], ";", "\n",
      RowBox[{"a", "=", "1.2"}], ";",
      RowBox[{"b", "=", "1.0"}], ";",
      RowBox[{"am", "=", "2.0"}], ";",
      RowBox[{"bm", "=", "2.0"}], ";", "\n",
      RowBox[{"gammaX", "=", "8.0"}], ";",
      RowBox[{"gammaR", "=", "8.0"}], ";",
      RowBox[{"sigma", "=", "0.01"}], ";", "\n",
      RowBox[{"lower", "=", RowBox[{"-", "10.0"}]}], ";",
      RowBox[{"upper", "=", "10.0"}], ";",
      RowBox[{"dt", "=", "0.001"}], ";",
      RowBox[{"tEnd", "=", "30.0"}], ";", "\n",
      RowBox[{
        RowBox[{"reference", "[", "t_", "]"}], ":=",
        RowBox[{
          RowBox[{"Sin", "[", RowBox[{"0.5", " ", "t"}], "]"}], "+",
          RowBox[{"0.5", " ", RowBox[{"Sin", "[", RowBox[{"1.5", " ", "t"}], "]"}]}]
        }]
      }], ";", "\n",
      RowBox[{
        RowBox[{"projectedRate", "[", RowBox[{"value_", ",", "raw_"}], "]"}], ":=",
        RowBox[{"Which", "[",
          RowBox[{
            RowBox[{RowBox[{"value", "<=", "lower"}], "&&", RowBox[{"raw", "<", "0"}]}],
            ",", "0", ",",
            RowBox[{RowBox[{"value", ">=", "upper"}], "&&", RowBox[{"raw", ">", "0"}]}],
            ",", "0", ",", "True", ",", "raw"
          }], "]"}]
      }], ";", "\n",
      RowBox[{"initial", "=", RowBox[{"{", RowBox[{"0.8", ",", "0.0", ",", "0.0", ",", "0.0"}], "}"}]}], ";", "\n",
      RowBox[{
        RowBox[{"rhs", "[", RowBox[{"state_", ",", "t_"}], "]"}], ":=",
        RowBox[{"Module", "[",
          RowBox[{
            RowBox[{"{",
              RowBox[{"x", ",", "xm", ",", "kx", ",", "kr", ",", "r", ",", "e", ",", "u", ",", "rawKx", ",", "rawKr"}],
            "}"}], ",",
            RowBox[{
              RowBox[{RowBox[{"{", RowBox[{"x", ",", "xm", ",", "kx", ",", "kr"}], "}"}], "=", "state"}], ";",
              RowBox[{"r", "=", RowBox[{"reference", "[", "t", "]"}]}], ";",
              RowBox[{"e", "=", RowBox[{"x", "-", "xm"}]}], ";",
              RowBox[{"u", "=", RowBox[{RowBox[{"kx", " ", "x"}], "+", RowBox[{"kr", " ", "r"}]}]}], ";",
              RowBox[{"rawKx", "=", RowBox[{RowBox[{"-", "gammaX", " ", "x", " ", "e"}], "-", RowBox[{"sigma", " ", "kx"}]}]}], ";",
              RowBox[{"rawKr", "=", RowBox[{RowBox[{"-", "gammaR", " ", "r", " ", "e"}], "-", RowBox[{"sigma", " ", "kr"}]}]}], ";",
              RowBox[{"{",
                RowBox[{
                  RowBox[{RowBox[{"-", "a", " ", "x"}], "+", RowBox[{"b", " ", "u"}]}], ",",
                  RowBox[{RowBox[{"-", "am", " ", "xm"}], "+", RowBox[{"bm", " ", "r"}]}], ",",
                  RowBox[{"projectedRate", "[", RowBox[{"kx", ",", "rawKx"}], "]"}], ",",
                  RowBox[{"projectedRate", "[", RowBox[{"kr", ",", "rawKr"}], "]"}]
                }], "}"}]
            }]
          }], "]"}]
      }], ";", "\n",
      RowBox[{
        RowBox[{"rk4", "[", RowBox[{"state_", ",", "t_"}], "]"}], ":=",
        RowBox[{"Module", "[",
          RowBox[{
            RowBox[{"{", RowBox[{"k1", ",", "k2", ",", "k3", ",", "k4", ",", "next"}], "}"}], ",",
            RowBox[{
              RowBox[{"k1", "=", RowBox[{"rhs", "[", RowBox[{"state", ",", "t"}], "]"}]}], ";",
              RowBox[{"k2", "=", RowBox[{"rhs", "[",
                RowBox[{RowBox[{"state", "+", RowBox[{"0.5", " ", "dt", " ", "k1"}]}], ",", RowBox[{"t", "+", RowBox[{"0.5", " ", "dt"}]}]}], "]"}]}], ";",
              RowBox[{"k3", "=", RowBox[{"rhs", "[",
                RowBox[{RowBox[{"state", "+", RowBox[{"0.5", " ", "dt", " ", "k2"}]}], ",", RowBox[{"t", "+", RowBox[{"0.5", " ", "dt"}]}]}], "]"}]}], ";",
              RowBox[{"k4", "=", RowBox[{"rhs", "[",
                RowBox[{RowBox[{"state", "+", RowBox[{"dt", " ", "k3"}]}], ",", RowBox[{"t", "+", "dt"}]}], "]"}]}], ";",
              RowBox[{"next", "=", RowBox[{"state", "+", RowBox[{RowBox[{"dt", "/", "6"}], " ",
                RowBox[{"(", RowBox[{"k1", "+", RowBox[{"2", " ", "k2"}], "+", RowBox[{"2", " ", "k3"}], "+", "k4"}], ")"}]}]}]}], ";",
              RowBox[{"ReplacePart", "[", RowBox[{"next", ",",
                RowBox[{"{",
                  RowBox[{
                    RowBox[{"3", "->", RowBox[{"Clip", "[", RowBox[{RowBox[{"next", "[", RowBox[{"[", "3", "]"}], "]"}], ",", RowBox[{"{", RowBox[{"lower", ",", "upper"}], "}"}]}], "]"}]}], ",",
                    RowBox[{"4", "->", RowBox[{"Clip", "[", RowBox[{RowBox[{"next", "[", RowBox[{"[", "4", "]"}], "]"}], ",", RowBox[{"{", RowBox[{"lower", ",", "upper"}], "}"}]}], "]"}]}]
                  }], "}"}]}], "]"}]
            }]
          }], "]"}]
      }], ";", "\n",
      RowBox[{"times", "=", RowBox[{"Range", "[", RowBox[{"0", ",", "tEnd", ",", "dt"}], "]"}]}], ";", "\n",
      RowBox[{"states", "=", RowBox[{"FoldList", "[", RowBox[{"rk4", ",", "initial", ",", RowBox[{"Most", "[", "times", "]"}]}], "]"}]}], ";", "\n",
      RowBox[{"xSeries", "=", RowBox[{"states", "[", RowBox[{"[", RowBox[{"All", ",", "1"}], "]"}], "]"}]}], ";",
      RowBox[{"xmSeries", "=", RowBox[{"states", "[", RowBox[{"[", RowBox[{"All", ",", "2"}], "]"}], "]"}]}], ";", "\n",
      RowBox[{"errorSeries", "=", RowBox[{"xSeries", "-", "xmSeries"}]}], ";", "\n",
      RowBox[{"idealKx", "=", RowBox[{RowBox[{"(", RowBox[{"a", "-", "am"}], ")"}], "/", "b"}]}], ";",
      RowBox[{"idealKr", "=", RowBox[{"bm", "/", "b"}]}], ";", "\n",
      RowBox[{"rms", "=", RowBox[{"Sqrt", "[", RowBox[{"Mean", "[",
        RowBox[{"Take", "[", RowBox[{RowBox[{"errorSeries", "^", "2"}], ",", RowBox[{"-", RowBox[{"Floor", "[", RowBox[{RowBox[{"Length", "[", "times", "]"}], "/", "2"}], "]"}]}]}], "]"}], "]"}], "]"}]}], ";", "\n",
      RowBox[{"Print", "[", RowBox[{"\"Ideal parameters: \"", ",", RowBox[{"{", RowBox[{"idealKx", ",", "idealKr"}], "}"}]}], "]"}], ";",
      RowBox[{"Print", "[", RowBox[{"\"Final estimates: \"", ",", RowBox[{"Last", "[", "states", "]"}], ",", "\"; RMS error = \"", ",", "rms"}], "]"}], ";", "\n",
      RowBox[{"ListLinePlot", "[",
        RowBox[{
          RowBox[{"{", RowBox[{RowBox[{"Transpose", "[", RowBox[{"{", RowBox[{"times", ",", "xSeries"}], "}"}], "]"}], ",",
          RowBox[{"Transpose", "[", RowBox[{"{", RowBox[{"times", ",", "xmSeries"}], "}"}], "]"}]}], "}"}], ",",
          RowBox[{"PlotLegends", "->", RowBox[{"{", RowBox[{"\"plant x\"", ",", "\"reference xm\""}], "}"}]}], ",",
          RowBox[{"Frame", "->", "True"}]
        }], "]"}]
    }]], "Input"]
},
WindowTitle -> "Chapter28_Lesson5",
Saveable -> True
]

20. Problems and Solutions

Problem 1 (Complete Scalar MRAC Proof): Consider \( \dot x=-ax+bu \), where \( b>0 \) is known and \( a \) is unknown. Let \( \dot x_m=-a_mx_m+b_mr \), with \( a_m>0 \), and use \( u=\widehat{k}_xx+\widehat{k}_rr \). Derive the ideal parameters and adaptive laws that ensure \( e=x-x_m\to 0 \).

Solution:

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

The error equation is

\[ \dot e = -a_me + b\widetilde{k}_x x + b\widetilde{k}_r r. \]

Choose

\[ V = \frac{1}{2}e^2 + \frac{b}{2\gamma_x}\widetilde{k}_x^2 + \frac{b}{2\gamma_r}\widetilde{k}_r^2. \]

Then

\[ \dot V = -a_me^2 + be\widetilde{k}_xx + be\widetilde{k}_rr + \frac{b}{\gamma_x}\widetilde{k}_x\dot{\widehat{k}}_x + \frac{b}{\gamma_r}\widetilde{k}_r\dot{\widehat{k}}_r. \]

Selecting \( \dot{\widehat{k}}_x=-\gamma_xxe \) and \( \dot{\widehat{k}}_r=-\gamma_rre \) gives \( \dot V=-a_me^2 \). All adaptive states are bounded, \( e\in\mathcal{L}_2 \), and bounded reference signals imply bounded \( \dot e \). Barbalat's lemma gives \( e(t)\to 0 \).

Problem 2 (Why Constant Commands Do Not Identify Both Gains): Suppose the reference command becomes a nonzero constant and the closed loop reaches steady state. Explain why perfect tracking need not imply \( \widehat{k}_x\to k_x^* \) and \( \widehat{k}_r\to k_r^* \).

Solution:

At steady state, \( x=x_m=\bar x \) and \( r=\bar r \). Zero tracking error only imposes the single scalar relation

\[ \widetilde{k}_x\bar x+\widetilde{k}_r\bar r=0. \]

This is one equation in two unknown parameter errors. The regressor \( [x\;\;r]^{\top} \) has rank one along the steady trajectory, so infinitely many parameter pairs are indistinguishable. A sufficiently rich time-varying command is needed to identify both.

Problem 3 (Ultimate Bound Under Leakage): For the scalar system in Problem 1, replace the update laws by \( \dot{\widehat{k}}_i=-\gamma_i\psi_i e-\sigma\widehat{k}_i \), where \( \psi_x=x \) and \( \psi_r=r \). Show why the origin of the complete error system is generally replaced by an ultimate bound.

Solution:

The leakage contribution to the derivative is

\[ -\frac{b\sigma}{\gamma_i} \widetilde{k}_i\widehat{k}_i = -\frac{b\sigma}{\gamma_i}\widetilde{k}_i^2 - \frac{b\sigma}{\gamma_i}\widetilde{k}_ik_i^*. \]

Using Young's inequality,

\[ -\widetilde{k}_ik_i^* \leq \frac{1}{2}\widetilde{k}_i^2 + \frac{1}{2}(k_i^*)^2. \]

Therefore

\[ \dot V \leq -a_me^2 - \sum_i \frac{b\sigma}{2\gamma_i}\widetilde{k}_i^2 + \sum_i \frac{b\sigma}{2\gamma_i}(k_i^*)^2. \]

The derivative is negative outside a compact set whose size depends on \( \sigma \) and the ideal parameters. This proves uniform ultimate boundedness rather than exact convergence to the complete error origin.

Problem 4 (DREM Scalar Error Monotonicity): Given \( \mathcal{Y}_i=\Delta\theta_i^* \) and \( \dot{\widehat{\theta}}_i =-\gamma_i\Delta(\Delta\widehat{\theta}_i-\mathcal{Y}_i) \), prove that \( |\widetilde{\theta}_i(t)| \) is nonincreasing.

Solution:

\[ \dot{\widetilde{\theta}}_i = -\gamma_i\Delta^2\widetilde{\theta}_i. \]

For \( W_i=\frac{1}{2}\widetilde{\theta}_i^2 \),

\[ \dot W_i = -\gamma_i\Delta^2\widetilde{\theta}_i^2 \leq 0. \]

Thus \( |\widetilde{\theta}_i| \) cannot increase. Furthermore, the explicit solution shows convergence when \( \int_0^{\infty}\Delta^2(\tau)d\tau=\infty \).

Problem 5 (Robust Adaptive Barrier Constraint): Suppose \( \dot{\mathbf{x}} =\mathbf{f}(\mathbf{x})+\mathbf{g}(\mathbf{x})u+ \mathbf{Y}(\mathbf{x})\boldsymbol{\theta}^* \) and the estimator satisfies \( \|\widetilde{\boldsymbol{\theta}}\|\leq\bar\theta \). Derive a conservative barrier inequality.

Solution:

The uncertain contribution to \( \dot h \) is \( \nabla h^{\top}\mathbf{Y}\boldsymbol{\theta}^* \). Write it as the estimated term plus parameter error and bound the latter with Cauchy-Schwarz:

\[ \nabla h^{\top}\mathbf{Y}\boldsymbol{\theta}^* \geq \nabla h^{\top}\mathbf{Y}\widehat{\boldsymbol{\theta}} - \|\nabla h^{\top}\mathbf{Y}\|\bar\theta. \]

A sufficient robust barrier constraint is therefore

\[ L_fh+L_gh\,u + \nabla h^{\top}\mathbf{Y}\widehat{\boldsymbol{\theta}} - \|\nabla h^{\top}\mathbf{Y}\|\bar\theta + \alpha(h) \geq 0. \]

If the parameter-error bound is valid and the inequality remains feasible, the uncertainty margin protects forward invariance.

Problem 6 (Research-Claim Audit): A paper reports lower RMS tracking error after replacing a conventional regressor by a deep neural network. List the minimum additional evidence required before concluding that the method advances adaptive-control theory.

Solution:

The paper should specify the uncertainty class and approximation domain; derive the exact closed-loop error dynamics; state and prove boundedness, convergence, robustness, or safety claims; quantify approximation error; define excitation or identifiability conditions; compare against fixed and robust baselines; test noise, delay, saturation, sampling, and unmodeled dynamics; report computational timing; and release enough details for reproduction. Improved simulation RMS alone does not provide a stability or generalization guarantee.

21. Summary

Adaptive control is best understood as a disciplined interaction among parameterization, feedback, online estimation, and stability analysis. Matching conditions create a tractable error model; Lyapunov design establishes boundedness and tracking; excitation and identifiability determine whether parameters converge; and projection, normalization, leakage, filtering, and saturation logic address practical robustness.

Current research seeks to learn from finite data, relax persistent excitation, integrate DREM and memory-based estimators, combine adaptive control with certified learning, enforce safety constraints, preserve feasibility in adaptive MPC, and extend rigorous guarantees to output-feedback, MIMO, networked, and computationally constrained systems. The enduring requirement is unchanged: every new learning mechanism must be connected to a verifiable closed-loop claim.

22. References

  1. Rohrs, C.E., Valavani, L., Athans, M., & Stein, G. (1985). Robustness of continuous-time adaptive control algorithms in the presence of unmodeled dynamics. IEEE Transactions on Automatic Control, 30(9), 881–889.
  2. Narendra, K.S., & Annaswamy, A.M. (1986). Robust adaptive control in the presence of bounded disturbances. IEEE Transactions on Automatic Control, 31(4), 306–315.
  3. Boyd, S., & Sastry, S.S. (1986). Necessary and sufficient conditions for parameter convergence in adaptive control. Automatica, 22(6), 629–639.
  4. Narendra, K.S., & Annaswamy, A.M. (1987). Persistent excitation in adaptive systems. International Journal of Control, 45(1), 127–160.
  5. Slotine, J.-J.E., & Li, W. (1987). On the adaptive control of robot manipulators. The International Journal of Robotics Research, 6(3), 49–59.
  6. Cao, C., & Hovakimyan, N. (2008). Design and analysis of a novel \( \mathcal{L}_1 \) adaptive control architecture with guaranteed transient performance. IEEE Transactions on Automatic Control, 53(2), 586–591.
  7. Chowdhary, G., Yucelen, T., Mühlegg, M., & Johnson, E.N. (2013). Concurrent learning adaptive control of linear systems with exponentially convergent bounds. International Journal of Adaptive Control and Signal Processing, 27(4), 280–301.
  8. Ames, A.D., Xu, X., Grizzle, J.W., & Tabuada, P. (2017). Control barrier function based quadratic programs for safety critical systems. IEEE Transactions on Automatic Control, 62(8), 3861–3876.
  9. Pan, Y., & Yu, H. (2018). Composite learning robot control with guaranteed parameter convergence. Automatica, 89, 398–406.
  10. Ortega, R., Praly, L., Aranovskiy, S., Yi, B., & Zhang, W. (2018). On dynamic regressor extension and mixing parameter estimators: Two Luenberger observers interpretations. Automatica, 95, 548–551.
  11. Loría, A., Panteley, E., & Maghenem, M. (2019). Strict Lyapunov functions for model reference adaptive control: Application to Lagrangian systems. IEEE Transactions on Automatic Control, 64(7), 3040–3045.
  12. Taylor, A.J., & Ames, A.D. (2020). Adaptive safety with control barrier functions. Proceedings of the American Control Conference, 1399–1405.
  13. Annaswamy, A.M., & Fradkov, A.L. (2021). A historical perspective of adaptive control and learning. Annual Reviews in Control, 52, 18–41.
  14. Annaswamy, A.M. (2023). Adaptive control and intersections with reinforcement learning. Annual Review of Control, Robotics, and Autonomous Systems, 6, 65–93.
Support CaaT Academy

Help keep these engineering tutorials free and growing

If these lessons, examples, and project pages help you, a small donation supports the continued creation and improvement of free control, robotics, software, and engineering education resources.

Created and maintained by Abolfazl Mohammadijoo.