引入

Riccati 方程是求解最优控制问题的关键。此前，我们都是用变分法与极大值原理理解的 Riccati 方程。事实上， Riccati 方程还可以用一种更直观的方法推导得来——不需要哈密顿量，不需要拉格朗日乘子，而是基于大家都知道的动态规划算法。

离散系统的动态规划

回顾动态规划算法

首先，我们回顾一下动态规划中最经典的最优路径问题。

例如有一块 $n \times n$ 的棋盘，棋子在棋盘的左上角。每一步棋子可以选择向右或向下走，直到抵达右下角的终点。此时，给棋盘上相邻格点间的每一段路径都设置一个代价，而我们的目标是找出那条是总路径所有代价之和最小的路径。

用动态规划的思想，我们知道我们可以从终点的格点开始，嵌套循环遍历求出每个格点到终点的最小代价。对每个格点，只需要比较右方格点的最小代价与向右路径代价的和，与下方格点的最小代价与向下路径代价的和，然后选出其中最小作为当前格点的最小代价，同时在格点数据结构的属性中记录代价最小的路径即可。

不难知道，整个路径的“最小代价”之和就是指标函数。为了模拟整个移动过程的代价累计，我们不妨设移动过程中的每一个状态都有一级最小代价，指标则是所有最小代价的和。

现在我们来写成公式形式。设走完全程需要 N 步，离散系统为：

$$\boldsymbol{x}_{k+1} = f(\boldsymbol{x}_k, \boldsymbol{u}_k)$$

对每个状态，定义每个状态的最小代价

$$V(\boldsymbol{x}_k) = \min_{\boldsymbol{u}_k}(l(\boldsymbol{x}_k, \boldsymbol{u}_k) + V(\boldsymbol{x}_{k+1})) = \min_{\boldsymbol{u}_k} [l(\boldsymbol{x}_k, \boldsymbol{u}_k) + V(f(\boldsymbol{x}_k, \boldsymbol{u}_k))], \qquad k \in [0, N-1]$$

同时我们有边界条件

$$V(\boldsymbol{\boldsymbol{x}}_N) = l_f$$

其中，$l(\boldsymbol{x}_k, \boldsymbol{u}_k)$ 表示路径的代价，$l_f$ 表示终端状态的代价，其中的参数 $\boldsymbol{u}_k$ 是输入，即向右还是向下。我们的目标是求出最优的输入使得指标最小。此时指标就等于：

$$J = V(\boldsymbol{x}_0) = l_f + \sum_{i=0}^{N-1} \min_{\boldsymbol{u}_i} l(\boldsymbol{x}_i,\boldsymbol{u}_i)$$

最优控制为：

$$\boldsymbol{u}^*_k(\boldsymbol{x}) = \mathrm{arg}\min_{\boldsymbol{u}_k} [l(\boldsymbol{x}_k, \boldsymbol{u}_k) + V(f(\boldsymbol{x}_k, \boldsymbol{u}_k))]$$

求解

设系统为线性系统

$$\boldsymbol{x}_{k+1} = \boldsymbol{A}_k \boldsymbol{x}_k + \boldsymbol{B}_k \boldsymbol{u}_k$$

我们不妨将性能指标定义为线性二次型形式，即

$$J = \boldsymbol{x}_N^T \boldsymbol{S} \boldsymbol{x}_N + \sum_{k=0}^{N-1} (\boldsymbol{x}^T_k \boldsymbol{Q}_k \boldsymbol{x}_k + \boldsymbol{u}^T_k \boldsymbol{R}_k \boldsymbol{u}_k)$$

则每一步的指标可以表示为 $l(\boldsymbol{x}_k, \boldsymbol{u}_k) = \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{u}^T \boldsymbol{R} \boldsymbol{u}$。对每一步的最优指标可以写出贝尔曼方程：

$$\begin{aligned} V(\boldsymbol{x}_k) &= \min_{\boldsymbol{u}_k} \left[ \boldsymbol{x}^T_k \boldsymbol{Q}_k \boldsymbol{x}_k + \boldsymbol{u}^T_k \boldsymbol{R}_k \boldsymbol{u}_k + V(\boldsymbol{x}_{k+1}) \right] \\ &= \min_{\boldsymbol{u}_k} \left[ \boldsymbol{x}^T_k \boldsymbol{Q}_k \boldsymbol{x}_k + \boldsymbol{u}^T_k \boldsymbol{R}_k \boldsymbol{u}_k + V(\boldsymbol{A}_k \boldsymbol{x}_k + \boldsymbol{B}_k \boldsymbol{u}_k) \right] \end{aligned}$$

边界条件为

$$V(\boldsymbol{x}_N) = \boldsymbol{x}_N^T \boldsymbol{S} \boldsymbol{x}_N$$

由于系统是线性系统，指标为二次型，我们可以定性得到 $V(\boldsymbol{x}_k)$ 也是二次型。我们不妨设 $V(\boldsymbol{x}_k) = \boldsymbol{x}^T_k \boldsymbol{P}_k \boldsymbol{x}_k$。则贝尔曼方程可以表示为：

$$V(\boldsymbol{x}_k) = \min_{\boldsymbol{u}_k} \left[ \boldsymbol{x}^T_k \boldsymbol{Q}_k \boldsymbol{x}_k + \boldsymbol{u}^T_k \boldsymbol{R}_k \boldsymbol{u}_k + (\boldsymbol{A}_k \boldsymbol{x}_k + \boldsymbol{B}_k \boldsymbol{u}_k)^T \boldsymbol{P}_{k+1} (\boldsymbol{A}_k \boldsymbol{x}_k + \boldsymbol{B}_k \boldsymbol{u}_k) \right]$$

令偏导等于零，

$$\frac{\partial V(\boldsymbol{x}_k)}{\partial \boldsymbol{u}^*_k} = 2\boldsymbol{R}_k \boldsymbol{u}_k + 2\boldsymbol{B}^T_k \boldsymbol{P}_{k+1} (\boldsymbol{A}_k \boldsymbol{x}_k + \boldsymbol{B}_k \boldsymbol{u}_k) = 0$$

解得

$$\boldsymbol{u}^*_k = - (\boldsymbol{R}_k + \boldsymbol{B}^T_k \boldsymbol{P}_{k+1} \boldsymbol{B}_k)^{-1} \boldsymbol{B}^T_k \boldsymbol{P}_{k+1} \boldsymbol{A}_k \boldsymbol{x}_k$$

二阶偏导为

$$\frac{\partial^2 V(\boldsymbol{x}_k)}{\partial \boldsymbol{u}_k^{*2}} = 2\boldsymbol{R}_k + 2\boldsymbol{B}^T_k \boldsymbol{P}_{k+1} \boldsymbol{B}_k$$

由 P 的定义可知 P 半正定，又由 R 的正定性，$\frac{\partial^2 V(\boldsymbol{x}_k)}{\partial \boldsymbol{u}_k^2} \succ 0$，故 $\boldsymbol{u}^*_k$ 是最优控制。不妨设 $\boldsymbol{K}_k = (\boldsymbol{R}_k + \boldsymbol{B}^T_k \boldsymbol{P}_{k+1} \boldsymbol{B}_k)^{-1} \boldsymbol{B}^T_k \boldsymbol{P}_{k+1} \boldsymbol{A}_k$，则 $\boldsymbol{u}^*_k = -\boldsymbol{K}_k \boldsymbol{x}_k$。接下来为了方便，将省略下标 k。将最优控制代回到贝尔曼方程，并改写 $V(\boldsymbol{x}_k)$：

$$\begin{aligned} \boldsymbol{x}^T \boldsymbol{P}_k \boldsymbol{x} &= \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{x}^T \boldsymbol{K}^T \boldsymbol{R} \boldsymbol{K} \boldsymbol{x} + \boldsymbol{x}^T (\boldsymbol{A} - \boldsymbol{B} \boldsymbol{K})^T \boldsymbol{P}_{k+1} (\boldsymbol{A} - \boldsymbol{B} \boldsymbol{K}) \boldsymbol{x} \\ \boldsymbol{P}_k &= \boldsymbol{Q} + \boldsymbol{K}^T \boldsymbol{R} \boldsymbol{K} + (\boldsymbol{A} - \boldsymbol{B} \boldsymbol{K})^T \boldsymbol{P}_{k+1} (\boldsymbol{A} - \boldsymbol{B} \boldsymbol{K}) \\ &= \boldsymbol{Q} + \boldsymbol{K}^T \boldsymbol{R} \boldsymbol{K} + \boldsymbol{A}^T \boldsymbol{P}_{k+1} \boldsymbol{A} + \boldsymbol{K}^T \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{B} \boldsymbol{K} - \boldsymbol{K}^T \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{A} - \boldsymbol{A}^T \boldsymbol{P}_{k+1} \boldsymbol{B} \boldsymbol{K} \\ &= \boldsymbol{Q} + \boldsymbol{A}^T \boldsymbol{P}_{k+1} \boldsymbol{A} + \boldsymbol{K}^T (\boldsymbol{R} + \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{B}) \boldsymbol{K} - \boldsymbol{K}^T \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{A} - \boldsymbol{A}^T \boldsymbol{P}_{k+1} \boldsymbol{B} \boldsymbol{K} \end{aligned}$$

为了化简两个非二次型形式的负项，我们需要回到 K 的定义：

$$\begin{aligned} \boldsymbol{K} &= (\boldsymbol{R} + \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{B}_k)^{-1} \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{A} \\ \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{A} &= (\boldsymbol{R} + \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{B}) \boldsymbol{K} \\ \boldsymbol{K}^T \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{A} &= \boldsymbol{K}^T (\boldsymbol{R} + \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{B}) \boldsymbol{K} \end{aligned}$$

将上式转置可得 $\boldsymbol{A}^T \boldsymbol{P}_{k+1} \boldsymbol{B} \boldsymbol{K} = \boldsymbol{K}^T (\boldsymbol{R} + \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{B})^T \boldsymbol{K} = \boldsymbol{K}^T (\boldsymbol{R} + \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{B}) \boldsymbol{K} = \boldsymbol{K}^T \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{A}$，故

$$\boldsymbol{P}_k = \boldsymbol{Q} + \boldsymbol{A}^T \boldsymbol{P}_{k+1} \boldsymbol{A} - \boldsymbol{K}^T \boldsymbol{B}^T \boldsymbol{P}_{k+1} \boldsymbol{A}$$

带入 $\boldsymbol{K}$ 的表达式，得

$$\boldsymbol{P}_k = \boldsymbol{Q}_k + \boldsymbol{A}_k^T \boldsymbol{P}_{k+1} \boldsymbol{A}_k - \boldsymbol{A}_k^T \boldsymbol{P}_{k+1} \boldsymbol{B}_k (\boldsymbol{R}_k + \boldsymbol{B}^T_k \boldsymbol{P}_{k+1} \boldsymbol{B}_k)^{-1} \boldsymbol{B}_k^T \boldsymbol{P}_{k+1} \boldsymbol{A}_k$$

这就是离散时间代数 Riccati 方程。负反馈增益由 $\boldsymbol{K}_k = (\boldsymbol{R}_k + \boldsymbol{B}^T_k \boldsymbol{P}_{k+1} \boldsymbol{B}_k)^{-1} \boldsymbol{B}^T_k \boldsymbol{P}_{k+1} \boldsymbol{A}_k$ 得到。

连续系统的动态规划

设连续系统

$$\dot{\boldsymbol{x}} = \boldsymbol{A}(t) \boldsymbol{x}(t) + \boldsymbol{B}(t) \boldsymbol{u}(t)$$

其中 $\boldsymbol{u}(t) \in \boldsymbol{U}$。这里我们使用线性二次型指标，性能指标为

$$J = \boldsymbol{x}^T(T) \boldsymbol{S} \boldsymbol{x}(T) + \int^T_{t_0} [\boldsymbol{x}^T(t) \boldsymbol{Q}(t) \boldsymbol{x}(t) + \boldsymbol{u}^T(t) \boldsymbol{R}(t) \boldsymbol{u}(t)] dt$$

则最小指标

$$\begin{aligned} V(\boldsymbol{x}(t), t) &= \min_{\boldsymbol{u}(t) \in \boldsymbol{U}} \left\{ \boldsymbol{x}^T(T) \boldsymbol{S} \boldsymbol{x}(T) + \int^T_{t} ( \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{u}^T \boldsymbol{R} \boldsymbol{u} ) dt \right\} \\ &= \min_{\boldsymbol{u}(t) \in \boldsymbol{U}} \left\{ \int^{t+\Delta t}_{t} ( \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{u}^T \boldsymbol{R} \boldsymbol{u} ) dt + \boldsymbol{x}^T(T) \boldsymbol{S} \boldsymbol{x}(T) + \int^T_{t+\Delta t} ( \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{u}^T \boldsymbol{R} \boldsymbol{u} ) dt \right\} \\ &= \min_{\boldsymbol{u}(t) \in \boldsymbol{U}} \left\{ \int^{t+\Delta t}_{t} ( \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{u}^T \boldsymbol{R} \boldsymbol{u} ) dt + V(\boldsymbol{x}(t + \Delta t), t + \Delta t) \right\} \end{aligned}$$

容易知道，连续系统无法直接求出迭代式，需要使用微分方程求解。其中 $V(\boldsymbol{x}(t + \Delta t))$ 使用泰勒公式展开

$$\begin{aligned} V(\boldsymbol{x}(t + \Delta t), t + \Delta t) &= V(\boldsymbol{x}(t), t) + \frac{dV}{dt} \Delta t + o(\Delta t)^2 \\ &= V(\boldsymbol{x}(t), t) + (\frac{\partial V}{\partial \boldsymbol{x}})^T \frac{d \boldsymbol{x}}{dt} \Delta t + \frac{\partial V}{\partial t} \Delta t + o(\Delta t)^2 \end{aligned}$$

两边 $V(\boldsymbol{x}(t))$ 抵消，有

$$0 = \min_{\boldsymbol{u}(t) \in \boldsymbol{U}} \int^{t+\Delta t}_{t} ( \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{u}^T \boldsymbol{R} \boldsymbol{u} ) dt + (\frac{\partial V}{\partial \boldsymbol{x}})^T \frac{d \boldsymbol{x}}{dt} \Delta t + \frac{\partial V}{\partial t} \Delta t + o(\Delta t)^2$$

现在令 $\Delta t \rightarrow 0$，对 t 求导有：

$$\begin{aligned} 0 &= \min_{\boldsymbol{u}(t) \in \boldsymbol{U}} \left[ ( \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{u}^T \boldsymbol{R} \boldsymbol{u} ) + (\frac{\partial V}{\partial \boldsymbol{x}})^T \frac{d \boldsymbol{x}}{dt} + \frac{\partial V}{\partial t} \right] \\ \frac{\partial V}{\partial t} &= - \min_{\boldsymbol{u}(t) \in \boldsymbol{U}} \left[ ( \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{u}^T \boldsymbol{R} \boldsymbol{u} ) + (\frac{\partial V}{\partial \boldsymbol{x}})^T (\boldsymbol{A}\boldsymbol{x} + \boldsymbol{B}\boldsymbol{u}) \right] \end{aligned}$$

我们不妨继续设 $V(\boldsymbol{x}(t), t) = \boldsymbol{x}^T(t) \boldsymbol{P}(t) \boldsymbol{x}(t)$，其中 $\boldsymbol{P}(t)$ 是对称矩阵。则 $\frac{\partial V}{\partial t} = \boldsymbol{x}^T \dot{\boldsymbol{P}} \boldsymbol{x}$，$\frac{\partial V}{\partial \boldsymbol{x}} = 2\boldsymbol{P}\boldsymbol{x}$。上式变为

$$\boldsymbol{x}^T \dot{\boldsymbol{P}} \boldsymbol{x} = - \min_{\boldsymbol{u}(t) \in \boldsymbol{U}} \left[ ( \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{u}^T \boldsymbol{R} \boldsymbol{u} ) + 2\boldsymbol{x}^T \boldsymbol{P}^T (\boldsymbol{A}\boldsymbol{x} + \boldsymbol{B}\boldsymbol{u}) \right]$$

令右式对u的偏导为零，

$$\begin{aligned} 2\boldsymbol{R}\boldsymbol{u}^* + 2\boldsymbol{B}^T \boldsymbol{P} \boldsymbol{x} = 0 \\ \boldsymbol{u}^* = - \boldsymbol{R}^{-1} \boldsymbol{B}^T \boldsymbol{P} \boldsymbol{x} \end{aligned}$$

设负反馈增益 $\boldsymbol{K}(t) = \boldsymbol{R}^{-1}(t) \boldsymbol{B}^T(t) \boldsymbol{P}(t)$，则 $\boldsymbol{u}^*(t) = - \boldsymbol{K}(t) \boldsymbol{x}(t)$。将最优控制代回，得到关于 $\boldsymbol{P}(t)$ 的微分方程：

$$\begin{aligned} -\boldsymbol{x}^T \dot{\boldsymbol{P}} \boldsymbol{x} &= \boldsymbol{x}^T \boldsymbol{Q} \boldsymbol{x} + \boldsymbol{x}^T \boldsymbol{K}^T \boldsymbol{R} \boldsymbol{K} \boldsymbol{x} + 2\boldsymbol{x}^T \boldsymbol{P}^T (\boldsymbol{A} - \boldsymbol{B}\boldsymbol{K}) \boldsymbol{x} \\ -\dot{\boldsymbol{P}} &= \boldsymbol{Q} + \boldsymbol{K}^T \boldsymbol{R} \boldsymbol{K} + 2 \boldsymbol{P}^T (\boldsymbol{A} - \boldsymbol{B}\boldsymbol{K}) \\ &= \boldsymbol{Q} + \boldsymbol{P}^T\boldsymbol{B}\boldsymbol{R}^{-1}\boldsymbol{B}^T\boldsymbol{P} + 2\boldsymbol{P}^T\boldsymbol{A} - 2\boldsymbol{P}^T\boldsymbol{B}\boldsymbol{R}^{-1}\boldsymbol{B}^T\boldsymbol{P} \\ &= \boldsymbol{Q} + 2\boldsymbol{P}^T\boldsymbol{A} - \boldsymbol{P}^T\boldsymbol{B}\boldsymbol{R}^{-1}\boldsymbol{B}^T\boldsymbol{P} \end{aligned}$$

即

$$\dot{\boldsymbol{P}} + 2\boldsymbol{P}\boldsymbol{A} - \boldsymbol{P}\boldsymbol{B}\boldsymbol{R}^{-1}\boldsymbol{B}^T\boldsymbol{P} + \boldsymbol{Q} = 0$$

这就是 $\boldsymbol{P}^T = \boldsymbol{P}$ 时的连续时间代数 Riccati 微分方程，负反馈增益为 $\boldsymbol{K}(t) = \boldsymbol{R}^{-1}(t) \boldsymbol{B}^T(t) \boldsymbol{P}(t)$。