Lesson 1 established the central difficulty: for a plant whose input channel contains an unknown nonzero coefficient, an ordinary feedback gain can accidentally inject energy instead of dissipating it. Consider the scalar plant

\[ \dot{x}(t)=a x(t)+b u(t), \qquad b\neq 0, \]

where \(a\) may be unknown and the sign of \(b\) is unavailable. If one chooses the conventional feedback \(u=-k x\) with \(k>0\), then

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

The feedback is stabilizing only when the product \(bk\) has the correct sign and sufficient magnitude. A Nussbaum design avoids committing permanently to one sign. It replaces the fixed gain by an oscillatory gain whose signed integral explores arbitrarily large positive and negative average values.

flowchart TD
  R["State or tracking error x"] --> A["Adaptive state: kappa_dot = x^2"]
  A --> N["Oscillatory gain N(kappa)"]
  N --> C["Control u = N(kappa) x"]
  C --> P["Plant with unknown constant input sign"]
  P --> R
  A --> I["kappa grows only while error energy persists"]
  I --> N
        

The adaptive variable \(\kappa\) is not, in general, an estimate of \(b\) or of \(\operatorname{sgn}(b)\). Its role is to index a sequence of alternating feedback directions and increasing gain magnitudes. The stability proof uses the integral history of this search, rather than convergence of a parameter estimate.

A continuous function \(N:[0,\infty)\rightarrow\mathbb{R}\) is called a classical Nussbaum-type function when

\[ \limsup_{\zeta\rightarrow\infty} \frac{1}{\zeta}\int_{0}^{\zeta}N(s)\,ds=+\infty \]

and simultaneously

\[ \liminf_{\zeta\rightarrow\infty} \frac{1}{\zeta}\int_{0}^{\zeta}N(s)\,ds=-\infty. \]

Define the normalized primitive

\[ \mathcal{A}_{N}(\zeta) := \frac{1}{\zeta}\int_{0}^{\zeta}N(s)\,ds, \qquad \zeta>0. \]

The definition requires more than sign changes. It requires the average accumulated gain to possess subsequences that diverge to both \(+\infty\) and \(-\infty\). Therefore, a bounded sinusoid such as \(N(\zeta)=\sin\zeta\) is not a Nussbaum function:

\[ \frac{1}{\zeta}\int_{0}^{\zeta}\sin s\,ds = \frac{1-\cos\zeta}{\zeta} \rightarrow 0. \]

Similarly, \(N(\zeta)=\zeta\sin\zeta\) is not sufficient, because

\[ \frac{1}{\zeta}\int_{0}^{\zeta}s\sin s\,ds = -\cos\zeta+\frac{\sin\zeta}{\zeta}, \]

whose limsup and liminf are finite. The oscillation amplitude must grow rapidly enough that successive favorable and unfavorable integral excursions become unbounded after normalization by \(\zeta\).

Frequently used examples include

\[ N_1(\kappa)=\kappa^2\cos\kappa, \qquad N_2(\kappa)=\kappa^2\sin\kappa, \]

\[ N_3(\kappa)=e^{\kappa^2}\cos\left(\frac{\pi\kappa}{2}\right), \qquad N_4(\kappa)=e^{\kappa^2}\sin\kappa. \]

The exponential forms produce especially dominant integral excursions, but they can overflow numerically and generate impractically large transient controls. Polynomial-amplitude choices are often preferable for teaching and numerical experiments, although the selected function must match the theorem used in the intended plant class.

Proposition 1 (scaling and sign reversal). If \(N\) is a Nussbaum function and \(c\neq 0\) is constant, then \(\widetilde{N}(\kappa)=cN(\kappa)\) is also a Nussbaum function.

Proof. Its normalized primitive is

\[ \mathcal{A}_{\widetilde{N}}(\zeta) = c\,\mathcal{A}_{N}(\zeta). \]

For \(c>0\), the positive and negative divergent subsequences are preserved. For \(c<0\), they are interchanged. In either case, the limsup is \(+\infty\) and the liminf is \(-\infty\).

Proposition 2 (bounded additive perturbation). If \(N\) is a Nussbaum function and \(q\) is bounded, then \(\widetilde{N}=N+q\) remains Nussbaum.

Proof. If \(|q(s)|\le M\), then

\[ \left| \frac{1}{\zeta}\int_{0}^{\zeta}q(s)\,ds \right| \le M. \]

Adding a bounded quantity cannot eliminate subsequences diverging to positive and negative infinity.

Integration by parts gives

\[ \int_{0}^{\zeta}s^2\cos s\,ds = \zeta^2\sin\zeta +2\zeta\cos\zeta -2\sin\zeta. \]

Hence

\[ \mathcal{A}_{N}(\zeta) = \zeta\sin\zeta +2\cos\zeta -\frac{2\sin\zeta}{\zeta}. \]

Select the sequence \(\zeta_n=\frac{\pi}{2}+2\pi n\). Since \(\sin\zeta_n=1\) and \(\cos\zeta_n=0\),

\[ \mathcal{A}_{N}(\zeta_n) = \zeta_n-\frac{2}{\zeta_n} \rightarrow +\infty. \]

Select instead \(\eta_n=\frac{3\pi}{2}+2\pi n\). Then \(\sin\eta_n=-1\) and \(\cos\eta_n=0\), so

\[ \mathcal{A}_{N}(\eta_n) = -\eta_n+\frac{2}{\eta_n} \rightarrow -\infty. \]

Both defining conditions hold; therefore \(N(\kappa)=\kappa^2\cos\kappa\) is a Nussbaum-type function.

The following lemma captures the role played by a Nussbaum function in many Lyapunov proofs.

Lemma (constant unknown control direction). Let \(V(t)\ge 0\) and \(\kappa(t)\ge 0\) be continuous, with \(\dot{\kappa}(t)\ge 0\). Suppose that for constants \(c_0,c_1\in\mathbb{R}\) and a constant \(b\neq 0\),

\[ V(t) \le c_0+ \int_{0}^{t} \left[bN(\kappa(\tau))+c_1\right]\dot{\kappa}(\tau)\,d\tau. \]

If \(N\) is Nussbaum, then \(V(t)\) and \(\kappa(t)\) are bounded.

Proof. Because \(\kappa\) is nondecreasing,

\[ V(t) \le C+ \int_{\kappa(0)}^{\kappa(t)} \left[bN(s)+c_1\right]ds \]

for a finite constant \(C\). Assume, for contradiction, that \(\kappa(t)\) is unbounded. Dividing the integral contribution by \(\kappa(t)\), its dominant average is

\[ b\, \frac{1}{\kappa(t)} \int_{0}^{\kappa(t)}N(s)\,ds+c_1. \]

If \(b>0\), the Nussbaum liminf supplies an unbounded sequence on which this average approaches \(-\infty\). If \(b<0\), the Nussbaum limsup supplies such a sequence after multiplication by \(b\). Along either sequence, the right-hand side becomes negative with arbitrarily large magnitude, contradicting \(V(t)\ge 0\). Thus \(\kappa\) is bounded. Substitution into the same inequality then establishes boundedness of \(V\).

flowchart TD
  S["Assume kappa becomes unbounded"] --> B["Use change of variable from time to kappa"]
  B --> A["Average term contains b times integral of N"]
  A --> Q{"Unknown b sign"}
  Q -->|"b positive"| L["Use negative-diverging subsequence"]
  Q -->|"b negative"| U["Use positive-diverging subsequence"]
  L --> C["Upper bound for nonnegative V becomes negative"]
  U --> C
  C --> X["Contradiction: kappa must remain bounded"]
        

Apply the controller

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

to

\[ \dot{x}=ax+bu, \qquad b\neq 0. \]

Let \(V=\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}. \]

Integrating gives

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

The lemma establishes boundedness of \(V\) and \(\kappa\). Because

\[ \kappa(t)-\kappa(0) = \int_{0}^{t}x^2(\tau)\,d\tau, \]

boundedness of \(\kappa\) implies \(x\in L_2\). Also, \(x\), \(\kappa\), \(N(\kappa)\), and \(u=N(\kappa)x\) are bounded, so \(\dot{x}\) is bounded. Therefore \(\frac{d}{dt}x^2=2x\dot{x}\) is bounded. By Barbalat's lemma,

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

The same implemented law works for both signs of \(b\). No explicit sign estimator appears in the controller.

The update \(\dot{\kappa}=x^2\) has three structural properties:

Monotonicity: \(\dot{\kappa}\ge 0\), enabling the change of variables from time to \(\kappa\) in the stability proof.

Error-driven exploration: when \(|x|\) is large, \(\kappa\) moves rapidly across alternating regions of \(N(\kappa)\). When the state is small, the exploration slows.

Energy accounting: \(\kappa(t)-\kappa(0)=\int_0^t x^2d\tau\). Thus boundedness of \(\kappa\) immediately provides finite error energy.

More general designs use

\[ \dot{\kappa} = \gamma\,\psi^2(t), \qquad \gamma>0, \]

where \(\psi\) may be a filtered tracking error, backstepping surface, output error, or another Lyapunov-related signal. The controller commonly has the structure

\[ u=N(\kappa)\,\alpha(t), \]

with \(\alpha(t)\) chosen so that the Lyapunov derivative can be organized into a Nussbaum-integral inequality.

The classical definition is sufficient for a constant unknown nonzero coefficient \(b\). Additional care is required when the input coefficient is time-varying:

\[ \dot{x}=f(x,t)+b(t)u. \]

Even when \(b(t)\) remains bounded away from zero and retains one unknown sign, a general classical Nussbaum function does not automatically satisfy every time-varying-coefficient boundedness lemma. Enhanced classes compare the positive and negative accumulated areas. Define

\[ N^{+}(s)=\max\{0,N(s)\}, \qquad N^{-}(s)=\max\{0,-N(s)\}. \]

A stronger type-B behavior requires both positive and negative normalized areas to diverge and, along suitable subsequences, each area to dominate the opposite area by an arbitrarily large ratio:

\[ \lim_{\zeta\rightarrow\infty} \frac{1}{\zeta}\int_{0}^{\zeta}N^{+}(s)\,ds = \lim_{\zeta\rightarrow\infty} \frac{1}{\zeta}\int_{0}^{\zeta}N^{-}(s)\,ds = +\infty, \]

\[ \limsup_{\zeta\rightarrow\infty} \frac{\int_{0}^{\zeta}N^{+}(s)\,ds} {\int_{0}^{\zeta}N^{-}(s)\,ds} = \limsup_{\zeta\rightarrow\infty} \frac{\int_{0}^{\zeta}N^{-}(s)\,ds} {\int_{0}^{\zeta}N^{+}(s)\,ds} = +\infty. \]

This lesson uses the classical constant-sign theorem. Lesson 3 develops the corresponding stability lemmas in greater depth. A coefficient that actually crosses zero or changes sign arbitrarily is a different and generally harder problem; classical constant-direction arguments cannot simply be reused.

Choose the theorem before the function. Verify whether the plant has a constant unknown sign, a time-varying coefficient with fixed sign, multiple unknown directions, or possible sign reversals.

Prefer analytically verifiable functions. A finite numerical plot of the normalized primitive is a diagnostic, not a proof of its asymptotic limsup and liminf.

Control numerical growth. Exponential functions such as \(e^{\kappa^2}\cos\kappa\) can overflow before the theoretical mechanism is visible. Polynomial growth is safer for instructional simulations.

Use sufficiently accurate integration. The controller may produce rapid transients when the gain changes sign. A coarse fixed step can create a false instability or miss a gain reversal. Adaptive solvers or a small RK4 step are appropriate for the scalar example.

Do not silently saturate in a proof-oriented experiment. Saturating \(u\) changes the closed-loop equations and invalidates the exact Lyapunov identity. Saturation is essential in practice, but it requires a modified analysis, addressed later in the course.

Expect transient excursions. During intervals in which the instantaneous direction is wrong, the state may grow. The Nussbaum argument guarantees that persistent growth forces \(\kappa\) to move into a later interval with a more favorable accumulated gain. It does not promise monotonic decay of the state.

Chapter16_Lesson2.py

This implementation uses NumPy, Matplotlib, and a from-scratch RK4 integrator. It simulates the same controller for \(b=+2\) and \(b=-2\), writes CSV files, and plots the normalized primitive.

"""
Chapter16_Lesson2.py

Nussbaum-gain demonstration for the scalar uncertain plant

    x_dot = a*x + b*u,

where b is constant, nonzero, and its sign is unknown to the controller.

Controller:
    u         = N(kappa)*x
    kappa_dot = x^2
    N(kappa)  = kappa^2*cos(kappa)

The same controller is simulated for positive and negative b.
"""

from __future__ import annotations

import csv
import math
from dataclasses import dataclass
from pathlib import Path

import matplotlib.pyplot as plt
import numpy as np


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


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


def normalized_primitive(zeta: float) -> float:
    """
    Return (1/zeta) * integral_0^zeta s^2 cos(s) ds.

    The exact primitive is
        zeta^2 sin(zeta) + 2 zeta cos(zeta) - 2 sin(zeta).
    """
    if abs(zeta) < 1.0e-12:
        return 0.0
    primitive = (
        zeta * zeta * math.sin(zeta)
        + 2.0 * zeta * math.cos(zeta)
        - 2.0 * math.sin(zeta)
    )
    return primitive / zeta


def rhs(state: np.ndarray, a: float, b: float) -> np.ndarray:
    """Closed-loop right-hand side for [x, kappa]."""
    x, kappa = float(state[0]), float(state[1])
    control = nussbaum(kappa) * x
    return np.array([a * x + b * control, x * x], dtype=float)


def rk4_step(state: np.ndarray, step: float, a: float, b: float) -> np.ndarray:
    """One fixed-step fourth-order Runge-Kutta update."""
    k1 = rhs(state, a, b)
    k2 = rhs(state + 0.5 * step * k1, a, b)
    k3 = rhs(state + 0.5 * step * k2, a, b)
    k4 = rhs(state + step * k3, a, b)
    return state + (step / 6.0) * (k1 + 2.0 * k2 + 2.0 * k3 + k4)


def simulate(b: float, cfg: SimulationConfig) -> dict[str, np.ndarray]:
    """Simulate one unknown-sign case and return all recorded signals."""
    count = int(round(cfg.final_time / cfg.step)) + 1
    time = np.linspace(0.0, cfg.final_time, count)
    state = np.array([cfg.x0, cfg.kappa0], dtype=float)

    x = np.empty(count)
    kappa = np.empty(count)
    gain = np.empty(count)
    control = np.empty(count)
    energy = np.empty(count)

    for index, _ in enumerate(time):
        x[index] = state[0]
        kappa[index] = state[1]
        gain[index] = nussbaum(state[1])
        control[index] = gain[index] * state[0]
        energy[index] = 0.5 * state[0] * state[0]

        if index + 1 < count:
            state = rk4_step(state, cfg.step, cfg.a, b)
            if not np.all(np.isfinite(state)):
                raise FloatingPointError(
                    "The numerical state became non-finite. Reduce the step size."
                )

    return {
        "time": time,
        "x": x,
        "kappa": kappa,
        "N": gain,
        "u": control,
        "V": energy,
        "b": np.full(count, b),
    }


def write_csv(path: Path, result: dict[str, np.ndarray]) -> None:
    """Write one simulation result without requiring pandas."""
    columns = ("time", "b", "x", "kappa", "N", "u", "V")
    with path.open("w", newline="", encoding="utf-8") as stream:
        writer = csv.writer(stream)
        writer.writerow(columns)
        for row in zip(*(result[name] for name in columns)):
            writer.writerow(row)


def main() -> None:
    cfg = SimulationConfig()
    results = {
        "positive_b": simulate(+2.0, cfg),
        "negative_b": simulate(-2.0, cfg),
    }

    output_directory = Path(__file__).resolve().parent
    for label, result in results.items():
        csv_path = output_directory / f"Chapter16_Lesson2_{label}.csv"
        write_csv(csv_path, result)
        print(
            f"{label}: final x={result['x'][-1]: .6e}, "
            f"final kappa={result['kappa'][-1]: .6f}, "
            f"max |x|={np.max(np.abs(result['x'])): .6f}"
        )

    zeta = np.linspace(0.01, 45.0, 6000)
    average = np.array([normalized_primitive(value) for value in zeta])
    plt.figure()
    plt.plot(zeta, average)
    plt.xlabel("zeta")
    plt.ylabel("(1/zeta) integral_0^zeta N(s) ds")
    plt.title("Nussbaum normalized primitive")
    plt.grid(True)
    plt.tight_layout()
    plt.savefig(output_directory / "Chapter16_Lesson2_integral_test.png", dpi=180)
    plt.close()

    plt.figure()
    for label, result in results.items():
        plt.plot(result["time"], result["x"], label=label)
    plt.xlabel("time")
    plt.ylabel("x")
    plt.title("Closed-loop state for both signs of b")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.savefig(output_directory / "Chapter16_Lesson2_state.png", dpi=180)
    plt.close()

    plt.figure()
    for label, result in results.items():
        plt.plot(result["time"], result["kappa"], label=label)
    plt.xlabel("time")
    plt.ylabel("kappa")
    plt.title("Adaptive Nussbaum variable")
    plt.grid(True)
    plt.legend()
    plt.tight_layout()
    plt.savefig(output_directory / "Chapter16_Lesson2_kappa.png", dpi=180)
    plt.close()


if __name__ == "__main__":
    main()

Chapter16_Lesson2.cpp

The C++17 program is dependency-free and writes CSV data for subsequent plotting or verification.

/*
Chapter16_Lesson2.cpp

Dependency-free RK4 simulation of a Nussbaum controller for
    x_dot = a*x + b*u,
with unknown constant sign of b.

Controller:
    u = kappa^2*cos(kappa)*x
    kappa_dot = x^2
*/

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

struct State {
    double x;
    double kappa;
};

struct Derivative {
    double xDot;
    double kappaDot;
};

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

double normalizedPrimitive(const double zeta) {
    if (std::abs(zeta) < 1.0e-12) {
        return 0.0;
    }
    const double primitive =
        zeta * zeta * std::sin(zeta)
        + 2.0 * zeta * std::cos(zeta)
        - 2.0 * std::sin(zeta);
    return primitive / zeta;
}

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

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

State rk4Step(const State& state, const double step, const double a, const double b) {
    const Derivative k1 = rhs(state, a, b);
    const Derivative k2 = rhs(addScaled(state, k1, 0.5 * step), a, b);
    const Derivative k3 = rhs(addScaled(state, k2, 0.5 * step), a, b);
    const Derivative k4 = rhs(addScaled(state, k3, step), a, b);

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

void simulate(const double b, const std::string& label) {
    constexpr double a = 1.0;
    constexpr double step = 1.0e-3;
    constexpr double finalTime = 12.0;
    State state{1.5, 0.0};

    std::ofstream output("Chapter16_Lesson2_" + label + ".csv");
    if (!output) {
        throw std::runtime_error("Could not create the CSV output file.");
    }

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

    double maximumAbsX = std::abs(state.x);
    const int steps = static_cast<int>(std::llround(finalTime / step));

    for (int index = 0; index <= steps; ++index) {
        const double time = index * step;
        const double gain = nussbaum(state.kappa);
        const double control = gain * state.x;
        const double energy = 0.5 * state.x * state.x;

        output << time << ',' << b << ',' << state.x << ','
               << state.kappa << ',' << gain << ','
               << control << ',' << energy << '\n';

        maximumAbsX = std::max(maximumAbsX, std::abs(state.x));

        if (index < steps) {
            state = rk4Step(state, step, a, b);
            if (!std::isfinite(state.x) || !std::isfinite(state.kappa)) {
                throw std::runtime_error(
                    "Numerical overflow detected; reduce the integration step."
                );
            }
        }
    }

    std::cout << label
              << ": final x=" << state.x
              << ", final kappa=" << state.kappa
              << ", max |x|=" << maximumAbsX << '\n';
}

int main() {
    try {
        simulate(+2.0, "positive_b");
        simulate(-2.0, "negative_b");

        std::ofstream diagnostic("Chapter16_Lesson2_integral_test.csv");
        diagnostic << "zeta,normalized_primitive\n";
        diagnostic << std::setprecision(16);
        for (double zeta = 0.01; zeta <= 45.0; zeta += 0.01) {
            diagnostic << zeta << ',' << normalizedPrimitive(zeta) << '\n';
        }
    } catch (const std::exception& error) {
        std::cerr << "Error: " << error.what() << '\n';
        return 1;
    }
    return 0;
}

Chapter16_Lesson2.java

The Java implementation uses records for the state and derivative and requires Java 17 or later.

/*
Chapter16_Lesson2.java

Dependency-free RK4 simulation of a Nussbaum controller for
    x_dot = a*x + b*u,
where the constant input gain b may be positive or negative.
*/

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_Lesson2 {
    private static final double A = 1.0;
    private static final double STEP = 1.0e-3;
    private static final double FINAL_TIME = 12.0;

    private Chapter16_Lesson2() {
    }

    private record State(double x, double kappa) {
    }

    private record Derivative(double xDot, double kappaDot) {
    }

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

    private static double normalizedPrimitive(double zeta) {
        if (Math.abs(zeta) < 1.0e-12) {
            return 0.0;
        }
        double primitive =
            zeta * zeta * Math.sin(zeta)
            + 2.0 * zeta * Math.cos(zeta)
            - 2.0 * Math.sin(zeta);
        return primitive / zeta;
    }

    private static Derivative rhs(State state, double b) {
        double control = nussbaum(state.kappa()) * state.x();
        return new Derivative(
            A * state.x() + b * control,
            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 b) {
        Derivative k1 = rhs(state, b);
        Derivative k2 = rhs(addScaled(state, k1, 0.5 * STEP), b);
        Derivative k3 = rhs(addScaled(state, k2, 0.5 * STEP), b);
        Derivative k4 = rhs(addScaled(state, k3, STEP), b);

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

    private static void simulate(double b, String label) throws IOException {
        State state = new State(1.5, 0.0);
        int steps = (int) Math.round(FINAL_TIME / STEP);
        double maximumAbsX = Math.abs(state.x());

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

            for (int index = 0; index <= steps; index++) {
                double time = index * STEP;
                double gain = nussbaum(state.kappa());
                double control = gain * state.x();
                double energy = 0.5 * state.x() * state.x();

                writer.write(String.format(
                    Locale.ROOT,
                    "%.12g,%.12g,%.12g,%.12g,%.12g,%.12g,%.12g%n",
                    time, b, state.x(), state.kappa(), gain, control, energy
                ));

                maximumAbsX = Math.max(maximumAbsX, Math.abs(state.x()));

                if (index < steps) {
                    state = rk4Step(state, b);
                    if (!Double.isFinite(state.x())
                        || !Double.isFinite(state.kappa())) {
                        throw new ArithmeticException(
                            "Numerical overflow detected; reduce STEP."
                        );
                    }
                }
            }
        }

        System.out.printf(
            Locale.ROOT,
            "%s: final x=% .6e, final kappa=% .6f, max |x|=% .6f%n",
            label, state.x(), state.kappa(), maximumAbsX
        );
    }

    private static void writeIntegralDiagnostic() throws IOException {
        Path path = Path.of("Chapter16_Lesson2_integral_test.csv");
        try (BufferedWriter writer = Files.newBufferedWriter(
            path,
            StandardCharsets.UTF_8
        )) {
            writer.write("zeta,normalized_primitive");
            writer.newLine();
            for (double zeta = 0.01; zeta <= 45.0; zeta += 0.01) {
                writer.write(String.format(
                    Locale.ROOT,
                    "%.12g,%.12g%n",
                    zeta,
                    normalizedPrimitive(zeta)
                ));
            }
        }
    }

    public static void main(String[] args) {
        try {
            simulate(+2.0, "positive_b");
            simulate(-2.0, "negative_b");
            writeIntegralDiagnostic();
        } catch (IOException | ArithmeticException exception) {
            System.err.println("Error: " + exception.getMessage());
            System.exit(1);
        }
    }
}

Chapter16_Lesson2.m

The script uses ode45 for the numerical study. Its final function can also construct a Simulink block diagram programmatically; set BUILD_SIMULINK_MODEL = true when Simulink is installed.

% Chapter16_Lesson2.m
%
% Nussbaum-gain demonstration for
%     x_dot = a*x + b*u,
% where b is constant and nonzero but its sign is not supplied to the
% controller.
%
% Controller:
%     u         = N(kappa)*x
%     kappa_dot = x^2
%     N(kappa)  = kappa^2*cos(kappa)

clear;
clc;
close all;

a = 1.0;
x0 = 1.5;
kappa0 = 0.0;
tspan = [0.0, 12.0];
bCases = [2.0, -2.0];
labels = ["positive b", "negative b"];

results = cell(numel(bCases), 1);

for caseIndex = 1:numel(bCases)
    b = bCases(caseIndex);
    initialState = [x0; kappa0];

    options = odeset( ...
        'RelTol', 1.0e-8, ...
        'AbsTol', 1.0e-10, ...
        'MaxStep', 1.0e-2);

    [time, state] = ode45( ...
        @(timeValue, stateValue) closedLoopRhs( ...
            timeValue, stateValue, a, b), ...
        tspan, initialState, options);

    x = state(:, 1);
    kappa = state(:, 2);
    N = arrayfun(@nussbaum, kappa);
    u = N .* x;
    V = 0.5 .* x.^2;

    results{caseIndex} = table( ...
        time, repmat(b, size(time)), x, kappa, N, u, V, ...
        'VariableNames', {'time', 'b', 'x', 'kappa', 'N', 'u', 'V'});

    fileLabel = strrep(labels(caseIndex), " ", "_");
    writetable( ...
        results{caseIndex}, ...
        "Chapter16_Lesson2_" + fileLabel + ".csv");

    fprintf( ...
        '%s: final x=% .6e, final kappa=% .6f, max |x|=% .6f\n', ...
        labels(caseIndex), x(end), kappa(end), max(abs(x)));
end

zeta = linspace(0.01, 45.0, 6000).';
averagePrimitive = arrayfun(@normalizedPrimitive, zeta);

figure;
plot(zeta, averagePrimitive, 'LineWidth', 1.2);
grid on;
xlabel('zeta');
ylabel('(1/zeta) integral_0^zeta N(s) ds');
title('Nussbaum normalized primitive');

figure;
hold on;
for caseIndex = 1:numel(results)
    plot( ...
        results{caseIndex}.time, ...
        results{caseIndex}.x, ...
        'DisplayName', labels(caseIndex));
end
grid on;
xlabel('time');
ylabel('x');
title('Closed-loop state for both signs of b');
legend('Location', 'best');
hold off;

figure;
hold on;
for caseIndex = 1:numel(results)
    plot( ...
        results{caseIndex}.time, ...
        results{caseIndex}.kappa, ...
        'DisplayName', labels(caseIndex));
end
grid on;
xlabel('time');
ylabel('kappa');
title('Adaptive Nussbaum variable');
legend('Location', 'best');
hold off;

% Set this flag to true when Simulink is installed and a block-diagram
% realization is desired.
BUILD_SIMULINK_MODEL = false;
if BUILD_SIMULINK_MODEL
    buildNussbaumSimulinkModel();
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 value = normalizedPrimitive(zeta)
    if abs(zeta) < 1.0e-12
        value = 0.0;
        return;
    end
    primitive = ...
        zeta.^2 .* sin(zeta) ...
        + 2.0 .* zeta .* cos(zeta) ...
        - 2.0 .* sin(zeta);
    value = primitive ./ zeta;
end

function buildNussbaumSimulinkModel()
    model = 'Chapter16_Lesson2_Simulink';
    if bdIsLoaded(model)
        close_system(model, 0);
    end
    if isfile(model + ".slx")
        delete(model + ".slx");
    end

    assignin('base', 'a', 1.0);
    assignin('base', 'b', -2.0);

    new_system(model);
    open_system(model);

    add_block( ...
        'simulink/Continuous/Integrator', model + "/x", ...
        'InitialCondition', '1.5', ...
        'Position', [520 110 550 140]);
    add_block( ...
        'simulink/Continuous/Integrator', model + "/kappa", ...
        'InitialCondition', '0', ...
        'Position', [260 280 290 310]);
    add_block( ...
        'simulink/Math Operations/Math Function', model + "/x_squared", ...
        'Operator', 'square', ...
        'Position', [120 275 165 315]);
    add_block( ...
        'simulink/Math Operations/Math Function', model + "/kappa_squared", ...
        'Operator', 'square', ...
        'Position', [330 245 380 285]);
    add_block( ...
        'simulink/Math Operations/Trigonometric Function', model + "/cos_kappa", ...
        'Operator', 'cos', ...
        'Position', [330 325 380 365]);
    add_block( ...
        'simulink/Math Operations/Product', model + "/Nussbaum_gain", ...
        'Position', [420 265 455 325]);
    add_block( ...
        'simulink/Math Operations/Product', model + "/control_u", ...
        'Position', [510 260 545 320]);
    add_block( ...
        'simulink/Math Operations/Gain', model + "/a_gain", ...
        'Gain', 'a', ...
        'Position', [330 75 380 105]);
    add_block( ...
        'simulink/Math Operations/Gain', model + "/b_gain", ...
        'Gain', 'b', ...
        'Position', [600 260 650 290]);
    add_block( ...
        'simulink/Math Operations/Sum', model + "/plant_sum", ...
        'Inputs', '++', ...
        'Position', [450 105 475 155]);
    add_block( ...
        'simulink/Signal Routing/Mux', model + "/scope_mux", ...
        'Inputs', '2', ...
        'Position', [690 115 695 185]);
    add_block( ...
        'simulink/Sinks/Scope', model + "/scope", ...
        'Position', [750 125 790 175]);

    add_line(model, 'x/1', 'x_squared/1', 'autorouting', 'on');
    add_line(model, 'x_squared/1', 'kappa/1', 'autorouting', 'on');
    add_line(model, 'kappa/1', 'kappa_squared/1', 'autorouting', 'on');
    add_line(model, 'kappa/1', 'cos_kappa/1', 'autorouting', 'on');
    add_line(model, 'kappa_squared/1', 'Nussbaum_gain/1', 'autorouting', 'on');
    add_line(model, 'cos_kappa/1', 'Nussbaum_gain/2', 'autorouting', 'on');
    add_line(model, 'Nussbaum_gain/1', 'control_u/1', 'autorouting', 'on');
    add_line(model, 'x/1', 'control_u/2', 'autorouting', 'on');
    add_line(model, 'control_u/1', 'b_gain/1', 'autorouting', 'on');
    add_line(model, 'x/1', 'a_gain/1', 'autorouting', 'on');
    add_line(model, 'a_gain/1', 'plant_sum/1', 'autorouting', 'on');
    add_line(model, 'b_gain/1', 'plant_sum/2', 'autorouting', 'on');
    add_line(model, 'plant_sum/1', 'x/1', 'autorouting', 'on');
    add_line(model, 'x/1', 'scope_mux/1', 'autorouting', 'on');
    add_line(model, 'kappa/1', 'scope_mux/2', 'autorouting', 'on');
    add_line(model, 'scope_mux/1', 'scope/1', 'autorouting', 'on');

    set_param(model, 'StopTime', '12');
    save_system(model);
    fprintf('Created %s.slx with b=-2 in the base workspace.\n', model);
end

Chapter16_Lesson2.nb

The notebook solves both control-direction cases with NDSolveValue, plots the state and adaptive variable, and exports sampled trajectories.

Notebook[{
  Cell["Chapter16_Lesson2.nb", "Title"], Cell["Nussbaum-Type Functions: Concept and Role", "Subtitle"],
  Cell["ClearAll[nussbaum, normalizedPrimitive, rhs, simulate];
nussbaum[kappa_?NumericQ] := kappa^2 Cos[kappa];
normalizedPrimitive[zeta_?NumericQ] := If[
  Abs[zeta] < 10^-12, 0.0,
  (zeta^2 Sin[zeta] + 2 zeta Cos[zeta] - 2 Sin[zeta])/zeta
];
rhs[a_?NumericQ, b_?NumericQ][{x_?NumericQ, kappa_?NumericQ}] := {
  a x + b nussbaum[kappa] x,
  x^2
};
simulate[b_?NumericQ, tFinal_: 12.0] := Module[
  {a = 1.0, x0 = 1.5, kappa0 = 0.0, solution},
  solution = NDSolveValue[
    {x'[t] == a x[t] + b nussbaum[kappa[t]] x[t],
     kappa'[t] == x[t]^2, x[0] == x0, kappa[0] == kappa0},
    {x, kappa}, {t, 0, tFinal},
    Method -> {\"StiffnessSwitching\"}, MaxStepFraction -> 1/2000
  ];
  <|\"b\" -> b, \"x\" -> solution[[1]], \"kappa\" -> solution[[2]]|>
];", "Input"],
  Cell["positiveCase = simulate[2.0];
negativeCase = simulate[-2.0];
statePlot = Plot[
  {positiveCase[\"x\"][t], negativeCase[\"x\"][t]}, {t, 0, 12},
  PlotLegends -> {\"b = +2\", \"b = -2\"}, AxesLabel -> {\"time\", \"x\"},
  PlotLabel -> \"Closed-loop state for both control directions\",
  PlotRange -> All
];
kappaPlot = Plot[
  {positiveCase[\"kappa\"][t], negativeCase[\"kappa\"][t]}, {t, 0, 12},
  PlotLegends -> {\"b = +2\", \"b = -2\"}, AxesLabel -> {\"time\", \"kappa\"},
  PlotLabel -> \"Adaptive Nussbaum variable\", PlotRange -> All
];
integralPlot = Plot[
  normalizedPrimitive[zeta], {zeta, 0.01, 45},
  AxesLabel -> {\"zeta\", \"(1/zeta) integral N\"},
  PlotLabel -> \"Nussbaum normalized primitive\", PlotRange -> All
];
Column[{statePlot, kappaPlot, integralPlot}]", "Input"],
  Cell["sampleTimes = Subdivide[0.0, 12.0, 12000]; exportCase[result_, label_String] := Module[
  {rows},
  rows = Table[
    With[
      {xValue = result[\"x\"][time],
       kappaValue = result[\"kappa\"][time], bValue = result[\"b\"]},
      {time, bValue, xValue, kappaValue, nussbaum[kappaValue],
       nussbaum[kappaValue] xValue, 0.5 xValue^2}
    ],
    {time, sampleTimes}
  ];
  Export[
    \"Chapter16_Lesson2_\" <> label <> \".csv\",
    Prepend[rows, {\"time\", \"b\", \"x\", \"kappa\", \"N\", \"u\", \"V\"}]
  ]
];
exportCase[positiveCase, \"positive_b\"]; exportCase[negativeCase, \"negative_b\"];", "Input"]
}, WindowTitle -> "Chapter16_Lesson2.nb"
]        
      

Problem 1 (Rejecting a false candidate): Determine whether \(N(\kappa)=\kappa\cos\kappa\) is a Nussbaum function.

Solution: Integration gives

\[ \int_{0}^{\zeta}s\cos s\,ds = \zeta\sin\zeta+\cos\zeta-1. \]

Therefore

\[ \mathcal{A}_{N}(\zeta) = \sin\zeta+\frac{\cos\zeta-1}{\zeta}. \]

Its limsup is \(1\) and its liminf is \(-1\), not positive and negative infinity. It is not a Nussbaum function.

Problem 2 (Phase shift): Prove that \(N_{\phi}(\kappa)=\kappa^2\cos(\kappa+\phi)\) is Nussbaum for every fixed \(\phi\in\mathbb{R}\).

Solution: Expand

\[ N_{\phi}(\kappa) = \cos\phi\,\kappa^2\cos\kappa - \sin\phi\,\kappa^2\sin\kappa. \]

More directly, choose sequences on which \(\sin(\zeta+\phi)=1\) and \(\sin(\zeta+\phi)=-1\). Integration by parts shows the normalized primitive has the leading term \(\zeta\sin(\zeta+\phi)\), while all remaining terms are bounded. Thus its limsup and liminf are \(+\infty\) and \(-\infty\).

Problem 3 (Unknown-sign scalar plant): For \(\dot{x}=ax+bu\), apply \(u=N(\kappa)x\) and \(\dot{\kappa}=x^2\). Show that finite-time escape is impossible.

Solution: The identity

\[ \frac{1}{2}x^2(t) = \frac{1}{2}x^2(0) + \int_{\kappa(0)}^{\kappa(t)} \left[a+bN(s)\right]ds \]

has the form required by the Nussbaum boundedness lemma. Hence \(x\) and \(\kappa\) are bounded on every maximal solution interval. Since the closed-loop vector field is locally Lipschitz and the state remains bounded, the continuation theorem extends the solution for all \(t\ge 0\); finite-time escape is impossible.

Problem 4 (Why monotonicity matters): Explain why the update \(\dot{\kappa}=x\) does not directly support the same proof.

Solution: The sign of \(\dot{\kappa}\) would change with \(x\), so \(\kappa(t)\) would not be monotone. The substitution \(\int N(\kappa)\dot{\kappa}\,dt=\int N(s)\,ds\) remains locally valid, but an unbounded-time trajectory could repeatedly move forward and backward through the same gain intervals. The standard contradiction based on an unbounded monotone upper limit would no longer follow. Choosing a nonnegative update such as \(x^2\) removes this ambiguity.

Problem 5 (Parameter estimate or search variable?): Suppose the scalar simulation converges with \(\kappa(t)\rightarrow 6.04\). Can \(6.04\) be interpreted as an estimate of \(b\)?

Solution: No. The update law contains no prediction error or regression model that would identify \(b\). Different initial conditions, plant coefficients, solver settings, and Nussbaum functions can yield different terminal values of \(\kappa\) while stabilizing the same plant. Its terminal value records where the sign-search process stopped after the error energy became small.

Problem 6 (Time-varying coefficient): Let \(b(t)\in[1,2]\). Why is it unsafe to assert that every classical Nussbaum function automatically yields the same lemma as for a constant \(b\)?

Solution: With constant \(b\), it can be taken outside the integral:

\[ \int_{0}^{t}bN(\kappa)\dot{\kappa}\,d\tau = b\int_{\kappa(0)}^{\kappa(t)}N(s)\,ds. \]

For \(b(t)\), the coefficient weights different positive and negative gain intervals unequally and cannot generally be removed from the integral. The classical cancellation geometry may be distorted. Stronger area-ratio properties, explicit coefficient bounds, or a specialized theorem are needed.

A Nussbaum-type function is an oscillatory gain whose normalized primitive has unbounded positive and negative excursions. Coupled with a monotone error-energy variable, this property converts an unknown constant control direction into a contradiction argument that bounds the adaptive variable and the Lyapunov function. The mechanism is not sign estimation: it is an error-driven search across alternating gain directions. Classical Nussbaum functions address constant unknown signs; time-varying, multivariable, or sign-switching coefficients require stronger definitions and lemmas.