Skip to content

Commit 7a21cbb

Browse files
feat: add model reference adaptive control (#223)
* add mrac * Apply suggestions from code review Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com> --------- Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>
1 parent 44c7da9 commit 7a21cbb

12 files changed

Lines changed: 507 additions & 494 deletions

File tree

‎README.md‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -26,7 +26,7 @@ Refer to the documentation to quickly integrate and utilize the library's signal
2626
| [Regularization](doc/regularization/README.md) | L1 (Lasso), L2 (Ridge) |
2727
| [Math](doc/math/README.md) | CORDIC, Quaternion, MatrixNorms, Step Response Metrics, MatrixExponential |
2828
| [Solvers](doc/solvers/README.md) | Gaussian Elimination, Levinson-Durbin, Durand-Kerner, Cholesky, DARE, Runge-Kutta ODE Integrators (RK4 + Dormand-Prince), Spectral Radius & Discrete Stability Margin, QR Decomposition (Householder / Givens), LU Decomposition with Partial Pivoting, Singular Value Decomposition (Golub-Kahan) |
29-
| [Nonlinear Control](doc/nonlinear_control/README.md) | Feedback Linearization, Backstepping Control |
29+
| [Nonlinear Control](doc/nonlinear_control/README.md) | Feedback Linearization, Backstepping Control, Model Reference Adaptive Control (MRAC) |
3030
| [Robust Control](doc/robust_control/README.md) | Active Disturbance Rejection Control (ADRC + ESO), Sliding Mode Control (SMC), Disturbance Observer (DOB), H∞ State-Feedback Control |
3131
| [Performance Optimization](doc/performance-optimization/README.md) | Compiler optimizations, SIMD |
3232

‎ROADMAP.md‎

Lines changed: 0 additions & 291 deletions
Large diffs are not rendered by default.
Lines changed: 130 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,130 @@
1+
# Model Reference Adaptive Control
2+
3+
## Overview & Motivation
4+
5+
Real hardware changes over time: motor resistance drifts with temperature, payload mass varies, actuators age. A fixed controller tuned at the factory cannot maintain performance across this variability. Model Reference Adaptive Control (MRAC) addresses this by continuously adjusting its own gains online, using only the tracking error between the real plant and a designer-specified reference model, until the plant's response matches the desired ideal. No knowledge of the exact plant parameters is needed — only the structure (the equations of motion) and the sign of the input gain.
6+
7+
This makes MRAC the canonical direct-adaptive scheme for the "known structure, unknown numbers" class of problems that appears throughout embedded control: a family of units sharing the same firmware but differing in parameter values, or a single unit whose parameters change slowly during operation.
8+
9+
## Mathematical Theory
10+
11+
### Problem Setup
12+
13+
Given a plant of the form
14+
15+
$$\dot{x} = a\,x + b\,u$$
16+
17+
where $a$ and $b$ are unknown constants with $b \neq 0$ and $\text{sgn}(b)$ known, the goal is to choose the control $u$ so that $x(t) \to x_m(t)$ as $t \to \infty$, where $x_m$ satisfies the reference model
18+
19+
$$\dot{x}_m = -a_m\,x_m + b_m\,r, \quad a_m > 0.$$
20+
21+
### Matching Conditions
22+
23+
The ideal control law that would make the plant identical to the reference model is
24+
25+
$$u^* = \theta_x^*\,x + \theta_r^*\,r$$
26+
27+
where
28+
29+
$$\theta_x^* = \frac{-(a_m + a)}{b}, \qquad \theta_r^* = \frac{b_m}{b}.$$
30+
31+
Since $a$ and $b$ are unknown, $\theta_x^*$ and $\theta_r^*$ cannot be computed directly; instead they are estimated online.
32+
33+
### Tracking Error Dynamics
34+
35+
Let $\hat{\theta}_x$ and $\hat{\theta}_r$ denote the current parameter estimates and define the parameter errors $\tilde{\theta} = \hat{\theta} - \theta^*$. The tracking error $e = x - x_m$ satisfies
36+
37+
$$\dot{e} = -a_m\,e + b\,(\tilde{\theta}_x\,x + \tilde{\theta}_r\,r).$$
38+
39+
### MIT Rule
40+
41+
The MIT rule minimises $J = \tfrac{1}{2}e^2$ by gradient descent on the parameter space:
42+
43+
$$\dot{\hat{\theta}}_x = -\gamma\,\text{sgn}(b)\,e\,x, \qquad \dot{\hat{\theta}}_r = -\gamma\,\text{sgn}(b)\,e\,r.$$
44+
45+
This is a first-order gradient update. It converges when $\gamma$ is small relative to the signal levels and the reference is persistently exciting. It has no global stability proof.
46+
47+
### Lyapunov Redesign
48+
49+
The Lyapunov method chooses the same update law but derives it from the Lyapunov function candidate
50+
51+
$$V(e, \tilde{\theta}) = \frac{e^2}{2} + \frac{|b|}{2\gamma}\left(\tilde{\theta}_x^2 + \tilde{\theta}_r^2\right).$$
52+
53+
Computing $\dot{V}$ and choosing the parameter update laws to make $\dot{V} \leq -a_m\,e^2 \leq 0$ yields precisely the same gradient form. The Lyapunov construction guarantees that $e$ and $\tilde{\theta}$ remain bounded for all $t \geq 0$ and that $e(t) \to 0$, even under large initial errors — a global stability guarantee that the pure MIT rule lacks.
54+
55+
### Discrete-Time Implementation
56+
57+
The continuous-time update is approximated by forward Euler with step $\Delta t$:
58+
59+
$$x_m[k+1] = x_m[k] + (A_m\,x_m[k] + B_m\,r[k])\,\Delta t$$
60+
61+
$$e[k] = x[k] - x_m[k+1]$$
62+
63+
$$u[k] = \hat{\theta}_x[k]\,x[k] + \hat{\theta}_r[k]\,r[k]$$
64+
65+
$$\hat{\theta}_x[k+1] = \hat{\theta}_x[k] - \gamma\,\text{sgn}(b)\,e[k]\,x[k]\,\Delta t$$
66+
67+
$$\hat{\theta}_r[k+1] = \hat{\theta}_r[k] - \gamma\,\text{sgn}(b)\,e[k]\,r[k]\,\Delta t$$
68+
69+
The same structure extends to the multi-input multi-output case using outer products of the error and regressor vectors, yielding matrix parameter estimates.
70+
71+
## Complexity Analysis
72+
73+
| Operation | Time | Space | Notes |
74+
|----------------|----------------------|----------------------------------------|------------------------------------------------------------------------|
75+
| ComputeControl | $O(n^2 + nm)$ | $O(1)$ working registers | Reference model step + two outer products + two matrix-vector products |
76+
| Memory | $O(n^2 + nm)$ static | Parameter matrices and reference state | All fixed-size; no heap allocation |
77+
78+
Here $n$ = StateSize and $m$ = InputSize. The dominant cost is the outer product update of the parameter matrices at each step.
79+
80+
## Step-by-Step Walkthrough
81+
82+
**First-order scalar example** ($n = m = 1$, $a = -1.5$, $b = 2$, $a_m = 1$, $b_m = 1$, $\gamma = 1$, $\Delta t = 0.01$):
83+
84+
1. At step $k=0$: $x = 2$, $r = 1$, $x_m = 0$, $\hat{\theta}_x = 0$, $\hat{\theta}_r = 0$.
85+
2. Advance reference model: $x_m \leftarrow 0 + (-1 \cdot 0 + 1 \cdot 1) \cdot 0.01 = 0.01$.
86+
3. Tracking error: $e = 2 - 0.01 = 1.99$.
87+
4. Control output: $u = 0 \cdot 2 + 0 \cdot 1 = 0$ (initial parameters are zero).
88+
5. Update $\hat{\theta}_x$: $0 - 1 \cdot (+1) \cdot 1.99 \cdot 2 \cdot 0.01 = -0.0398$.
89+
6. Update $\hat{\theta}_r$: $0 - 1 \cdot (+1) \cdot 1.99 \cdot 1 \cdot 0.01 = -0.0199$.
90+
7. At the next step the control becomes $u = -0.0398 \cdot x - 0.0199 \cdot r$, already driving the plant toward the reference model. After many steps, $\hat{\theta}_x \to \theta_x^* = 0.25$ and $\hat{\theta}_r \to \theta_r^* = 0.5$.
91+
92+
## Pitfalls & Edge Cases
93+
94+
- **Adaptation gain too large**: with the MIT rule, $\gamma$ large relative to signal power causes $\hat{\theta}$ to overshoot and the closed-loop to go unstable. The Lyapunov law is more forgiving but still requires reasonable $\gamma$.
95+
- **Wrong sign of $b$**: if $\text{sgn}(b)$ is set incorrectly, adaptation drives parameters in the wrong direction and the error grows. This is a hard fault — the algorithm is designed around knowing the plant's sign.
96+
- **No persistent excitation**: if $r$ is constant or bandlimited, $\hat{\theta}$ converges to some values that achieve tracking but not necessarily to $\theta^*$. The plant still tracks the reference model; only parameter identification fails.
97+
- **Parameter drift**: in the presence of noise or unmodelled disturbances, $\hat{\theta}$ drifts even when $e$ is small. Standard remedies are $\sigma$-modification ($\dot{\hat{\theta}} = -\gamma\,e\,\phi - \sigma\,\hat{\theta}$) and parameter projection onto a compact set.
98+
- **Euler discretisation error**: the Euler step introduces $O(\Delta t)$ error in the reference model trajectory and the parameter update. Small $\Delta t$ is needed for accuracy; a higher-order integrator (RK4) can be used for the reference model when $\Delta t$ is large.
99+
100+
## Variants & Generalizations
101+
102+
- **Indirect MRAC**: first identifies the plant parameters (via Recursive Least Squares) and then recomputes the matching gains; separates identification from control but requires a plant model structure assumption.
103+
- **$\sigma$-modification**: adds a leakage term $-\sigma\,\hat{\theta}$ to the parameter update, bounding parameter drift at the cost of a small steady-state bias.
104+
- **e-modification**: replaces $\sigma$ with $\sigma\,|e|\,\hat{\theta}$, so leakage is active only when the error is large and vanishes at steady state.
105+
- **Parameter projection**: constrains $\hat{\theta}$ to a known compact set, preventing unbounded drift without degrading tracking.
106+
- **Adaptive backstepping**: embeds MRAC-style parameter adaptation inside the Backstepping recursive design for strict-feedback plants with unknown parameters.
107+
- **MRAC with reference model order > 1**: the same gradient law extends to any LTI reference model by tracking the full state vector $x_m \in \mathbb{R}^n$.
108+
109+
## Applications
110+
111+
- Motor drives with varying load inertia or resistance: MRAC maintains speed/position bandwidth despite parameter changes.
112+
- Aerospace: autopilot adaptation to changing mass, fuel consumption, or aerodynamic coefficients.
113+
- Robotic manipulators: payload-varying adaptive torque controllers.
114+
- Power electronics: adaptive current controllers for converters with uncertain filter inductance.
115+
- Embedded firmware for a product family: a single binary adapts to unit-to-unit hardware variation at startup.
116+
117+
## Connections to Other Algorithms
118+
119+
- **Backstepping Control**: MRAC handles unknown parameters in a fixed-structure plant; adaptive backstepping merges both, providing recursive stability proofs for uncertain strict-feedback systems.
120+
- **Feedback Linearization**: cancels nonlinearities via an explicit model; MRAC is model-free in the parameter sense and is more robust when the model is uncertain.
121+
- **Recursive Least Squares (RLS)**: the indirect-MRAC identification step; direct MRAC avoids explicit parameter estimation by updating control gains rather than plant parameters.
122+
- **LQR**: optimal fixed-gain control for known linear plants; MRAC extends the design to plants with unknown parameters at the cost of the adaptation transient.
123+
- **math::LinearTimeInvariant**: provides the reference model state-space container; the MRAC controller holds a reference to an LTI instance to step the reference trajectory.
124+
125+
## References & Further Reading
126+
127+
- K. J. Åström, B. Wittenmark, *Adaptive Control*, 2nd ed., Addison-Wesley, 1995.
128+
- K. S. Narendra, A. M. Annaswamy, *Stable Adaptive Systems*, Prentice-Hall, 1989; Dover reprint 2005.
129+
- S. Sastry, M. Bodson, *Adaptive Control: Stability, Convergence and Robustness*, Prentice-Hall, 1989.
130+
- P. A. Ioannou, J. Sun, *Robust Adaptive Control*, Prentice-Hall, 1996 (free PDF available from the authors).

‎doc/nonlinear_control/README.md‎

Lines changed: 3 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -6,5 +6,6 @@ Algorithms for nonlinear control design: controllers that exploit or cancel plan
66

77
| Algorithm | Description |
88
|----------------------------------------------------|---------------------------------------------------------------------------------------------------------------------------------------------------------------------------------|
9-
| [Backstepping Control](BacksteppingControl.md) | Lyapunov-based recursive design for strict-feedback nonlinear cascades that stabilises each integrator stage in sequence, yielding a provably stable controller by construction |
10-
| [Feedback Linearization](FeedbackLinearization.md) | Cancels a control-affine plant's known nonlinear dynamics via an inner control law, leaving decoupled integrator chains that a simple outer PD/LQR loop drives |
9+
| [Backstepping Control](BacksteppingControl.md) | Lyapunov-based recursive design for strict-feedback nonlinear cascades that stabilises each integrator stage in sequence, yielding a provably stable controller by construction |
10+
| [Feedback Linearization](FeedbackLinearization.md) | Cancels a control-affine plant's known nonlinear dynamics via an inner control law, leaving decoupled integrator chains that a simple outer PD/LQR loop drives |
11+
| [Model Reference Adaptive Control](ModelReferenceAdaptiveControl.md) | Online gradient-descent adaptation of control gains to match plant response to a stable reference model, with Lyapunov stability guarantee for bounded error under unknown plant parameters |

‎numerical/nonlinear_control/CMakeLists.txt‎

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -13,11 +13,13 @@ target_link_libraries(numerical.nonlinear_control ${NUMERICAL_VISIBILITY}
1313
target_sources(numerical.nonlinear_control PRIVATE
1414
BacksteppingControl.hpp
1515
FeedbackLinearization.hpp
16+
ModelReferenceAdaptiveControl.hpp
1617
)
1718

1819
numerical_add_coverage_sources(numerical.nonlinear_control
1920
BacksteppingControl.cpp
2021
FeedbackLinearization.cpp
22+
ModelReferenceAdaptiveControl.cpp
2123
)
2224

2325
add_subdirectory(test)
Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,6 @@
1+
#include "numerical/nonlinear_control/ModelReferenceAdaptiveControl.hpp"
2+
3+
namespace nonlinear_control
4+
{
5+
template class ModelReferenceAdaptiveControl<float, 1, 1>;
6+
}
Lines changed: 113 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,113 @@
1+
#pragma once
2+
3+
#if defined(__GNUC__) || defined(__clang__)
4+
#pragma GCC optimize("O3", "fast-math")
5+
#endif
6+
7+
#include "numerical/math/CompilerOptimizations.hpp"
8+
#include "numerical/math/LinearTimeInvariant.hpp"
9+
#include "numerical/math/Matrix.hpp"
10+
#include <cstddef>
11+
#include <type_traits>
12+
13+
namespace nonlinear_control
14+
{
15+
enum class AdaptationLaw
16+
{
17+
MitRule,
18+
Lyapunov
19+
};
20+
21+
template<typename T, std::size_t StateSize, std::size_t InputSize>
22+
class ModelReferenceAdaptiveControl
23+
{
24+
static_assert(std::is_floating_point_v<T>,
25+
"ModelReferenceAdaptiveControl supports floating-point types");
26+
static_assert(StateSize > 0, "ModelReferenceAdaptiveControl requires StateSize > 0");
27+
static_assert(InputSize > 0, "ModelReferenceAdaptiveControl requires InputSize > 0");
28+
29+
public:
30+
using StateVector = math::Vector<T, StateSize>;
31+
using InputVector = math::Vector<T, InputSize>;
32+
using FeedbackMatrix = math::Matrix<T, InputSize, StateSize>;
33+
using FeedforwardMatrix = math::Matrix<T, InputSize, InputSize>;
34+
using ReferenceModel = math::LinearTimeInvariant<T, StateSize, InputSize, StateSize>;
35+
36+
ModelReferenceAdaptiveControl(const ReferenceModel& referenceModel,
37+
T gamma, T signB, AdaptationLaw law);
38+
39+
OPTIMIZE_FOR_SPEED InputVector ComputeControl(
40+
const StateVector& x, const InputVector& r, T dt);
41+
42+
void Reset();
43+
44+
[[nodiscard]] const StateVector& GetReferenceState() const
45+
{
46+
return xm;
47+
}
48+
49+
[[nodiscard]] const FeedbackMatrix& GetThetaX() const
50+
{
51+
return thetaX;
52+
}
53+
54+
[[nodiscard]] const FeedforwardMatrix& GetThetaR() const
55+
{
56+
return thetaR;
57+
}
58+
59+
private:
60+
const ReferenceModel& reference;
61+
StateVector xm{};
62+
FeedbackMatrix thetaX{};
63+
FeedforwardMatrix thetaR{};
64+
T gamma;
65+
T signB;
66+
AdaptationLaw law;
67+
};
68+
69+
template<typename T, std::size_t StateSize, std::size_t InputSize>
70+
ModelReferenceAdaptiveControl<T, StateSize, InputSize>::ModelReferenceAdaptiveControl(
71+
const ReferenceModel& referenceModel, T gamma, T signB, AdaptationLaw law)
72+
: reference{ referenceModel }
73+
, gamma{ gamma }
74+
, signB{ signB }
75+
, law{ law }
76+
{}
77+
78+
template<typename T, std::size_t StateSize, std::size_t InputSize>
79+
OPTIMIZE_FOR_SPEED typename ModelReferenceAdaptiveControl<T, StateSize, InputSize>::InputVector
80+
ModelReferenceAdaptiveControl<T, StateSize, InputSize>::ComputeControl(
81+
const StateVector& x, const InputVector& r, T dt)
82+
{
83+
xm = xm + (reference.A * xm + reference.B * r) * dt;
84+
85+
const StateVector e{ x - xm };
86+
87+
const InputVector u{ thetaX * x + thetaR * r };
88+
89+
const T scale{ gamma * signB * dt };
90+
for (std::size_t i = 0; i < InputSize; ++i)
91+
{
92+
for (std::size_t j = 0; j < StateSize; ++j)
93+
thetaX.at(i, j) -= scale * e.at(i, 0) * x.at(j, 0);
94+
95+
for (std::size_t j = 0; j < InputSize; ++j)
96+
thetaR.at(i, j) -= scale * e.at(i, 0) * r.at(j, 0);
97+
}
98+
99+
return u;
100+
}
101+
102+
template<typename T, std::size_t StateSize, std::size_t InputSize>
103+
void ModelReferenceAdaptiveControl<T, StateSize, InputSize>::Reset()
104+
{
105+
xm = StateVector{};
106+
thetaX = FeedbackMatrix{};
107+
thetaR = FeedforwardMatrix{};
108+
}
109+
110+
#ifdef NUMERICAL_TOOLBOX_COVERAGE_BUILD
111+
extern template class ModelReferenceAdaptiveControl<float, 1, 1>;
112+
#endif
113+
}

‎numerical/nonlinear_control/test/CMakeLists.txt‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -10,4 +10,5 @@ target_link_libraries(numerical.nonlinear_control_test PUBLIC
1010
target_sources(numerical.nonlinear_control_test PRIVATE
1111
TestBacksteppingControl.cpp
1212
TestFeedbackLinearization.cpp
13+
TestModelReferenceAdaptiveControl.cpp
1314
)

0 commit comments

Comments
 (0)