After completing this lesson, students should be able to:

• recognize a nonlinear system in strict-feedback form;
• identify where unknown constant parameters enter each subsystem;
• separate known nonlinearities, linearly parameterized uncertainty, and unmodeled dynamics;
• construct the first backstepping error coordinate and virtual control;
• derive a tuning function that cancels parameter-error cross terms;
• verify a complete Lyapunov proof for a second-order uncertain system; and
• simulate the resulting adaptive closed loop in several programming environments.

Students are assumed to know Lyapunov stability, nonlinear state feedback, ordinary backstepping for known systems, and standard numerical integration. No later adaptive-backstepping concepts, such as command filters or dynamic surface control, are required.

Let \( ar{x}_i=[x_1,\ldots,x_i]^\mathsf{T} \). A single-input nonlinear system is in strict-feedback form when it can be written as

\[ \begin{aligned} \dot x_1 &= f_1(x_1)+g_1(x_1)x_2+ \boldsymbol{\phi}_1^\mathsf{T}(x_1)\boldsymbol{\theta},\\ \dot x_2 &= f_2(\bar{x}_2)+g_2(\bar{x}_2)x_3+ \boldsymbol{\phi}_2^\mathsf{T}(\bar{x}_2)\boldsymbol{\theta},\\ &\hspace{1.2cm}\vdots\\ \dot x_{n-1} &= f_{n-1}(\bar{x}_{n-1}) +g_{n-1}(\bar{x}_{n-1})x_n +\boldsymbol{\phi}_{n-1}^\mathsf{T}(\bar{x}_{n-1}) \boldsymbol{\theta},\\ \dot x_n &= f_n(x)+g_n(x)u+ \boldsymbol{\phi}_n^\mathsf{T}(x)\boldsymbol{\theta},\\ y &= x_1 . \end{aligned} \]

The state \(x_{i+1}\) enters the \(i\)th equation affinely and acts as a virtual control during the recursive design. The actual control input appears only in the final equation. The triangular dependence is essential: the \(i\)th subsystem depends only on \(\bar{x}_i\), the next state, and the common unknown parameter vector.

flowchart TD
  X1["Subsystem 1: x1_dot"] -->|"virtual input x2"| X2["Subsystem 2: x2_dot"]
  X2 -->|"virtual input x3"| X3["Intermediate subsystems"]
  X3 -->|"virtual input xn"| XN["Final subsystem: xn_dot"]
  U["actual input u"] --> XN
  TH["unknown constant parameter vector theta"] --> X1
  TH --> X2
  TH --> X3
  TH --> XN
  X1 --> Y["measured output y = x1"]
        

The elementary adaptive-backstepping construction in this chapter uses the following assumptions.

1. Known smooth structure: \(f_i\), \(g_i\), and \(\boldsymbol{\phi}_i\) are known and sufficiently differentiable.
2. Known control directions: each virtual-control coefficient has known sign and does not vanish:

   \[ 0 < g_{i,\min} \le |g_i(\bar{x}_i)| \le g_{i,\max}, \qquad i=1,\ldots,n . \]

3. Constant parametric uncertainty:

   \[ \boldsymbol{\theta}\in\mathbb{R}^p, \qquad \dot{\boldsymbol{\theta}}=0, \qquad \boldsymbol{\theta}\text{ is unknown}. \]

4. Linear parameterization: unknown quantities appear as known regressors multiplied by the same unknown parameter vector.
5. Full-state availability: all states required by the state-feedback design are measured or otherwise available.
6. Regular reference: the desired output \(r(t)\) and the derivatives required by the recursion are bounded and known.

Unknown control direction, unmodeled dynamics, actuator saturation, and output-feedback reconstruction are separate complications treated in other chapters. They must not be silently absorbed into \(\boldsymbol{\theta}\).

Adaptive backstepping exploits uncertainty that is nonlinear in the state but linear in unknown constants:

\[ \Delta_i(\bar{x}_i) =\boldsymbol{\phi}_i^\mathsf{T}(\bar{x}_i) \boldsymbol{\theta} =\sum_{j=1}^p \phi_{ij}(\bar{x}_i)\theta_j . \]

The functions \(\phi_{ij}\) may be highly nonlinear. Linearity is required only with respect to the unknown coefficients. For example,

\[ \theta_1x_1^3+\theta_2\sin(x_1x_2) = \begin{bmatrix}x_1^3 & \sin(x_1x_2)\end{bmatrix} \begin{bmatrix}\theta_1\\\theta_2\end{bmatrix} . \]

For output tracking, define the first error coordinate

\[ z_1=x_1-r . \]

A stabilizing function \(\alpha_1\) is designed as the desired value of \(x_2\). The second error coordinate is

\[ z_2=x_2-\alpha_1 . \]

Recursively,

\[ z_i=x_i-\alpha_{i-1}, \qquad i=2,\ldots,n, \qquad \alpha_0=r . \]

In nonadaptive backstepping, each \(\alpha_i\) is chosen to remove known nonlinearities and create a negative quadratic term. Adaptive backstepping adds a parameter estimate \(\hat{\boldsymbol{\theta}}\) and constructs the control and update law together so that every parameter-error cross term is canceled or dominated.

\[ \tilde{\boldsymbol{\theta}} =\hat{\boldsymbol{\theta}}-\boldsymbol{\theta} . \]

This sign convention is used throughout the lesson.

Consider the first subsystem

\[ \dot x_1=f_1(x_1)+g_1(x_1)x_2+ \boldsymbol{\phi}_1^\mathsf{T}(x_1)\boldsymbol{\theta} . \]

Because \(x_2=z_2+\alpha_1\),

\[ \dot z_1=f_1+g_1(z_2+\alpha_1) +\boldsymbol{\phi}_1^\mathsf{T}\boldsymbol{\theta} -\dot r . \]

A certainty-equivalent stabilizing function is

\[ \alpha_1= \frac{1}{g_1} \left[ -f_1+\dot r-k_1z_1 -\boldsymbol{\phi}_1^\mathsf{T} \hat{\boldsymbol{\theta}} \right], \qquad k_1>0 . \]

Substitution gives the exact error equation

\[ \dot z_1=-k_1z_1+g_1z_2 -\boldsymbol{\phi}_1^\mathsf{T} \tilde{\boldsymbol{\theta}} . \]

With a first-step Lyapunov candidate

\[ V_1=\frac{1}{2}z_1^2+ \frac{1}{2} \tilde{\boldsymbol{\theta}}^\mathsf{T} \boldsymbol{\Gamma}^{-1} \tilde{\boldsymbol{\theta}}, \qquad \boldsymbol{\Gamma}=\boldsymbol{\Gamma}^\mathsf{T}>0 , \]

its derivative is

\[ \dot V_1= -k_1z_1^2+g_1z_1z_2+ \tilde{\boldsymbol{\theta}}^\mathsf{T} \left[ \boldsymbol{\Gamma}^{-1} \dot{\hat{\boldsymbol{\theta}}} -\boldsymbol{\phi}_1z_1 \right]. \]

If the design ended here, the update law \(\dot{\hat{\boldsymbol{\theta}}} =\boldsymbol{\Gamma}\boldsymbol{\phi}_1z_1\) would cancel the parameter-error term. In a higher-order system, however, later coordinate transformations generate additional parameter-error terms. Therefore the final update law is postponed and assembled using tuning functions.

\[ \boldsymbol{\tau}_1 =\boldsymbol{\phi}_1z_1 \]

is the first tuning function. Lesson 3 develops the complete recursion.

Consider the uncertain strict-feedback plant

\[ \dot x_1=x_2+\theta x_1, \qquad \dot x_2=u, \qquad y=x_1 , \]

where \(\theta\) is an unknown constant. Let \(r,\dot r,\ddot r\) be bounded and available. Define

\[ z_1=x_1-r, \qquad \alpha_1=-k_1z_1-\hat\theta x_1+\dot r, \qquad z_2=x_2-\alpha_1 . \]

Then

\[ \dot z_1=-k_1z_1+z_2-\tilde\theta x_1, \qquad \tilde\theta=\hat\theta-\theta . \]

Since \(\alpha_1\) depends on \(x_1\) and \(\hat\theta\), define

\[ \alpha_{1,x} =\frac{\partial\alpha_1}{\partial x_1} =-k_1-\hat\theta . \]

Using \(\theta=\hat\theta-\tilde\theta\), the derivative of the virtual control can be decomposed into known and parameter-error parts:

\[ \dot\alpha_1= \underbrace{ \alpha_{1,x}(x_2+\hat\theta x_1) +k_1\dot r+\ddot r }_{\dot\alpha_{1,\mathrm{known}}} -\alpha_{1,x}x_1\tilde\theta -x_1\dot{\hat\theta} . \]

Therefore,

\[ \dot z_2= u-\dot\alpha_{1,\mathrm{known}} +\alpha_{1,x}x_1\tilde\theta +x_1\dot{\hat\theta} . \]

The second tuning function is selected as

\[ \tau_2=x_1z_1-\alpha_{1,x}x_1z_2, \qquad \dot{\hat\theta}=\gamma\tau_2, \qquad \gamma>0 . \]

The actual control is

\[ u= \dot\alpha_{1,\mathrm{known}} -z_1-k_2z_2-x_1\dot{\hat\theta}, \qquad k_2>0 . \]

Stability Theorem

For any positive \(k_1,k_2,\gamma\), the preceding controller and update law make \(z_1,z_2,\tilde\theta\) bounded and guarantee \(z_1(t)→0\) and \(z_2(t)→0\).

Proof

Choose

\[ V=\frac{1}{2}z_1^2+ \frac{1}{2}z_2^2+ \frac{1}{2\gamma}\tilde\theta^2 . \]

Using the two error equations,

\[ \begin{aligned} \dot V &= z_1(-k_1z_1+z_2-\tilde\theta x_1) +z_2\left( u-\dot\alpha_{1,\mathrm{known}} +\alpha_{1,x}x_1\tilde\theta +x_1\dot{\hat\theta} \right) +\frac{1}{\gamma}\tilde\theta\dot{\hat\theta} . \end{aligned} \]

Insert the control law:

\[ \begin{aligned} \dot V &=-k_1z_1^2-k_2z_2^2\\ &\quad+ \tilde\theta \left( -x_1z_1+\alpha_{1,x}x_1z_2+ \frac{1}{\gamma}\dot{\hat\theta} \right). \end{aligned} \]

The update law \(\dot{\hat\theta}= \gamma(x_1z_1-\alpha_{1,x}x_1z_2)\) cancels the bracketed term exactly, so

\[ \boxed{\dot V=-k_1z_1^2-k_2z_2^2\le0} . \]

Hence \(V(t)\le V(0)\), so \(z_1,z_2,\tilde\theta\in\mathcal{L}_\infty\). Boundedness of the reference and smooth closed-loop expressions implies bounded derivatives of the tracking coordinates. Moreover,

\[ \int_0^\infty \left(k_1z_1^2+k_2z_2^2\right)dt \le V(0), \]

so \(z_1,z_2\in\mathcal{L}_2\). Barbalat's lemma then gives \(z_1(t)→0\) and \(z_2(t)→0\). Parameter convergence is not implied without an excitation condition. This completes the proof.

flowchart TD
  A["Write plant in strict-feedback form"] --> B["Define z1 = x1 - reference"]
  B --> C["Treat x2 as virtual control"]
  C --> D["Choose alpha1 using parameter estimates"]
  D --> E["Expose parameter-error cross terms"]
  E --> F["Define first tuning function"]
  F --> G["Set z2 = x2 - alpha1"]
  G --> H["Differentiate alpha1 and collect new cross terms"]
  H --> I["Extend tuning function"]
  I --> J["Choose actual control and final update law"]
  J --> K["Prove negative-semidefinite Lyapunov derivative"]
        

The key point is that adaptation is not appended after designing a nonadaptive controller. The virtual controls, actual control, tuning functions, and Lyapunov function are constructed as one recursive object.

The identity \(\dot V=-k_1z_1^2-k_2z_2^2\) establishes tracking without requiring the estimate to converge to the true parameter. The estimate is bounded because it is part of the Lyapunov function, but convergence of \(\hat\theta\) requires sufficient information in the regressor.

\[ z_1,z_2→0 \quad\not\Rightarrow\quad \tilde\theta→0 . \]

In the numerical example, the sinusoidal reference and transient state trajectory usually provide enough excitation for the scalar estimate to approach the true value. That simulation outcome should not be mistaken for a general parameter-convergence theorem.

Increasing \(k_1\) and \(k_2\) increases nominal error damping. Increasing \(\gamma\) accelerates parameter adaptation but may amplify measurement noise and numerical stiffness in practical implementations.

This implementation uses NumPy for arrays, SciPy for adaptive-step integration, and Matplotlib for visualization.

Chapter17_Lesson1.py

"""
Chapter17_Lesson1.py
Adaptive backstepping simulation for
    x1_dot = x2 + theta*x1
    x2_dot = u
with sinusoidal reference tracking.

Dependencies:
    numpy, scipy, matplotlib
Install:
    python -m pip install numpy scipy matplotlib
"""

from __future__ import annotations

import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt


THETA_TRUE = 1.2
K1 = 2.0
K2 = 2.5
GAMMA = 3.0


def reference(t: float) -> tuple[float, float, float]:
    """Return r(t), r_dot(t), and r_ddot(t)."""
    r = 0.5 * np.sin(0.7 * t)
    r_dot = 0.35 * np.cos(0.7 * t)
    r_ddot = -0.245 * np.sin(0.7 * t)
    return float(r), float(r_dot), float(r_ddot)


def closed_loop(t: float, state: np.ndarray) -> np.ndarray:
    """Closed-loop plant and parameter-estimator dynamics."""
    x1, x2, theta_hat = state
    r, r_dot, r_ddot = reference(t)

    phi = x1
    z1 = x1 - r
    alpha1 = -K1 * z1 - theta_hat * phi + r_dot
    z2 = x2 - alpha1

    # Partial derivative of alpha1 with respect to x1 for phi(x1)=x1.
    alpha1_x1 = -K1 - theta_hat

    # Tuning function and adaptive law.
    tau2 = phi * z1 - alpha1_x1 * phi * z2
    theta_hat_dot = GAMMA * tau2

    # Known part of alpha1_dot after replacing theta by theta_hat.
    alpha1_dot_known = (
        alpha1_x1 * (x2 + theta_hat * phi)
        + K1 * r_dot
        + r_ddot
    )

    # Actual control input.
    u = (
        alpha1_dot_known
        - z1
        - K2 * z2
        - phi * theta_hat_dot
    )

    x1_dot = x2 + THETA_TRUE * phi
    x2_dot = u
    return np.array([x1_dot, x2_dot, theta_hat_dot], dtype=float)


def main() -> None:
    initial_state = np.array([1.2, -0.5, 0.0], dtype=float)
    t_eval = np.linspace(0.0, 20.0, 4001)

    solution = solve_ivp(
        closed_loop,
        (t_eval[0], t_eval[-1]),
        initial_state,
        t_eval=t_eval,
        method="RK45",
        rtol=1e-9,
        atol=1e-11,
        max_step=0.02,
    )
    if not solution.success:
        raise RuntimeError(solution.message)

    t = solution.t
    x1, x2, theta_hat = solution.y
    r = 0.5 * np.sin(0.7 * t)
    tracking_error = x1 - r

    print(f"Final tracking error: {tracking_error[-1]: .6e}")
    print(f"Final parameter estimate: {theta_hat[-1]: .6f}")
    print(f"True parameter: {THETA_TRUE: .6f}")

    plt.figure()
    plt.plot(t, x1, label="x1")
    plt.plot(t, r, "--", label="reference")
    plt.xlabel("Time (s)")
    plt.ylabel("Output")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()

    plt.figure()
    plt.plot(t, tracking_error)
    plt.xlabel("Time (s)")
    plt.ylabel("Tracking error z1")
    plt.grid(True)
    plt.tight_layout()

    plt.figure()
    plt.plot(t, theta_hat, label="theta_hat")
    plt.axhline(THETA_TRUE, linestyle="--", label="theta_true")
    plt.xlabel("Time (s)")
    plt.ylabel("Parameter")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()

    plt.show()


if __name__ == "__main__":
    main()

The C++ version uses only the ISO C++ standard library and implements fourth-order Runge-Kutta integration from scratch. It writes a CSV file suitable for plotting with Python, MATLAB, GNUplot, or spreadsheet tools.

Chapter17_Lesson1.cpp

/*
Chapter17_Lesson1.cpp
From-scratch RK4 simulation of a second-order adaptive-backstepping loop.

Build:
    g++ -std=c++17 -O2 Chapter17_Lesson1.cpp -o Chapter17_Lesson1
Run:
    ./Chapter17_Lesson1
*/

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

namespace {

constexpr double theta_true = 1.2;
constexpr double k1 = 2.0;
constexpr double k2 = 2.5;
constexpr double gamma_gain = 3.0;

using State = std::array<double, 3>;

struct Reference {
    double r;
    double r_dot;
    double r_ddot;
};

Reference reference(double t) {
    return {
        0.5 * std::sin(0.7 * t),
        0.35 * std::cos(0.7 * t),
        -0.245 * std::sin(0.7 * t)
    };
}

State dynamics(double t, const State& state) {
    const double x1 = state[0];
    const double x2 = state[1];
    const double theta_hat = state[2];
    const Reference ref = reference(t);

    const double phi = x1;
    const double z1 = x1 - ref.r;
    const double alpha1 = -k1 * z1 - theta_hat * phi + ref.r_dot;
    const double z2 = x2 - alpha1;
    const double alpha1_x1 = -k1 - theta_hat;

    const double tau2 = phi * z1 - alpha1_x1 * phi * z2;
    const double theta_hat_dot = gamma_gain * tau2;

    const double alpha1_dot_known =
        alpha1_x1 * (x2 + theta_hat * phi)
        + k1 * ref.r_dot
        + ref.r_ddot;

    const double u =
        alpha1_dot_known
        - z1
        - k2 * z2
        - phi * theta_hat_dot;

    const double x1_dot = x2 + theta_true * phi;
    const double x2_dot = u;
    return {x1_dot, x2_dot, theta_hat_dot};
}

State add_scaled(const State& a, const State& b, double scale) {
    return {
        a[0] + scale * b[0],
        a[1] + scale * b[1],
        a[2] + scale * b[2]
    };
}

State rk4_step(double t, const State& state, double h) {
    const State k_1 = dynamics(t, state);
    const State k_2 = dynamics(t + 0.5 * h, add_scaled(state, k_1, 0.5 * h));
    const State k_3 = dynamics(t + 0.5 * h, add_scaled(state, k_2, 0.5 * h));
    const State k_4 = dynamics(t + h, add_scaled(state, k_3, h));

    return {
        state[0] + h * (k_1[0] + 2.0 * k_2[0] + 2.0 * k_3[0] + k_4[0]) / 6.0,
        state[1] + h * (k_1[1] + 2.0 * k_2[1] + 2.0 * k_3[1] + k_4[1]) / 6.0,
        state[2] + h * (k_1[2] + 2.0 * k_2[2] + 2.0 * k_3[2] + k_4[2]) / 6.0
    };
}

}  // namespace

int main() {
    const double h = 0.002;
    const double final_time = 20.0;
    State state{1.2, -0.5, 0.0};

    std::ofstream csv("Chapter17_Lesson1_cpp.csv");
    if (!csv) {
        throw std::runtime_error("Could not open output CSV file.");
    }
    csv << "t,x1,x2,theta_hat,reference,tracking_error\n";
    csv << std::setprecision(12);

    for (double t = 0.0; t <= final_time + 0.5 * h; t += h) {
        const Reference ref = reference(t);
        csv << t << ','
            << state[0] << ','
            << state[1] << ','
            << state[2] << ','
            << ref.r << ','
            << state[0] - ref.r << '\n';

        state = rk4_step(t, state, h);
    }

    const Reference final_ref = reference(final_time);
    std::cout << std::setprecision(8)
              << "Final tracking error: " << state[0] - final_ref.r << '\n'
              << "Final parameter estimate: " << state[2] << '\n'
              << "True parameter: " << theta_true << '\n'
              << "Data written to Chapter17_Lesson1_cpp.csv\n";
    return 0;
}

The Java implementation uses the standard JDK only. It mirrors the C++ RK4 simulation and writes the complete state and tracking history to CSV.

Chapter17_Lesson1.java

/*
Chapter17_Lesson1.java
From-scratch RK4 simulation of a second-order adaptive-backstepping loop.

Build:
    javac Chapter17_Lesson1.java
Run:
    java Chapter17_Lesson1
*/

import java.io.BufferedWriter;
import java.io.FileWriter;
import java.io.IOException;
import java.util.Locale;

public final class Chapter17_Lesson1 {
    private static final double THETA_TRUE = 1.2;
    private static final double K1 = 2.0;
    private static final double K2 = 2.5;
    private static final double GAMMA = 3.0;

    private Chapter17_Lesson1() {
    }

    private static double[] reference(double t) {
        return new double[] {
            0.5 * Math.sin(0.7 * t),
            0.35 * Math.cos(0.7 * t),
            -0.245 * Math.sin(0.7 * t)
        };
    }

    private static double[] dynamics(double t, double[] state) {
        double x1 = state[0];
        double x2 = state[1];
        double thetaHat = state[2];
        double[] ref = reference(t);

        double phi = x1;
        double z1 = x1 - ref[0];
        double alpha1 = -K1 * z1 - thetaHat * phi + ref[1];
        double z2 = x2 - alpha1;
        double alpha1X1 = -K1 - thetaHat;

        double tau2 = phi * z1 - alpha1X1 * phi * z2;
        double thetaHatDot = GAMMA * tau2;

        double alpha1DotKnown =
            alpha1X1 * (x2 + thetaHat * phi)
            + K1 * ref[1]
            + ref[2];

        double u =
            alpha1DotKnown
            - z1
            - K2 * z2
            - phi * thetaHatDot;

        double x1Dot = x2 + THETA_TRUE * phi;
        double x2Dot = u;
        return new double[] {x1Dot, x2Dot, thetaHatDot};
    }

    private static double[] addScaled(double[] a, double[] b, double scale) {
        return new double[] {
            a[0] + scale * b[0],
            a[1] + scale * b[1],
            a[2] + scale * b[2]
        };
    }

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

        return new double[] {
            state[0] + h * (k1[0] + 2.0 * k2[0] + 2.0 * k3[0] + k4[0]) / 6.0,
            state[1] + h * (k1[1] + 2.0 * k2[1] + 2.0 * k3[1] + k4[1]) / 6.0,
            state[2] + h * (k1[2] + 2.0 * k2[2] + 2.0 * k3[2] + k4[2]) / 6.0
        };
    }

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

        double h = 0.002;
        double finalTime = 20.0;
        double[] state = {1.2, -0.5, 0.0};

        try (BufferedWriter writer =
                 new BufferedWriter(new FileWriter("Chapter17_Lesson1_java.csv"))) {
            writer.write("t,x1,x2,theta_hat,reference,tracking_error");
            writer.newLine();

            for (double t = 0.0; t <= finalTime + 0.5 * h; t += h) {
                double[] ref = reference(t);
                writer.write(String.format(
                    Locale.US,
                    "%.12f,%.12f,%.12f,%.12f,%.12f,%.12f",
                    t, state[0], state[1], state[2], ref[0], state[0] - ref[0]
                ));
                writer.newLine();
                state = rk4Step(t, state, h);
            }
        } catch (IOException exception) {
            System.err.println("I/O error: " + exception.getMessage());
            System.exit(1);
        }

        double[] finalRef = reference(finalTime);
        System.out.printf(Locale.US, "Final tracking error: %.8e%n",
                          state[0] - finalRef[0]);
        System.out.printf(Locale.US, "Final parameter estimate: %.8f%n", state[2]);
        System.out.printf(Locale.US, "True parameter: %.8f%n", THETA_TRUE);
        System.out.println("Data written to Chapter17_Lesson1_java.csv");
    }
}

The first file uses ode45. The second file creates a Simulink model programmatically with a Clock block, a vector Integrator, a MATLAB Function block containing the adaptive closed-loop dynamics, and a workspace sink.

Chapter17_Lesson1.m

% Chapter17_Lesson1.m
% Adaptive backstepping simulation using ode45.
%
% Required product:
%   MATLAB
% Optional product:
%   Simulink (see Chapter17_Lesson1_Simulink.m)

clear; clc; close all;

thetaTrue = 1.2;
k1 = 2.0;
k2 = 2.5;
gamma = 3.0;

initialState = [1.2; -0.5; 0.0];
timeSpan = [0, 20];

options = odeset( ...
    'RelTol', 1e-9, ...
    'AbsTol', 1e-11, ...
    'MaxStep', 0.02);

[t, state] = ode45(@closedLoop, timeSpan, initialState, options);

x1 = state(:, 1);
x2 = state(:, 2);
thetaHat = state(:, 3);
referenceSignal = 0.5 * sin(0.7 * t);
trackingError = x1 - referenceSignal;

fprintf('Final tracking error: %.6e\n', trackingError(end));
fprintf('Final parameter estimate: %.6f\n', thetaHat(end));
fprintf('True parameter: %.6f\n', thetaTrue);

figure;
plot(t, x1, 'LineWidth', 1.4); hold on;
plot(t, referenceSignal, '--', 'LineWidth', 1.4);
grid on;
xlabel('Time (s)');
ylabel('Output');
legend('x_1', 'reference', 'Location', 'best');
title('Adaptive-Backstepping Tracking');

figure;
plot(t, trackingError, 'LineWidth', 1.4);
grid on;
xlabel('Time (s)');
ylabel('z_1');
title('Tracking Error');

figure;
plot(t, thetaHat, 'LineWidth', 1.4); hold on;
yline(thetaTrue, '--', 'LineWidth', 1.4);
grid on;
xlabel('Time (s)');
ylabel('Parameter');
legend('\hat{\theta}', '\theta', 'Location', 'best');
title('Parameter Estimate');

results = table(t, x1, x2, thetaHat, referenceSignal, trackingError);
writetable(results, 'Chapter17_Lesson1_matlab.csv');

function derivative = closedLoop(t, state)
    thetaTrue = 1.2;
    k1 = 2.0;
    k2 = 2.5;
    gamma = 3.0;

    x1 = state(1);
    x2 = state(2);
    thetaHat = state(3);

    r = 0.5 * sin(0.7 * t);
    rDot = 0.35 * cos(0.7 * t);
    rDDot = -0.245 * sin(0.7 * t);

    phi = x1;
    z1 = x1 - r;
    alpha1 = -k1 * z1 - thetaHat * phi + rDot;
    z2 = x2 - alpha1;
    alpha1X1 = -k1 - thetaHat;

    tau2 = phi * z1 - alpha1X1 * phi * z2;
    thetaHatDot = gamma * tau2;

    alpha1DotKnown = ...
        alpha1X1 * (x2 + thetaHat * phi) ...
        + k1 * rDot ...
        + rDDot;

    u = ...
        alpha1DotKnown ...
        - z1 ...
        - k2 * z2 ...
        - phi * thetaHatDot;

    x1Dot = x2 + thetaTrue * phi;
    x2Dot = u;

    derivative = [x1Dot; x2Dot; thetaHatDot];
end

Chapter17_Lesson1_Simulink.m

% Chapter17_Lesson1_Simulink.m
% Programmatically construct and simulate a Simulink model for the
% second-order adaptive-backstepping example.
%
% Required products:
%   MATLAB, Simulink

clear; clc;

model = 'Chapter17_Lesson1_Simulink_Model';
if bdIsLoaded(model)
    close_system(model, 0);
end
if isfile([model, '.slx'])
    delete([model, '.slx']);
end

new_system(model);
open_system(model);

add_block('simulink/Sources/Clock', ...
    [model, '/Clock'], ...
    'Position', [40, 80, 70, 100]);

add_block('simulink/Continuous/Integrator', ...
    [model, '/State Integrator'], ...
    'InitialCondition', '[1.2; -0.5; 0.0]', ...
    'Position', [410, 65, 445, 115]);

add_block('simulink/User-Defined Functions/MATLAB Function', ...
    [model, '/Adaptive Closed Loop'], ...
    'Position', [160, 45, 330, 135]);

add_block('simulink/Sinks/To Workspace', ...
    [model, '/State Data'], ...
    'VariableName', 'stateData', ...
    'SaveFormat', 'Structure With Time', ...
    'Position', [520, 70, 610, 110]);

add_line(model, 'Clock/1', 'Adaptive Closed Loop/1', 'autorouting', 'on');
add_line(model, 'State Integrator/1', ...
    'Adaptive Closed Loop/2', 'autorouting', 'on');
add_line(model, 'Adaptive Closed Loop/1', ...
    'State Integrator/1', 'autorouting', 'on');
add_line(model, 'State Integrator/1', ...
    'State Data/1', 'autorouting', 'on');

root = sfroot;
chart = root.find('-isa', 'Stateflow.EMChart', ...
                  'Path', [model, '/Adaptive Closed Loop']);

chart.Script = sprintf([ ...
'function derivative = fcn(t, state)\n' ...
'%%#codegen\n' ...
'thetaTrue = 1.2;\n' ...
'k1 = 2.0;\n' ...
'k2 = 2.5;\n' ...
'gamma = 3.0;\n' ...
'x1 = state(1);\n' ...
'x2 = state(2);\n' ...
'thetaHat = state(3);\n' ...
'r = 0.5*sin(0.7*t);\n' ...
'rDot = 0.35*cos(0.7*t);\n' ...
'rDDot = -0.245*sin(0.7*t);\n' ...
'phi = x1;\n' ...
'z1 = x1-r;\n' ...
'alpha1 = -k1*z1-thetaHat*phi+rDot;\n' ...
'z2 = x2-alpha1;\n' ...
'alpha1X1 = -k1-thetaHat;\n' ...
'tau2 = phi*z1-alpha1X1*phi*z2;\n' ...
'thetaHatDot = gamma*tau2;\n' ...
'alpha1DotKnown = alpha1X1*(x2+thetaHat*phi)+k1*rDot+rDDot;\n' ...
'u = alpha1DotKnown-z1-k2*z2-phi*thetaHatDot;\n' ...
'x1Dot = x2+thetaTrue*phi;\n' ...
'x2Dot = u;\n' ...
'derivative = [x1Dot; x2Dot; thetaHatDot];\n' ...
'end\n']);

set_param(model, ...
    'StopTime', '20', ...
    'Solver', 'ode45', ...
    'MaxStep', '0.02');

save_system(model);
simulationOutput = sim(model);

time = simulationOutput.stateData.time;
state = simulationOutput.stateData.signals.values;
referenceSignal = 0.5 * sin(0.7 * time);

figure;
plot(time, state(:, 1), 'LineWidth', 1.4); hold on;
plot(time, referenceSignal, '--', 'LineWidth', 1.4);
grid on;
xlabel('Time (s)');
ylabel('Output');
legend('x_1', 'reference', 'Location', 'best');
title('Simulink Adaptive-Backstepping Tracking');

disp(['Saved model: ', model, '.slx']);

The notebook uses NDSolveValue and produces tracking, error, and parameter-estimate plots.

Chapter17_Lesson1.nb

        Notebook[{
  Cell[
    "Chapter 17, Lesson 1: Strict-Feedback Systems with Parametric Uncertainty",
    "Title"
  ],
  Cell[
    "Adaptive backstepping simulation for a second-order uncertain strict-feedback system.",
    "Text"
  ],
  Cell[
"ClearAll[\"Global`*\"];

thetaTrue = 1.2; k1 = 2.0; k2 = 2.5; gamma = 3.0;

r[t_] := 0.5 Sin[0.7 t];
rDot[t_] := 0.35 Cos[0.7 t];
rDDot[t_] := -0.245 Sin[0.7 t];

phi[t_] := x1[t];
z1[t_] := x1[t] - r[t];
alpha1[t_] := -k1 z1[t] - thetaHat[t] phi[t] + rDot[t];
z2[t_] := x2[t] - alpha1[t];
alpha1X1[t_] := -k1 - thetaHat[t];

tau2[t_] := phi[t] z1[t] - alpha1X1[t] phi[t] z2[t];
thetaHatDot[t_] := gamma tau2[t];

alpha1DotKnown[t_] :=
  alpha1X1[t] (x2[t] + thetaHat[t] phi[t]) +
  k1 rDot[t] + rDDot[t];

control[t_] :=
  alpha1DotKnown[t] - z1[t] - k2 z2[t] -
  phi[t] thetaHatDot[t];

solution = NDSolveValue[
  {
    x1'[t] == x2[t] + thetaTrue phi[t],
    x2'[t] == control[t],
    thetaHat'[t] == thetaHatDot[t],
    x1[0] == 1.2, x2[0] == -0.5, thetaHat[0] == 0.0
  },
  {x1, x2, thetaHat}, {t, 0, 20},
  Method -> {
    \"TimeIntegration\" -> {
      \"ExplicitRungeKutta\", \"DifferenceOrder\" -> 5
    }
  }
];

x1Solution[t_] := solution[[1]][t];
x2Solution[t_] := solution[[2]][t];
thetaHatSolution[t_] := solution[[3]][t];

trackingPlot = Plot[
  {x1Solution[t], r[t]}, {t, 0, 20},
  PlotLegends -> {\"x1\", \"reference\"},
  AxesLabel -> {\"Time (s)\", \"Output\"},
  PlotRange -> All, GridLines -> Automatic
];

errorPlot = Plot[
  x1Solution[t] - r[t], {t, 0, 20},
  AxesLabel -> {\"Time (s)\", \"z1\"},
  PlotRange -> All, GridLines -> Automatic
];

parameterPlot = Plot[
  {thetaHatSolution[t], thetaTrue}, {t, 0, 20},
  PlotLegends -> {\"thetaHat\", \"thetaTrue\"},
  AxesLabel -> {\"Time (s)\", \"Parameter\"},
  PlotRange -> All, GridLines -> Automatic
];

Print[trackingPlot]; Print[errorPlot]; Print[parameterPlot];
Print[\"Final tracking error = \", N[x1Solution[20] - r[20], 8]];
Print[\"Final parameter estimate = \", N[thetaHatSolution[20], 8]];
",
    "Input"
  ]
},
  WindowSize -> {1200, 800},
  StyleDefinitions -> "Default.nb"
]

Use the common settings \(\theta=1.2\), \(k_1=2\), \(k_2=2.5\), \(\gamma=3\), and

\[ r(t)=0.5\sin(0.7t). \]

The initial conditions are \(x_1(0)=1.2\), \(x_2(0)=-0.5\), and \(\hat\theta(0)=0\). A correct implementation should show the following:

• the output approaches the sinusoidal reference;
• the tracking coordinate \(z_1\) approaches zero;
• the virtual-control error \(z_2\) approaches zero;
• all simulated states remain bounded; and
• for this particular experiment, the scalar estimate approaches 1.2.

A useful implementation check is to evaluate the theoretical Lyapunov derivative from logged data:

\[ \dot V_{\mathrm{theory}} =-k_1z_1^2-k_2z_2^2 . \]

Finite-difference estimates of \(\dot V\) may show small positive spikes because of discretization error, but the analytical expression is nonpositive.

1. Confusing strict feedback with pure feedback: in strict feedback, the next state enters affinely. A nonlinear occurrence such as \(\sin(x_2)\) in the first equation generally creates a pure-feedback problem.
2. Using a separate estimate at every step: this causes overparameterization when the same physical parameter appears in several equations. Tuning functions retain a single estimate vector.
3. Canceling an unknown quantity directly: virtual controls may use \(\hat\theta\), not \(\theta\).
4. Ignoring derivatives of parameter estimates: differentiating \(\alpha_i\) produces terms containing \(\dot{\hat\theta}\); these must appear consistently in the actual control.
5. Claiming parameter convergence from tracking convergence: asymptotic tracking alone does not establish \(\hat\theta→\theta\).

Problem 1 (Structural Classification): Determine whether each system is in strict-feedback form:

\[ \begin{aligned} \text{(a)}\quad &\dot x_1=x_1^2+x_2,\qquad \dot x_2=\theta\sin x_1+u,\\ \text{(b)}\quad &\dot x_1=x_1+\sin x_2,\qquad \dot x_2=u,\\ \text{(c)}\quad &\dot x_1=x_2+\theta_1x_1^3,\qquad \dot x_2=x_3+\theta_2x_1x_2,\qquad \dot x_3=u . \end{aligned} \]

Solution: Systems (a) and (c) are strict-feedback systems. In (a), \(x_2\) enters the first equation affinely, and the actual input enters the second equation affinely. System (c) has the required lower-triangular chain. System (b) is not in standard strict-feedback form because \(x_2\) enters through \(\sin x_2\); it is a pure-feedback-type structure.

Problem 2 (First Virtual Control): For

\[ \dot x_1=f_1(x_1)+g_1(x_1)x_2+ \phi_1(x_1)\theta , \qquad z_1=x_1-r, \]

derive a virtual control that gives \(\dot z_1=-k_1z_1+g_1z_2-\phi_1\tilde\theta\).

Solution: Set \(z_2=x_2-\alpha_1\). Then

\[ \dot z_1=f_1+g_1(z_2+\alpha_1)+\phi_1\theta-\dot r . \]

Choose

\[ \alpha_1= \frac{-f_1+\dot r-k_1z_1-\phi_1\hat\theta}{g_1} . \]

Since \(\theta-\hat\theta=-\tilde\theta\), substitution yields the required equation.

Problem 3 (Lyapunov Cancellation): For the second-order example, start from

\[ \dot V= -k_1z_1^2-k_2z_2^2+ \tilde\theta \left( -x_1z_1+\alpha_{1,x}x_1z_2+ \frac{1}{\gamma}\dot{\hat\theta} \right). \]

Find the update law that removes the indefinite term.

Solution: Require the bracket to be zero:

\[ \frac{1}{\gamma}\dot{\hat\theta} =x_1z_1-\alpha_{1,x}x_1z_2 . \]

Thus

\[ \dot{\hat\theta} =\gamma\left(x_1z_1-\alpha_{1,x}x_1z_2\right), \]

and the derivative reduces to \(\dot V=-k_1z_1^2-k_2z_2^2\).

Problem 4 (Tracking Does Not Imply Identification): Explain why \(z_1,z_2→0\) does not necessarily imply \(\tilde\theta→0\).

Solution: The Lyapunov derivative contains only the tracking coordinates:

\[ \dot V=-k_1z_1^2-k_2z_2^2 . \]

It guarantees square integrability and convergence of \(z_1,z_2\), but it does not contain a negative term proportional to \(\tilde\theta^2\). Once the tracking errors become zero, the adaptation signal may also become zero for many different parameter estimates. Additional excitation is needed to distinguish the true parameter from these alternatives.

Problem 5 (Vector Parameter Extension): Replace the scalar uncertainty in the first equation by

\[ \boldsymbol{\phi}^\mathsf{T}(x_1) \boldsymbol{\theta} . \]

State the natural vector generalization of the first tuning function and adaptive gain.

Solution: Let \(\boldsymbol{\Gamma}= \boldsymbol{\Gamma}^\mathsf{T}>0\). The first tuning function and tentative update law are

\[ \boldsymbol{\tau}_1= \boldsymbol{\phi}(x_1)z_1, \qquad \dot{\hat{\boldsymbol{\theta}}} =\boldsymbol{\Gamma}\boldsymbol{\tau}_1 . \]

In a complete higher-order design, \(\boldsymbol{\tau}_1\) is extended at subsequent steps before the final update law is applied.

Problem 6 (Implementation Exercise): Modify any supplied program to use \(r(t)=0.4\cos(0.5t)\). Give the required derivatives.

Solution:

\[ r(t)=0.4\cos(0.5t), \qquad \dot r(t)=-0.2\sin(0.5t), \qquad \ddot r(t)=-0.1\cos(0.5t). \]

Replace the three reference expressions consistently. Replacing only \(r\) while leaving its derivatives unchanged breaks the analytical cancellation used in the control law.

Strict-feedback systems possess a lower-triangular structure in which each next state can be treated as a virtual control. When uncertainty is linearly parameterized by unknown constants, adaptive backstepping combines recursive coordinate transformations, certainty-equivalent virtual controls, tuning functions, a parameter update law, and an augmented Lyapunov function. The complete second-order derivation showed how differentiating a virtual control creates new parameter-error terms and how the tuning function cancels them exactly. The next lesson develops virtual controls and parameter dependence more systematically for general strict-feedback chains.