1. Learning Objectives and Scope

After completing this lesson, students should be able to:

distinguish direct MRAC from an indirect STR using signal flow and parameter meaning;
derive a first-order Lyapunov-based MRAC law and identify its assumptions;
derive a certainty-equivalent STR control law from an estimated discrete-time plant;
explain why plant-estimation convergence and closed-loop tracking are different claims;
analyze sensitivity to weak excitation, small estimated input gain, noise, saturation, and unmodeled dynamics;
select between STR and MRAC using model structure, performance specifications, and verification needs.

The comparison uses concepts already developed in earlier chapters: linearly parameterized models, Lyapunov functions, projection, normalization, persistent excitation, and direct versus indirect adaptation. Detailed recursive least squares is intentionally deferred to Chapter 12; the STR implementation here uses a normalized-gradient identifier so that the architectural comparison remains the focus.

2. A Common Control Objective

Consider the uncertain first-order plant \( \dot y(t)=-a y(t)+b u(t) \), where \(a\) and \(b\) are unknown constants, while the sign of \(b\) is known and positive. The desired response is generated by the stable reference model

\[ \dot y_m(t)=-a_m y_m(t)+b_m r(t), \qquad a_m>0. \]

Define the tracking error \(e(t)=y(t)-y_m(t)\). Both MRAC and STR attempt to make this error small, but they adapt different quantities:

\[ \begin{aligned} \text{MRAC:}\quad & \hat{\boldsymbol{k}}(t) \text{ is updated directly from tracking information},\\ \text{STR:}\quad & \hat{\boldsymbol{\theta}}(t) \text{ identifies the plant, and } \hat{\boldsymbol{k}}(t)=\mathcal{D}\!\left( \hat{\boldsymbol{\theta}}(t)\right). \end{aligned} \]

Here \(\mathcal{D}\) denotes the controller-design map: pole placement, minimum variance, linear-quadratic design, PID tuning, or another synthesis rule applied to the current plant estimate.

3. Architectural Distinction

flowchart TD
  R1["Reference r"] --> RM["Reference model"]
  RM --> E["Tracking error"]
  E --> AL["MRAC adaptation law"]
  AL --> CK["Controller gains"]
  R1 --> MC["MRAC controller"]
  CK --> MC
  MC --> P1["Plant"]
  P1 --> E

  R2["Reference r"] --> CD["Controller design map"]
  ID["Online plant identifier"] --> CD
  CD --> SC["STR controller"]
  R2 --> SC
  SC --> P2["Plant"]
  P2 --> ID

In direct MRAC, the adapted parameters normally have the meaning of controller gains. In an indirect STR, the adapted parameters normally have the meaning of plant coefficients. The latter are passed through a controller-design map before they affect the actuator command.

This distinction creates different failure channels. MRAC can suffer from an unsuitable matching structure or a fragile tracking-error adaptation loop. STR can suffer from both poor identification and a badly conditioned redesign map, even when each block appears reasonable in isolation.

4. Direct MRAC: Derivation and Lyapunov Argument

Choose the direct adaptive control law

\[ u(t)=\hat k_r(t)r(t)-\hat k_y(t)y(t). \]

If the plant parameters were known, exact model matching would require

\[ a+b k_y^=a_m, \qquad b k_r^=b_m, \]

and therefore

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

Define controller-parameter errors \(\tilde k_y=\hat k_y-k_y^*\) and \(\tilde k_r=\hat k_r-k_r^*\). Substitution gives

\[ \dot e=-a_m e-b\tilde k_y y+b\tilde k_r r. \]

For known positive input-gain sign, consider

\[ V=\frac{1}{2}e^2+\frac{b}{2\gamma_y}\tilde k_y^2+ \frac{b}{2\gamma_r}\tilde k_r^2, \qquad \gamma_y>0,\;\gamma_r>0. \]

Its derivative is

\[ \dot V=e\left(-a_m e-b\tilde k_y y+b\tilde k_r r\right) +\frac{b}{\gamma_y}\tilde k_y\dot{\hat k}_y +\frac{b}{\gamma_r}\tilde k_r\dot{\hat k}_r. \]

Choose the adaptive laws

\[ \dot{\hat k}_y=\gamma_y e y, \qquad \dot{\hat k}_r=-\gamma_r e r. \]

The cross terms cancel, yielding

\[ \dot V=-a_m e^2\le 0. \]

Thus \(e\), \(\tilde k_y\), and \(\tilde k_r\) are bounded. Under the usual bounded-signal conditions, \(\dot e\) is bounded and Barbalat's lemma gives \(e(t)\) tending to zero. This conclusion is a tracking result; it does not imply that the controller gains converge to their ideal values unless sufficient excitation is present.

Primary MRAC strength: the adaptation law is constructed together with a stability proof for the combined tracking and parameter dynamics. Primary limitation: the proof depends on the assumed plant class, matching conditions, known control direction, and absence of neglected dynamics that invalidate the error model.

5. Indirect STR: Identification Followed by Controller Redesign

Sample the first-order plant and represent it by the one-step model

\[ y_{k+1}=A y_k+B u_k, \qquad \boldsymbol{\phi}_k=\begin{bmatrix}y_k & u_k\end{bmatrix}^T, \qquad \boldsymbol{\theta}=\begin{bmatrix}A & B\end{bmatrix}^T. \]

The prediction and prediction error are

\[ \hat y_{k+1|k}=\hat{\boldsymbol{\theta}}_k^T\boldsymbol{\phi}_k, \qquad \varepsilon_{k+1}=y_{k+1}- \hat{\boldsymbol{\theta}}_k^T\boldsymbol{\phi}_k. \]

A normalized-gradient identifier with projection is

\[ \hat{\boldsymbol{\theta}}_{k+1}= \operatorname{Proj}_\Theta\!\left[ \hat{\boldsymbol{\theta}}_k+ \frac{\gamma\boldsymbol{\phi}_k\varepsilon_{k+1}} {\delta+\boldsymbol{\phi}_k^T\boldsymbol{\phi}_k} \right], \qquad \gamma>0,\;\delta>0. \]

Let the desired sampled closed-loop relation be \(y_{k+1}=A_m y_k+B_m r_k\). The certainty-equivalent control law is obtained by replacing the unknown coefficients with their estimates:

\[ u_k=\frac{(A_m-\hat A_k)y_k+B_m r_k}{\hat B_k}. \]

This formula shows both the attraction and danger of an STR. Once a plant model is available, a familiar controller-design equation can be reused. However, the map becomes singular when the estimated input gain approaches zero. Projection, a known gain sign, a lower gain bound, and control saturation are therefore not cosmetic implementation details.

5.1 Error Decomposition

Substitution into the true plant gives

\[ y_{k+1}=A_m y_k+B_m r_k+ \left(A-\hat A_k\right)y_k+ \left(B-\hat B_k\right)u_k. \]

The departure from the target model is therefore driven by prediction error components filtered through closed-loop signals. Unlike the simple MRAC derivation above, an estimator convergence result alone does not automatically establish stable transients for every intermediate controller generated by \(\mathcal{D}(\hat{\boldsymbol{\theta}})\).

5.2 Controller-Map Sensitivity

For the scalar pole-placement law,

\[ \frac{\partial u_k}{\partial \hat A_k}=-\frac{y_k}{\hat B_k}, \qquad \frac{\partial u_k}{\partial \hat B_k}=- \frac{(A_m-\hat A_k)y_k+B_m r_k}{\hat B_k^2}. \]

The second derivative scales with \(1/\hat B_k^2\). Consequently, a modest input-gain estimation error can cause a large control perturbation when \(\hat B_k\) is small. This is a characteristic STR drawback that is absent from controller parameterizations which do not divide by an estimated plant coefficient.

6. Advantages and Drawbacks: Structured Comparison

Criterion	Self-Tuning Regulator (STR)	Model Reference Adaptive Control (MRAC)
Adapted quantity	Plant-model parameters, then controller gains through a design map.	Controller parameters directly, usually from tracking error.
Interpretability	High when identified coefficients correspond to physical or transfer-function parameters.	Controller gains may be less useful as explicit plant estimates.
Use of classical synthesis	Strong advantage: estimated models can feed pole placement, minimum variance, LQ, or PID redesign.	Controller structure and adaptive law are usually derived together for a prescribed reference model.
Stability proof	Often conditional and modular: identification, admissible estimates, and redesign must jointly preserve stability.	Lyapunov design can directly certify boundedness and tracking for the assumed plant class.
Excitation requirements	Weak excitation degrades the plant estimate and can produce unreliable redesign.	Tracking can converge without full parameter convergence, but gain identification still requires excitation.
Transient behavior	May show bursts when estimated models cross poorly conditioned regions or redesign changes rapidly.	May show aggressive transients when adaptation gains are large or reference-model demands are fast.
Noise sensitivity	Prediction-error identification is directly affected by measurement noise and regressor conditioning.	Tracking-error adaptation can integrate noise into gain motion and parameter drift.
Model changes	Useful when explicit model estimates are needed for supervision, diagnosis, or controller redesign.	Useful when the main requirement is direct tracking rather than online plant interpretation.
Structural assumptions	Requires a valid identifiable model order and a safe controller design for all admissible estimates.	Requires model matching, suitable relative degree, and known control direction in standard forms.
Implementation complexity	More modular but usually more components: estimator, projection, redesign, validation, and fallback logic.	Can be compact for low-order matched plants, although robust modifications add complexity.

7. Stability, Convergence, and Robustness Implications

7.1 Frozen-Time Stability Is Not Enough

Suppose every fixed estimate in an admissible set produces a stable controller. This does not by itself prove stability of the time-varying adaptive loop. If the estimate moves rapidly, the closed-loop system is time varying, and repeated switching among individually stable dynamics can still produce large transients or instability.

\[ x_{k+1}=A_{cl}\!\left(\hat{\boldsymbol{\theta}}_k\right)x_k, \qquad \rho\!\left(A_{cl}(\vartheta)\right)<1 \text{ for each fixed }\vartheta \not\Rightarrow \text{uniform stability for arbitrary }\hat{\boldsymbol{\theta}}_k. \]

A rigorous STR analysis therefore needs more than nominal pole placement: bounds on parameter motion, an invariant admissible set, common Lyapunov arguments, dwell-time reasoning, or a specific adaptive-control theorem.

7.2 Tracking Convergence Is Not Parameter Convergence

In the MRAC proof, \(\dot V=-a_m e^2\) ensures an integrability property for the tracking error, but it does not force \(\tilde{\boldsymbol{k}}\) to zero. Similarly, an STR can regulate the output so effectively that its regressor loses excitation, leaving several plant models indistinguishable in closed loop. The persistent-excitation condition introduced in Chapter 10 is

\[ \sum_{i=k}^{k+N-1}\boldsymbol{\phi}_i \boldsymbol{\phi}_i^T\succeq \alpha I, \qquad N\in\mathbb{N},\;\alpha>0. \]

Without this richness condition, bounded tracking may coexist with nonunique or drifting parameter estimates.

7.3 Unmodeled Dynamics and Disturbances

Both architectures can react to a modeling error as though it were a correctable parameter error. In MRAC this may drive controller gains; in STR it may bias the plant identifier and repeatedly redesign the controller around a fictitious low-order model. The robust modifications from Chapters 8 and 9—leakage, dead zones, normalization, and projection—are therefore relevant to both methods.

\[ \begin{aligned} \text{measured mismatch} &= \text{parametric mismatch} + \text{unmodeled dynamics} \\ &\quad + \text{disturbance} + \text{measurement noise}. \end{aligned} \]

An adaptive law generally cannot separate these terms without additional structural information or filtering.

7.4 Saturation Changes the Mathematics

If the commanded input is \(u_c\) but the applied input is \(u=\operatorname{sat}(u_c)\), the nominal error model contains an additional term:

\[ \dot e=-a_m e-b\tilde k_y y+b\tilde k_r r+ b\left[\operatorname{sat}(u_c)-u_c\right]. \]

The simple cancellation in the MRAC Lyapunov proof is no longer exact. In an STR, saturation also corrupts identification if the estimator uses the commanded input instead of the actual actuator input. Practical implementations must log and use the applied input.

8. Selection Logic: When Each Architecture Is More Natural

flowchart TD
  S["Start with uncertainty and \nperformance requirements"] --> M["Need explicit online \nplant estimates?"]
  M -->|"yes"| C["Can every admissible estimate \nproduce a safe controller?"]
  C -->|"yes"| STR["Prefer STR with projection \nand redesign validation"]
  C -->|"no"| SAFE["Restrict model set or \nuse a supervisory fallback"]
  M -->|"no"| R["Does a matched reference-model \nstructure exist?"]
  R -->|"yes"| MRAC["Prefer MRAC with \nLyapunov-based adaptation"]
  R -->|"no"| ALT["Reconsider parameterization or\n another adaptive architecture"]
  STR --> V["Verify excitation, \nsaturation, noise, \nand rate limits"]
  MRAC --> V
  SAFE --> V

Choose STR when explicit online identification has independent value, the controller-design map is well conditioned over a certified parameter set, and the application benefits from modular redesign. Choose MRAC when a reliable matching model exists and direct tracking-error adaptation allows a transparent Lyapunov analysis. In safety-critical systems, either method normally requires projection, monitoring, actuator limits, and a fixed fallback controller.

9. Comparative Numerical Experiment

The implementations below simulate two copies of the same plant. The continuous coefficients change halfway through the experiment:

\[ (a,b)=\begin{cases} (1.2,1.0), & k<600,\\ (0.7,0.65), & k\ge 600. \end{cases} \]

The sampling interval is \(T_s=0.02\) seconds. The reference-model sampled coefficients are

\[ A_m=e^{-a_mT_s}, \qquad B_m=1-A_m, \qquad a_m=1.8. \]

The STR uses normalized-gradient estimation, projection, a small probing signal, certainty-equivalent pole placement, and actuator saturation. The MRAC uses the Lyapunov-gradient laws derived earlier, discretized by forward Euler, with gain projection and the same actuator limit.

This experiment is pedagogical rather than a theorem-level benchmark. Saturation and numerical discretization mean that the exact continuous Lyapunov proof does not apply without additional analysis. Also, a lower tracking error in one tuning does not establish general superiority; it only illustrates the different transient mechanisms.

10. Python Implementation

Required libraries are numpy and matplotlib. The standard library supplies data classes and mathematical functions.

`Chapter11_Lesson4.py`

"""Chapter 11, Lesson 4: STR versus MRAC.

Pedagogical comparison on the same uncertain first-order plant.  The STR uses
an online normalized-gradient identifier and certainty-equivalent pole
placement.  The MRAC directly adjusts controller gains from tracking error.
"""

from __future__ import annotations

import math
from dataclasses import dataclass

import matplotlib.pyplot as plt
import numpy as np


@dataclass(frozen=True)
class Config:
    sample_time: float = 0.02
    samples: int = 1200
    model_rate: float = 1.8
    control_limit: float = 8.0


def reference_signal(samples: int) -> np.ndarray:
    k = np.arange(samples + 1)
    return np.where((k // 200) % 2 == 0, 1.0, -0.5)


def plant_parameters(k: int, samples: int) -> tuple[float, float]:
    """Return continuous-time coefficients in y_dot = -a*y + b*u."""
    return (1.2, 1.0) if k < samples // 2 else (0.7, 0.65)


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


def simulate() -> dict[str, np.ndarray]:
    cfg = Config()
    ts, n = cfg.sample_time, cfg.samples
    r = reference_signal(n)

    # Exact sampled reference model: y_m[k+1] = A_m*y_m[k] + B_m*r[k].
    a_m_d = math.exp(-cfg.model_rate * ts)
    b_m_d = 1.0 - a_m_d

    y_model = np.zeros(n + 1)
    y_mrac = np.zeros(n + 1)
    y_str = np.zeros(n + 1)
    u_mrac = np.zeros(n)
    u_str = np.zeros(n)

    # Direct MRAC gains for u = k_r*r - k_y*y.
    k_y_hat, k_r_hat = 0.0, 0.2
    gamma_y, gamma_r = 1.4, 1.8
    k_y_hist = np.zeros(n + 1)
    k_r_hist = np.zeros(n + 1)
    k_y_hist[0], k_r_hist[0] = k_y_hat, k_r_hat

    # STR estimate of sampled plant y[k+1] = A*y[k] + B*u[k].
    theta_hat = np.array([0.94, 0.018], dtype=float)
    theta_hist = np.zeros((n + 1, 2))
    theta_hist[0] = theta_hat
    identifier_gain, normalization = 0.55, 1.0e-4
    previous_phi = np.zeros(2)

    for k in range(n):
        a_true, b_true = plant_parameters(k, n)
        a_true_d = 1.0 - ts * a_true
        b_true_d = ts * b_true

        # Reference model shared by the two adaptive controllers.
        y_model[k + 1] = a_m_d * y_model[k] + b_m_d * r[k]

        # ---------------- Direct MRAC ----------------
        tracking_error = y_mrac[k] - y_model[k]
        k_y_hat += ts * gamma_y * tracking_error * y_mrac[k]
        k_r_hat += ts * (-gamma_r * tracking_error * r[k])
        k_y_hat = project(k_y_hat, -20.0, 20.0)
        k_r_hat = project(k_r_hat, -20.0, 20.0)

        u_mrac[k] = np.clip(
            k_r_hat * r[k] - k_y_hat * y_mrac[k],
            -cfg.control_limit,
            cfg.control_limit,
        )
        y_mrac[k + 1] = y_mrac[k] + ts * (
            -a_true * y_mrac[k] + b_true * u_mrac[k]
        )
        k_y_hist[k + 1], k_r_hist[k + 1] = k_y_hat, k_r_hat

        # ---------------- Indirect STR ----------------
        # At time k, y_str[k] is available.  It identifies the model from the
        # previous regressor [y_str[k-1], u_str[k-1]].
        if k > 0:
            prediction_error = y_str[k] - float(theta_hat @ previous_phi)
            theta_hat += (
                identifier_gain
                * previous_phi
                * prediction_error
                / (normalization + float(previous_phi @ previous_phi))
            )
            # Projection prevents a near-zero estimated input gain and limits
            # unreasonable sampled-pole estimates.
            theta_hat[0] = project(theta_hat[0], 0.0, 1.20)
            theta_hat[1] = project(theta_hat[1], 0.004, 0.060)

        a_hat_d, b_hat_d = theta_hat
        probing_signal = 0.025 * math.sin(0.11 * k)
        u_ce = ((a_m_d - a_hat_d) * y_str[k] + b_m_d * r[k]) / b_hat_d
        u_str[k] = np.clip(
            u_ce + probing_signal, -cfg.control_limit, cfg.control_limit
        )
        y_str[k + 1] = a_true_d * y_str[k] + b_true_d * u_str[k]
        previous_phi = np.array([y_str[k], u_str[k]], dtype=float)
        theta_hist[k + 1] = theta_hat

    time = ts * np.arange(n + 1)
    return {
        "time": time,
        "reference": r,
        "model": y_model,
        "mrac": y_mrac,
        "str": y_str,
        "u_mrac": u_mrac,
        "u_str": u_str,
        "k_y": k_y_hist,
        "k_r": k_r_hist,
        "theta": theta_hist,
    }


def main() -> None:
    data = simulate()
    rmse_mrac = float(np.sqrt(np.mean((data["mrac"] - data["model"]) ** 2)))
    rmse_str = float(np.sqrt(np.mean((data["str"] - data["model"]) ** 2)))
    print(f"MRAC tracking RMSE: {rmse_mrac:.6f}")
    print(f"STR tracking RMSE:  {rmse_str:.6f}")
    print("The plant parameters change halfway through the run.")

    plt.figure(figsize=(10, 5))
    plt.plot(data["time"], data["model"], label="reference model")
    plt.plot(data["time"], data["mrac"], label="MRAC output")
    plt.plot(data["time"], data["str"], label="STR output")
    plt.xlabel("time [s]")
    plt.ylabel("output")
    plt.title("STR versus MRAC on a changing first-order plant")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()

    plt.figure(figsize=(10, 5))
    plt.plot(data["time"][:-1], data["u_mrac"], label="MRAC control")
    plt.plot(data["time"][:-1], data["u_str"], label="STR control")
    plt.xlabel("time [s]")
    plt.ylabel("control input")
    plt.title("Adaptive control effort")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.show()


if __name__ == "__main__":
    main()

11. C++ Implementation

This version uses only the C++17 standard library. It exports a CSV file that can be plotted with Python, MATLAB, a spreadsheet, or another tool.

`Chapter11_Lesson4.cpp`

// Chapter 11, Lesson 4: STR versus MRAC.
// Build: g++ -std=c++17 -O2 Chapter11_Lesson4.cpp -o Chapter11_Lesson4

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

namespace {
constexpr double kTs = 0.02;
constexpr int kSamples = 1200;
constexpr double kModelRate = 1.8;
constexpr double kControlLimit = 8.0;

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

double reference(int k) {
    return ((k / 200) % 2 == 0) ? 1.0 : -0.5;
}

void plantParameters(int k, double& a, double& b) {
    if (k < kSamples / 2) {
        a = 1.2;
        b = 1.0;
    } else {
        a = 0.7;
        b = 0.65;
    }
}
}  // namespace

int main() {
    const double am = std::exp(-kModelRate * kTs);
    const double bm = 1.0 - am;

    std::vector<double> model(kSamples + 1, 0.0);
    std::vector<double> mrac(kSamples + 1, 0.0);
    std::vector<double> str(kSamples + 1, 0.0);
    std::vector<double> uMrac(kSamples, 0.0);
    std::vector<double> uStr(kSamples, 0.0);

    // Direct MRAC parameters for u = kr*r - ky*y.
    double kyHat = 0.0;
    double krHat = 0.2;
    const double gammaY = 1.4;
    const double gammaR = 1.8;

    // STR estimate for y[k+1] = A*y[k] + B*u[k].
    double aHat = 0.94;
    double bHat = 0.018;
    const double identifierGain = 0.55;
    const double normalization = 1.0e-4;
    double previousY = 0.0;
    double previousU = 0.0;

    double squaredErrorMrac = 0.0;
    double squaredErrorStr = 0.0;

    std::ofstream csv("Chapter11_Lesson4_results.csv");
    if (!csv) {
        throw std::runtime_error("Cannot create Chapter11_Lesson4_results.csv");
    }
    csv << "time,reference_model,mrac_output,str_output,mrac_control,str_control,A_hat,B_hat\n";
    csv << std::setprecision(10);

    for (int k = 0; k < kSamples; ++k) {
        const double r = reference(k);
        double aTrue = 0.0;
        double bTrue = 0.0;
        plantParameters(k, aTrue, bTrue);
        const double aTrueD = 1.0 - kTs * aTrue;
        const double bTrueD = kTs * bTrue;

        model[k + 1] = am * model[k] + bm * r;

        // ---------------- Direct MRAC ----------------
        const double e = mrac[k] - model[k];
        kyHat = clip(kyHat + kTs * gammaY * e * mrac[k], -20.0, 20.0);
        krHat = clip(krHat - kTs * gammaR * e * r, -20.0, 20.0);
        uMrac[k] = clip(krHat * r - kyHat * mrac[k],
                        -kControlLimit, kControlLimit);
        mrac[k + 1] = mrac[k] + kTs * (-aTrue * mrac[k] + bTrue * uMrac[k]);

        // ---------------- Indirect STR ----------------
        if (k > 0) {
            const double prediction = aHat * previousY + bHat * previousU;
            const double epsilon = str[k] - prediction;
            const double denominator = normalization
                                     + previousY * previousY
                                     + previousU * previousU;
            aHat += identifierGain * previousY * epsilon / denominator;
            bHat += identifierGain * previousU * epsilon / denominator;
            aHat = clip(aHat, 0.0, 1.20);
            bHat = clip(bHat, 0.004, 0.060);
        }

        const double probingSignal = 0.025 * std::sin(0.11 * k);
        const double certaintyEquivalent =
            ((am - aHat) * str[k] + bm * r) / bHat;
        uStr[k] = clip(certaintyEquivalent + probingSignal,
                       -kControlLimit, kControlLimit);
        str[k + 1] = aTrueD * str[k] + bTrueD * uStr[k];
        previousY = str[k];
        previousU = uStr[k];

        squaredErrorMrac += std::pow(mrac[k] - model[k], 2.0);
        squaredErrorStr += std::pow(str[k] - model[k], 2.0);

        csv << k * kTs << ',' << model[k] << ',' << mrac[k] << ',' << str[k]
            << ',' << uMrac[k] << ',' << uStr[k] << ',' << aHat << ',' << bHat
            << '\n';
    }

    const double rmseMrac = std::sqrt(squaredErrorMrac / kSamples);
    const double rmseStr = std::sqrt(squaredErrorStr / kSamples);
    std::cout << std::fixed << std::setprecision(6)
              << "MRAC tracking RMSE: " << rmseMrac << '\n'
              << "STR tracking RMSE:  " << rmseStr << '\n'
              << "Results written to Chapter11_Lesson4_results.csv\n";
    return 0;
}

12. Java Implementation

This implementation uses the Java standard library and writes the same comparison data to a CSV file.

`Chapter11_Lesson4.java`

// Chapter 11, Lesson 4: STR versus MRAC.
// Build and run: javac Chapter11_Lesson4.java && java Chapter11_Lesson4

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

public final class Chapter11_Lesson4 {
    private static final double TS = 0.02;
    private static final int SAMPLES = 1200;
    private static final double MODEL_RATE = 1.8;
    private static final double CONTROL_LIMIT = 8.0;

    private Chapter11_Lesson4() {
    }

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

    private static double reference(int k) {
        return ((k / 200) % 2 == 0) ? 1.0 : -0.5;
    }

    private static double[] plantParameters(int k) {
        return (k < SAMPLES / 2)
                ? new double[] {1.2, 1.0}
                : new double[] {0.7, 0.65};
    }

    public static void main(String[] args) throws IOException {
        Locale.setDefault(Locale.US);
        final double modelA = Math.exp(-MODEL_RATE * TS);
        final double modelB = 1.0 - modelA;

        double[] model = new double[SAMPLES + 1];
        double[] mrac = new double[SAMPLES + 1];
        double[] str = new double[SAMPLES + 1];
        double[] uMrac = new double[SAMPLES];
        double[] uStr = new double[SAMPLES];

        // Direct MRAC parameters for u = kr*r - ky*y.
        double kyHat = 0.0;
        double krHat = 0.2;
        final double gammaY = 1.4;
        final double gammaR = 1.8;

        // STR estimate for y[k+1] = A*y[k] + B*u[k].
        double aHat = 0.94;
        double bHat = 0.018;
        final double identifierGain = 0.55;
        final double normalization = 1.0e-4;
        double previousY = 0.0;
        double previousU = 0.0;

        double squaredErrorMrac = 0.0;
        double squaredErrorStr = 0.0;

        Path output = Path.of("Chapter11_Lesson4_results.csv");
        try (BufferedWriter writer = Files.newBufferedWriter(output)) {
            writer.write("time,reference_model,mrac_output,str_output,"
                    + "mrac_control,str_control,A_hat,B_hat\n");

            for (int k = 0; k < SAMPLES; k++) {
                double r = reference(k);
                double[] plant = plantParameters(k);
                double aTrue = plant[0];
                double bTrue = plant[1];
                double aTrueD = 1.0 - TS * aTrue;
                double bTrueD = TS * bTrue;

                model[k + 1] = modelA * model[k] + modelB * r;

                // ---------------- Direct MRAC ----------------
                double trackingError = mrac[k] - model[k];
                kyHat = clip(
                        kyHat + TS * gammaY * trackingError * mrac[k],
                        -20.0,
                        20.0);
                krHat = clip(
                        krHat - TS * gammaR * trackingError * r,
                        -20.0,
                        20.0);
                uMrac[k] = clip(
                        krHat * r - kyHat * mrac[k],
                        -CONTROL_LIMIT,
                        CONTROL_LIMIT);
                mrac[k + 1] = mrac[k]
                        + TS * (-aTrue * mrac[k] + bTrue * uMrac[k]);

                // ---------------- Indirect STR ----------------
                if (k > 0) {
                    double prediction = aHat * previousY + bHat * previousU;
                    double predictionError = str[k] - prediction;
                    double denominator = normalization
                            + previousY * previousY
                            + previousU * previousU;
                    aHat += identifierGain * previousY
                            * predictionError / denominator;
                    bHat += identifierGain * previousU
                            * predictionError / denominator;
                    aHat = clip(aHat, 0.0, 1.20);
                    bHat = clip(bHat, 0.004, 0.060);
                }

                double probingSignal = 0.025 * Math.sin(0.11 * k);
                double certaintyEquivalent =
                        ((modelA - aHat) * str[k] + modelB * r) / bHat;
                uStr[k] = clip(
                        certaintyEquivalent + probingSignal,
                        -CONTROL_LIMIT,
                        CONTROL_LIMIT);
                str[k + 1] = aTrueD * str[k] + bTrueD * uStr[k];
                previousY = str[k];
                previousU = uStr[k];

                squaredErrorMrac += Math.pow(mrac[k] - model[k], 2.0);
                squaredErrorStr += Math.pow(str[k] - model[k], 2.0);

                writer.write(String.format(
                        Locale.US,
                        "%.10f,%.10f,%.10f,%.10f,%.10f,%.10f,%.10f,%.10f%n",
                        k * TS,
                        model[k],
                        mrac[k],
                        str[k],
                        uMrac[k],
                        uStr[k],
                        aHat,
                        bHat));
            }
        }

        double rmseMrac = Math.sqrt(squaredErrorMrac / SAMPLES);
        double rmseStr = Math.sqrt(squaredErrorStr / SAMPLES);
        System.out.printf("MRAC tracking RMSE: %.6f%n", rmseMrac);
        System.out.printf("STR tracking RMSE:  %.6f%n", rmseStr);
        System.out.println("Results written to " + output.toAbsolutePath());
    }
}

13. MATLAB Implementation

The MATLAB script simulates the two adaptive loops, plots outputs, control actions, and parameter histories, and exports the results.

`Chapter11_Lesson4.m`

% Chapter 11, Lesson 4: STR versus MRAC.
% Pedagogical comparison on the same changing first-order plant.

clear; clc; close all;

Ts = 0.02;
N = 1200;
modelRate = 1.8;
controlLimit = 8.0;

k = 0:N;
time = k * Ts;
r = ones(1, N + 1);
r(mod(floor(k / 200), 2) == 1) = -0.5;

Am = exp(-modelRate * Ts);
Bm = 1 - Am;

yModel = zeros(1, N + 1);
yMRAC = zeros(1, N + 1);
ySTR = zeros(1, N + 1);
uMRAC = zeros(1, N);
uSTR = zeros(1, N);

% Direct MRAC: u = krHat*r - kyHat*y.
kyHat = 0.0;
krHat = 0.2;
gammaY = 1.4;
gammaR = 1.8;
kyHistory = zeros(1, N + 1);
krHistory = zeros(1, N + 1);
kyHistory(1) = kyHat;
krHistory(1) = krHat;

% Indirect STR model: y(k+1) = A*y(k) + B*u(k).
thetaHat = [0.94; 0.018];
thetaHistory = zeros(2, N + 1);
thetaHistory(:, 1) = thetaHat;
identifierGain = 0.55;
normalization = 1.0e-4;
previousPhi = [0; 0];

for index = 1:N
    sample = index - 1;
    if sample < N / 2
        aTrue = 1.2;
        bTrue = 1.0;
    else
        aTrue = 0.7;
        bTrue = 0.65;
    end
    Atrue = 1 - Ts * aTrue;
    Btrue = Ts * bTrue;

    yModel(index + 1) = Am * yModel(index) + Bm * r(index);

    % ---------------- Direct MRAC ----------------
    trackingError = yMRAC(index) - yModel(index);
    kyHat = clipValue(kyHat + Ts * gammaY * trackingError * yMRAC(index), ...
                      -20, 20);
    krHat = clipValue(krHat - Ts * gammaR * trackingError * r(index), ...
                      -20, 20);
    uMRAC(index) = clipValue(krHat * r(index) - kyHat * yMRAC(index), ...
                             -controlLimit, controlLimit);
    yMRAC(index + 1) = yMRAC(index) + Ts * ...
        (-aTrue * yMRAC(index) + bTrue * uMRAC(index));
    kyHistory(index + 1) = kyHat;
    krHistory(index + 1) = krHat;

    % ---------------- Indirect STR ----------------
    if index > 1
        predictionError = ySTR(index) - thetaHat.' * previousPhi;
        thetaHat = thetaHat + identifierGain * previousPhi * ...
            predictionError / (normalization + previousPhi.' * previousPhi);
        thetaHat(1) = clipValue(thetaHat(1), 0.0, 1.20);
        thetaHat(2) = clipValue(thetaHat(2), 0.004, 0.060);
    end

    probingSignal = 0.025 * sin(0.11 * sample);
    certaintyEquivalent = ((Am - thetaHat(1)) * ySTR(index) ...
                          + Bm * r(index)) / thetaHat(2);
    uSTR(index) = clipValue(certaintyEquivalent + probingSignal, ...
                            -controlLimit, controlLimit);
    ySTR(index + 1) = Atrue * ySTR(index) + Btrue * uSTR(index);
    previousPhi = [ySTR(index); uSTR(index)];
    thetaHistory(:, index + 1) = thetaHat;
end

rmseMRAC = sqrt(mean((yMRAC - yModel).^2));
rmseSTR = sqrt(mean((ySTR - yModel).^2));
fprintf('MRAC tracking RMSE: %.6f\n', rmseMRAC);
fprintf('STR tracking RMSE:  %.6f\n', rmseSTR);

figure('Name', 'STR versus MRAC outputs');
plot(time, yModel, 'LineWidth', 1.5); hold on;
plot(time, yMRAC, 'LineWidth', 1.2);
plot(time, ySTR, 'LineWidth', 1.2);
grid on;
xlabel('Time [s]');
ylabel('Output');
title('STR versus MRAC on a changing first-order plant');
legend('Reference model', 'MRAC output', 'STR output', ...
       'Location', 'best');

figure('Name', 'Adaptive control effort');
plot(time(1:end-1), uMRAC, 'LineWidth', 1.2); hold on;
plot(time(1:end-1), uSTR, 'LineWidth', 1.2);
grid on;
xlabel('Time [s]');
ylabel('Control input');
title('Adaptive control effort');
legend('MRAC control', 'STR control', 'Location', 'best');

figure('Name', 'Adaptive parameters');
tiledlayout(2, 1);
nexttile;
plot(time, kyHistory, time, krHistory, 'LineWidth', 1.2);
grid on;
ylabel('MRAC gains');
legend('k_y estimate', 'k_r estimate', 'Location', 'best');
nexttile;
plot(time, thetaHistory(1, :), time, thetaHistory(2, :), ...
     'LineWidth', 1.2);
grid on;
xlabel('Time [s]');
ylabel('STR estimates');
legend('A estimate', 'B estimate', 'Location', 'best');

results = table(time.', r.', yModel.', yMRAC.', ySTR.', ...
    [uMRAC, NaN].', [uSTR, NaN].', ...
    'VariableNames', {'Time', 'Reference', 'Model', 'MRAC', 'STR', ...
                      'ControlMRAC', 'ControlSTR'});
writetable(results, 'Chapter11_Lesson4_results.csv');

function value = clipValue(value, lowerBound, upperBound)
    value = min(max(value, lowerBound), upperBound);
end

14. Simulink Architecture Builder

The following MATLAB function programmatically creates and saves a Simulink model that exposes the structural difference between the two adaptive architectures. The detailed numerical laws remain in Chapter11_Lesson4.m; students can place those laws inside MATLAB Function blocks as a follow-on laboratory task.

`Chapter11_Lesson4_Simulink.m`

function Chapter11_Lesson4_Simulink
% Chapter 11, Lesson 4: programmatically create a Simulink architecture.
% The generated model emphasizes the structural distinction:
%   MRAC: tracking error directly updates controller gains.
%   STR: identifier estimates plant parameters, then a controller-design
%        block computes the controller gains.
%
% The numerical adaptive laws are implemented in Chapter11_Lesson4.m.

modelName = 'Chapter11_Lesson4_STR_vs_MRAC_Model';
if bdIsLoaded(modelName)
    close_system(modelName, 0);
end
if exist([modelName '.slx'], 'file')
    delete([modelName '.slx']);
end

new_system(modelName);
open_system(modelName);
set_param(modelName, 'StopTime', '24');

add_block('simulink/Sources/Step', [modelName '/Reference'], ...
    'Time', '1', 'Before', '0', 'After', '1', ...
    'Position', [30 125 60 155]);

add_block('simulink/Ports & Subsystems/Subsystem', ...
    [modelName '/Reference Model'], ...
    'Position', [120 55 270 115]);
add_block('simulink/Ports & Subsystems/Subsystem', ...
    [modelName '/Direct MRAC'], ...
    'Position', [120 155 270 225]);
add_block('simulink/Ports & Subsystems/Subsystem', ...
    [modelName '/Plant for MRAC'], ...
    'Position', [335 155 475 225]);

add_block('simulink/Ports & Subsystems/Subsystem', ...
    [modelName '/Online Identifier'], ...
    'Position', [120 300 270 365]);
add_block('simulink/Ports & Subsystems/Subsystem', ...
    [modelName '/Controller Design'], ...
    'Position', [335 300 475 365]);
add_block('simulink/Ports & Subsystems/Subsystem', ...
    [modelName '/Plant for STR'], ...
    'Position', [540 300 680 365]);

add_block('simulink/Signal Routing/Mux', [modelName '/Output Mux'], ...
    'Inputs', '3', 'Position', [745 125 750 335]);
add_block('simulink/Sinks/Scope', [modelName '/Comparison Scope'], ...
    'Position', [810 205 845 240]);

add_line(modelName, 'Reference/1', 'Reference Model/1', 'autorouting', 'on');
add_line(modelName, 'Reference/1', 'Direct MRAC/1', 'autorouting', 'on');
add_line(modelName, 'Direct MRAC/1', 'Plant for MRAC/1', 'autorouting', 'on');
add_line(modelName, 'Reference/1', 'Controller Design/1', 'autorouting', 'on');
add_line(modelName, 'Online Identifier/1', 'Controller Design/2', ...
    'autorouting', 'on');
add_line(modelName, 'Controller Design/1', 'Plant for STR/1', ...
    'autorouting', 'on');
add_line(modelName, 'Reference Model/1', 'Output Mux/1', 'autorouting', 'on');
add_line(modelName, 'Plant for MRAC/1', 'Output Mux/2', 'autorouting', 'on');
add_line(modelName, 'Plant for STR/1', 'Output Mux/3', 'autorouting', 'on');
add_line(modelName, 'Output Mux/1', 'Comparison Scope/1', 'autorouting', 'on');

set_param([modelName '/Reference Model'], 'AttributesFormatString', ...
    'Desired closed-loop dynamics');
set_param([modelName '/Direct MRAC'], 'AttributesFormatString', ...
    'e updates controller gains directly');
set_param([modelName '/Online Identifier'], 'AttributesFormatString', ...
    'prediction error updates plant estimates');
set_param([modelName '/Controller Design'], 'AttributesFormatString', ...
    'estimated plant maps to controller gains');

Simulink.Annotation(modelName, ...
    'MRAC path: reference model -> tracking error -> gain adaptation');
annotation = Simulink.Annotation(modelName, ...
    'STR path: measured data -> identifier -> controller redesign');
annotation.Position = [110 405 560 435];

save_system(modelName);
fprintf('Created %s.slx\n', modelName);
end

15. Wolfram Mathematica Implementation

The downloadable notebook contains the following Wolfram Language input cell and exports a CSV file after simulation.

`Chapter11_Lesson4.nb`

ClearAll["Global`*"];

ts = 0.02;
n = 1200;
modelRate = 1.8;
controlLimit = 8.0;
modelA = Exp[-modelRate ts];
modelB = 1.0 - modelA;

reference[k_Integer] := If[EvenQ[Quotient[k, 200]], 1.0, -0.5];
plantParameters[k_Integer] := If[k < n/2, {1.2, 1.0}, {0.7, 0.65}];
clip[value_, lower_, upper_] := Clip[value, {lower, upper}];

yModel = ConstantArray[0.0, n + 1];
yMRAC = ConstantArray[0.0, n + 1];
ySTR = ConstantArray[0.0, n + 1];
uMRAC = ConstantArray[0.0, n];
uSTR = ConstantArray[0.0, n];

(* Direct MRAC parameters for u = krHat r - kyHat y. *)
kyHat = 0.0;
krHat = 0.2;
gammaY = 1.4;
gammaR = 1.8;

(* STR estimate for y(k+1) = A y(k) + B u(k). *)
thetaHat = {0.94, 0.018};
identifierGain = 0.55;
normalization = 1.0*^-4;
previousPhi = {0.0, 0.0};

Do[
  r = reference[k];
  {aTrue, bTrue} = plantParameters[k];
  aTrueD = 1.0 - ts aTrue;
  bTrueD = ts bTrue;

  yModel[[k + 2]] = modelA yModel[[k + 1]] + modelB r;

  (* ---------------- Direct MRAC ---------------- *)
  trackingError = yMRAC[[k + 1]] - yModel[[k + 1]];
  kyHat = clip[
    kyHat + ts gammaY trackingError yMRAC[[k + 1]], -20.0, 20.0];
  krHat = clip[krHat - ts gammaR trackingError r, -20.0, 20.0];
  uMRAC[[k + 1]] = clip[
    krHat r - kyHat yMRAC[[k + 1]], -controlLimit, controlLimit];
  yMRAC[[k + 2]] = yMRAC[[k + 1]] +
    ts (-aTrue yMRAC[[k + 1]] + bTrue uMRAC[[k + 1]]);

  (* ---------------- Indirect STR ---------------- *)
  If[k > 0,
    predictionError = ySTR[[k + 1]] - thetaHat.previousPhi;
    thetaHat = thetaHat + identifierGain previousPhi predictionError/
      (normalization + previousPhi.previousPhi);
    thetaHat[[1]] = clip[thetaHat[[1]], 0.0, 1.20];
    thetaHat[[2]] = clip[thetaHat[[2]], 0.004, 0.060];
  ];

  probingSignal = 0.025 Sin[0.11 k];
  certaintyEquivalent =
    ((modelA - thetaHat[[1]]) ySTR[[k + 1]] + modelB r)/thetaHat[[2]];
  uSTR[[k + 1]] = clip[
    certaintyEquivalent + probingSignal, -controlLimit, controlLimit];
  ySTR[[k + 2]] = aTrueD ySTR[[k + 1]] + bTrueD uSTR[[k + 1]];
  previousPhi = {ySTR[[k + 1]], uSTR[[k + 1]]};
, {k, 0, n - 1}];

time = Range[0, n] ts;
rmseMRAC = Sqrt[Mean[(yMRAC - yModel)^2]];
rmseSTR = Sqrt[Mean[(ySTR - yModel)^2]];
Print["MRAC tracking RMSE: ", NumberForm[rmseMRAC, {8, 6}]];
Print["STR tracking RMSE:  ", NumberForm[rmseSTR, {8, 6}]];

outputPlot = ListLinePlot[
  {Transpose[{time, yModel}], Transpose[{time, yMRAC}],
   Transpose[{time, ySTR}]},
  PlotLegends -> {"Reference model", "MRAC output", "STR output"},
  Frame -> True,
  FrameLabel -> {"Time [s]", "Output"},
  PlotLabel -> "STR versus MRAC on a changing first-order plant",
  ImageSize -> Large
];

controlPlot = ListLinePlot[
  {Transpose[{Most[time], uMRAC}], Transpose[{Most[time], uSTR}]},
  PlotLegends -> {"MRAC control", "STR control"},
  Frame -> True,
  FrameLabel -> {"Time [s]", "Control input"},
  PlotLabel -> "Adaptive control effort",
  ImageSize -> Large
];

Export["Chapter11_Lesson4_results.csv",
  Prepend[
    Transpose[{time, Table[reference[k], {k, 0, n}], yModel, yMRAC, ySTR}],
    {"time", "reference", "reference_model", "mrac_output", "str_output"}
  ]
];

Column[{outputPlot, controlPlot}]

16. Interpreting the Simulation Correctly

Do not rank architectures from one RMSE value. The STR and MRAC tunings have different gains, excitation mechanisms, and transient objectives.
Inspect control effort. Small output error obtained by repeated saturation or large high-frequency control is not necessarily acceptable performance.
Inspect parameter histories. Tracking recovery after the parameter change can occur before estimates settle, and estimates may stop converging when the reference becomes insufficiently rich.
Use applied actuator input for identification. Replacing it with an unsaturated command creates a false regression model.
Test model mismatch. Add a neglected pole, delay, noise, bias, and disturbance before drawing implementation conclusions.

16.1 Relevant Software Ecosystems

Python: NumPy for arrays, SciPy for integration and signal processing, Matplotlib for visualization, and the control package for linear-system utilities.
C++: Eigen can simplify vector and matrix estimators; Boost.Odeint can support continuous-time simulation.
Java: EJML or Apache Commons Math can support higher-order matrix adaptive laws.
MATLAB/Simulink: Control System Toolbox supports model analysis; Simulink and Simulink Control Design support block-level experiments and adaptive-control workflows.
Wolfram Mathematica: NDSolve, RecurrenceTable, symbolic algebra, and optimization functions can support analysis and reproducible derivations.

17. Problems and Solutions

Problem 1 — MRAC Matching and Adaptive Law: For \(\dot y=-a y+b u\), reference model \(\dot y_m=-a_m y_m+b_m r\), and controller \(u=\hat k_r r-\hat k_y y\), derive the ideal gains and a Lyapunov adaptive law for positive \(b\).

Solution: Matching coefficients gives

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

The error equation is

\[ \dot e=-a_m e-b\tilde k_y y+b\tilde k_r r. \]

With

\[ V=\frac{1}{2}e^2+\frac{b}{2\gamma_y}\tilde k_y^2+ \frac{b}{2\gamma_r}\tilde k_r^2, \]

select

\[ \dot{\hat k}_y=\gamma_y e y, \qquad \dot{\hat k}_r=-\gamma_r e r, \]

which yields \(\dot V=-a_m e^2\). The result establishes bounded adaptive signals and asymptotic tracking under the standard boundedness conditions, but not parameter convergence without excitation.

Problem 2 — STR Gain-Singularity Sensitivity: For \(u=((A_m-\hat A)y+B_m r)/\hat B\), compare the control sensitivity at \(\hat B=0.1\) and \(\hat B=0.01\), assuming the numerator equals \(0.2\).

Solution: The derivative with respect to the estimated input gain is

\[ \frac{\partial u}{\partial \hat B}=-\frac{0.2}{\hat B^2}. \]

Therefore,

\[ \left.\frac{\partial u}{\partial \hat B}\right|_{\hat B=0.1}=-20, \qquad \left.\frac{\partial u}{\partial \hat B}\right|_{\hat B=0.01}=-2000. \]

Reducing \(\hat B\) by a factor of ten increases the sensitivity by a factor of one hundred. A lower projection bound on \(|\hat B|\) is therefore essential.

Problem 3 — One Normalized-Gradient STR Update: Let \(\hat{\boldsymbol{\theta}}_k=[0.90,\;0.10]^T\), \(\boldsymbol{\phi}_k=[1.0,\;0.5]^T\), measured \(y_{k+1}=1.02\), \(\gamma=0.5\), and \(\delta=0.01\). Compute the unprojected update.

Solution: First compute

\[ \hat y_{k+1|k}=0.90(1.0)+0.10(0.5)=0.95, \qquad \varepsilon_{k+1}=1.02-0.95=0.07. \]

The denominator is

\[ \delta+\boldsymbol{\phi}_k^T\boldsymbol{\phi}_k =0.01+1^2+0.5^2=1.26. \]

Hence

\[ \hat{\boldsymbol{\theta}}_{k+1}= \begin{bmatrix}0.90\\0.10\end{bmatrix}+ \frac{0.5(0.07)}{1.26} \begin{bmatrix}1.0\\0.5\end{bmatrix} \approx \begin{bmatrix}0.92778\\0.11389\end{bmatrix}. \]

Problem 4 — Loss of Excitation: An adaptive loop reaches nearly perfect constant-reference regulation, so that \(y_k\) and \(u_k\) become almost constant. Explain the implications for STR and MRAC.

Solution:

In STR, the regressor vectors become nearly collinear. The information matrix loses rank, so several combinations of \(A\) and \(B\) can explain the same steady data. Plant estimates need not converge uniquely.
In MRAC, the tracking error can approach zero and adaptation can stop, even though the controller gains are not equal to their ideal values.
The control objective may still be met. Parameter convergence is a stronger requirement than output tracking and requires richer signals.

Problem 5 — Architecture Selection: A thermal process has slowly changing gain and time constant. Operators require online model estimates for diagnostics, and a certified PI tuning rule is valid over a bounded parameter region. A projection operator and a fixed fallback PI controller are available. Which architecture is more natural, and what checks remain necessary?

Solution:

An indirect STR is more natural because the estimated model has explicit diagnostic value and can feed the certified PI redesign map. The design must still verify excitation, estimator bias under disturbances, safe behavior throughout the projected parameter region, redesign rate limits, actuator saturation, bumpless gain updates, and fallback activation. The existence of a stable PI controller for each frozen estimate is not alone a proof of stability for the complete time-varying adaptive loop.

18. Summary

MRAC and STR solve related but structurally different adaptive-control problems. Direct MRAC adjusts controller parameters from tracking error and can admit a compact Lyapunov proof when matching assumptions hold. An indirect STR estimates a plant model and repeatedly applies a controller synthesis rule, offering modularity and interpretable plant estimates at the cost of an additional identification-redesign coupling.

MRAC is not automatically robust merely because a Lyapunov function exists, and STR is not automatically stable merely because every frozen estimated model yields stable poles. Both require careful treatment of excitation, projection, normalization, noise, saturation, unmodeled dynamics, and implementation rate. The appropriate choice follows from the uncertainty structure and verification requirements, not from a universal ranking.

19. References

Åström, K.J., & Wittenmark, B. (1973). On self-tuning regulators. Automatica, 9(2), 185–199.
Åström, K.J., Borisson, U., Ljung, L., & Wittenmark, B. (1977). Theory and applications of self-tuning regulators. Automatica, 13(5), 457–476.
Parks, P.C. (1966). Liapunov redesign of model reference adaptive control systems. IEEE Transactions on Automatic Control, 11(3), 362–367.
Monopoli, R.V. (1974). Model reference adaptive control with an augmented error signal. IEEE Transactions on Automatic Control, 19(5), 474–484.
Narendra, K.S., & Valavani, L.S. (1978). Stable adaptive controller design—direct control. IEEE Transactions on Automatic Control, 23(4), 570–583.
Morse, A.S. (1980). Global stability of parameter-adaptive control systems. IEEE Transactions on Automatic Control, 25(3), 433–439.
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.
Middleton, R.H., Goodwin, G.C., Hill, D.J., & Mayne, D.Q. (1988). Design issues in adaptive control. IEEE Transactions on Automatic Control, 33(1), 50–58.

Chapter 11: Indirect Adaptive Control and Self-Tuning Regulators (STR)

Lesson 4: STR vs MRAC – Advantages and Drawbacks

1. Learning Objectives and Scope

2. A Common Control Objective

\[ \dot y_m(t)=-a_m y_m(t)+b_m r(t), \qquad a_m>0. \]

3. Architectural Distinction

4. Direct MRAC: Derivation and Lyapunov Argument

\[ u(t)=\hat k_r(t)r(t)-\hat k_y(t)y(t). \]

\[ a+b k_y^*=a_m, \qquad b k_r^*=b_m, \]

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

\[ \dot e=-a_m e-b\tilde k_y y+b\tilde k_r r. \]

\[ V=\frac{1}{2}e^2+\frac{b}{2\gamma_y}\tilde k_y^2+ \frac{b}{2\gamma_r}\tilde k_r^2, \qquad \gamma_y>0,\;\gamma_r>0. \]

\[ \dot V=e\left(-a_m e-b\tilde k_y y+b\tilde k_r r\right) +\frac{b}{\gamma_y}\tilde k_y\dot{\hat k}_y +\frac{b}{\gamma_r}\tilde k_r\dot{\hat k}_r. \]

\[ \dot{\hat k}_y=\gamma_y e y, \qquad \dot{\hat k}_r=-\gamma_r e r. \]

\[ \dot V=-a_m e^2\le 0. \]

5. Indirect STR: Identification Followed by Controller Redesign

\[ y_{k+1}=A y_k+B u_k, \qquad \boldsymbol{\phi}_k=\begin{bmatrix}y_k & u_k\end{bmatrix}^T, \qquad \boldsymbol{\theta}=\begin{bmatrix}A & B\end{bmatrix}^T. \]

\[ \hat y_{k+1|k}=\hat{\boldsymbol{\theta}}_k^T\boldsymbol{\phi}_k, \qquad \varepsilon_{k+1}=y_{k+1}- \hat{\boldsymbol{\theta}}_k^T\boldsymbol{\phi}_k. \]

\[ \hat{\boldsymbol{\theta}}_{k+1}= \operatorname{Proj}_\Theta\!\left[ \hat{\boldsymbol{\theta}}_k+ \frac{\gamma\boldsymbol{\phi}_k\varepsilon_{k+1}} {\delta+\boldsymbol{\phi}_k^T\boldsymbol{\phi}_k} \right], \qquad \gamma>0,\;\delta>0. \]

\[ u_k=\frac{(A_m-\hat A_k)y_k+B_m r_k}{\hat B_k}. \]

5.1 Error Decomposition

\[ y_{k+1}=A_m y_k+B_m r_k+ \left(A-\hat A_k\right)y_k+ \left(B-\hat B_k\right)u_k. \]

5.2 Controller-Map Sensitivity

\[ \frac{\partial u_k}{\partial \hat A_k}=-\frac{y_k}{\hat B_k}, \qquad \frac{\partial u_k}{\partial \hat B_k}=- \frac{(A_m-\hat A_k)y_k+B_m r_k}{\hat B_k^2}. \]

6. Advantages and Drawbacks: Structured Comparison

7. Stability, Convergence, and Robustness Implications

7.1 Frozen-Time Stability Is Not Enough

\[ x_{k+1}=A_{cl}\!\left(\hat{\boldsymbol{\theta}}_k\right)x_k, \qquad \rho\!\left(A_{cl}(\vartheta)\right)<1 \text{ for each fixed }\vartheta \not\Rightarrow \text{uniform stability for arbitrary }\hat{\boldsymbol{\theta}}_k. \]

7.2 Tracking Convergence Is Not Parameter Convergence

\[ \sum_{i=k}^{k+N-1}\boldsymbol{\phi}_i \boldsymbol{\phi}_i^T\succeq \alpha I, \qquad N\in\mathbb{N},\;\alpha>0. \]

7.3 Unmodeled Dynamics and Disturbances

\[ \begin{aligned} \text{measured mismatch} &= \text{parametric mismatch} + \text{unmodeled dynamics} \\ &\quad + \text{disturbance} + \text{measurement noise}. \end{aligned} \]

7.4 Saturation Changes the Mathematics

\[ \dot e=-a_m e-b\tilde k_y y+b\tilde k_r r+ b\left[\operatorname{sat}(u_c)-u_c\right]. \]

8. Selection Logic: When Each Architecture Is More Natural

9. Comparative Numerical Experiment

\[ (a,b)=\begin{cases} (1.2,1.0), & k<600,\\ (0.7,0.65), & k\ge 600. \end{cases} \]

\[ A_m=e^{-a_mT_s}, \qquad B_m=1-A_m, \qquad a_m=1.8. \]

10. Python Implementation

Chapter11_Lesson4.py

11. C++ Implementation

Chapter11_Lesson4.cpp

12. Java Implementation

Chapter11_Lesson4.java

13. MATLAB Implementation

Chapter11_Lesson4.m

14. Simulink Architecture Builder

Chapter11_Lesson4_Simulink.m

15. Wolfram Mathematica Implementation

Chapter11_Lesson4.nb

16. Interpreting the Simulation Correctly

16.1 Relevant Software Ecosystems

17. Problems and Solutions

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

\[ \dot e=-a_m e-b\tilde k_y y+b\tilde k_r r. \]

\[ V=\frac{1}{2}e^2+\frac{b}{2\gamma_y}\tilde k_y^2+ \frac{b}{2\gamma_r}\tilde k_r^2, \]

\[ \dot{\hat k}_y=\gamma_y e y, \qquad \dot{\hat k}_r=-\gamma_r e r, \]

\[ \frac{\partial u}{\partial \hat B}=-\frac{0.2}{\hat B^2}. \]

\[ \left.\frac{\partial u}{\partial \hat B}\right|_{\hat B=0.1}=-20, \qquad \left.\frac{\partial u}{\partial \hat B}\right|_{\hat B=0.01}=-2000. \]

\[ \hat y_{k+1|k}=0.90(1.0)+0.10(0.5)=0.95, \qquad \varepsilon_{k+1}=1.02-0.95=0.07. \]

\[ \delta+\boldsymbol{\phi}_k^T\boldsymbol{\phi}_k =0.01+1^2+0.5^2=1.26. \]

\[ \hat{\boldsymbol{\theta}}_{k+1}= \begin{bmatrix}0.90\\0.10\end{bmatrix}+ \frac{0.5(0.07)}{1.26} \begin{bmatrix}1.0\\0.5\end{bmatrix} \approx \begin{bmatrix}0.92778\\0.11389\end{bmatrix}. \]

18. Summary

19. References

Previous

Next

Help keep these engineering tutorials free and growing

Similar Tutorials

\[ a+b k_y^=a_m, \qquad b k_r^=b_m, \]

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

`Chapter11_Lesson4.py`

`Chapter11_Lesson4.cpp`

`Chapter11_Lesson4.java`

`Chapter11_Lesson4.m`

`Chapter11_Lesson4_Simulink.m`

`Chapter11_Lesson4.nb`

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