After completing this lesson, students should be able to:

1. formulate a scalar adaptive-control problem with unknown high-frequency-gain sign;
2. verify that a candidate function satisfies the Nussbaum oscillation property;
3. derive a Nussbaum-gain controller and tuning law;
4. prove global boundedness and asymptotic state regulation;
5. reproduce the case study in Python, C++, Java, MATLAB/Simulink, and Wolfram Mathematica; and
6. distinguish a stability guarantee from acceptable transient performance.

The plant is deliberately minimal so every proof step can be inspected. The unknown-direction issue remains fundamental: a conventional fixed-sign feedback can stabilize one plant and destabilize another whose only difference is the actuator-gain sign.

Consider the scalar continuous-time plant

\[ \dot{x}(t)=a x(t)+b u(t), \qquad a\in\mathbb{R}, \qquad b\in\mathbb{R}\setminus\{0\}. \]

The constants \(a\) and \(b\) are unknown. In particular, the controller does not know whether \(b>0\) or \(b<0\). The objective is

\[ \lim_{t→\infty}x(t)=0 \]

while keeping every closed-loop signal bounded. Suppose one ignores the unknown sign and applies \(u=-kx\), where \(k>0\). Then

\[ \dot{x}=(a-bk)x. \]

If \(b>0\), sufficiently large positive \(k\) can make \(a-bk<0\). If \(b<0\), increasing the same gain produces \(a-bk=a+|b|k\), which is more unstable. Thus the uncertainty concerns feedback polarity, not merely gain magnitude.

A continuous function \(N:[0,\infty)→\mathbb{R}\) is a Nussbaum-type function when its averaged integral has unbounded positive and negative excursions:

\[ \limsup_{\zeta→\infty} \frac{1}{\zeta}\int_0^\zeta N(s)\,ds=+\infty, \qquad \liminf_{\zeta→\infty} \frac{1}{\zeta}\int_0^\zeta N(s)\,ds=-\infty. \]

The case study uses

\[ N(\zeta)=\zeta^2\cos\zeta. \]

Its integral can be evaluated exactly:

\[ \int_0^\zeta s^2\cos s\,ds = \zeta^2\sin\zeta+2\zeta\cos\zeta-2\sin\zeta. \]

Choose \(\zeta_n^+=\frac{\pi}{2}+2\pi n\) and \(\zeta_n^-=\frac{3\pi}{2}+2\pi n\). Along the first sequence, the normalized integral grows like \(+\zeta_n^+\); along the second, it grows like \(-\zeta_n^-\). Hence the required conditions hold.

The function does not estimate the sign of \(b\) as a discrete variable. Its oscillatory, increasing-amplitude gain explores both feedback polarities until the Lyapunov argument prevents further growth of the tuning variable.

Use the control and tuning laws

\[ u=N(\kappa)x, \qquad \dot{\kappa}=x^2, \qquad \kappa(0)\geq 0, \]

with

\[ N(\kappa)=\kappa^2\cos\kappa. \]

Since \(\dot{\kappa}=x^2\geq0\), the tuning variable is continuous and nondecreasing. It records accumulated state energy:

\[ \kappa(t)=\kappa(0)+\int_0^t x^2(\tau)\,d\tau. \]

flowchart TD
  X["Measured state x"] --> A["kappa_dot = x^2"]
  A --> K["Adaptive state kappa"]
  K --> N["N(kappa) = kappa^2 cos(kappa)"]
  N --> U["u = N(kappa) x"]
  X --> U
  U --> P["Plant: x_dot = a x + b u"]
  P --> X
        

Unlike an identifier-based design, there is no estimate \(\hat{a}\), no estimate \(\hat{b}\), and no sign classifier. This is a direct adaptive stabilizer.

Substitution gives

\[ \dot{x}=\left(a+bN(\kappa)\right)x, \qquad \dot{\kappa}=x^2. \]

Theorem. For every finite initial condition \((x(0),\kappa(0))\), every \(a\in\mathbb{R}\), and every constant \(b\neq0\), the controller yields a global solution satisfying

\[ x,\;\kappa,\;u\in\mathcal{L}_\infty, \qquad x\in\mathcal{L}_2, \qquad \lim_{t→\infty}x(t)=0. \]

Proof, Step 1: Lyapunov identity. Choose

\[ V(x)=\frac{1}{2}x^2. \]

Then

\[ \dot{V} = x\dot{x} = \left(a+bN(\kappa)\right)x^2 = \left(a+bN(\kappa)\right)\dot{\kappa}. \]

Integration gives

\[ V(t) = V(0) +a\left(\kappa(t)-\kappa(0)\right) +b\int_{\kappa(0)}^{\kappa(t)}N(s)\,ds. \]

Proof, Step 2: boundedness of \(\kappa\). Assume for contradiction that \(\kappa(t)\) is unbounded. Because it is continuous and nondecreasing, it reaches arbitrarily large values. Along a sequence for which \(b\int N\) has a sufficiently large negative averaged value, the Nussbaum term dominates the finite linear term \(a\kappa\). The identity would then imply

\[ V(t_j)<0, \]

contradicting \(V=\frac{1}{2}x^2\geq0\). Therefore \(\kappa\in\mathcal{L}_\infty\).

Proof, Step 3: boundedness of the remaining signals. Since \(\kappa\) is bounded and \(N\) is continuous, \(N(\kappa)\) is bounded. The Lyapunov identity then implies bounded \(V\) and \(x\). Consequently, \(u=N(\kappa)x\) and \(\dot{x}=(a+bN(\kappa))x\) are bounded. No finite-escape time is possible, so the solution is global.

Proof, Step 4: asymptotic regulation. The tuning law and boundedness of \(\kappa\) imply

\[ \int_0^\infty x^2(t)\,dt = \lim_{t→\infty} \left(\kappa(t)-\kappa(0)\right) <\infty. \]

Thus \(x\in\mathcal{L}_2\). Moreover, \(\frac{d}{dt}x^2=2x\dot{x}\) is bounded, so \(x^2\) is uniformly continuous. Barbalat's lemma gives

\[ x^2(t)→0, \qquad x(t)→0. \qquad \blacksquare \]

The implementations use

\[ a=0.8, \qquad x(0)=1.2, \qquad \kappa(0)=0, \]

and compare

\[ b=+2 \qquad\text{and}\qquad b=-2. \]

A fixed-step fourth-order Runge-Kutta method is used with \(h=10^{-3}\) s over 20 s. Representative results are:

Both plants converge, but the negative-gain case experiences a much larger search transient. Asymptotic stability does not imply small overshoot, modest control effort, or compatibility with actuator limits.

flowchart TD
  S["Set a, x0, kappa0, time step"] --> B["Run case b = +2"]
  S --> C["Run case b = -2"]
  B --> R["At each RK4 stage compute N(kappa), u, x_dot, kappa_dot"]
  C --> R
  R --> L["Log x, kappa, u, V"]
  L --> Q["Check boundedness and convergence"]
  Q --> T["Compare transient peaks for opposite control directions"]
        

Required libraries are NumPy and Matplotlib. The controller and RK4 integrator are implemented explicitly.

Chapter16_Lesson5.py

# Chapter16_Lesson5.py
# Case study: adaptive stabilization with unknown input-gain sign.
#
# Plant:
#     x_dot = a*x + b*u
# Controller:
#     u       = N(kappa)*x
#     kappa_dot = x^2
#     N(kappa) = kappa^2*cos(kappa)
#
# The same controller is simulated for b > 0 and b < 0.

from __future__ import annotations

import math
from dataclasses import dataclass

import matplotlib.pyplot as plt
import numpy as np


@dataclass(frozen=True)
class SimulationConfig:
    a: float = 0.8
    x0: float = 1.2
    kappa0: float = 0.0
    final_time: float = 20.0
    step: float = 1.0e-3


def nussbaum(kappa: float) -> float:
    """Polynomial Nussbaum function N(kappa) = kappa^2 cos(kappa)."""
    return kappa * kappa * math.cos(kappa)


def closed_loop_rhs(x: float, kappa: float, a: float, b: float) -> tuple[float, float]:
    """Return (x_dot, kappa_dot) for the closed-loop system."""
    u = nussbaum(kappa) * x
    x_dot = a * x + b * u
    kappa_dot = x * x
    return x_dot, kappa_dot


def rk4_step(x: float, kappa: float, dt: float, a: float, b: float) -> tuple[float, float]:
    """One classical fourth-order Runge-Kutta step."""
    k1x, k1k = closed_loop_rhs(x, kappa, a, b)
    k2x, k2k = closed_loop_rhs(
        x + 0.5 * dt * k1x,
        kappa + 0.5 * dt * k1k,
        a,
        b,
    )
    k3x, k3k = closed_loop_rhs(
        x + 0.5 * dt * k2x,
        kappa + 0.5 * dt * k2k,
        a,
        b,
    )
    k4x, k4k = closed_loop_rhs(
        x + dt * k3x,
        kappa + dt * k3k,
        a,
        b,
    )

    x_next = x + (dt / 6.0) * (k1x + 2.0 * k2x + 2.0 * k3x + k4x)
    kappa_next = kappa + (dt / 6.0) * (k1k + 2.0 * k2k + 2.0 * k3k + k4k)
    return x_next, kappa_next


def simulate(b: float, cfg: SimulationConfig) -> dict[str, np.ndarray]:
    """Simulate one unknown-control-direction case."""
    count = int(round(cfg.final_time / cfg.step)) + 1
    time = np.linspace(0.0, cfg.final_time, count)
    x = np.empty(count)
    kappa = np.empty(count)
    u = np.empty(count)
    lyapunov = np.empty(count)

    x[0] = cfg.x0
    kappa[0] = cfg.kappa0
    u[0] = nussbaum(kappa[0]) * x[0]
    lyapunov[0] = 0.5 * x[0] * x[0]

    for i in range(count - 1):
        x[i + 1], kappa[i + 1] = rk4_step(
            x[i], kappa[i], cfg.step, cfg.a, b
        )
        u[i + 1] = nussbaum(kappa[i + 1]) * x[i + 1]
        lyapunov[i + 1] = 0.5 * x[i + 1] * x[i + 1]

        if not (
            math.isfinite(x[i + 1])
            and math.isfinite(kappa[i + 1])
            and math.isfinite(u[i + 1])
        ):
            raise FloatingPointError(
                "Non-finite numerical value encountered. Reduce the integration step."
            )

    return {
        "time": time,
        "x": x,
        "kappa": kappa,
        "u": u,
        "V": lyapunov,
        "b": np.full(count, b),
    }


def print_summary(result: dict[str, np.ndarray], a: float) -> None:
    """Print key transient and terminal quantities."""
    b = float(result["b"][0])
    x = result["x"]
    kappa = result["kappa"]
    u = result["u"]

    print(
        f"a={a:+.3f}, b={b:+.3f}: "
        f"x(T)={x[-1]:+.6e}, "
        f"kappa(T)={kappa[-1]:.6f}, "
        f"max|x|={np.max(np.abs(x)):.6f}, "
        f"max|u|={np.max(np.abs(u)):.6f}"
    )


def plot_results(results: list[dict[str, np.ndarray]]) -> None:
    """Plot state, adaptive variable, control input, and Lyapunov function."""
    figure, axes = plt.subplots(4, 1, figsize=(10, 11), sharex=True)

    for result in results:
        b = float(result["b"][0])
        label = f"b = {b:+g}"
        axes[0].plot(result["time"], result["x"], label=label)
        axes[1].plot(result["time"], result["kappa"], label=label)
        axes[2].plot(result["time"], result["u"], label=label)
        axes[3].semilogy(
            result["time"],
            np.maximum(result["V"], np.finfo(float).tiny),
            label=label,
        )

    axes[0].set_ylabel("x(t)")
    axes[1].set_ylabel("kappa(t)")
    axes[2].set_ylabel("u(t)")
    axes[3].set_ylabel("V(t)")
    axes[3].set_xlabel("Time [s]")

    for axis in axes:
        axis.grid(True)
        axis.legend()

    figure.suptitle("Nussbaum-Gain Adaptive Stabilization for Unknown Sign of b")
    figure.tight_layout()
    figure.savefig("Chapter16_Lesson5_results.png", dpi=180)
    plt.show()


def main() -> None:
    cfg = SimulationConfig()
    results = [simulate(+2.0, cfg), simulate(-2.0, cfg)]

    for result in results:
        print_summary(result, cfg.a)

    plot_results(results)


if __name__ == "__main__":
    main()

This implementation uses only the C++17 standard library and writes CSV files for subsequent plotting.

Chapter16_Lesson5.cpp

// Chapter16_Lesson5.cpp
// Case study: Nussbaum-gain adaptive stabilization with unknown sign of b.
//
// Build:
//   g++ -std=c++17 -O2 Chapter16_Lesson5.cpp -o Chapter16_Lesson5
//
// Output:
//   Chapter16_Lesson5_b_pos.csv
//   Chapter16_Lesson5_b_neg.csv

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

struct State {
    double x;
    double kappa;
};

struct Derivative {
    double x_dot;
    double kappa_dot;
};

double nussbaum(double kappa) {
    return kappa * kappa * std::cos(kappa);
}

Derivative rhs(const State& state, double a, double b) {
    const double u = nussbaum(state.kappa) * state.x;
    return {
        a * state.x + b * u,
        state.x * state.x
    };
}

State add_scaled(const State& state, const Derivative& derivative, double scale) {
    return {
        state.x + scale * derivative.x_dot,
        state.kappa + scale * derivative.kappa_dot
    };
}

State rk4_step(const State& state, double dt, double a, double b) {
    const Derivative k1 = rhs(state, a, b);
    const Derivative k2 = rhs(add_scaled(state, k1, 0.5 * dt), a, b);
    const Derivative k3 = rhs(add_scaled(state, k2, 0.5 * dt), a, b);
    const Derivative k4 = rhs(add_scaled(state, k3, dt), a, b);

    return {
        state.x + (dt / 6.0) *
            (k1.x_dot + 2.0 * k2.x_dot + 2.0 * k3.x_dot + k4.x_dot),
        state.kappa + (dt / 6.0) *
            (k1.kappa_dot + 2.0 * k2.kappa_dot +
             2.0 * k3.kappa_dot + k4.kappa_dot)
    };
}

void simulate(
    double a,
    double b,
    double x0,
    double kappa0,
    double final_time,
    double dt,
    const std::string& output_path
) {
    std::ofstream output(output_path);
    if (!output) {
        throw std::runtime_error("Cannot open output file: " + output_path);
    }

    output << "time,x,kappa,u,V\n";
    output << std::setprecision(16);

    State state{x0, kappa0};
    const int steps = static_cast<int>(std::llround(final_time / dt));

    double max_abs_x = std::abs(state.x);
    double max_abs_u = std::abs(nussbaum(state.kappa) * state.x);

    for (int i = 0; i <= steps; ++i) {
        const double time = i * dt;
        const double u = nussbaum(state.kappa) * state.x;
        const double V = 0.5 * state.x * state.x;

        output << time << ','
               << state.x << ','
               << state.kappa << ','
               << u << ','
               << V << '\n';

        max_abs_x = std::max(max_abs_x, std::abs(state.x));
        max_abs_u = std::max(max_abs_u, std::abs(u));

        if (i < steps) {
            state = rk4_step(state, dt, a, b);
            if (!std::isfinite(state.x) || !std::isfinite(state.kappa)) {
                throw std::runtime_error(
                    "Non-finite numerical value encountered. Reduce dt."
                );
            }
        }
    }

    std::cout << std::showpos << std::fixed << std::setprecision(3)
              << "a=" << a << ", b=" << b
              << std::noshowpos << std::setprecision(6)
              << ": x(T)=" << state.x
              << ", kappa(T)=" << state.kappa
              << ", max|x|=" << max_abs_x
              << ", max|u|=" << max_abs_u
              << '\n';
}

int main() {
    try {
        constexpr double a = 0.8;
        constexpr double x0 = 1.2;
        constexpr double kappa0 = 0.0;
        constexpr double final_time = 20.0;
        constexpr double dt = 1.0e-3;

        simulate(
            a, +2.0, x0, kappa0, final_time, dt,
            "Chapter16_Lesson5_b_pos.csv"
        );
        simulate(
            a, -2.0, x0, kappa0, final_time, dt,
            "Chapter16_Lesson5_b_neg.csv"
        );
    } catch (const std::exception& error) {
        std::cerr << "Error: " << error.what() << '\n';
        return 1;
    }

    return 0;
}

The Java version uses the standard library, an immutable state class, and CSV output compatible with spreadsheet or plotting software.

Chapter16_Lesson5.java

// Chapter16_Lesson5.java
// Case study: Nussbaum-gain adaptive stabilization with unknown sign of b.
//
// Build and run:
//   javac Chapter16_Lesson5.java
//   java Chapter16_Lesson5
//
// Output:
//   Chapter16_Lesson5_b_pos.csv
//   Chapter16_Lesson5_b_neg.csv

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

public final class Chapter16_Lesson5 {
    private Chapter16_Lesson5() {
        // Utility class.
    }

    private static final class State {
        final double x;
        final double kappa;

        State(double x, double kappa) {
            this.x = x;
            this.kappa = kappa;
        }
    }

    private static final class Derivative {
        final double xDot;
        final double kappaDot;

        Derivative(double xDot, double kappaDot) {
            this.xDot = xDot;
            this.kappaDot = kappaDot;
        }
    }

    private static double nussbaum(double kappa) {
        return kappa * kappa * Math.cos(kappa);
    }

    private static Derivative rhs(State state, double a, double b) {
        double u = nussbaum(state.kappa) * state.x;
        return new Derivative(
            a * state.x + b * u,
            state.x * state.x
        );
    }

    private static State addScaled(
        State state,
        Derivative derivative,
        double scale
    ) {
        return new State(
            state.x + scale * derivative.xDot,
            state.kappa + scale * derivative.kappaDot
        );
    }

    private static State rk4Step(State state, double dt, double a, double b) {
        Derivative k1 = rhs(state, a, b);
        Derivative k2 = rhs(addScaled(state, k1, 0.5 * dt), a, b);
        Derivative k3 = rhs(addScaled(state, k2, 0.5 * dt), a, b);
        Derivative k4 = rhs(addScaled(state, k3, dt), a, b);

        return new State(
            state.x + (dt / 6.0) * (
                k1.xDot + 2.0 * k2.xDot + 2.0 * k3.xDot + k4.xDot
            ),
            state.kappa + (dt / 6.0) * (
                k1.kappaDot
                    + 2.0 * k2.kappaDot
                    + 2.0 * k3.kappaDot
                    + k4.kappaDot
            )
        );
    }

    private static void simulate(
        double a,
        double b,
        double x0,
        double kappa0,
        double finalTime,
        double dt,
        String outputPath
    ) throws IOException {
        int steps = (int) Math.round(finalTime / dt);
        State state = new State(x0, kappa0);

        double maxAbsX = Math.abs(state.x);
        double maxAbsU = Math.abs(nussbaum(state.kappa) * state.x);

        try (
            BufferedWriter writer = Files.newBufferedWriter(
                Path.of(outputPath),
                StandardCharsets.UTF_8
            )
        ) {
            writer.write("time,x,kappa,u,V");
            writer.newLine();

            for (int i = 0; i <= steps; i++) {
                double time = i * dt;
                double u = nussbaum(state.kappa) * state.x;
                double V = 0.5 * state.x * state.x;

                writer.write(
                    String.format(
                        Locale.US,
                        "%.16g,%.16g,%.16g,%.16g,%.16g%n",
                        time,
                        state.x,
                        state.kappa,
                        u,
                        V
                    )
                );

                maxAbsX = Math.max(maxAbsX, Math.abs(state.x));
                maxAbsU = Math.max(maxAbsU, Math.abs(u));

                if (i < steps) {
                    state = rk4Step(state, dt, a, b);
                    if (
                        !Double.isFinite(state.x)
                            || !Double.isFinite(state.kappa)
                    ) {
                        throw new ArithmeticException(
                            "Non-finite value encountered. Reduce dt."
                        );
                    }
                }
            }
        }

        System.out.printf(
            Locale.US,
            "a=%+.3f, b=%+.3f: x(T)=%+.6e, "
                + "kappa(T)=%.6f, max|x|=%.6f, max|u|=%.6f%n",
            a,
            b,
            state.x,
            state.kappa,
            maxAbsX,
            maxAbsU
        );
    }

    public static void main(String[] args) {
        final double a = 0.8;
        final double x0 = 1.2;
        final double kappa0 = 0.0;
        final double finalTime = 20.0;
        final double dt = 1.0e-3;

        try {
            simulate(
                a,
                +2.0,
                x0,
                kappa0,
                finalTime,
                dt,
                "Chapter16_Lesson5_b_pos.csv"
            );
            simulate(
                a,
                -2.0,
                x0,
                kappa0,
                finalTime,
                dt,
                "Chapter16_Lesson5_b_neg.csv"
            );
        } catch (IOException | ArithmeticException error) {
            System.err.println("Error: " + error.getMessage());
            System.exit(1);
        }
    }
}

The MATLAB script contains a from-scratch RK4 simulation. It also contains an optional function that programmatically creates a Simulink model. Set BUILD_SIMULINK_MODEL = true when Simulink is installed.

Chapter16_Lesson5.m

% Chapter16_Lesson5.m
% Case study: adaptive stabilization with unknown input-gain sign.
%
% Plant:
%   x_dot = a*x + b*u
%
% Controller:
%   u         = N(kappa)*x
%   kappa_dot = x^2
%   N(kappa)  = kappa^2*cos(kappa)
%
% The script first runs a from-scratch RK4 simulation for b=+2 and b=-2.
% Set BUILD_SIMULINK_MODEL=true to generate a continuous-time Simulink
% block diagram for one selected value of b.

clear;
clc;
close all;

a = 0.8;
x0 = 1.2;
kappa0 = 0.0;
finalTime = 20.0;
dt = 1.0e-3;
bCases = [2.0, -2.0];

results = cell(size(bCases));

for index = 1:numel(bCases)
    b = bCases(index);
    results{index} = simulateCase(a, b, x0, kappa0, finalTime, dt);

    fprintf( ...
        ['a=%+.3f, b=%+.3f: x(T)=%+.6e, kappa(T)=%.6f, ' ...
         'max|x|=%.6f, max|u|=%.6f\n'], ...
        a, ...
        b, ...
        results{index}.x(end), ...
        results{index}.kappa(end), ...
        max(abs(results{index}.x)), ...
        max(abs(results{index}.u)) ...
    );

    if b > 0
        fileName = 'Chapter16_Lesson5_b_pos.csv';
    else
        fileName = 'Chapter16_Lesson5_b_neg.csv';
    end

    tableOutput = table( ...
        results{index}.time, ...
        results{index}.x, ...
        results{index}.kappa, ...
        results{index}.u, ...
        results{index}.V, ...
        'VariableNames', {'time', 'x', 'kappa', 'u', 'V'} ...
    );
    writetable(tableOutput, fileName);
end

figure('Name', 'Unknown input-gain sign case study');

subplot(4, 1, 1);
hold on;
for index = 1:numel(results)
    plot(results{index}.time, results{index}.x, 'LineWidth', 1.2);
end
grid on;
ylabel('x(t)');
legend('b=+2', 'b=-2');

subplot(4, 1, 2);
hold on;
for index = 1:numel(results)
    plot(results{index}.time, results{index}.kappa, 'LineWidth', 1.2);
end
grid on;
ylabel('\kappa(t)');

subplot(4, 1, 3);
hold on;
for index = 1:numel(results)
    plot(results{index}.time, results{index}.u, 'LineWidth', 1.2);
end
grid on;
ylabel('u(t)');

subplot(4, 1, 4);
hold on;
for index = 1:numel(results)
    semilogy( ...
        results{index}.time, ...
        max(results{index}.V, realmin), ...
        'LineWidth', ...
        1.2 ...
    );
end
grid on;
ylabel('V(t)');
xlabel('Time [s]');

exportgraphics(gcf, 'Chapter16_Lesson5_results.png', 'Resolution', 180);

% Optional Simulink model generation.
BUILD_SIMULINK_MODEL = false;
if BUILD_SIMULINK_MODEL
    selectedB = -2.0;
    buildSimulinkModel( ...
        'Chapter16_Lesson5_Simulink', ...
        a, ...
        selectedB, ...
        x0, ...
        kappa0, ...
        finalTime ...
    );
end

function result = simulateCase(a, b, x0, kappa0, finalTime, dt)
    count = round(finalTime / dt) + 1;
    time = linspace(0.0, finalTime, count).';

    x = zeros(count, 1);
    kappa = zeros(count, 1);
    u = zeros(count, 1);
    V = zeros(count, 1);

    x(1) = x0;
    kappa(1) = kappa0;
    u(1) = nussbaum(kappa(1)) * x(1);
    V(1) = 0.5 * x(1)^2;

    for index = 1:(count - 1)
        current = [x(index); kappa(index)];
        next = rk4Step(current, dt, a, b);

        x(index + 1) = next(1);
        kappa(index + 1) = next(2);
        u(index + 1) = nussbaum(kappa(index + 1)) * x(index + 1);
        V(index + 1) = 0.5 * x(index + 1)^2;

        if any(~isfinite(next))
            error('Non-finite value encountered. Reduce the RK4 step.');
        end
    end

    result = struct( ...
        'time', time, ...
        'x', x, ...
        'kappa', kappa, ...
        'u', u, ...
        'V', V, ...
        'b', b ...
    );
end

function next = rk4Step(state, dt, a, b)
    k1 = closedLoopRhs(state, a, b);
    k2 = closedLoopRhs(state + 0.5 * dt * k1, a, b);
    k3 = closedLoopRhs(state + 0.5 * dt * k2, a, b);
    k4 = closedLoopRhs(state + dt * k3, a, b);

    next = state + (dt / 6.0) * (k1 + 2.0 * k2 + 2.0 * k3 + k4);
end

function derivative = closedLoopRhs(state, a, b)
    x = state(1);
    kappa = state(2);

    u = nussbaum(kappa) * x;
    derivative = [
        a * x + b * u;
        x^2
    ];
end

function value = nussbaum(kappa)
    value = kappa.^2 .* cos(kappa);
end

function buildSimulinkModel(modelName, a, b, x0, kappa0, finalTime)
    if ~license('test', 'Simulink')
        error('Simulink is not available in this MATLAB installation.');
    end

    if bdIsLoaded(modelName)
        close_system(modelName, 0);
    end

    new_system(modelName);
    open_system(modelName);

    add_block( ...
        'simulink/Continuous/Integrator', ...
        [modelName '/x_integrator'], ...
        'InitialCondition', num2str(x0), ...
        'Position', [520, 100, 550, 130] ...
    );
    add_block( ...
        'simulink/Continuous/Integrator', ...
        [modelName '/kappa_integrator'], ...
        'InitialCondition', num2str(kappa0), ...
        'Position', [520, 280, 550, 310] ...
    );

    add_block( ...
        'simulink/Math Operations/Gain', ...
        [modelName '/a_gain'], ...
        'Gain', num2str(a), ...
        'Position', [170, 70, 230, 110] ...
    );
    add_block( ...
        'simulink/Math Operations/Gain', ...
        [modelName '/b_gain'], ...
        'Gain', num2str(b), ...
        'Position', [360, 145, 420, 185] ...
    );
    add_block( ...
        'simulink/Math Operations/Sum', ...
        [modelName '/x_dot_sum'], ...
        'Inputs', '++', ...
        'Position', [455, 95, 480, 155] ...
    );

    add_block( ...
        'simulink/Math Operations/Product', ...
        [modelName '/kappa_squared'], ...
        'Inputs', '**', ...
        'Position', [170, 245, 215, 285] ...
    );
    add_block( ...
        'simulink/Math Operations/Trigonometric Function', ...
        [modelName '/cos_kappa'], ...
        'Operator', 'cos', ...
        'Position', [170, 315, 230, 345] ...
    );
    add_block( ...
        'simulink/Math Operations/Product', ...
        [modelName '/nussbaum_product'], ...
        'Inputs', '**', ...
        'Position', [270, 260, 320, 310] ...
    );
    add_block( ...
        'simulink/Math Operations/Product', ...
        [modelName '/control_product'], ...
        'Inputs', '**', ...
        'Position', [285, 145, 330, 195] ...
    );
    add_block( ...
        'simulink/Math Operations/Product', ...
        [modelName '/x_squared'], ...
        'Inputs', '**', ...
        'Position', [360, 280, 410, 320] ...
    );

    add_block( ...
        'simulink/Sinks/To Workspace', ...
        [modelName '/x_output'], ...
        'VariableName', 'x_simulink', ...
        'SaveFormat', 'Structure With Time', ...
        'Position', [610, 95, 710, 125] ...
    );
    add_block( ...
        'simulink/Sinks/To Workspace', ...
        [modelName '/kappa_output'], ...
        'VariableName', 'kappa_simulink', ...
        'SaveFormat', 'Structure With Time', ...
        'Position', [610, 275, 710, 305] ...
    );

    add_line(modelName, 'x_integrator/1', 'a_gain/1');
    add_line(modelName, 'a_gain/1', 'x_dot_sum/1');
    add_line(modelName, 'x_dot_sum/1', 'x_integrator/1');
    add_line(modelName, 'x_integrator/1', 'x_output/1');

    add_line(modelName, 'kappa_integrator/1', 'kappa_squared/1');
    add_line(modelName, 'kappa_integrator/1', 'kappa_squared/2');
    add_line(modelName, 'kappa_integrator/1', 'cos_kappa/1');
    add_line(modelName, 'kappa_squared/1', 'nussbaum_product/1');
    add_line(modelName, 'cos_kappa/1', 'nussbaum_product/2');

    add_line(modelName, 'nussbaum_product/1', 'control_product/1');
    add_line(modelName, 'x_integrator/1', 'control_product/2');
    add_line(modelName, 'control_product/1', 'b_gain/1');
    add_line(modelName, 'b_gain/1', 'x_dot_sum/2');

    add_line(modelName, 'x_integrator/1', 'x_squared/1');
    add_line(modelName, 'x_integrator/1', 'x_squared/2');
    add_line(modelName, 'x_squared/1', 'kappa_integrator/1');
    add_line(modelName, 'kappa_integrator/1', 'kappa_output/1');

    set_param( ...
        modelName, ...
        'StopTime', num2str(finalTime), ...
        'Solver', 'ode45' ...
    );

    save_system(modelName, [modelName '.slx']);
    fprintf('Generated %s.slx\n', modelName);
end

The notebook contains one executable input cell. Evaluating it runs the same RK4 case study and generates state, adaptive-gain, control, and Lyapunov plots.

Chapter16_Lesson5.nb

Notebook[{
  Cell["Chapter 16, Lesson 5: Unknown Input-Gain Sign", "Title"],
  Cell["Evaluate the input cell to run the RK4 simulation for b=+2 and b=-2.", "Text"],
  Cell[BoxData[
    RowBox[{"ToExpression", "[", "\"ClearAll[\"Global`*\"];
a = 0.8; x0 = 1.2; kappa0 = 0.0; finalTime = 20.0; dt = 0.001;
bCases = {2.0, -2.0};
nussbaum[kappa_?NumericQ] := kappa^2 Cos[kappa];
rhs[{x_?NumericQ, kappa_?NumericQ}, b_?NumericQ] := {
  a x + b nussbaum[kappa] x, x^2
};
rk4Step[state_List, b_?NumericQ] := Module[{k1, k2, k3, k4},
  k1 = rhs[state, b];
  k2 = rhs[state + 0.5 dt k1, b];
  k3 = rhs[state + 0.5 dt k2, b];
  k4 = rhs[state + dt k3, b];
  state + (dt/6.0) (k1 + 2.0 k2 + 2.0 k3 + k4)
];
simulate[b_?NumericQ] := Module[
  {steps, states, times, xValues, kappaValues, uValues, vValues},
  steps = Round[finalTime/dt];
  states = NestList[rk4Step[#, b] &, {x0, kappa0}, steps];
  times = Range[0, steps] dt; xValues = states[[All, 1]];
  kappaValues = states[[All, 2]];
  uValues = MapThread[nussbaum[#2] #1 &, {xValues, kappaValues}];
  vValues = 0.5 xValues^2;
  <|\"b\" -> b, \"Time\" -> times, \"x\" -> xValues,
    \"Kappa\" -> kappaValues, \"u\" -> uValues, \"V\" -> vValues|>
];
results = simulate /@ bCases;
Do[Print[Row[{
  \"a=\", NumberForm[a, {4, 3}], \", b=\", NumberForm[result[\"b\"], {4, 1}],
  \": x(T)=\", ScientificForm[Last[result[\"x\"]], 6],
  \", kappa(T)=\", NumberForm[Last[result[\"Kappa\"]], {8, 6}],
  \", max|x|=\", NumberForm[Max[Abs[result[\"x\"]]], {8, 6}],
  \", max|u|=\", NumberForm[Max[Abs[result[\"u\"]]], {10, 6}]
}]], {result, results}];
statePlot = ListLinePlot[
  Table[Transpose[{result[\"Time\"], result[\"x\"]}], {result, results}],
  PlotLegends -> {\"b = +2\", \"b = -2\"}, Frame -> True,
  FrameLabel -> {\"Time [s]\", \"x(t)\"}, PlotRange -> All, ImageSize -> Large];
kappaPlot = ListLinePlot[
  Table[Transpose[{result[\"Time\"], result[\"Kappa\"]}], {result, results}],
  PlotLegends -> {\"b = +2\", \"b = -2\"}, Frame -> True,
  FrameLabel -> {\"Time [s]\", \"kappa(t)\"}, PlotRange -> All, ImageSize -> Large];
controlPlot = ListLinePlot[
  Table[Transpose[{result[\"Time\"], result[\"u\"]}], {result, results}],
  PlotLegends -> {\"b = +2\", \"b = -2\"}, Frame -> True,
  FrameLabel -> {\"Time [s]\", \"u(t)\"}, PlotRange -> All, ImageSize -> Large];
lyapunovPlot = ListLogPlot[
  Table[Transpose[{result[\"Time\"],
    Map[Max[#, $MinMachineNumber] &, result[\"V\"]]}], {result, results}],
  PlotLegends -> {\"b = +2\", \"b = -2\"}, Frame -> True,
  FrameLabel -> {\"Time [s]\", \"V(t)\"}, PlotRange -> All, ImageSize -> Large];
GraphicsGrid[{ {statePlot}, {kappaPlot}, {controlPlot}, {lyapunovPlot} }]
\"", "]"}]
  ], "Input"]
},
WindowSize -> {1100, 850},
StyleDefinitions -> "Default.nb"]        
      

12.1 Why \(\kappa\) stops changing

The tuning state increases only while \(x\neq0\). Once the gain enters a region that produces stabilizing feedback for the actual sign of \(b\), the state decays, \(x^2\) becomes integrable, and \(\kappa\) approaches a finite constant. That constant is not an estimate of \(b\).

12.2 Why the final \(\kappa\) differs for opposite signs

The effective coefficient is \(a+bN(\kappa)\). For \(b=+2\), stabilization requires a sufficiently negative value of \(N(\kappa)\). For \(b=-2\), it requires a sufficiently positive value.

12.3 No parameter-convergence claim

The theorem establishes state regulation and boundedness, not convergence of a plant-parameter estimate. There is no persistent excitation requirement because the design is not recovering \(a\) or \(b\).

12.4 Numerical stiffness and step size

During the search phase, \(|N(\kappa)|\) can grow rapidly. A coarse explicit integration step can miss these dynamics or create numerical instability even when the continuous-time theorem is valid.

The result depends on constant nonzero \(b\), exact continuous-time implementation, and an unsaturated actuator. Practical complications include:

1. Large transient control: the polarity search can create high peaks in \(u\).
2. Actuator saturation: saturation destroys the exact Lyapunov identity.
3. Measurement noise: \(\dot{\kappa}=x^2\) remains positive under zero-mean noise and may drift.
4. Time-varying or sign-changing gain: the constant-gain proof does not transfer automatically.
5. Sampling: a digital implementation requires separate analysis.

Practical implementations may use supervisory logic, state and input limits, normalization, dead zones, switching restrictions, or alternative designs with explicit transient bounds. Each modification creates a new closed-loop system requiring its own proof.

Problem 1 (Verification of a Nussbaum function): Prove directly that \(N(\zeta)=\zeta^2\cos\zeta\) satisfies the Nussbaum conditions.

Solution: Integration by parts gives

\[ F(\zeta) = \int_0^\zeta s^2\cos s\,ds = \zeta^2\sin\zeta+2\zeta\cos\zeta-2\sin\zeta. \]

At \(\zeta_n^+=\frac{\pi}{2}+2\pi n\),

\[ \frac{F(\zeta_n^+)}{\zeta_n^+} = \zeta_n^+ - \frac{2}{\zeta_n^+} →+\infty. \]

At \(\zeta_n^-=\frac{3\pi}{2}+2\pi n\),

\[ \frac{F(\zeta_n^-)}{\zeta_n^-} = -\zeta_n^- + \frac{2}{\zeta_n^-} →-\infty. \]

Problem 2 (Why ordinary high gain is insufficient): For \(a=0.8\), determine whether \(u=-5x\) stabilizes the plant for \(b=2\) and \(b=-2\).

Solution: The closed-loop pole is

\[ \lambda_{cl}=a-bk. \]

For \(b=2\), \(\lambda_{cl}=0.8-10=-9.2<0\). For \(b=-2\), \(\lambda_{cl}=0.8+10=10.8>0\). The same controller stabilizes one plant and destabilizes the other.

Problem 3 (Energy interpretation): Show that bounded \(\kappa(t)\) implies \(x\in\mathcal{L}_2\).

Solution:

\[ \int_0^T x^2(t)\,dt = \kappa(T)-\kappa(0). \]

If \(\kappa(T)\leq\bar{\kappa}\), then

\[ \int_0^\infty x^2(t)\,dt \leq \bar{\kappa}-\kappa(0) <\infty. \]

Problem 4 (Initial motion): For \(a=0.8\), \(x(0)=1.2\), and \(\kappa(0)=0\), compute \(u(0)\), \(\dot{x}(0)\), and \(\dot{\kappa}(0)\).

Solution: Since \(N(0)=0\),

\[ u(0)=0, \qquad \dot{x}(0)=0.8(1.2)=0.96, \qquad \dot{\kappa}(0)=1.2^2=1.44. \]

These values are identical for both signs of \(b\) because the initial control is zero.

Problem 5 (Actuator saturation): Suppose \(u_{sat}=\operatorname{sat}_{u_{max} }(N(\kappa)x)\). Identify the proof step that fails.

Solution: With saturation,

\[ \dot{V}=ax^2+b x u_{sat}, \]

which cannot generally be rewritten as \((a+bN(\kappa))\dot{\kappa}\). The exact integral identity is lost, so the contradiction argument proving bounded \(\kappa\) is unavailable without a new design and proof.

The scalar plant \(\dot{x}=ax+bu\) cannot be stabilized by a fixed feedback polarity when the sign of \(b\) is unknown. The Nussbaum-gain law \(u=N(\kappa)x\), \(\dot{\kappa}=x^2\) replaces explicit sign identification with an oscillatory adaptive search. A Lyapunov identity proves boundedness of the tuning state, square integrability of the plant state, and asymptotic regulation. The simulations also expose the principal practical weakness: potentially large state and control transients.