In the previous lessons, persistent excitation was introduced as a sufficient condition for parameter convergence, and the distinction between state convergence and parameter convergence was established. The central practical fact is now stated more sharply: a controller can achieve small tracking error without uniquely recovering the physical parameter vector.

Consider a linear regression embedded in an adaptive loop:

\[ y(t)=\boldsymbol{\phi}^{T}(t)\boldsymbol{\theta}^{\ast}+d(t), \qquad \hat{y}(t)=\boldsymbol{\phi}^{T}(t)\hat{\boldsymbol{\theta} }(t) \]

where \( \boldsymbol{\theta}^{\ast}\in\mathbb{R}^{p} \) is the unknown constant parameter vector, \( \boldsymbol{\phi}(t)\in\mathbb{R}^{p} \) is the measured regressor, and \( d(t) \) collects noise, disturbances, and unmodeled effects. Define

\[ \tilde{\boldsymbol{\theta} }(t)= \hat{\boldsymbol{\theta} }(t)-\boldsymbol{\theta}^{\ast}, \qquad \varepsilon(t)=\hat{y}(t)-y(t) =\boldsymbol{\phi}^{T}(t)\tilde{\boldsymbol{\theta} }(t)-d(t). \]

If \( d(t)=0 \), then \( \varepsilon(t)=0 \) only requires \( \boldsymbol{\phi}^{T}(t)\tilde{\boldsymbol{\theta} }(t)=0 \). This scalar condition does not generally imply \( \tilde{\boldsymbol{\theta} }(t)=\mathbf{0} \) when several unknown parameters are present.

flowchart TD
  A["Adaptive closed loop"] --> B["Tracking error becomes small"]
  B --> C["Regulation suppresses state and input variation"]
  C --> D["Regressor loses independent directions"]
  D --> E["Information matrix becomes singular or ill-conditioned"]
  E --> F["Many parameter vectors explain nearly the same data"]
  F --> G["Good control, incomplete identification"]
        

A normalized gradient estimator, consistent with the normalization ideas of Chapter 9, is

\[ \dot{\hat{\boldsymbol{\theta} } }(t) =-\frac{\boldsymbol{\Gamma}\boldsymbol{\phi}(t)} {1+\boldsymbol{\phi}^{T}(t)\boldsymbol{\phi}(t)}\, \varepsilon(t), \qquad \boldsymbol{\Gamma}=\boldsymbol{\Gamma}^{T}\gt 0. \]

In the ideal noise-free case, substitution of \( \varepsilon=\boldsymbol{\phi}^{T}\tilde{\boldsymbol{\theta} } \) gives

\[ \dot{\tilde{\boldsymbol{\theta} } }(t) =-\frac{\boldsymbol{\Gamma}\boldsymbol{\phi}(t) \boldsymbol{\phi}^{T}(t)} {1+\boldsymbol{\phi}^{T}(t)\boldsymbol{\phi}(t)} \tilde{\boldsymbol{\theta} }(t). \]

With the Lyapunov function \( V=\tfrac{1}{2}\tilde{\boldsymbol{\theta} }^{T} \boldsymbol{\Gamma}^{-1}\tilde{\boldsymbol{\theta} } \),

\[ \dot{V}(t)= -\frac{\left(\boldsymbol{\phi}^{T}(t) \tilde{\boldsymbol{\theta} }(t)\right)^{2} } {1+\boldsymbol{\phi}^{T}(t)\boldsymbol{\phi}(t)}\le 0. \]

Therefore, the parameter error is bounded and the prediction error can approach zero. However, this derivative is only negative semidefinite. Every parameter-error vector orthogonal to the currently active regressor lies in an instantaneous zero-decrease direction.

Uniform persistent excitation requires constants \( T\gt 0 \) and \( \alpha\gt 0 \) such that

\[ \int_{t}^{t+T}\boldsymbol{\phi}(\tau) \boldsymbol{\phi}^{T}(\tau)\,d\tau \ge \alpha\mathbf{I}_{p}, \qquad \text{for every } t\ge 0. \]

This condition demands repeated information in every direction of the \( p \)-dimensional parameter space. Boundedness of \( \boldsymbol{\phi} \) alone is insufficient; large signals may still be collinear and therefore uninformative about some parameter combinations.

Over a finite interval \( [t_0,t_1] \), define the information or excitation Gramian

\[ \mathbf{G}(t_0,t_1)= \int_{t_0}^{t_1}\boldsymbol{\phi}(\tau) \boldsymbol{\phi}^{T}(\tau)\,d\tau. \]

For any fixed direction \( \mathbf{v}\in\mathbb{R}^{p} \),

\[ \mathbf{v}^{T}\mathbf{G}\mathbf{v} =\int_{t_0}^{t_1}\left(\boldsymbol{\phi}^{T}(\tau) \mathbf{v}\right)^{2}d\tau. \]

Hence \( \mathbf{v}^{T}\mathbf{G}\mathbf{v}=0 \) if and only if the data are completely insensitive to parameter changes in the direction \( \mathbf{v} \) almost everywhere on the interval. The null space \( \mathcal{N}(\mathbf{G}) \) is precisely the set of unidentifiable directions for that data record.

\[ \operatorname{rank}(\mathbf{G})=r\lt p \quad\Longrightarrow\quad p-r \text{ independent parameter directions are not identified}. \]

flowchart TD
  T["True parameter vector"] --> P["Project onto excited subspace"]
  P --> O["Observable parameter combinations"]
  T --> N["Project onto Gramian null space"]
  N --> U["Unobservable parameter combinations"]
  O --> Y["Influence measured output"]
  U --> S["Can change without changing measured data"]
        

If \( \lambda_{1}\ge\cdots\ge\lambda_{p}\ge 0 \) are the eigenvalues of \( \mathbf{G} \), then \( \lambda_{p} \) measures the weakest excited direction. A useful finite-data diagnostic is

\[ \kappa(\mathbf{G})= \frac{\lambda_{\max}(\mathbf{G})} {\lambda_{\min}(\mathbf{G})}. \]

A very large condition number indicates practical non-identifiability: the model is theoretically full rank over the chosen record, but one or more directions are so weakly excited that noise or modeling error can dominate their estimates.

Suppose the noise-free data satisfy \( y(t)=\boldsymbol{\phi}^{T}(t)\boldsymbol{\theta}^{\ast} \) on \( [t_0,t_1] \), and let \( \mathbf{v}\in\mathcal{N}(\mathbf{G}) \) with \( \mathbf{v}\ne\mathbf{0} \). Then

\[ 0=\mathbf{v}^{T}\mathbf{G}\mathbf{v} =\int_{t_0}^{t_1} \left(\boldsymbol{\phi}^{T}(\tau)\mathbf{v}\right)^{2}d\tau. \]

Because the integrand is nonnegative, it follows that \( \boldsymbol{\phi}^{T}(t)\mathbf{v}=0 \) almost everywhere. For any scalar \( c \), define \( \boldsymbol{\theta}_{c}=\boldsymbol{\theta}^{\ast}+c\mathbf{v} \). Then

\[ \boldsymbol{\phi}^{T}(t)\boldsymbol{\theta}_{c} =\boldsymbol{\phi}^{T}(t)\boldsymbol{\theta}^{\ast} +c\boldsymbol{\phi}^{T}(t)\mathbf{v} =y(t). \]

Thus infinitely many distinct parameter vectors generate exactly the same measured trajectory. No estimator can select the true vector from these data alone without additional excitation, prior information, constraints, or a different model structure.

Let \( \boldsymbol{\phi}(t)=[s(t),\;2s(t)]^{T} \). Then

\[ y(t)=s(t)\theta_{1}^{\ast}+2s(t)\theta_{2}^{\ast} =s(t)\left(\theta_{1}^{\ast}+2\theta_{2}^{\ast}\right). \]

Only the scalar combination \( \eta^{\ast}=\theta_{1}^{\ast}+2\theta_{2}^{\ast} \) affects the data. The finite-window Gramian is

\[ \mathbf{G}= \left(\int_{t_0}^{t_1}s^{2}(\tau)d\tau\right) \begin{bmatrix}1&2\\2&4\end{bmatrix}. \]

Its eigenvalues are \( 5\int s^{2}d\tau \) and \( 0 \). A null-space vector is \( [-2,\;1]^{T} \), because

\[ [1,\;2]\begin{bmatrix}-2\\1\end{bmatrix}=0. \]

Therefore, all vectors satisfying \( \theta_{1}+2\theta_{2}=\eta^{\ast} \) are indistinguishable. An adaptive law may converge to one point on this line, but that point depends on initialization, gain geometry, normalization, leakage, projection bounds, noise, and transient data—not solely on the physical parameters.

Persistent excitation is especially difficult in feedback systems because the controller is designed to suppress error and variation. When the reference becomes constant and disturbances are rejected, state and input trajectories may approach steady values. Regressors constructed from those signals then lose richness.

For a regressor generated by a stable closed-loop system, suppose \( \boldsymbol{\phi}(t)→\boldsymbol{\phi}_{\infty} \). Over a window of length \( T \),

\[ \int_{t}^{t+T}\boldsymbol{\phi}(\tau) \boldsymbol{\phi}^{T}(\tau)d\tau → T\boldsymbol{\phi}_{\infty} \boldsymbol{\phi}_{\infty}^{T}. \]

The limiting matrix has rank at most one, regardless of the number of unknown parameters. Thus a high-quality regulator can naturally remove the very signal variation required for complete identification.

This creates a dual-purpose conflict. Control seeks to minimize a performance measure such as

\[ J_{c}=\int_{0}^{T_f} \left(e^{T}(t)\mathbf{Q}e(t)+u^{T}(t)\mathbf{R}u(t)\right)dt, \]

whereas identification benefits from increasing a measure such as \( \lambda_{\min}(\mathbf{G}) \) or \( \log\det(\mathbf{G}+\delta\mathbf{I}) \). Aggressive probing may improve information while worsening tracking, actuator usage, energy consumption, comfort, wear, or safety margins.

1. Insufficient reference richness. Constant commands, slowly varying trajectories, or a single sinusoid may excite fewer independent modes than the parameter dimension requires.

2. Correlated regressors. Different regressor components may be nearly proportional in closed loop. The Gramian is then full rank numerically but badly conditioned.

3. Overparameterization. A model may contain more adjustable coefficients than the input-output behavior can uniquely determine. Redundant basis functions create equivalent parameterizations.

4. Structural mismatch. Unmodeled dynamics, delays, nonlinearities, hysteresis, or disturbances can be absorbed into parameter estimates. The estimate may become an effective coefficient rather than a physical parameter.

5. Measurement noise. Weakly excited directions have low signal-to-noise ratio. Increasing adaptation gain may accelerate apparent learning but also amplifies parameter fluctuations.

6. Actuator saturation and safety constraints. Inputs that would generate adequate excitation may be prohibited. Saturation also changes the assumed regression relation.

7. Finite operating time. Persistent excitation is an infinite-time uniform condition. A finite experiment may contain useful excitation only briefly or may not repeat it sufficiently.

8. Time-varying parameters. If the true parameter drifts, exact convergence to a fixed vector is not a meaningful objective. The estimator must trade tracking of variation against noise rejection.

9. Limited sensing. Regressors formed from unmeasured or poorly reconstructed states may not preserve the excitation present in the physical plant.

10. Robust modifications. Leakage, dead zones, and projection improve boundedness and robustness, but they generally bias the steady parameter estimate or stop adaptation when errors become small.

With additive disturbance \( d(t) \), the normalized gradient error dynamics become

\[ \dot{\tilde{\boldsymbol{\theta} } } =-\frac{\boldsymbol{\Gamma}\boldsymbol{\phi} \boldsymbol{\phi}^{T} } {1+\boldsymbol{\phi}^{T}\boldsymbol{\phi} } \tilde{\boldsymbol{\theta} } +\frac{\boldsymbol{\Gamma}\boldsymbol{\phi} } {1+\boldsymbol{\phi}^{T}\boldsymbol{\phi} }d. \]

For the same Lyapunov function,

\[ \dot{V}=-\frac{z^{2} }{1+\|\boldsymbol{\phi}\|^{2} } +\frac{z d}{1+\|\boldsymbol{\phi}\|^{2} }, \qquad z=\boldsymbol{\phi}^{T}\tilde{\boldsymbol{\theta} }. \]

Applying \( ab\le\tfrac{1}{2}a^{2}+\tfrac{1}{2}b^{2} \) yields

\[ \dot{V}\le -\frac{1}{2}\frac{z^{2} }{1+\|\boldsymbol{\phi}\|^{2} } +\frac{1}{2}\frac{d^{2} }{1+\|\boldsymbol{\phi}\|^{2} }. \]

The disturbance injects energy into the estimator. Even under good excitation, the parameter error generally approaches a neighborhood whose size depends on disturbance magnitude, excitation strength, adaptation gain, and model mismatch. In weak directions, the neighborhood can be large because the restoring information is small.

A local least-squares interpretation gives the same warning. If \( \mathbf{G} \) is nonsingular and measurement noise has variance \( \sigma^{2} \), the estimate covariance scales approximately as

\[ \operatorname{Cov}(\hat{\boldsymbol{\theta} }) \approx \sigma^{2}\mathbf{G}^{-1}. \]

Thus the variance along an eigenvector associated with \( \lambda_i \) is proportional to \( 1/\lambda_i \). A small positive eigenvalue can be almost as problematic as an exact zero in practice.

Parameter convergence should not be inferred from tracking error alone. A practical implementation should separately monitor:

Prediction residual: \( \varepsilon(t)=\hat{y}(t)-y(t) \).

Parameter increment: \( \|\dot{\hat{\boldsymbol{\theta} } }(t)\| \) or its sampled equivalent.

Windowed information matrix:

\[ \mathbf{G}_{T}(t)= \int_{t-T}^{t}\boldsymbol{\phi}(\tau) \boldsymbol{\phi}^{T}(\tau)d\tau. \]

Weakest excitation level: \( \lambda_{\min}(\mathbf{G}_{T}) \).

Condition number: \( \kappa(\mathbf{G}_{T}) \).

Constraint activity: projection boundaries, actuator saturation, dead zones, and adaptation freezes.

Model validity indicators: residual bias, residual correlation, unexplained frequency content, and systematic dependence on operating point.

A small prediction residual together with a small \( \lambda_{\min}(\mathbf{G}_{T}) \) means the data fit is good but full identification is not supported. Parameter estimates should then be reported as nonunique or weakly identifiable rather than as verified physical constants.

Use task-relevant parameter combinations. If the control law depends only on a lower-dimensional combination, estimating that combination may be sufficient and better conditioned.

Reduce the parameterization. Remove redundant basis functions and fix parameters known accurately from design data or offline tests.

Design safe probing signals. Small multisine, chirp, or bounded dither signals may improve information, but their amplitudes and frequencies must respect actuator, state, and performance constraints.

Use multiple operating regimes. Data collected under different references, loads, or equilibrium points can excite complementary parameter directions.

Separate control and calibration phases. A commissioning experiment can prioritize identification, after which online adaptation is used mainly for bounded compensation.

Retain uncertainty information. Report confidence regions, excitation eigenvalues, and sensitivity rather than only one parameter vector.

Prefer performance guarantees when identification is unnecessary. In many adaptive-control problems, boundedness and tracking are the actual objectives. Exact physical identification should not be imposed as an unstated requirement.

More advanced methods can relax classical persistent-excitation demands by exploiting stored data or alternative regressions. Those approaches belong to later topics; they do not remove the fundamental need for enough independent information somewhere in the experiment.

The implementations below use the same true parameter vector \( \boldsymbol{\theta}^{\ast}=[1.5,-0.8]^{T} \) and the same normalized gradient estimator. Two regressors are compared:

\[ \boldsymbol{\phi}_{\mathrm{nonPE} }(t)= \begin{bmatrix}\sin t\\2\sin t\end{bmatrix}, \qquad \boldsymbol{\phi}_{\mathrm{PE} }(t)= \begin{bmatrix}\sin t\\\cos(\sqrt{2}t)\end{bmatrix}. \]

The first Gramian has rank one. The second becomes well conditioned over a sufficiently long interval because its components are not permanently collinear. Each program reports final parameter error, prediction-error RMS, Gramian eigenvalues, and condition number.

Chapter10_Lesson5.py

Libraries: NumPy for numerical linear algebra and Matplotlib for parameter-history plots.

# Chapter10_Lesson5.py
# Demonstrates why good tracking/prediction does not imply full parameter identification.
# Requirements: numpy, matplotlib

import numpy as np
import matplotlib.pyplot as plt

DT = 0.002
T_END = 60.0
GAMMA = 1.2
THETA_TRUE = np.array([1.5, -0.8], dtype=float)


def regressor(t: float, pe_case: bool) -> np.ndarray:
    """Return a rank-deficient or persistently exciting two-parameter regressor."""
    if pe_case:
        return np.array([np.sin(t), np.cos(np.sqrt(2.0) * t)])
    s = np.sin(t)
    return np.array([s, 2.0 * s])


def symmetric_eigenvalues_2x2(matrix: np.ndarray) -> tuple[float, float]:
    values = np.linalg.eigvalsh(matrix)
    return float(values[0]), float(values[1])


def simulate(pe_case: bool) -> dict[str, np.ndarray | float]:
    time = np.arange(0.0, T_END + DT, DT)
    theta_hat = np.zeros(2)
    gramian = np.zeros((2, 2))

    theta_history = np.zeros((time.size, 2))
    prediction_error = np.zeros(time.size)

    for k, t in enumerate(time):
        phi = regressor(t, pe_case)
        deterministic_noise = 0.02 * np.sin(13.7 * t) + 0.015 * np.cos(7.3 * t)
        y = float(phi @ THETA_TRUE + deterministic_noise)
        error = y - float(phi @ theta_hat)

        # Normalized gradient update (normalization was introduced in Chapter 9).
        normalization = 1.0 + float(phi @ phi)
        theta_dot = GAMMA * phi * error / normalization
        theta_hat = theta_hat + DT * theta_dot

        gramian = gramian + DT * np.outer(phi, phi)
        theta_history[k] = theta_hat
        prediction_error[k] = error

    lambda_min, lambda_max = symmetric_eigenvalues_2x2(gramian)
    condition_number = np.inf if lambda_min <= 1e-12 else lambda_max / lambda_min
    parameter_error = float(np.linalg.norm(theta_hat - THETA_TRUE))
    rms_prediction_error = float(np.sqrt(np.mean(prediction_error**2)))

    return {
        "time": time,
        "theta_history": theta_history,
        "prediction_error": prediction_error,
        "gramian": gramian,
        "lambda_min": lambda_min,
        "lambda_max": lambda_max,
        "condition_number": condition_number,
        "parameter_error": parameter_error,
        "rms_prediction_error": rms_prediction_error,
    }


def report(name: str, result: dict[str, np.ndarray | float]) -> None:
    theta_final = np.asarray(result["theta_history"])[-1]
    print(f"\n{name}")
    print(f"  final theta_hat        = {theta_final}")
    print(f"  true theta             = {THETA_TRUE}")
    print(f"  ||theta_hat-theta||    = {result['parameter_error']:.6f}")
    print(f"  RMS prediction error   = {result['rms_prediction_error']:.6f}")
    print(f"  Gramian eigenvalues    = ({result['lambda_min']:.6e}, {result['lambda_max']:.6e})")
    print(f"  Gramian condition no.  = {result['condition_number']}")


def main() -> None:
    non_pe = simulate(pe_case=False)
    pe = simulate(pe_case=True)

    report("Rank-deficient regressor: phi=[sin(t), 2 sin(t)]", non_pe)
    report("Full-rank exciting regressor: phi=[sin(t), cos(sqrt(2)t)]", pe)

    time = np.asarray(non_pe["time"])
    non_pe_theta = np.asarray(non_pe["theta_history"])
    pe_theta = np.asarray(pe["theta_history"])

    plt.figure(figsize=(10, 5))
    plt.plot(time, non_pe_theta[:, 0], label="theta1 estimate (non-PE)")
    plt.plot(time, non_pe_theta[:, 1], label="theta2 estimate (non-PE)")
    plt.axhline(THETA_TRUE[0], linestyle="--", label="true theta1")
    plt.axhline(THETA_TRUE[1], linestyle=":", label="true theta2")
    plt.xlabel("Time [s]")
    plt.ylabel("Parameter value")
    plt.title("Rank-deficient data: output fit can improve while parameters remain wrong")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()

    plt.figure(figsize=(10, 5))
    plt.plot(time, pe_theta[:, 0], label="theta1 estimate (PE)")
    plt.plot(time, pe_theta[:, 1], label="theta2 estimate (PE)")
    plt.axhline(THETA_TRUE[0], linestyle="--", label="true theta1")
    plt.axhline(THETA_TRUE[1], linestyle=":", label="true theta2")
    plt.xlabel("Time [s]")
    plt.ylabel("Parameter value")
    plt.title("Full-rank excitation: both parameter directions are learned")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.show()


if __name__ == "__main__":
    main()

Chapter10_Lesson5.cpp

Implementation: C++17 standard library, explicit two-by-two Gramian and eigenvalue calculation, and CSV output for plotting in another tool.

// Chapter10_Lesson5.cpp
// C++17 implementation of the rank-deficient versus full-rank excitation example.

#include <array>
#include <cmath>
#include <fstream>
#include <iomanip>
#include <iostream>
#include <limits>
#include <stdexcept>
#include <string>

struct Result {
    std::array<double, 2> theta_hat{};
    std::array<std::array<double, 2>, 2> gramian{};
    double lambda_min{};
    double lambda_max{};
    double condition_number{};
    double parameter_error{};
    double rms_prediction_error{};
};

constexpr double kDt = 0.002;
constexpr double kEndTime = 60.0;
constexpr double kGamma = 1.2;
constexpr std::array<double, 2> kThetaTrue{1.5, -0.8};
constexpr double kSqrt2 = 1.4142135623730950488;

std::array<double, 2> regressor(double t, bool pe_case) {
    if (pe_case) {
        return {std::sin(t), std::cos(kSqrt2 * t)};
    }
    const double s = std::sin(t);
    return {s, 2.0 * s};
}

std::array<double, 2> eigenvaluesSymmetric2x2(
    const std::array<std::array<double, 2>, 2>& matrix) {
    const double a = matrix[0][0];
    const double b = matrix[0][1];
    const double d = matrix[1][1];
    const double trace = a + d;
    const double discriminant = std::sqrt((a - d) * (a - d) + 4.0 * b * b);
    return {0.5 * (trace - discriminant), 0.5 * (trace + discriminant)};
}

Result simulate(bool pe_case, const std::string& csv_name) {
    Result result;
    std::array<double, 2> theta_hat{0.0, 0.0};
    std::array<std::array<double, 2>, 2> gramian{ { {0.0, 0.0}, {0.0, 0.0} } };

    std::ofstream csv(csv_name);
    if (!csv) {
        throw std::runtime_error("Cannot open output CSV: " + csv_name);
    }
    csv << "time,theta1_hat,theta2_hat,prediction_error\n";

    double squared_error_sum = 0.0;
    std::size_t samples = 0;
    const int steps = static_cast<int>(std::round(kEndTime / kDt));

    for (int k = 0; k <= steps; ++k) {
        const double t = k * kDt;
        const auto phi = regressor(t, pe_case);
        const double noise = 0.02 * std::sin(13.7 * t) + 0.015 * std::cos(7.3 * t);
        const double y = phi[0] * kThetaTrue[0] + phi[1] * kThetaTrue[1] + noise;
        const double y_hat = phi[0] * theta_hat[0] + phi[1] * theta_hat[1];
        const double error = y - y_hat;
        const double normalization = 1.0 + phi[0] * phi[0] + phi[1] * phi[1];

        theta_hat[0] += kDt * kGamma * phi[0] * error / normalization;
        theta_hat[1] += kDt * kGamma * phi[1] * error / normalization;

        gramian[0][0] += kDt * phi[0] * phi[0];
        gramian[0][1] += kDt * phi[0] * phi[1];
        gramian[1][0] += kDt * phi[1] * phi[0];
        gramian[1][1] += kDt * phi[1] * phi[1];

        squared_error_sum += error * error;
        ++samples;

        if (k % 50 == 0) {
            csv << std::setprecision(12) << t << ',' << theta_hat[0] << ','
                << theta_hat[1] << ',' << error << '\n';
        }
    }

    const auto eigenvalues = eigenvaluesSymmetric2x2(gramian);
    const double e0 = theta_hat[0] - kThetaTrue[0];
    const double e1 = theta_hat[1] - kThetaTrue[1];

    result.theta_hat = theta_hat;
    result.gramian = gramian;
    result.lambda_min = eigenvalues[0];
    result.lambda_max = eigenvalues[1];
    result.condition_number = eigenvalues[0] <= 1e-12
        ? std::numeric_limits<double>::infinity()
        : eigenvalues[1] / eigenvalues[0];
    result.parameter_error = std::sqrt(e0 * e0 + e1 * e1);
    result.rms_prediction_error = std::sqrt(squared_error_sum / static_cast<double>(samples));
    return result;
}

void report(const std::string& name, const Result& result) {
    std::cout << "\n" << name << '\n'
              << "  final theta_hat       = [" << result.theta_hat[0] << ", "
              << result.theta_hat[1] << "]\n"
              << "  parameter error norm  = " << result.parameter_error << '\n'
              << "  RMS prediction error  = " << result.rms_prediction_error << '\n'
              << "  Gramian eigenvalues   = [" << result.lambda_min << ", "
              << result.lambda_max << "]\n"
              << "  condition number      = " << result.condition_number << '\n';
}

int main() {
    try {
        const Result non_pe = simulate(false, "Chapter10_Lesson5_non_pe_cpp.csv");
        const Result pe = simulate(true, "Chapter10_Lesson5_pe_cpp.csv");
        report("Rank-deficient regressor", non_pe);
        report("Full-rank exciting regressor", pe);
    } catch (const std::exception& exception) {
        std::cerr << "Error: " << exception.what() << '\n';
        return 1;
    }
    return 0;
}

Chapter10_Lesson5.java

Implementation: Java standard library, immutable result record, explicit numerical integration, and CSV output.

// Chapter10_Lesson5.java
// Standard-library Java implementation of the excitation/identifiability example.

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 Chapter10_Lesson5 {
    private static final double DT = 0.002;
    private static final double END_TIME = 60.0;
    private static final double GAMMA = 1.2;
    private static final double[] THETA_TRUE = {1.5, -0.8};

    private Chapter10_Lesson5() {
    }

    private record Result(
            double[] thetaHat,
            double[][] gramian,
            double lambdaMin,
            double lambdaMax,
            double conditionNumber,
            double parameterError,
            double rmsPredictionError) {
    }

    private static double[] regressor(double t, boolean peCase) {
        if (peCase) {
            return new double[]{Math.sin(t), Math.cos(Math.sqrt(2.0) * t)};
        }
        double s = Math.sin(t);
        return new double[]{s, 2.0 * s};
    }

    private static double[] eigenvaluesSymmetric2x2(double[][] matrix) {
        double a = matrix[0][0];
        double b = matrix[0][1];
        double d = matrix[1][1];
        double trace = a + d;
        double discriminant = Math.sqrt((a - d) * (a - d) + 4.0 * b * b);
        return new double[]{0.5 * (trace - discriminant), 0.5 * (trace + discriminant)};
    }

    private static Result simulate(boolean peCase, Path csvPath) throws IOException {
        double[] thetaHat = {0.0, 0.0};
        double[][] gramian = { {0.0, 0.0}, {0.0, 0.0} };
        double squaredErrorSum = 0.0;
        int samples = 0;
        int steps = (int) Math.round(END_TIME / DT);

        try (BufferedWriter writer = Files.newBufferedWriter(csvPath)) {
            writer.write("time,theta1_hat,theta2_hat,prediction_error\n");

            for (int k = 0; k <= steps; k++) {
                double t = k * DT;
                double[] phi = regressor(t, peCase);
                double noise = 0.02 * Math.sin(13.7 * t) + 0.015 * Math.cos(7.3 * t);
                double y = phi[0] * THETA_TRUE[0] + phi[1] * THETA_TRUE[1] + noise;
                double yHat = phi[0] * thetaHat[0] + phi[1] * thetaHat[1];
                double error = y - yHat;
                double normalization = 1.0 + phi[0] * phi[0] + phi[1] * phi[1];

                thetaHat[0] += DT * GAMMA * phi[0] * error / normalization;
                thetaHat[1] += DT * GAMMA * phi[1] * error / normalization;

                gramian[0][0] += DT * phi[0] * phi[0];
                gramian[0][1] += DT * phi[0] * phi[1];
                gramian[1][0] += DT * phi[1] * phi[0];
                gramian[1][1] += DT * phi[1] * phi[1];

                squaredErrorSum += error * error;
                samples++;

                if (k % 50 == 0) {
                    writer.write(String.format(
                            Locale.US,
                            "%.12f,%.12f,%.12f,%.12f%n",
                            t, thetaHat[0], thetaHat[1], error));
                }
            }
        }

        double[] eigenvalues = eigenvaluesSymmetric2x2(gramian);
        double conditionNumber = eigenvalues[0] <= 1.0e-12
                ? Double.POSITIVE_INFINITY
                : eigenvalues[1] / eigenvalues[0];
        double e0 = thetaHat[0] - THETA_TRUE[0];
        double e1 = thetaHat[1] - THETA_TRUE[1];
        double parameterError = Math.hypot(e0, e1);
        double rmsPredictionError = Math.sqrt(squaredErrorSum / samples);

        return new Result(
                thetaHat.clone(), gramian, eigenvalues[0], eigenvalues[1],
                conditionNumber, parameterError, rmsPredictionError);
    }

    private static void report(String name, Result result) {
        System.out.println("\n" + name);
        System.out.printf(Locale.US, "  final theta_hat       = [%.8f, %.8f]%n",
                result.thetaHat()[0], result.thetaHat()[1]);
        System.out.printf(Locale.US, "  parameter error norm  = %.8f%n", result.parameterError());
        System.out.printf(Locale.US, "  RMS prediction error  = %.8f%n", result.rmsPredictionError());
        System.out.printf(Locale.US, "  Gramian eigenvalues   = [%.8e, %.8e]%n",
                result.lambdaMin(), result.lambdaMax());
        System.out.println("  condition number      = " + result.conditionNumber());
    }

    public static void main(String[] args) {
        try {
            Result nonPe = simulate(false, Path.of("Chapter10_Lesson5_non_pe_java.csv"));
            Result pe = simulate(true, Path.of("Chapter10_Lesson5_pe_java.csv"));
            report("Rank-deficient regressor", nonPe);
            report("Full-rank exciting regressor", pe);
        } catch (IOException exception) {
            System.err.println("I/O error: " + exception.getMessage());
            System.exit(1);
        }
    }
}

Chapter10_Lesson5.m

Implementation: base MATLAB. The same script can be translated into a Simulink model using Integrator blocks for the two parameter states, a MATLAB Function block for the regressor, and Matrix Multiply blocks for the normalized gradient law.

% Chapter10_Lesson5.m
% Rank-deficient versus full-rank excitation in a two-parameter gradient estimator.
% This script uses base MATLAB only.

clear; clc; close all;

dt = 0.002;
Tend = 60.0;
gamma = 1.2;
thetaTrue = [1.5; -0.8];

nonPE = runCase(false, dt, Tend, gamma, thetaTrue);
PE = runCase(true, dt, Tend, gamma, thetaTrue);

printReport('Rank-deficient regressor', nonPE, thetaTrue);
printReport('Full-rank exciting regressor', PE, thetaTrue);

figure;
plot(nonPE.t, nonPE.thetaHat(1, :), 'LineWidth', 1.2); hold on;
plot(nonPE.t, nonPE.thetaHat(2, :), 'LineWidth', 1.2);
yline(thetaTrue(1), '--', 'true theta_1');
yline(thetaTrue(2), ':', 'true theta_2');
grid on;
xlabel('Time [s]'); ylabel('Parameter value');
title('Rank-deficient data: prediction can improve while parameters remain wrong');
legend('theta_1 estimate', 'theta_2 estimate', 'Location', 'best');

figure;
plot(PE.t, PE.thetaHat(1, :), 'LineWidth', 1.2); hold on;
plot(PE.t, PE.thetaHat(2, :), 'LineWidth', 1.2);
yline(thetaTrue(1), '--', 'true theta_1');
yline(thetaTrue(2), ':', 'true theta_2');
grid on;
xlabel('Time [s]'); ylabel('Parameter value');
title('Full-rank excitation: both parameter directions are learned');
legend('theta_1 estimate', 'theta_2 estimate', 'Location', 'best');

writematrix([nonPE.t(:), nonPE.thetaHat.', nonPE.predictionError(:)], ...
    'Chapter10_Lesson5_non_pe_matlab.csv');
writematrix([PE.t(:), PE.thetaHat.', PE.predictionError(:)], ...
    'Chapter10_Lesson5_pe_matlab.csv');

function result = runCase(peCase, dt, Tend, gamma, thetaTrue)
    t = 0:dt:Tend;
    thetaHat = zeros(2, numel(t));
    predictionError = zeros(1, numel(t));
    gramian = zeros(2, 2);

    estimate = [0; 0];
    for k = 1:numel(t)
        tk = t(k);
        if peCase
            phi = [sin(tk); cos(sqrt(2) * tk)];
        else
            phi = [sin(tk); 2 * sin(tk)];
        end

        noise = 0.02 * sin(13.7 * tk) + 0.015 * cos(7.3 * tk);
        y = phi.' * thetaTrue + noise;
        error = y - phi.' * estimate;
        normalization = 1 + phi.' * phi;
        estimateDot = gamma * phi * error / normalization;
        estimate = estimate + dt * estimateDot;

        gramian = gramian + dt * (phi * phi.');
        thetaHat(:, k) = estimate;
        predictionError(k) = error;
    end

    eigenvalues = sort(eig(gramian));
    if eigenvalues(1) <= 1e-12
        conditionNumber = Inf;
    else
        conditionNumber = eigenvalues(2) / eigenvalues(1);
    end

    result.t = t;
    result.thetaHat = thetaHat;
    result.predictionError = predictionError;
    result.gramian = gramian;
    result.lambdaMin = eigenvalues(1);
    result.lambdaMax = eigenvalues(2);
    result.conditionNumber = conditionNumber;
    result.parameterError = norm(estimate - thetaTrue);
    result.rmsPredictionError = sqrt(mean(predictionError.^2));
end

function printReport(name, result, thetaTrue)
    finalEstimate = result.thetaHat(:, end);
    fprintf('\n%s\n', name);
    fprintf('  final thetaHat       = [%.8f, %.8f]\n', finalEstimate(1), finalEstimate(2));
    fprintf('  true theta           = [%.8f, %.8f]\n', thetaTrue(1), thetaTrue(2));
    fprintf('  parameter error norm = %.8f\n', result.parameterError);
    fprintf('  RMS prediction error = %.8f\n', result.rmsPredictionError);
    fprintf('  Gramian eigenvalues  = [%.8e, %.8e]\n', result.lambdaMin, result.lambdaMax);
    fprintf('  condition number     = %.8e\n', result.conditionNumber);
end

Chapter10_Lesson5_Simulink.m

Implementation: programmatic Simulink model construction. The function creates either a full-rank exciting model or a rank-deficient model using standard source, product, sum, unit-delay, and workspace logging blocks.

function Chapter10_Lesson5_Simulink(peCase)
% Chapter10_Lesson5_Simulink.m
% Programmatically builds a discrete Simulink model for the normalized
% gradient estimator used in Chapter 10, Lesson 5.
%
% Usage:
%   Chapter10_Lesson5_Simulink(true)   % full-rank exciting regressor
%   Chapter10_Lesson5_Simulink(false)  % rank-deficient regressor
%
% Requires Simulink.

if nargin == 0
    peCase = true;
end

if ~license('test', 'Simulink')
    error('Simulink is required to build this model.');
end

dt = 0.002;
Tend = 60.0;
gamma = 1.2;
thetaTrue = [1.5; -0.8];

t = (0:dt:Tend).';
phi1 = sin(t);
if peCase
    phi2 = cos(sqrt(2) * t);
    caseName = 'PE';
else
    phi2 = 2 * sin(t);
    caseName = 'NonPE';
end
noise = 0.02 * sin(13.7 * t) + 0.015 * cos(7.3 * t);
y = phi1 * thetaTrue(1) + phi2 * thetaTrue(2) + noise;

assignin('base', 'Chapter10_Lesson5_phi1_ts', timeseries(phi1, t));
assignin('base', 'Chapter10_Lesson5_phi2_ts', timeseries(phi2, t));
assignin('base', 'Chapter10_Lesson5_y_ts', timeseries(y, t));
assignin('base', 'Chapter10_Lesson5_dt', dt);
assignin('base', 'Chapter10_Lesson5_gamma', gamma);

model = ['Chapter10_Lesson5_' caseName '_Model'];
if bdIsLoaded(model)
    close_system(model, 0);
end
new_system(model);
open_system(model);

% Sources and constants.
add_block('simulink/Sources/From Workspace', [model '/phi1'], ...
    'VariableName', 'Chapter10_Lesson5_phi1_ts', ...
    'Position', [30 50 150 80]);
add_block('simulink/Sources/From Workspace', [model '/phi2'], ...
    'VariableName', 'Chapter10_Lesson5_phi2_ts', ...
    'Position', [30 130 150 160]);
add_block('simulink/Sources/From Workspace', [model '/measured_y'], ...
    'VariableName', 'Chapter10_Lesson5_y_ts', ...
    'Position', [30 230 150 260]);
add_block('simulink/Sources/Constant', [model '/gamma'], ...
    'Value', 'Chapter10_Lesson5_gamma', ...
    'Position', [250 330 310 360]);
add_block('simulink/Sources/Constant', [model '/one'], ...
    'Value', '1', ...
    'Position', [250 430 310 460]);
add_block('simulink/Sources/Constant', [model '/dt'], ...
    'Value', 'Chapter10_Lesson5_dt', ...
    'Position', [760 330 820 360]);

% Parameter memories break the feedback loop and implement forward Euler.
add_block('simulink/Discrete/Unit Delay', [model '/theta1_memory'], ...
    'InitialCondition', '0', 'SampleTime', 'Chapter10_Lesson5_dt', ...
    'Position', [980 70 1040 100]);
add_block('simulink/Discrete/Unit Delay', [model '/theta2_memory'], ...
    'InitialCondition', '0', 'SampleTime', 'Chapter10_Lesson5_dt', ...
    'Position', [980 170 1040 200]);

% Predicted output y_hat = phi1*theta1 + phi2*theta2.
add_block('simulink/Math Operations/Product', [model '/phi1_theta1'], ...
    'Inputs', '**', 'Position', [250 50 310 90]);
add_block('simulink/Math Operations/Product', [model '/phi2_theta2'], ...
    'Inputs', '**', 'Position', [250 130 310 170]);
add_block('simulink/Math Operations/Sum', [model '/predicted_y'], ...
    'Inputs', '++', 'Position', [380 90 410 140]);
add_block('simulink/Math Operations/Sum', [model '/prediction_error'], ...
    'Inputs', '+-', 'Position', [480 210 510 260]);

% Normalization denominator 1 + phi1^2 + phi2^2.
add_block('simulink/Math Operations/Product', [model '/phi1_squared'], ...
    'Inputs', '**', 'Position', [250 500 310 540]);
add_block('simulink/Math Operations/Product', [model '/phi2_squared'], ...
    'Inputs', '**', 'Position', [250 570 310 610]);
add_block('simulink/Math Operations/Sum', [model '/normalization'], ...
    'Inputs', '+++', 'Position', [410 490 440 560]);

% theta_dot_i = gamma*phi_i*error/normalization.
add_block('simulink/Math Operations/Product', [model '/theta1_numerator'], ...
    'Inputs', '***', 'Position', [600 40 660 90]);
add_block('simulink/Math Operations/Product', [model '/theta2_numerator'], ...
    'Inputs', '***', 'Position', [600 140 660 190]);
add_block('simulink/Math Operations/Product', [model '/theta1_dot'], ...
    'Inputs', '*/', 'Position', [730 50 790 90]);
add_block('simulink/Math Operations/Product', [model '/theta2_dot'], ...
    'Inputs', '*/', 'Position', [730 150 790 190]);

% Forward-Euler increments and state updates.
add_block('simulink/Math Operations/Product', [model '/theta1_increment'], ...
    'Inputs', '**', 'Position', [850 50 910 90]);
add_block('simulink/Math Operations/Product', [model '/theta2_increment'], ...
    'Inputs', '**', 'Position', [850 150 910 190]);
add_block('simulink/Math Operations/Sum', [model '/theta1_next'], ...
    'Inputs', '++', 'Position', [940 40 970 100]);
add_block('simulink/Math Operations/Sum', [model '/theta2_next'], ...
    'Inputs', '++', 'Position', [940 140 970 200]);

% Logged outputs.
add_block('simulink/Sinks/To Workspace', [model '/theta1_log'], ...
    'VariableName', 'Chapter10_Lesson5_theta1_sim', ...
    'SaveFormat', 'Structure With Time', ...
    'Position', [1120 60 1240 90]);
add_block('simulink/Sinks/To Workspace', [model '/theta2_log'], ...
    'VariableName', 'Chapter10_Lesson5_theta2_sim', ...
    'SaveFormat', 'Structure With Time', ...
    'Position', [1120 160 1240 190]);
add_block('simulink/Sinks/To Workspace', [model '/error_log'], ...
    'VariableName', 'Chapter10_Lesson5_error_sim', ...
    'SaveFormat', 'Structure With Time', ...
    'Position', [590 230 710 260]);

% Connections for prediction.
add_line(model, 'phi1/1', 'phi1_theta1/1', 'autorouting', 'on');
add_line(model, 'theta1_memory/1', 'phi1_theta1/2', 'autorouting', 'on');
add_line(model, 'phi2/1', 'phi2_theta2/1', 'autorouting', 'on');
add_line(model, 'theta2_memory/1', 'phi2_theta2/2', 'autorouting', 'on');
add_line(model, 'phi1_theta1/1', 'predicted_y/1', 'autorouting', 'on');
add_line(model, 'phi2_theta2/1', 'predicted_y/2', 'autorouting', 'on');
add_line(model, 'measured_y/1', 'prediction_error/1', 'autorouting', 'on');
add_line(model, 'predicted_y/1', 'prediction_error/2', 'autorouting', 'on');

% Connections for normalization.
add_line(model, 'phi1/1', 'phi1_squared/1', 'autorouting', 'on');
add_line(model, 'phi1/1', 'phi1_squared/2', 'autorouting', 'on');
add_line(model, 'phi2/1', 'phi2_squared/1', 'autorouting', 'on');
add_line(model, 'phi2/1', 'phi2_squared/2', 'autorouting', 'on');
add_line(model, 'one/1', 'normalization/1', 'autorouting', 'on');
add_line(model, 'phi1_squared/1', 'normalization/2', 'autorouting', 'on');
add_line(model, 'phi2_squared/1', 'normalization/3', 'autorouting', 'on');

% Connections for theta1 update.
add_line(model, 'gamma/1', 'theta1_numerator/1', 'autorouting', 'on');
add_line(model, 'phi1/1', 'theta1_numerator/2', 'autorouting', 'on');
add_line(model, 'prediction_error/1', 'theta1_numerator/3', 'autorouting', 'on');
add_line(model, 'theta1_numerator/1', 'theta1_dot/1', 'autorouting', 'on');
add_line(model, 'normalization/1', 'theta1_dot/2', 'autorouting', 'on');
add_line(model, 'theta1_dot/1', 'theta1_increment/1', 'autorouting', 'on');
add_line(model, 'dt/1', 'theta1_increment/2', 'autorouting', 'on');
add_line(model, 'theta1_memory/1', 'theta1_next/1', 'autorouting', 'on');
add_line(model, 'theta1_increment/1', 'theta1_next/2', 'autorouting', 'on');
add_line(model, 'theta1_next/1', 'theta1_memory/1', 'autorouting', 'on');

% Connections for theta2 update.
add_line(model, 'gamma/1', 'theta2_numerator/1', 'autorouting', 'on');
add_line(model, 'phi2/1', 'theta2_numerator/2', 'autorouting', 'on');
add_line(model, 'prediction_error/1', 'theta2_numerator/3', 'autorouting', 'on');
add_line(model, 'theta2_numerator/1', 'theta2_dot/1', 'autorouting', 'on');
add_line(model, 'normalization/1', 'theta2_dot/2', 'autorouting', 'on');
add_line(model, 'theta2_dot/1', 'theta2_increment/1', 'autorouting', 'on');
add_line(model, 'dt/1', 'theta2_increment/2', 'autorouting', 'on');
add_line(model, 'theta2_memory/1', 'theta2_next/1', 'autorouting', 'on');
add_line(model, 'theta2_increment/1', 'theta2_next/2', 'autorouting', 'on');
add_line(model, 'theta2_next/1', 'theta2_memory/1', 'autorouting', 'on');

% Logging branches.
add_line(model, 'theta1_memory/1', 'theta1_log/1', 'autorouting', 'on');
add_line(model, 'theta2_memory/1', 'theta2_log/1', 'autorouting', 'on');
add_line(model, 'prediction_error/1', 'error_log/1', 'autorouting', 'on');

set_param(model, ...
    'StopTime', num2str(Tend), ...
    'SolverType', 'Fixed-step', ...
    'Solver', 'FixedStepDiscrete', ...
    'FixedStep', num2str(dt));

save_system(model);
fprintf('Created %s.slx. Run it with: sim(''%s'')\n', model, model);
end

Chapter10_Lesson5.nb

Implementation: Wolfram Mathematica notebook expression using Module, Eigenvalues, and ListLinePlot.

Notebook[{
  Cell["Chapter 10, Lesson 5: Practical Limits of Full Parameter Identification", "Title"],
  Cell["Chapter10_Lesson5.nb", "Subtitle"],
  Cell[BoxData[
    RowBox[{"ClearAll", "[", "\"Global`*\"", "]"}]], "Input"],
  Cell[BoxData[
    RowBox[{
      RowBox[{"dt", "=", "0.002"}], ";",
      RowBox[{"tEnd", "=", "60.0"}], ";",
      RowBox[{"gamma", "=", "1.2"}], ";",
      RowBox[{"thetaTrue", "=", RowBox[{"{", RowBox[{"1.5", ",", RowBox[{"-", "0.8"}]}], "}"}]}], ";"
    }]], "Input"],
  Cell[BoxData[
    RowBox[{
      RowBox[{"regressor", "[", RowBox[{"t_", ",", "peCase_"}], "]"}], ":=",
      RowBox[{"If", "[", RowBox[{"peCase", ",",
        RowBox[{"{", RowBox[{RowBox[{"Sin", "[", "t", "]"}], ",", RowBox[{"Cos", "[", RowBox[{RowBox[{"Sqrt", "[", "2", "]"}], " ", "t"}], "]"}]}], "}"}], ",",
        RowBox[{"{", RowBox[{RowBox[{"Sin", "[", "t", "]"}], ",", RowBox[{"2", " ", RowBox[{"Sin", "[", "t", "]"}]}]}], "}"}]}], "]"}]
    }]], "Input"],
  Cell[BoxData[
    RowBox[{
      RowBox[{"runCase", "[", "peCase_", "]"}], ":=",
      RowBox[{"Module", "[",
        RowBox[{
          RowBox[{"{", RowBox[{"times", ",", "estimate", ",", "gramian", ",", "history", ",", "errors", ",", "phi", ",", "noise", ",", "y", ",", "error", ",", "normalization", ",", "eigenvalues"}], "}"}], ",",
          RowBox[{
            RowBox[{"times", "=", RowBox[{"Range", "[", RowBox[{"0.0", ",", "tEnd", ",", "dt"}], "]"}]}], ";",
            RowBox[{"estimate", "=", RowBox[{"{", RowBox[{"0.0", ",", "0.0"}], "}"}]}], ";",
            RowBox[{"gramian", "=", RowBox[{"ConstantArray", "[", RowBox[{"0.0", ",", RowBox[{"{", RowBox[{"2", ",", "2"}], "}"}]}], "]"}]}], ";",
            RowBox[{"history", "=", RowBox[{"{", "}"}]}], ";",
            RowBox[{"errors", "=", RowBox[{"{", "}"}]}], ";",
            RowBox[{"Do", "[",
              RowBox[{
                RowBox[{
                  RowBox[{"phi", "=", RowBox[{"regressor", "[", RowBox[{"t", ",", "peCase"}], "]"}]}], ";",
                  RowBox[{"noise", "=", RowBox[{RowBox[{"0.02", " ", RowBox[{"Sin", "[", RowBox[{"13.7", " ", "t"}], "]"}]}], "+", RowBox[{"0.015", " ", RowBox[{"Cos", "[", RowBox[{"7.3", " ", "t"}], "]"}]}]}]}], ";",
                  RowBox[{"y", "=", RowBox[{RowBox[{"phi", ".", "thetaTrue"}], "+", "noise"}]}], ";",
                  RowBox[{"error", "=", RowBox[{"y", "-", RowBox[{"phi", ".", "estimate"}]}]}], ";",
                  RowBox[{"normalization", "=", RowBox[{"1.0", "+", RowBox[{"phi", ".", "phi"}]}]}], ";",
                  RowBox[{"estimate", "=", RowBox[{"estimate", "+", RowBox[{"dt", " ", "gamma", " ", "phi", " ", "error", "/", "normalization"}]}]}], ";",
                  RowBox[{"gramian", "=", RowBox[{"gramian", "+", RowBox[{"dt", " ", RowBox[{"Outer", "[", RowBox[{"Times", ",", "phi", ",", "phi"}], "]"}]}]}]}], ";",
                  RowBox[{"history", "=", RowBox[{"Append", "[", RowBox[{"history", ",", "estimate"}], "]"}]}], ";",
                  RowBox[{"errors", "=", RowBox[{"Append", "[", RowBox[{"errors", ",", "error"}], "]"}]}]
                }], ",", RowBox[{"{", RowBox[{"t", ",", "times"}], "}"}]
              }], "]"}], ";",
            RowBox[{"eigenvalues", "=", RowBox[{"Sort", "[", RowBox[{"Eigenvalues", "[", "gramian", "]"}], "]"}]}], ";",
            RowBox[{"Association", "[",
              RowBox[{
                RowBox[{"\"Times\"", "->", "times"}], ",",
                RowBox[{"\"History\"", "->", "history"}], ",",
                RowBox[{"\"Errors\"", "->", "errors"}], ",",
                RowBox[{"\"Gramian\"", "->", "gramian"}], ",",
                RowBox[{"\"Eigenvalues\"", "->", "eigenvalues"}], ",",
                RowBox[{"\"ParameterError\"", "->", RowBox[{"Norm", "[", RowBox[{"estimate", "-", "thetaTrue"}], "]"}]}], ",",
                RowBox[{"\"RMSError\"", "->", RowBox[{"Sqrt", "[", RowBox[{"Mean", "[", RowBox[{"errors", "^", "2"}], "]"}], "]"}]}]
              }], "]"}]
          }]
        }], "]"}]
    }]], "Input"],
  Cell[BoxData[
    RowBox[{
      RowBox[{"nonPE", "=", RowBox[{"runCase", "[", "False", "]"}]}], ";",
      RowBox[{"pe", "=", RowBox[{"runCase", "[", "True", "]"}]}], ";",
      RowBox[{"Print", "[", RowBox[{"\"Non-PE eigenvalues: \"", ",", RowBox[{"nonPE", "[", "\"Eigenvalues\"", "]"}], ",", "\"; parameter error: \"", ",", RowBox[{"nonPE", "[", "\"ParameterError\"", "]"}]}], "]"}], ";",
      RowBox[{"Print", "[", RowBox[{"\"PE eigenvalues: \"", ",", RowBox[{"pe", "[", "\"Eigenvalues\"", "]"}], ",", "\"; parameter error: \"", ",", RowBox[{"pe", "[", "\"ParameterError\"", "]"}]}], "]"}]
    }]], "Input"],
  Cell[BoxData[
    RowBox[{"ListLinePlot", "[",
      RowBox[{
        RowBox[{"Transpose", "[", RowBox[{"pe", "[", "\"History\"", "]"}], "]"}], ",",
        RowBox[{"PlotLegends", "->", RowBox[{"{", RowBox[{"\"theta1 estimate\"", ",", "\"theta2 estimate\""}], "}"}]}], ",",
        RowBox[{"PlotLabel", "->", "\"Full-rank excitation\""}], ",",
        RowBox[{"AxesLabel", "->", RowBox[{"{", RowBox[{"\"sample\"", ",", "\"parameter value\""}], "}"}]}]
      }], "]"}]], "Input"]
}, WindowSize -> {1100, 760}, StyleDefinitions -> "Default.nb"]

Problem 1 (Null-space ambiguity): Let \( \boldsymbol{\phi}(t)=[q(t),-3q(t)]^{T} \) and \( y(t)=\boldsymbol{\phi}^{T}(t)\boldsymbol{\theta}^{\ast} \). Find the identifiable parameter combination and a basis for the unidentifiable subspace.

Solution: The output is

\[ y(t)=q(t)(\theta_{1}^{\ast}-3\theta_{2}^{\ast}). \]

Therefore, only \( \eta^{\ast}=\theta_{1}^{\ast}-3\theta_{2}^{\ast} \) is identifiable. A vector \( \mathbf{v} \) is unidentifiable if \( [1,-3]\mathbf{v}=0 \). Taking \( \mathbf{v}=[3,1]^{T} \) gives a basis. Every vector \( \boldsymbol{\theta}^{\ast}+c[3,1]^{T} \) produces the same output.

Problem 2 (Gramian eigenvalues): For \( \boldsymbol{\phi}(t)=[\sin t,2\sin t]^{T} \) on \( [0,2\pi] \), compute the Gramian and its eigenvalues.

Solution: Since \( \int_{0}^{2\pi}\sin^{2}t\,dt=\pi \),

\[ \mathbf{G}=\pi \begin{bmatrix}1&2\\2&4\end{bmatrix}. \]

The base matrix has trace \( 5 \) and determinant \( 0 \), so its eigenvalues are \( 5 \) and \( 0 \). Thus the Gramian eigenvalues are \( 5\pi \) and \( 0 \); PE fails for a two-dimensional parameter vector.

Problem 3 (Steady regressor limitation): Suppose \( \boldsymbol{\phi}(t)→\boldsymbol{\phi}_{\infty}\ne\mathbf{0} \) for a model with \( p=4 \) unknown parameters. Show why the limiting windowed Gramian cannot be positive definite.

Solution: For fixed window length \( T \),

\[ \mathbf{G}_{T}(t)→ T\boldsymbol{\phi}_{\infty}\boldsymbol{\phi}_{\infty}^{T}. \]

An outer product of one nonzero vector has rank one. A positive-definite four-by-four matrix must have rank four. Therefore, at least three parameter directions remain unexcited asymptotically.

Problem 4 (Conditioning and noise): A two-parameter information matrix has eigenvalues \( 100 \) and \( 0.01 \). If noise variance is \( \sigma^{2}=0.04 \), estimate the parameter variance along each eigenvector using \( \operatorname{Cov}(\hat{\boldsymbol{\theta} })\approx \sigma^{2}\mathbf{G}^{-1} \).

Solution: Along an eigenvector with eigenvalue \( \lambda_i \), the variance is approximately \( \sigma^{2}/\lambda_i \). Hence

\[ \operatorname{var}_{1}\approx\frac{0.04}{100}=0.0004, \qquad \operatorname{var}_{2}\approx\frac{0.04}{0.01}=4. \]

The weak direction has ten thousand times more variance. Although the matrix is mathematically nonsingular, the second parameter direction is practically unreliable.

Problem 5 (Tracking versus identification): An adaptive regulator reaches \( e(t)→0 \) for a constant command, while its parameter estimates settle at values different from laboratory measurements. Is this result necessarily a failure of the controller?

Solution: No. Tracking convergence verifies the control objective but not full parameter convergence. With a constant command, the closed-loop regressor may become constant or low rank. The estimates may lie in an equivalence set that generates the same closed-loop behavior. The controller has failed only if exact identification was an explicit and justified requirement. The engineer should inspect the windowed Gramian, residuals, model structure, projection activity, and sensitivity before assigning physical meaning to the estimates.

Problem 6 (Safe excitation design): A plant permits an added probing input \( u_p(t)=a_1\sin(\omega_1 t)+a_2\sin(\omega_2 t) \). Explain how amplitudes and frequencies should be selected conceptually.

Solution: The frequencies should generate independent regressor directions over the relevant model bandwidth and should avoid frequencies that are strongly attenuated or nearly redundant. Amplitudes must be large enough to raise the smallest Gramian eigenvalue above the noise floor, yet small enough to satisfy actuator, state, comfort, wear, and safety limits. The resulting experiment should be evaluated using both control-performance metrics and excitation metrics; maximizing signal amplitude alone is not a valid design rule.

Full parameter identification is difficult because adaptive control data are generated inside a feedback loop whose purpose is usually to suppress variation. Prediction or tracking convergence constrains only the parameter combinations visible through the regressor. The excitation Gramian exposes the geometry: its null space contains completely unidentifiable directions, and its small eigenvalues contain weakly identifiable directions that are highly sensitive to noise and modeling error. Closed-loop regulation, correlated regressors, overparameterization, finite data, constraints, unmodeled dynamics, and robust modifications all weaken identification. Engineers should therefore distinguish control performance from physical parameter recovery, monitor excitation explicitly, reduce unnecessary parameterization, and report uncertainty rather than treating every converged estimate as a verified plant constant.