iLQR & DDP
A derivation of differential dynamic programming (DDP) and its Gauss-Newton cousin iterative linear quadratic regulator (iLQR), the workhorse algorithms for trajectory optimization. We start from the LQR backward pass, regularize and invert the action-value Hessian with Cholesky, line-search along the dynamics manifold, and handle inequality constraints with an outer augmented Lagrangian loop.
Preamble
iLQR is used everywhere in robotics today. I first came across it around 2015 when I started my "AV career" and it was the hot new thing. Coming from a controls background where LQR was standard, I was honestly surprised that something this simple (once you really understand LQR) could solve such complex problems inside tight real-time budgets where much more advanced algorithms would fail (mostly due to timing requirements). It felt real fishy.
Part of why it works so well is structural. Classical solvers treat the dynamics as equality constraints that have to be enforced. iLQR doesn't. It bakes the dynamics into its forward rollout (thank you LQR!), so they're just part of how the algorithm searches.
I spent most of the following decade working closely with it. The first few years in controls and decision-making, where you mostly need to understand its limitations. And since 2018 in trajectory planning where understanding its inner workings is often needed.
I never actually implemented it from scratch outside some quick Python prototyping when I first came across it, and I wanted to take the time to write it up properly, especially the augmented Lagrangian part. Mostly for fun, and maybe to mark an end of an chapter in my work life.
Problem setup
We consider a discrete-time optimal control problem over a horizon of $N$ stages with state $x_t \in \mathbb{R}^n$ and control $u_t \in \mathbb{R}^m$:
$$\begin{aligned} \min_{\{u_t\}} \quad & \sum_{t=0}^{N-1} \ell_t(x_t, u_t) \;+\; \ell_f(x_N) \\ \text{s.t.} \quad & x_{t+1} = f_t(x_t, u_t), \quad t = 0, \dots, N-1, \\ & x_0 = x_{\text{init}}. \end{aligned}$$The dynamics $f_t$ are nonlinear and assumed twice differentiable. Stage costs $\ell_t$ and the terminal cost $\ell_f$ are scalar-valued but not necessarily convex. Since the state is fully determined by $x_0$ and the control sequence, the real decision variables are just the controls $\{u_0, \dots, u_{N-1}\}$.
We write Jacobians as $A_t = \partial f_t / \partial x$ and $B_t = \partial f_t / \partial u$, evaluated at the current nominal point. Subscripts on costs denote partial derivatives: $\ell_x$, $\ell_{xx}$, $\ell_{ux}$, and so on.
Bellman recursion & the LQR base case
Define the optimal cost-to-go (or value function) from stage $t$ onward as
$$V_t(x) = \min_{u_t, \dots, u_{N-1}} \;\; \ell_f(x_N) + \sum_{\tau = t}^{N-1} \ell_\tau(x_\tau, u_\tau), \quad x_t = x.$$Bellman's principle of optimality gives the backward recursion
$$V_t(x) = \min_{u} \;\Big\{ \underbrace{\ell_t(x, u) + V_{t+1}\big(f_t(x, u)\big)}_{\displaystyle Q_t(x,\,u)} \Big\}, \qquad V_N(x) = \ell_f(x).$$The inner expression $Q_t(x, u)$ is the action-value function. The minimizer over $u$ defines the optimal policy $u_t^\star(x)$.
In general $V_t$ has no closed form. But in the linear-quadratic case (linear dynamics $f_t(x, u) = A_t x + B_t u$ and convex quadratic costs), the value function stays quadratic, $V_t(x) = \tfrac{1}{2} x^\top P_t x + p_t^\top x + s_t$, and the optimal policy is affine, $u_t^\star(x) = K_t x + k_t$. The matrices $P_t, K_t$ come from the discrete-time Riccati recursion. A full derivation of LQR is given by Tedrake, Underactuated Robotics §8.
DDP backward pass
Fix a nominal trajectory $\bar\tau = (\bar x_0, \bar u_0, \dots, \bar u_{N-1}, \bar x_N)$, obtained by rolling out some initial control sequence. Consider a state perturbation $\delta x$ and a control perturbation $\delta u$ at stage $t$. Expand the action-value to second order around $(\bar x_t, \bar u_t)$:
$$Q_t(\bar x_t + \delta x,\, \bar u_t + \delta u) \;\approx\; Q_t^{(0)} + \begin{pmatrix} Q_x \\ Q_u \end{pmatrix}^{\!\top}\!\! \begin{pmatrix} \delta x \\ \delta u \end{pmatrix} + \tfrac{1}{2} \begin{pmatrix} \delta x \\ \delta u \end{pmatrix}^{\!\top}\!\! \begin{pmatrix} Q_{xx} & Q_{xu} \\ Q_{ux} & Q_{uu} \end{pmatrix}\! \begin{pmatrix} \delta x \\ \delta u \end{pmatrix}.$$Suppose by induction that the next value function is also quadratic: $V_{t+1}(\bar x_{t+1} + \delta x') = V_{t+1}^{(0)} + V_x'^{\,\top} \delta x' + \tfrac{1}{2}\delta x'^{\,\top} V_{xx}' \,\delta x'$. Applying the chain rule to $V_{t+1}(f_t(x, u))$ and combining with the cost expansion gives:
$$\begin{aligned} Q_x &= \ell_x + A_t^\top V_x' \\ Q_u &= \ell_u + B_t^\top V_x' \\ Q_{xx} &= \ell_{xx} + A_t^\top V_{xx}' A_t \;+\; \boxed{V_x' \cdot f_{xx}} \\ Q_{ux} &= \ell_{ux} + B_t^\top V_{xx}' A_t \;+\; \boxed{V_x' \cdot f_{ux}} \\ Q_{uu} &= \ell_{uu} + B_t^\top V_{xx}' B_t \;+\; \boxed{V_x' \cdot f_{uu}} \end{aligned}$$A quick reminder on what's what in these formulas, since $\ell$ and $V$ play very different roles. The $\ell_x, \ell_u, \ell_{xx}, \ldots$ are derivatives of the immediate stage cost $\ell_t(x_t, u_t)$ at the nominal point. The primed $V_x', V_{xx}'$ are derivatives of the cost-to-go $V_{t+1}$ from the next stage onward, which already sums over all future stages.
The boxed terms are tensors built from the value gradient and the second derivatives of the dynamics. Concretely, $(V_x' \cdot f_{xx})_{ij} = \sum_k [V_x']_k \,[\partial^2 f_t^{(k)} / \partial x_i \partial x_j]$, and similarly for the other blocks.
$Q$ is quadratic in $\delta u$, so the unconstrained minimum comes from setting $\partial_{\delta u} Q = 0$:
$$Q_{uu}\,\delta u^\star + Q_{ux}\,\delta x + Q_u = 0 \quad\Longrightarrow\quad \boxed{\;\delta u^\star(\delta x) = K_t\,\delta x + k_t\;}$$where the feedback gain and feedforward correction are
$$K_t = -Q_{uu}^{-1} Q_{ux}, \qquad k_t = -Q_{uu}^{-1} Q_u.$$Substituting $\delta u^\star$ back into the quadratic model produces a quadratic in $\delta x$, which is exactly the value function $V_t$ at the current stage. After simplification (the cross terms cancel out):
$$\begin{aligned} V_x &= Q_x + K_t^\top Q_u \;=\; Q_x - Q_{xu} Q_{uu}^{-1} Q_u, \\ V_{xx} &= Q_{xx} + K_t^\top Q_{ux} \;=\; Q_{xx} - Q_{xu} Q_{uu}^{-1} Q_{ux}, \\ \Delta V &= \tfrac{1}{2}\, k_t^\top Q_u \;=\; -\tfrac{1}{2}\, Q_u^\top Q_{uu}^{-1} Q_u \;\le\; 0. \end{aligned}$$These are the value derivatives at stage $t$. On the next pass through the loop (stage $t-1$), they take the role of the primed quantities $V_x', V_{xx}'$ in the $Q$ formulas above. The constant offset $\Delta V$ accumulates the predicted cost reduction across stages, which we will use during the line search.
So the backward sweep runs like this:
- Initialize at the horizon with the terminal cost: $V_x \leftarrow \ell_{f,x}(\bar x_N)$ and $V_{xx} \leftarrow \ell_{f,xx}(\bar x_N)$.
- For $t = N{-}1, \dots, 0$: assemble $Q_x, Q_u, Q_{xx}, Q_{ux}, Q_{uu}$ from the current $V_x, V_{xx}$ (which are the value derivatives at stage $t{+}1$). Solve for $(K_t, k_t)$ and store them. Overwrite $V_x, V_{xx}$ with the formulas above, which are now the value derivatives at stage $t$, ready for the next iteration.
iLQR: dropping the dynamics Hessians
The boxed tensor terms in the $Q$ Hessians are what make DDP genuinely second-order in the dynamics. They're also expensive: computing $f_{xx}, f_{ux}, f_{uu}$ requires $n$ Hessian evaluations. iLQR drops them entirely:
$$\begin{aligned} Q_{xx}^{\text{iLQR}} &= \ell_{xx} + A_t^\top V_{xx}' A_t, \\ Q_{ux}^{\text{iLQR}} &= \ell_{ux} + B_t^\top V_{xx}' A_t, \\ Q_{uu}^{\text{iLQR}} &= \ell_{uu} + B_t^\top V_{xx}' B_t. \end{aligned}$$This is the Gauss-Newton approximation of the action-value Hessian. We treat the dynamics as locally linear and use a true quadratic on costs only. Three consequences:
- Each backward pass only needs $A_t, B_t$, not $f_{xx}, f_{ux}, f_{uu}$.
- If the original costs are convex (so $\ell_{xx}, \ell_{uu} \succeq 0$) and $V_{xx}' \succeq 0$, the iLQR $Q_{uu}$ is automatically PSD. DDP can fail this condition even on simple problems because $V_x' \cdot f_{uu}$ can become negative.
In practice iLQR is the default. DDP's extra terms rarely justify their cost unless dynamics curvature really matters for your problem.
Regularizing $Q_{uu}$
The backward pass requires inverting $Q_{uu}$. Worse, for $\delta u^\star = -Q_{uu}^{-1}(Q_u + Q_{ux}\delta x)$ to be a descent direction, we need $Q_{uu} \succ 0$. If $Q_{uu}$ is indefinite, the local quadratic model has unbounded negative directions and the "minimum" doesn't mean much. You could shoot off to infinity while minimizing.
To ensure positive definiteness we use a Levenberg-Marquardt style regularizer:
$$\tilde Q_{uu} = Q_{uu} + \mu I, \qquad \mu \ge 0.$$Increasing $\mu$ shifts every eigenvalue of $Q_{uu}$ up by $\mu$, pushing them towards positive values. As $\mu \to \infty$, the step $k_t = -\tilde Q_{uu}^{-1} Q_u \to -\tfrac{1}{\mu} Q_u$ becomes a tiny gradient step (safe but slow). As $\mu \to 0$, we recover the pure Newton step.
In practice we never form $\tilde Q_{uu}^{-1}$ explicitly. We use the Cholesky factorization $\tilde Q_{uu} = L L^\top$, which exists iff $\tilde Q_{uu}$ is positive definite. The same operation then does two jobs at once:
- PSD check. Cholesky succeeds iff the matrix is positive definite. We attempt the factorization at the current $\mu$, and if it fails (non-positive diagonal pivot), we increase $\mu$ and retry. No separate eigenvalue computation needed.
- Fast linear solves. Once we have $L$, $k_t$ and $K_t$ come from two triangular back-substitutions: $L y = -Q_u$ then $L^\top k_t = y$, and similarly $L Y = -Q_{ux}$ then $L^\top K_t = Y$.
Forward pass & line search
Given the gains $\{(K_t, k_t)\}$ from the backward pass, the forward pass rolls out a new trajectory. The natural step is
$$\delta u_t = K_t\, \delta x_t + k_t,$$but the local quadratic model is only trustworthy near the nominal trajectory. Taking the full step often overshoots, so we need a line-search parameter $\alpha \in (0, 1]$. The right thing to do is scale only the feedforward term:
$$\boxed{\; \delta u_t(\alpha) = K_t\,\delta x_t(\alpha) + \alpha\, k_t \;}$$where $\delta x_t(\alpha)$ is the state perturbation that the nonlinear rollout actually produces:
$$\begin{aligned} x_0^+ &= x_{\text{init}}, \\ \delta x_t &= x_t^+ - \bar x_t, \\ u_t^+ &= \bar u_t + K_t\,\delta x_t + \alpha\, k_t, \\ x_{t+1}^+ &= f_t(x_t^+,\, u_t^+). \end{aligned}$$Why scale only $k_t$?
The two terms play very different roles:
- $K_t \delta x_t$ is a feedback stabilizer. It corrects for the gap between the linearized model and the true nonlinear evolution. Without it, even a perfect open-loop step would drift as the linearization approximation gets worse from stage to stage.
- $\alpha k_t$ is the feedforward improvement, the open-loop correction we'd apply if the model were exact. This contains most of the uncertainty and shrinking it is considered safe.
Scaling both terms by $\alpha$ would weaken the feedback as we shrink $\alpha$, leaving the rollout unstable just when we need it stable. Keeping $K_t$ at full strength while reducing $k_t$ effectively does a line search along the manifold of dynamically feasible trajectories. As $\alpha \to 0$ the feedback policy drives $\delta x_t \to 0$ and we recover the nominal trajectory, so the worst case is "no progress" rather than "blow up".
Acceptance uses the standard sufficient-decrease pattern. The model predicts a cost reduction (a closed form falls out of the backward pass):
$$\Delta J_{\text{pred}}(\alpha) = -\Big(\alpha - \tfrac{\alpha^2}{2}\Big)\sum_{t=0}^{N-1} k_t^\top Q_u^{(t)} \;\;>\;0 \;\;\text{for}\;\; \alpha \in (0, 1].$$Let $J^+ = J(\tau^+(\alpha))$ be the actual rolled-out cost. The step is accepted if
$$\zeta \;=\; \frac{J(\bar\tau) - J^+}{\Delta J_{\text{pred}}(\alpha)} \;\ge\; \eta,$$typically with $\beta \in [10^{-4},\, 1/4]$. If $\zeta < \beta$, halve $\alpha$ and re-roll. If $\alpha$ falls below a floor, raise $\mu$ and re-run the backward pass.
Augmented Lagrangian for $g(x, u) \le 0$
Now add stage inequality constraints to the problem:
$$g_t(x_t, u_t) \;\le\; 0, \qquad t = 0, \dots, N-1.$$The high-level picture
AL handles each inequality $g_t \le 0$ by absorbing it into the stage cost as a smooth penalty term. We trade one constrained problem (minimize the trajectory cost subject to $g_t \le 0$) for a sequence of unconstrained ones, each of which is just a standard iLQR solve on an augmented cost. The multiplier $\lambda_t \ge 0$ and penalty $\rho > 0$ are updated between solves and steer the successive solutions toward feasibility.
Concretely, the original constrained problem
$$\min_{\{u_t\}} \sum_{t=0}^{N-1} \ell_t(x_t, u_t) + \ell_f(x_N) \quad\text{s.t.}\quad x_{t+1} = f_t(x_t, u_t),\;\; g_t(x_t, u_t) \le 0$$is replaced by a sequence of unconstrained (well, equality-only) subproblems
$$\min_{\{u_t\}} \sum_{t=0}^{N-1} \tilde\ell_t(x_t, u_t;\, \lambda_t, \rho) + \ell_f(x_N) \quad\text{s.t.}\quad x_{t+1} = f_t(x_t, u_t),$$indexed by an evolving $(\lambda_t, \rho)$. The dynamics constraint stays exactly where iLQR likes it (baked into the rollout); only the inequalities move into the cost.
The augmented stage cost
The standard derivation (introduce slack $s_t \ge 0$ via $g_t + s_t = 0$, minimize analytically over $s_t$, drop a constant) gives the per-stage augmented cost
$$\tilde\ell_t(x_t, u_t;\, \lambda_t, \rho) \;=\; \ell_t(x_t, u_t) \;+\; \frac{1}{2\rho}\!\left[\,\max\{0,\, \lambda_t + \rho\, g_t(x_t, u_t)\}^2 \;-\; \lambda_t^2\,\right].$$The first term is the original stage cost; the second is the smooth penalty that encodes the inequality. The $\max$ implements the active set. When $\lambda_t + \rho\, g_t \le 0$ (constraint inactive with slack) the penalty collapses to the constant $-\lambda_t^2/(2\rho)$ and contributes nothing to the gradient. When $\lambda_t + \rho\, g_t > 0$ (violation, or active with no slack) it becomes a smooth quadratic penalty on $g_t$.
Define the active multiplier estimate at the current point as
$$\hat\lambda_t \;\equiv\; \max\{0,\, \lambda_t + \rho\, g_t(x_t, u_t)\}.$$The AL derivatives then have a clean form, which we can include directly into the stage cost gradient and Hessian that feed the backward pass:
$$\nabla \tilde\ell_t = \nabla \ell_t + \hat\lambda_t\,\nabla g_t, \qquad \nabla^2 \tilde\ell_t = \nabla^2 \ell_t + \hat\lambda_t\,\nabla^2 g_t + \rho\, \mathbb{1}_{\{\hat\lambda_t > 0\}}\, \nabla g_t\, \nabla g_t^\top.$$All gradients and Hessians here are with respect to $(x_t, u_t)$. The indicator $\mathbb{1}_{\{\hat\lambda_t > 0\}}$ zeros out the rank-one penalty Hessian whenever this stage's constraint is inactive. The Hessian is discontinuous across the active-set boundary $g_t = -\lambda_t/\rho$, but it's PSD on each side, which is consistent with what the Cholesky regularizer expects.
The outer loop alternates two things:
- Inner solve. Run iLQR/DDP to (approximately) minimize the AL-augmented total cost $$J_{\text{AL}}(\{u_t\};\, \{\lambda_t\}, \rho) \;\equiv\; \sum_{t=0}^{N-1} \tilde\ell_t(x_t, u_t;\, \lambda_t, \rho) \;+\; \ell_f(x_N),$$ subject only to the dynamics. This scalar $J_{\text{AL}}$ is what the next subsection refers to.
- Outer update. Compare the worst constraint violation $v = \max_t \|\max\{0, g_t\}\|_\infty$ against the current feasibility tolerance $\eta$:
- If $v \le \eta$ (feasibility is on track): update multipliers $\lambda_t^+ = \max\{0,\, \lambda_t + \rho\, g_t(x_t^\star, u_t^\star)\}$, keep $\rho$, and tighten both tolerances. The multiplier update is the first-order estimate from the AL gradient at the inner optimum, and the $\max$ automatically enforces $\lambda \ge 0$.
- If $v > \eta$ (still too infeasible): keep multipliers, inflate $\rho \leftarrow 10\rho$, and reset tolerances.
The tolerance schedule
A subtle but important property of AL is that the inner solves don't need to be exact. A pure quadratic penalty would need $\rho \to \infty$ to drive $g \to 0$, which would leave the inner Hessian badly conditioned. AL avoids this. Once $\lambda$ gets close to the true dual variable, the penalty term $\rho \|g\|^2$ stays small at the optimum, so we can keep $\rho$ bounded. The downside is that there are now two things to watch, how well the inner solve converged and how feasible the iterate is, and we tighten both as the outer loop goes on.
This gives the algorithm two pairs of tolerances:
- Optimality tolerance $\omega$ (current) and $\omega^\star$ (final). At each outer iteration we ask the inner iLQR to drive $\|\nabla_u J_{\text{AL}}\| \le \omega$. Note that $J_{\text{AL}}$ is scalar but $\nabla_u J_{\text{AL}} \in \mathbb{R}^{mN}$ is the gradient with respect to the entire control trajectory $(u_0, \ldots, u_{N-1})$. Conveniently, $\partial J_{\text{AL}} / \partial u_t = Q_u^{(t)}$, the partial at stage $t$ already bundles in the cost-to-go from $t{+}1$ onward via the chain rule through the dynamics, which is exactly what the backward pass computed. So we just stack the $Q_u^{(t)}$'s as they fall out of the backward sweep; the norm is essentially free.
- Feasibility tolerance $\eta$ (current) and $\eta^\star$ (final). The threshold above which the outer step inflates $\rho$ instead of updating multipliers.
A typical schedule (Bertsekas-style, also in Nocedal & Wright ch. 17):
- Initial: $\omega_0 = 1/\rho_0$, $\eta_0 = 1/\rho_0^{0.1}$.
- Successful outer step ($v \le \eta$): $\omega \leftarrow \omega/\rho$, $\eta \leftarrow \eta/\rho^{0.9}$.
- Failed outer step ($v > \eta$, after $\rho$ is inflated): $\omega \leftarrow 1/\rho_{\text{new}}$, $\eta \leftarrow 1/\rho_{\text{new}}^{0.1}$.
Final convergence is declared when $v \le \eta^\star$ and $\|\nabla_u J_{\text{AL}}\| \le \omega^\star$ both hold.
Full AL iLQR algorithm
Putting it all together: an outer AL loop over multipliers and penalty wraps an inner iLQR loop over trajectory iterates. The inner loop alternates a regularized backward pass and a line-searched forward pass. The pseudocode below makes the bookkeeping explicit, including the tolerance schedule and iteration budgets.
$x_{\text{init}}$ : initial state $\{\bar u_t\}_{t=0}^{N-1}$ : warm-start control trajectory $\omega^\star,\, \eta^\star$ : final optimality / feasibility tolerances $\rho_0,\, \mu_0$ : initial AL penalty / LM regularizer $\beta \in [10^{-4}, \tfrac{1}{4}],\, \alpha_{\min}$ : line-search threshold and $\alpha$ floor $\texttt{max\_outer},\, \texttt{max\_inner}$ : iteration budgets
initialize
$\lambda_t \leftarrow 0$ for all $t$ : one multiplier per constraint per stage $\rho \leftarrow \rho_0,\; \mu \leftarrow \mu_0$ $\omega \leftarrow 1/\rho,\; \eta \leftarrow 1/\rho^{0.1}$ : initial inner tolerances roll out $\{\bar x_t\}$ from $\{\bar u_t\}$ through $f_t$
for $k = 1, \ldots, \texttt{max\_outer}$: // outer AL loop
for $j = 1, \ldots, \texttt{max\_inner}$: // inner iLQR loop, $(\lambda, \rho)$ fixed backward pass: $V_x \leftarrow \ell_{f,x}(\bar x_N),\; V_{xx} \leftarrow \ell_{f,xx}(\bar x_N)$ for $t = N{-}1, \ldots, 0$: assemble $Q_x, Q_u, Q_{xx}, Q_{ux}, Q_{uu}$ // from AL-augmented stage cost attempt Cholesky $L L^\top = Q_{uu} + \mu I$ if Cholesky fails: $\mu \leftarrow 10\mu$, restart backward pass solve $L L^\top k_t = -Q_u,\; L L^\top K_t = -Q_{ux}$ update $V_x, V_{xx}$, accumulate $\Delta J_{\text{pred}}$
forward pass: // line search on $\alpha$ $\alpha \leftarrow 1$ while $\alpha \ge \alpha_{\min}$: roll out $u_t^+ = \bar u_t + K_t(x_t^+ - \bar x_t) + \alpha\, k_t$ $\zeta \leftarrow (J(\bar\tau) - J(\tau^+)) / \Delta J_{\text{pred}}(\alpha)$ if $\zeta \ge \beta$: accept ($\bar\tau \leftarrow \tau^+,\, \mu \leftarrow \mu/2$), break else: $\alpha \leftarrow \alpha/2$ if $\alpha < \alpha_{\min}$ (line search collapsed): $\mu \leftarrow 10\mu$, retry inner iteration
if $\|\nabla_u J_{\text{AL}}\| \le \omega$: break inner loop // inner converged
outer step: // feasibility check, multiplier / penalty update $v \leftarrow \max_t \|\max\{0, g_t(\bar x_t, \bar u_t)\}\|_\infty$ // worst constraint violation if $v \le \eta^\star$ and $\|\nabla_u J_{\text{AL}}\| \le \omega^\star$: return $\bar\tau$ // fully converged
if $v \le \eta$: // feasibility on track, update $\lambda$, tighten tolerances $\lambda_t \leftarrow \max\{0,\, \lambda_t + \rho\, g_t(\bar x_t, \bar u_t)\}$ for all $t$ $\omega \leftarrow \omega/\rho,\; \eta \leftarrow \eta/\rho^{0.9}$ else: // not feasible enough, inflate penalty $\rho \leftarrow 10\rho$ $\omega \leftarrow 1/\rho,\; \eta \leftarrow 1/\rho^{0.1}$
A poor man's approach to constraint
The augmented Lagrangian above handles any inequality $g(x, u) \le 0$, but it is an outer loop with its own penalty schedule. For the very common case of plain box limits on the control, $u_{\min} \le u_t \le u_{\max}$, there is a much cheaper trick that lives entirely inside the backward pass. I used a version of this on an AV planner and it held up well in production.
Recall the backward pass builds the stage law by minimizing the second-order action-value over $\delta u$. Collecting the terms that depend on $\delta u$, it solves
$$\delta u^\star = \arg\min_{\delta u}\ \tfrac{1}{2}\,\delta u^\top Q_{uu}\,\delta u + \big(Q_u + Q_{ux}\,\delta x\big)^\top \delta u \;=\; K_t\,\delta x + k_t.$$Box limits on the control just add a constraint to that same minimization, turning it into a small QP at each stage:
$$\min_{\delta u}\ \tfrac{1}{2}\,\delta u^\top Q_{uu}\,\delta u + \big(Q_u + Q_{ux}\,\delta x\big)^\top \delta u \quad\text{s.t.}\quad u_{\min} \le \bar u_t + \delta u \le u_{\max}.$$If the control is at most two-dimensional (true for most vehicle models, think longitudinal jerk and a lateral input), the box has at most four corners and each input is in one of three states: free, at its lower bound, or at its upper bound. So you can read the active set straight off the sign of the gradient $Q_{uu}\,\delta u + Q_u$ at the box corners, with a handful of if statements. Once you know which bounds are active, you correct the gain and feedforward $(K_t, k_t)$ analytically from the active-set KKT system, so the feedback law itself respects the active limit, not just the open-loop term.
All other constraints have to be handled through the cost formulation, with some kind of barrier term. We often call those soft constraints. With careful design this approach works well in practice and is very efficient.
The catch
The active set is computed at the nominal point ($\delta x = 0$), so the corrected law is only feasible there. During the forward rollout the feedback term $K_t\,\delta x$ can push the realized control back outside the box. Clamping would be the obvious fix, but it breaks the feedback law and the cost-to-go bookkeeping, so the rollout is left unclamped. The line search only enforces a cost decrease, not feasibility. What keeps the controls in bounds is the convergence test: the solver refuses to accept a solution while any control is out of bounds, so it just keeps iterating. As the trajectory settles ($\delta x \to 0$, steps shrink) the feedback overshoot vanishes and the controls come back inside the box. So box feasibility is only guaranteed at convergence; intermediate iterates can violate it. In practice this is rarely a problem and the method works well, but it is good to know where the approximation is.
Word of caution / Reality check
iLQR is basically an iterative local quadratic approximation. It's not a global search and it won't magically solve every problem you throw at it. The cost function and the dynamic model matter a lot, and the solver will exploit anything you get wrong in either.
It's also very sensitive to warm start, especially on non-convex problems, which is most of the interesting ones. Obstacle avoidance is the obvious example. You really want a warm-start module that can both generate candidate seeds and evaluate them to find the best one. What "best" means here is most of the challenge.
You also want a strong feasibility/safety checker downstream. The trajectories the solver gives you aren't automatically feasible in any real sense. They satisfy the dynamics they were given, but the model is always an approximation and the cost is always a proxy.
The pipeline I've seen actually work in production looks something like this:
The good news is that the barebones solver (no augmented Lagrangian, no DDP, just iLQR with a reasonable regularizer and a working line search) is pretty quick to get running. I've shipped systems that intentionally skipped state constraints and used soft cost penalties instead. They were fine, and stayed fast. I also got to push iLQR to the limit of vehicle handling using a full dynamic model instead of the kinematic model that's standard in AV, and was genuinely surprised that it didn't blow up and actually performed really well. Someday I'd like to point it at a new dynamical system and see how that AV experience holds up there.
Getting help from ML
iLQR is a classical solver, but there are a few places where a learned model can quietly make it better. Three that I find interesting:
Estimating the dynamics residual. Split the model into a known nominal part and a small learned correction:
$$x_{k+1} = f_\text{nom}(x_k, u_k) + g_\theta(x_k, u_k).$$iLQR only ever touches the dynamics through the Jacobians it linearizes at each step, so the residual just adds its own Jacobians into the mix:
$$A_k = \frac{\partial f_\text{nom}}{\partial x} + \frac{\partial g_\theta}{\partial x}, \qquad B_k = \frac{\partial f_\text{nom}}{\partial u} + \frac{\partial g_\theta}{\partial u}.$$Nothing else in the backward pass changes. I think this is generally not needed for a kinematic model, where the equations are simple and trustworthy, but it could be useful for a highly dynamical system, since more dynamics usually means more uncertainty in the model.
Estimating the terminal cost. A learned terminal cost lets you shorten the horizon: instead of planning all the way out, you stop early and let the learned value carry the rest. It is a quite challenging problem for AV given how complex the cost landscape can be, but it could lead to big performance gains and naturally bakes long-term planning into the trajectory.
Picking a warm start. The warm-start module from above is a natural place for ML. You can learn to propose a good seed, or to score candidate seeds and pick the best one. This is something I have used in the past and it helped a lot on the non-convex problems.
References
- D. Q. Mayne and D. H. Jacobson, Differential Dynamic Programming, Elsevier, 1970. The original.
- Y. Tassa, N. Mansard, E. Todorov, Control-Limited Differential Dynamic Programming, ICRA 2014. The modern reference for regularization, line search, and box constraints.
- R. Tedrake, Underactuated Robotics §8. LQR backward pass derived from scratch.
- J. Nocedal and S. Wright, Numerical Optimization, ch. 17. Augmented Lagrangian and modified Cholesky.
- T. A. Howell, B. Jackson, Z. Manchester, ALTRO: A Fast Solver for Constrained Trajectory Optimization, IROS 2019. A clean implementation of the AL-iLQR pipeline.
Terminal reward
If you read all that congrats! Treat yourself by running the secret shell command ./ahhhhh.sh in the terminal of the homepage.