企业官网建设流程全解析

最近把一篇关于离散时间非线性系统自适应动态规划（ADP）的经典文章复现做完了，趁热打铁记录一下。这篇博客针对的是标题里那类“非线性值迭代自适应动态规划”研究，目标是用Matlab实现离散时间非线性系统的策略迭代/值迭代算法。如果你正在做最优控制方向，或者被论文里的HJB方程、递增迭代序列、容许控制这些概念弄得焦头烂额，这篇笔记应该能帮你少走不少弯路。我会从算法原理出发，把值迭代和策略迭代的数学逻辑讲清楚，然后给出一套可以直接抄作业的Matlab实现方案，最后分享几个我实测踩过的坑。

1. 项目背景与核心问题拆解

1.1 非线性最优控制为什么难：HJB方程这堵墙

先聊点背景。学过线性系统最优控制的都知道LQR，给定一个线性系统，Q和R一选，Riccati方程一解，最优状态反馈增益就出来了，整个过程非常标准。但现实世界的系统几乎都是非线性的，倒立摆、机械臂、电力系统、飞行器，哪个都不是简单的dx/dt = Ax + Bu。对非线性系统做最优控制，理论上的解法是把Hamilton-Jacobi-Bellman方程解出来，得到一个关于状态的值函数。问题是，HJB方程是一个非线性偏微分方程，除了极少数特殊系统，根本不存在解析解。

动态规划倒是能从数值上逼近最优解，但它需要对连续状态空间逐点求值，每增加一个状态维度，计算量就指数膨胀，这就是所谓的“维数灾难”。ADP的出现就是冲着这个问题去的。它不再逐点求值，而是用一个参数化的函数结构——多项式、RBF基函数或者神经网络——去逼近值函数，把偏微分方程的求解问题转化为参数估计问题。这是整个ADP家族的核心思想，也是这篇文章复现的基础。

1.2 值迭代和策略迭代到底在迭代什么

标题里同时出现了“值迭代”和“策略迭代”，乍看有点乱，实际是同一类算法的两个分支。

值迭代的出发点是：随便给一个初始值函数，然后反复套用Bellman方程更新它，每轮更新都隐式地改进一次策略。它的优点是无需初始稳定策略，条件放宽很多。缺点是收敛速度通常比较慢，尤其当初始值函数与真实最优值差得很远时，前几十轮迭代基本都在“热身”。

策略迭代则是从一条明确的策略出发，交替执行策略评估和策略改进。策略评估是在固定策略下求解Bellman方程，得到对应的值函数；策略改进是根据新值函数更新控制律。它的收敛速度通常更快，但代价是需要一个初始的稳定策略——也就是能镇定系统的控制律，这一步往往不是白给的。

下面这个表格直观对比了两者的差异，复现之前最好心里有个数。

说句实在话，论文里这两种方法经常放在一起讲，不是要你二选一，而是帮你理解同一个最优控制问题的两条求解路径。复现时先把值迭代跑通，再改造成策略迭代，理解起来会顺很多。

2. 算法原理：两大迭代机制的数学逻辑

2.1 系统建模与性能指标

回到基础。考虑一类离散时间非线性系统，状态方程写成：

$$x_{k+1} = f(x_k) + g(x_k)u_k$$

其中$x_k \in \mathbb{R}^n$是状态，$u_k \in \mathbb{R}^m$是控制输入。这里要求$g(x)$已知，这在很多机械系统里是成立的，因为执行器往往以线性方式注入控制量。

性能指标用无限时域折扣或非折扣代价函数：

$$J(x_0, u) = \sum_{k=0}^{\infty} (x_k^T Q x_k + u_k^T R u_k)$$

其中$Q$和$R$是对称正定矩阵。最优控制问题的目标就是找到控制序列 $u_k = h(x_k)$，让这个代价最小。根据Bellman最优性原理，最优值函数满足：

$$V^(x_k) = \min_{u_k} { x_k^T Q x_k + u_k^T R u_k + V^(x_{k+1}) }$$

这就是离散时间的HJB方程。它把无限时域的最优性问题转化成一个函数方程。但$V^(x)$未知，而且嵌套在$V^(x_{k+1})$里，没法直接解。ADP的思路是：既然直接解不行，那就一层一层迭代逼近它。

2.2 值迭代更新方程与收敛逻辑

值迭代的更新规则很直接。设第$i$轮的值函数为$V_i(x)$，首先根据当前值函数求一个“贪心”控制：

$$h_i(x_k) = \arg\min_{u_k} { x_k^T Q x_k + u_k^T R u_k + V_i(x_{k+1}) }$$

然后用这个控制去更新值函数：

$$V_{i+1}(x_k) = x_k^T Q x_k + h_i^T(x_k) R h_i(x_k) + V_i(x_{k+1})$$

从任意非负初始函数$V_0$出发，在比较温和的条件下，这个序列会单调收敛到最优值函数$V^*$。这就是为什么我前面说“不需要初始稳定策略”。实际仿真中初始$V_0$常直接设为零，或者用一个随机的小权重向量。这样做的风险是前期迭代中策略质量很差，状态轨迹可能发散，但由于只是值函数在更新，数值上并不会崩溃，只要采样区域覆盖合理就能收回来。

需要留意的是，理论上的收敛性建立在“精确更新”的前提下。一旦引入函数近似，每次更新的投影误差就会累积。所以复现时控制逼近误差很关键，这也是后面调参的重点。

2.3 策略迭代的评估-改进循环

策略迭代走的是另一条路。假设当前策略为$h_i$，第一步是策略评估：把策略固定下来，求解Bellman方程

$$V_i(x_k) = x_k^T Q x_k + h_i^T(x_k) R h_i(x_k) + V_i(x_{k+1})$$

这是一个关于$V_i$的函数方程，比HJB方程简单，因为它去掉了对$u$的极小化操作。第二步是策略改进，利用新值函数更新控制：

$$h_{i+1}(x_k) = \arg\min_{u_k} { x_k^T Q x_k + u_k^T R u_k + V_i(x_{k+1}) }$$

两个步骤交替，直到策略不再变化。这个流程的关键条件：初始策略$h_0$必须是容许控制，也就是能让系统在该策略下稳定，同时使代价函数有限。如何获得这样一个初始策略？工程上最省事的办法就是用LQR先算一个线性反馈作为起点，只要系统的非线性项不是太强，这个反馈通常能镇定系统。

策略迭代与值迭代在“策略改进”这一步用的是同一个极小化操作，核心区别在于值函数如何更新。前者需要先评估当前策略，后者直接递推值函数。

3. Matlab仿真环境搭建与模型参数确定

3.1 复现需要准备的工具与环境

复现这套算法不需要多新的环境。我用的版本是Matlab R2022b，实测在R2018b上也能跑通。核心依赖只有几个：基本矩阵运算、优化工具箱的fminbnd函数（用于策略改进中求极小值）、控制系统工具箱的dlqr函数（用于生成初始线性策略）。如果你连控制系统工具箱都没有，手写Riccati方程迭代也完全可行，后面我会提一句。

神经网络工具箱不是必须的。虽然不少论文用BP神经网络做值函数逼近，但我在复现时采用了“线性参数化基函数 + 最小二乘”的方案，也就是把值函数写成

$$V(x) = W^T \phi(x)$$

其中$\phi(x)$是选定的基函数向量，$W$是待定权重。这个方案的好处是权重更新可以直接用最小二乘求解，绕开了BP训练的收敛震荡问题。很多学生一上来就上工具箱，被trainlm的学习率、动量项折腾半天，其实完全没必要。对复现论文而言，线性参数化是最稳的起点。

3.2 仿真系统模型设计

我采用的仿真对象是一个近线性系统配非线性扰动，这种模型在很多ADP论文中出现过，既能体现非线性特征，又方便用LQR生成初始稳定策略。具体状态方程：

$$x_{k+1} = A x_k + B u_k + \Delta(x_k)$$

其中

$$A = \begin{bmatrix} 0.8 & 0.1 \ 0.05 & 0.9 \end{bmatrix}, \quad B = \begin{bmatrix} 0 \ 1 \end{bmatrix}, \quad \Delta(x) = \begin{bmatrix} 0.1 \sin(x_1) \ 0.05 x_2 \sin^2(x_1) \end{bmatrix}$$

性能指标中取$Q=I_{2\times 2}$，$R=1$。系统有两个状态，单控制输入。$\Delta(x)$是带有界扰动项，满足Lipschitz条件，保证系统在非线性扰动下仍然有良好行为。采样范围取$x_1, x_2 \in [-2, 2]$，网格步长0.2，一共401个状态点，这个规模在Matlab单次循环中大概几十毫秒能算完，调试时非常友好。

为什么选这个模型？关键因素是$B$向量使得系统关于控制输入是线性的，这样在策略改进步骤中，对每个状态点求最优$u$的问题就退化成一个单变量极值问题，fminbnd处理起来很稳定。如果系统在$u$上也是非线性，就需要更复杂的优化手段，复现难度会陡增。

4. 核心代码实现与关键环节拆解

4.1 基函数选取与值函数参数化

线性参数化方案里，基函数的选择直接决定逼近精度。我测试了两类基函数：多项式基函数和高斯RBF基函数。多项式形式简单，比如二次型加交叉项，但在状态空间跨度大时容易产生数值不稳定。RBF形式更稳，表达式为

$$\phi_j(x) = \exp \left( -\frac{|x - c_j|^2}{\sigma^2} \right)$$

其中$c_j$是中心点，$\sigma$是宽度参数。我在$[-2, 2]^2$区域均匀布了9个中心点（3x3网格），$\sigma$取1.2。这个配置下，即使状态跑到采样边界之外，RBF的输出也不会像多项式那样极速膨胀，安全性好得多。

权重个数上，我用了9个基函数加上一个常数项，总共10个参数。值函数写成矩阵形式$V_W(x) = W^T \phi(x)$，这个参数规模对当前二维系统来说完全够用。如果你复现的是更高维系统，经验是：基函数个数要大于状态维数的5到10倍，中心点覆盖采样区域，宽度不能太小否则容易造成“孤岛效应”。

4.2 值迭代MATLAB核心实现

先给出值迭代的主循环。初始化$V_0 = 0$，对应$W$全零，然后在每个采样状态点$X(:, j)$上求优化控制，再更新值函数。

% 采样状态点 [X1, X2] = meshgrid(-2:0.2:2, -2:0.2:2); X = [X1(:)'; X2(:)']; N = size(X, 2); % RBF基函数初始化 centers = [-1.5 -0.5 0.5 1.5 -1.5 -0.5 0.5 1.5 -0.5; -1.5 -1.5 -1.5 -1.5 0.5 0.5 0.5 0.5 1.5]; % 简化的9个中心 sigma = 1.2; phi = @(x) [1; exp(-sum((x - centers).^2, 1)' / sigma^2)]; % 值迭代参数 tol = 1e-5; maxIter = 500; errList = zeros(maxIter, 1); W = zeros(10, 1); for iter = 1:maxIter V_new = zeros(1, N); u_new = zeros(1, N); for j = 1:N x = X(:, j); % 策略改进：对当前状态求最优控制 fun = @(u) x'*Q*x + u'*R*u + W' * phi(f(x) + g(x)*u); u_opt = fminbnd(fun, -5, 5); u_new(j) = u_opt; V_new(j) = fun(u_opt); end % 最小二乘更新权重：W' * PhiPhi ≈ V_new PhiMat = zeros(10, N); for j = 1:N PhiMat(:, j) = phi(X(:, j)); end W_new = (PhiMat * PhiMat') \ (PhiMat * V_new'); % 收敛判据 err = max(abs(W_new - W)); errList(iter) = err; W = W_new; if err < tol fprintf('值迭代收敛于第 %d 轮，误差 %e\n', iter, err); break; end end

这里有两个关键细节需要展开说。第一，fminbnd优化时边界选$[-5, 5]$，这个边界要在仿真初期就测试好。边界太大浪费计算，太小可能漏掉最优控制点。4.2中第二处关键是权重更新那行：$W = (\Phi \Phi^T)^{-1} \Phi V_{new}^T$，这是典型的最小二乘正规方程。在实际运行中，如果$\Phi \Phi^T$接近奇异，需要用pinv代替反斜杠。我测试的9个RBF中心在网格上位置分布均匀，矩阵条件数在可接受范围内，直接用反斜杠没问题。

这段代码运行完，W就是值函数网络参数。接着用这个值函数做闭环仿真，验证控制效果。

4.3 策略迭代MATLAB核心实现

策略迭代的代码结构有些不同，核心在于策略评估和策略改进两阶段。初始策略我直接用dlqr基于线性部分$A, B$求$K_0$，保证系统在初始策略下稳定。策略评估阶段需要解关于值函数权重的线性方程组。

% 初始线性策略 K = dlqr(A, B, Q, R); h_fun = @(x, K) -K * x; % 策略迭代参数 tol = 1e-5; maxIterPI = 100; W_pi = zeros(10, 1); for iter = 1:maxIterPI % ---- 策略评估：固定当前策略，求 W 满足 % V(x) = x'Qx + u'Ru + V(f(x)+g(x)u) A_mat = zeros(N, 10); b_vec = zeros(N, 1); for j = 1:N x = X(:, j); u = h_fun(x, K); x_next = f(x) + g(x) * u; A_mat(j, :) = phi(x)' - phi(x_next)'; b_vec(j) = x'*Q*x + u'*R*u; end W_pi = (A_mat' * A_mat) \ (A_mat' * b_vec); % ---- 策略改进：更新控制增益（这里用神经网络的局部线性化近似） K_new = update_gain_from_value(X, W_pi, phi, f, g, Q, R, K); gain_change = norm(K_new - K); K = K_new; fprintf('策略迭代第 %d 轮，增益变化 %e\n', iter, gain_change); if gain_change < 1e-6 break; end end

策略改进部分的实现，我采用的思路是：在每个采样点计算出贪心控制$u^(x) = \arg\min_u { x^TQx + u^TRu + W^T \phi(f(x)+g(x)u) }$，得到一组$(x, u^)$数据，然后用最小二乘拟合一个线性反馈$u = -Kx$。这样做的好处是保持策略的线性结构，策略迭代过程非常稳定。如果你希望策略提高表达能力，可以用RBF网络去拟合任意非线性策略，但要注意避免过拟合。

4.4 主流程与闭环仿真验证

学习完成后，需要做闭环仿真验证。从初始状态$x_0 = [1.5; -1]$出发，用学到的控制律跑50个时间步，记录状态和控制轨迹。

% 闭环仿真 x0 = [1.5; -1]; Tmax = 50; x_traj = zeros(2, Tmax+1); u_traj = zeros(1, Tmax); x_traj(:, 1) = x0; for k = 1:Tmax x = x_traj(:, k); % 使用值迭代学到的策略（min_u）或策略迭代的线性反馈 fun = @(u) x'*Q*x + u'*R*u + W' * phi(f(x) + g(x)*u); u_k = fminbnd(fun, -5, 5); u_traj(k) = u_k; x_traj(:, k+1) = f(x) + g(x) * u_k; end % 绘图 figure; subplot(2, 1, 1); plot(0:Tmax, x_traj(1, :), 'b-', 'LineWidth', 1.5); hold on; plot(0:Tmax, x_traj(2, :), 'r--', 'LineWidth', 1.5); grid on; legend('x1', 'x2'); title('闭环状态轨迹'); subplot(2, 1, 2); plot(0:Tmax-1, u_traj, 'k-', 'LineWidth', 1.5); grid on; title('控制轨迹');

把闭环轨迹画出来，这是最能直接反映算法好坏的图。

5. 仿真结果分析与收敛性验证

5.1 判断算法是否真正收敛了

复现时最容易犯的错误是：迭代误差降下去了，就宣布收敛，结果闭环仿真一塌糊涂。判断收敛不能只看单一指标，至少要看三个信号。

第一个信号是相邻迭代的权重差。程序中我用$err = \max |W_{new} - W|$作为停止条件。这个值反映的是参数层面的稳定程度，权重变化小说明值函数逼近基本饱和。

第二个信号是值函数曲面形态。把学到的值函数画成三维曲面，最优点应该是一个光滑的“漏斗”形状，底部在原点附近，没有任何异常的尖峰或者凹陷。如果曲面出现大范围震荡，基本可以判定RBF宽度或中心点设置有问题。

第三个信号是闭环状态轨迹。稳定收敛的控制策略应当让状态快速收敛到原点，控制量不过大、不震荡。这三个信号同时满足，才能说算法真的收住了。我在实验中发现，有时权重误差已经很小，但状态轨迹仍然有缓慢漂移，原因是最优值函数在原点邻域内的曲率没有被充分刻画。解决办法是增加原点附近的采样点密度，或者在基函数里加二次型硬编码项。

5.2 三张必出的仿真图与判读方法

复现论文时，审稿人或者导师最关心的图一般有三类。第一类是收敛曲线图，横轴迭代轮次，纵轴权重误差或值函数误差，这条曲线应当单调下降然后趋于平缓。第二类是最终的闭环状态轨迹图，反映控制器对系统的镇定能力。第三类是控制输入图，检查控制量是否在合理幅度内。

如果做对比实验，建议把值迭代和策略迭代的收敛曲线叠在一张图里。策略迭代通常会在更少的轮次内收敛，但每一轮的计算量更大，因为策略评估需要解一次线性方程组。值迭代虽然轮次多，但单轮计算简单，整体时间未必差很多。把两种算法的运行时间和收敛轮次记录下来，对比会非常有说服力。

5.3 参数敏感性：哪些参数决定成败

参数选择对ADP算法的影响是决定性的。我整理了一个快速参考表，所有结论都来自实测。

这里面最容易出问题的参数是$\sigma$。$\sigma$太小时，每个RBF中心只能影响极小区域，基函数之间没有重叠，值函数曲面会出现一个个“鼓包”；$\sigma$太大时，所有基函数几乎一样，逼近能力直线下降。一个实用的调参方法是：先固定中心点，把$\sigma$从0.5开始逐渐增大，观察拟合误差的变化，找到误差平台期的起始点。

6. 常见问题与调参经验

6.1 值迭代收敛慢或不收敛怎么办

值迭代收敛慢往往是三个原因：初始值函数离最优解太远、采样范围过大导致无效点占比高、fminbnd的搜索边界太窄。我的经验是把采样区域分阶段考虑：先用较粗网格跑通流程，再用细网格精调。不要一上来就是密密麻麻的状态点，几百个点循环一次都要半天，排错效率太低了。

另一个常见问题是不收敛。如果你看到误差曲线在某个数值附近反复震荡，首先检查值函数更新里是不是忘了包含代价项$x^TQx$。这个错误非常隐蔽，因为程序不会报错，但迭代就是降不下去。其次检查最小二乘的矩阵条件数，条件数飙升时用pinv代替反斜杠，问题通常会立刻缓解。

6.2 策略迭代初始策略不稳定的问题

策略迭代对初始策略的稳定性极其敏感。我最初偷懒用全零控制作为初始策略，结果策略评估阶段解出的值函数直接爆炸，数值上出现NaN。后来老老实实先用LQR求线性反馈K0，仿真才顺利跑通。

如果你的系统非线性太强，LQR无法镇定系统，可以改用状态相关Riccati方程（SDRE）生成初始策略。SDRE逐点求解线性矩阵方程，非线性项被当作“冻结”状态进入矩阵，基本能处理大部分工程系统。再不行，可以先人为给一个阻尼项，比如$u = -Kx - \alpha x^3$，让系统先稳定下来再做迭代。

6.3 网络逼近过拟合与欠拟合的表现

线性参数化RBF网络也会过拟合或者欠拟合。欠拟合的表现是：值函数曲面过于平滑，细节特征丢失，闭环轨迹收敛很慢。过拟合的表现是：值函数在采样点上完全贴合，但采样点之间剧烈跳变，闭环仿真时控制量来回振荡。判断方法是随机抽出一部分状态点做验证集，比较验证误差和训练误差。

缓解过拟合的办法很简单：增加采样点密度、适当增大$\sigma$、减少RBF中心数量。缓解欠拟合则相反。很多时候，你只需要调整中心点的布局方式，把网格中心改成“中间密、四周疏”的非均匀分布，就能在不增加参数数量的前提下提升逼近效果。

6.4 Matlab实现路上的三个坑

最后说几个Matlab特有的坑，每一个我都亲自踩过。

第一个坑是关于fminbnd的：这个函数默认使用黄金分割搜索，要求目标函数在边界内是单峰的。我们的代价函数通常不是严格凸的，如果最优解落在边界上，fminbnd返回的结果可能不可靠。处理方式是把边界放宽，或者在搜索后检查解是否在边界上，如果是，扩大边界重新搜索。

第二个坑是反斜杠运算与伪逆的取舍。矩阵$\Phi\Phi^T$条件数过大时，反斜杠结果数值不稳定。我建议在复现阶段一律改成pinv，虽然稍微慢一点，但换来的是稳定的收敛行为。

第三个坑是要特别注意避免变量名覆盖内置函数。我曾经把存储状态网格的变量命名为data，结果无意中覆盖了某工具箱的数据类引用，程序到处报错。后来养成习惯，所有自定义变量都用带前缀的名字，比如xGrid、WCritic、KPolicy，从此再没犯过这种低级错误。

最后分享一个我个人的习惯。复现这类文章时，我会刻意把值迭代和策略迭代放在同一个脚本里，共享同一套基函数和采样点。这样做的好处是，对比实验时只需要切换一个开关，就能同时看到两种算法的收敛曲线和闭环效果。这套代码跑通之后，换任何系统模型只需要改f和g函数，基函数自动按新状态的取值范围重新生成中心。做ADP方向的研究生，这套代码基本能应付绝大多数论文复现任务。