1. 项目概述：为什么LMI工具箱是解决复杂优化问题的利器

在控制系统设计、鲁棒滤波、信号处理乃至金融工程等领域，我们常常会遇到一类“棘手”的优化问题：目标函数或约束条件中包含了矩阵不等式。这类问题，传统的基于梯度的优化算法（如fmincon）往往力不从心，要么难以求解，要么得到的解可靠性存疑。这时，线性矩阵不等式（Linear Matrix Inequality, LMI）理论及其在MATLAB中的实现——LMI工具箱，就成为了工程师和研究员手中的一把瑞士军刀。

简单来说，LMI就是一个关于矩阵变量的线性不等式，其形式通常为 F(x) = F0 + x1F1 + ... + xnFn < 0，其中“< 0”表示矩阵是负定的。许多复杂的系统分析与设计问题，如判断系统稳定性、求解线性二次型最优控制（LQR）的Riccati方程、设计H∞控制器等，都可以转化为一个或多个LMI的可行性问题或优化问题。MATLAB的LMI工具箱，正是为了高效、可靠地求解这类问题而生的。

然而，工具箱的强大也伴随着一定的入门门槛。其核心函数lmivar和lmiterm的语法相对独特，需要使用者以一种“声明式”的思维来构建问题，这与我们熟悉的命令式编程有所不同。网上很多资料要么过于理论化，要么示例零散，让初学者望而却步。本文将从一个实践者的角度，手把手带你拆解这两个核心函数，并通过一个完整的控制器设计案例，展示如何从问题描述到代码实现，最终求解出可靠的解。无论你是正在做课题的研究生，还是需要解决实际工程问题的工程师，这篇超详细的教程都将为你提供一条清晰的路径。

2. LMI问题构建的核心思维：从数学描述到工具箱语言

在动手写代码之前，我们必须先建立正确的思维模型。LMI工具箱的求解流程，类似于我们向一个“求解器”提交一份“问题说明书”。这份说明书需要明确两点：变量和约束。lmivar函数负责定义变量，lmiterm函数负责描述约束中的每一项。整个构建过程是“分而治之”的：先定义所有矩阵变量，再针对每一个LMI约束，逐一累加其构成项。

2.1 理解LMI的标准形式

LMI工具箱处理的问题一般形式为： 找到矩阵变量 X1, X2, ..., Xk，使得一组线性矩阵不等式成立： N1^T * L(X1, ..., Xk) * N1 < M1^T * R(X1, ..., Xk) * M1 N2^T * L(X1, ..., Xk) * N2 < M2^T * R(X1, ..., Xk) * M2 ... 其中 L(.) 和 R(.) 是矩阵变量的仿射函数。在大多数基础应用中，我们遇到的是更简单的形式：A^T * X * A - X + Q < 0或A*X + X*A^T + B*B^T < 0等。

工具箱要求我们将这些不等式，都以“小于零”的形式，写在等号的同一侧。例如，不等式P > 0（P正定）等价于-P < 0。不等式A^T*P + P*A + Q < 0则已经是在“小于零”的一侧了。

2.2 变量定义函数 lmivar 详解

lmivar用于定义LMI问题中的未知矩阵变量。其基本调用语法为：[var_id, n, s] = lmivar(type, struct)

• var_id: 输出参数，是该变量的句柄（ID号），在后续的lmiterm中通过此ID来引用该变量。
• n: 输出参数，当type=1时，表示变量的维度；当type=2或3时，表示该变量在最终决策变量向量中的块大小（可暂时忽略）。
• type: 输入参数，定义变量类型，是关键所在。
  • type = 1: 对称块对角矩阵变量。这是最常见的情况，用于定义对称矩阵（如Lyapunov函数中的P矩阵）。struct参数是一个 k×2 的矩阵，每一行描述一个对角块。例如：
    • struct = [n, 0]: 表示一个 n×n 的满对称块。
    • struct = [n, 1]: 表示一个 n×n 的标量矩阵（即单位矩阵乘以一个标量变量）。
    • struct = [n, -1]: 表示一个 n×n 的全零块（占位，无变量）。
  • type = 2: 长向量变量（即标准的决策向量 x）。struct = [m, n]表示一个 m×n 的矩形矩阵变量，它将被按列拉直成一个长向量。这种类型在将传统优化问题转化为LMI形式时有用。
  • type = 3: 其他结构。struct是一个与变量维度相同的矩阵，其中每个元素只能是 0, x, -x。0表示该位置固定为0；x表示该位置是一个新的标量决策变量；-x表示该位置是-1乘以一个已由x定义的变量（用于强制对称性或特殊结构）。这种类型最灵活，但也最复杂。

实操心得一：变量定义的策略对于初学者，我的建议是优先使用type=1来定义所有对称矩阵变量。除非你非常确定需要type=2或3，否则它们容易引入错误。在定义type=1时，想清楚你的矩阵是“完全自由的对称矩阵”（用[n, 0]）还是“单位矩阵的倍数”（用[n, 1]）。后者可以显著减少变量个数，加快求解速度，适用于某些特定结构（例如，你需要寻找一个公共的Lyapunov矩阵比例因子）。

2.3 项描述函数 lmiterm 详解

lmiterm用于向指定的LMI约束中添加一项。这是构建LMI最核心、也最容易出错的一步。其语法为：lmiterm(termID, A, B, flag)

• termID: 一个四元向量[p, q, x, y]，它指明了当前项添加到哪个位置。
  • p: LMI的编号。正数表示第p个“小于零”不等式（< 0），负数表示第p个“大于零”不等式（> 0，内部会转换为-LMI < 0）。我们通常全用正数，手动处理正负号更清晰。
  • q: 项在LMI中的块行和块列索引。通常设为0，表示该项作用于整个矩阵（1x1块）。只有在处理分块矩阵LMI时才会用到非零值。
  • x和y: 矩阵变量ID。x指定包含变量的矩阵，y指定其转置（如果涉及）。具体规则见下。
• A,B: 系数矩阵。它们与变量（或其转置）相乘。如果项中不包含变量，则A或B中有一个是标量常数（通常是1）。
• flag: 可选字符串，默认为's'（对称项）。当添加形如X*A或A*X的项时，如果希望自动添加其转置项以保持整体对称性，则使用's'。对于常数项或明显不对称的项，使用'n'（非对称）。

x和y的取值规则（重中之重！）：

• 如果项是常数矩阵（如Q），则x=0, y=0。此时A就是该常数矩阵（或标量）。
• 如果项是变量矩阵（如P），则x为该变量的ID（var_id），y=0。此时A通常是左乘的系数矩阵，若A=1则表示变量本身。
• 如果项是变量矩阵的转置（如P*A中的P如果写在左边，但实际是A'*P的一部分），则x=0, y为该变量的ID。此时B是右乘的系数矩阵。
• 如果项是两个变量相乘（较少见），则x和y分别为两个变量的ID。

核心技巧：拆解项与符号处理lmiterm一次只添加一个“项”。一个复杂的表达式需要拆成多个项。例如，表达式A'*P + P*A包含两个项：A'*P和P*A。我们需要分别添加它们。所有项都必须放在不等式的一侧。例如，对于 Lyapunov 不等式A'*P + P*A < -Q，我们需要将其改写为A'*P + P*A + Q < 0。那么，在代码中我们就需要添加三项：A'*P,P*A,Q。

实操心得二：lmiterm的“分项累加”思维我习惯像搭积木一样构建LMI。首先，在纸上或注释里把目标不等式整理成Sum(terms) < 0的形式。然后，对求和号里的每一项，单独写一个lmiterm。每写一项，都检查其termID是否指向正确的LMI编号，以及x,y赋值是否符合上述规则。一个非常有效的调试方法是：在构建完所有LMI后，用lmiedit命令打开GUI查看器，直观地检查你构建的矩阵是否正确。

3. 完整案例：连续系统状态反馈控制器设计与求解

现在，我们通过一个经典问题——连续系统状态反馈镇定——来串联整个流程。问题描述：给定一个线性时不变系统dx/dt = A*x + B*u，设计状态反馈控制律u = K*x，使得闭环系统dx/dt = (A+B*K)*x渐近稳定。这可以转化为一个LMI可行性问题。

步骤1：问题转化为LMI根据Lyapunov稳定性理论，寻找一个对称正定矩阵P和一个矩阵K，使得闭环系统满足：(A+B*K)' * P + P * (A+B*K) < 0这个不等式关于P和K不是线性的（因为包含了P*B*K和K'*B'*P）。为了将其线性化，我们引入一个变量替换。令Y = K * X，其中X = inv(P)。对上式左右同时乘以X（合同变换），并代入Y = K*X，可以得到一个关于X和Y的LMI：A*X + X*A' + B*Y + Y'*B' < 0同时，P > 0等价于X > 0。 因此，我们的LMI问题为： 寻找对称矩阵X和矩阵Y，使得：

1. A*X + X*A' + B*Y + Y'*B' < 0
2. X > 0(即-X < 0) 求解得到X和Y后，控制器增益为K = Y * inv(X)。

步骤2：MATLAB代码实现

% 案例：连续系统状态反馈控制器设计 clear; clc; % 1. 定义系统矩阵 (不稳定系统) A = [1 2; -1 0]; B = [0; 1]; n = size(A, 1); % 状态维度 m = size(B, 2); % 输入维度 % 2. 初始化LMI系统描述 setlmis([]); % 开始描述LMI系统 % 3. 定义矩阵变量 % X 是 n x n 对称正定矩阵 (type=1, 满块) [X_id, nX, sX] = lmivar(1, [n, 1]); % struct=[n,1] 表示标量矩阵？这里错了！应该是[n,0] % 更正：对于自由对称矩阵P，应使用 [n, 0] [X_id, nX, sX] = lmivar(1, [n, 0]); % Y 是 m x n 的矩形矩阵 (type=2) [Y_id, nY, sY] = lmivar(2, [m, n]); % 4. 定义第一个LMI: A*X + X*A' + B*Y + Y'*B' < 0 lmi_id = 1; % 第一个不等式编号 % 项1: A*X lmiterm([lmi_id, 1, 1, X_id], A, 1); % [1,1, X_id]: 添加到LMI1的(1,1)块，项为 A * X * 1 % 注意：这里 A 左乘 X，所以 A 是系数，放在第二个参数位置。1表示右乘单位阵。 % 项2: X*A' lmiterm([lmi_id, 1, 1, X_id], 1, A', 's'); % 使用's'标志，自动添加对称部分 (X*A') + (A*X') % 实际上，因为项1已经加了A*X，这里用's'会自动补上其转置(A*X)' = X'*A' = X*A' (X对称)。 % 更稳妥的写法是分别添加两项，但's'在这里更简洁且不易出错。 % 项3: B*Y lmiterm([lmi_id, 1, 1, Y_id], B, 1); % B * Y * 1 % 项4: Y'*B' lmiterm([lmi_id, 1, 1, Y_id], 1, B', 's'); % 使用's'自动补全 (B*Y) + (Y'*B') % 5. 定义第二个LMI: X > 0 (等价于 -X < 0) lmi_id = 2; lmiterm([lmi_id, 1, 1, X_id], -1, 1); % 添加项：-1 * X * 1 < 0 => X > 0 % 6. 完成LMI系统描述 lmisys = getlmis; % 7. 求解LMI可行性问题 [tmin, xfeas] = feasp(lmisys); % 8. 检查可行性并提取解 if tmin < 0 % tmin是可行性测度，小于0表示可行 disp('LMI可行！'); % 从解向量xfeas中提取矩阵变量 X_sol = dec2mat(lmisys, xfeas, X_id); Y_sol = dec2mat(lmisys, xfeas, Y_id); % 计算控制器增益 K = Y * X^{-1} K = Y_sol / X_sol; % 等价于 Y_sol * inv(X_sol) disp('求解得到的 X:'); disp(X_sol); disp('求解得到的 Y:'); disp(Y_sol); disp('状态反馈增益矩阵 K:'); disp(K); % 验证闭环系统稳定性 A_cl = A + B*K; eig_cl = eig(A_cl); disp('闭环系统特征值:'); disp(eig_cl); if all(real(eig_cl) < 0) disp('闭环系统稳定！'); else disp('警告：闭环系统不稳定，请检查求解结果。'); end else disp('未找到可行解。tmin = '); disp(tmin); end

步骤3：代码关键点解析

1. setlmis([])与getlmis: 这是构建LMI问题的固定框架。setlmis([])初始化一个空的LMI系统，之后所有的lmivar和lmiterm都在向这个系统添加内容。getlmis最终获取完整的LMI系统描述对象lmisys，用于后续求解。
2. feasp函数: 这是求解LMI可行性问题的核心求解器。它尝试找到一组决策变量，使得所有LMI约束都被满足。返回值tmin是最小化约束违背的量，tmin < 0意味着找到了可行解（所有LMI严格小于零）。xfeas是决策变量的解向量。
3. dec2mat函数: 求解器返回的解xfeas是一个向量。我们需要用dec2mat函数，根据之前定义的变量结构（lmisys），将其还原成我们熟悉的矩阵形式。参数依次为：LMI系统、解向量、变量ID。
4. 符号处理: 注意第二个LMIX > 0，我们是通过添加项-X到< 0的不等式中实现的。这是LMI工具箱的标准做法：将所有约束统一为“小于零”形式。

4. 高级应用与性能优化技巧

掌握了基础构建方法后，我们可以处理更复杂的问题，并优化求解过程。

4.1 处理多个LMI约束与分块矩阵

有时一个问题包含多个LMI，或者一个LMI本身是分块矩阵。lmiterm的termID中第二、三个参数p, q就派上了用场。

示例：带有性能约束的H∞状态反馈问题：设计u = Kx，使得闭环系统满足||Tzw||∞ < γ。这可以转化为以下LMI组（简化形式）： 寻找X > 0,Y, 使得：

[ A*X + X*A' + B*Y + Y'*B' Bw (Cz*X+Dzu*Y)' ] [ Bw' -γ*I Dzw' ] < 0 [ Cz*X+Dzu*Y Dzw -γ*I ]

以及X > 0。 这是一个2x2的分块矩阵LMI（实际上左上角块是另一个矩阵不等式）。

% ... 假设已定义 A, B, Bw, Cz, Dzu, Dzw, gamma ... setlmis([]); [X_id, ~, ~] = lmivar(1, [n, 0]); [Y_id, ~, ~] = lmivar(2, [m, n]); lmi_id = 1; % H∞性能LMI % 块(1,1): A*X + X*A' + B*Y + Y'*B' lmiterm([lmi_id, 1, 1, X_id], A, 1); lmiterm([lmi_id, 1, 1, X_id], 1, A', 's'); lmiterm([lmi_id, 1, 1, Y_id], B, 1, 's'); % 's' 会自动添加 B*Y + Y'*B' % 块(1,2): Bw lmiterm([lmi_id, 1, 2, 0], Bw); % 常数项，x=0,y=0 % 块(1,3): (Cz*X + Dzu*Y)' % 先添加 Cz*X lmiterm([lmi_id, 3, 1, X_id], Cz, 1); % 注意：这项在(3,1)位置，是(1,3)的转置 % 再添加 Dzu*Y lmiterm([lmi_id, 3, 1, Y_id], Dzu, 1); % 块(2,1): Bw' (是块(1,2)的转置，由工具箱自动保证对称性，通常只需定义上三角或下三角) % 因为我们用's'或完整定义，这里通常只需定义(1,2)，对称部分会自动处理。但为清晰，也可定义： % lmiterm([lmi_id, 2, 1, 0], Bw'); % 块(2,2): -gamma*I lmiterm([lmi_id, 2, 2, 0], -gamma, eye(size(Bw,2))); % 块(2,3): Dzw' lmiterm([lmi_id, 2, 3, 0], Dzw'); % 块(3,1): (已定义) % 块(3,2): Dzw lmiterm([lmi_id, 3, 2, 0], Dzw); % 块(3,3): -gamma*I lmiterm([lmi_id, 3, 3, 0], -gamma, eye(size(Cz,1))); % 第二个LMI: X > 0 lmi_id = 2; lmiterm([lmi_id, 1, 1, X_id], -1, 1); % ... 后续求解步骤同上

在定义分块矩阵时，关键是理清每个块的位置(p, q)。工具箱不强制要求定义对称位置，但为了清晰和避免遗漏，建议至少定义完上三角部分，并利用's'标志处理对称项。

4.2 优化问题求解：mincx 与 gevp

除了可行性问题feasp，LMI工具箱还能解决两类优化问题：

1. 线性目标最小化 (mincx)：最小化c' * x，满足LMI约束。其中c是用户定义的向量，x是决策变量向量。这常用于最小化线性组合的变量，例如最小化矩阵的迹trace(X)。
   • 关键：需要用defcx函数来定义目标函数中的系数向量c。defcx(lmisys, k, X_id, i, j)用于设置与变量X_id的(i,j)位置元素相关的系数。对于最小化迹，通常是对角线元素系数设为1。

   % 在 getlmis 之后，定义目标为最小化 trace(X) n = decnbr(lmisys); % 获取决策变量数量 c = zeros(n, 1); for j = 1:n [Xj, ~, ~] = defcx(lmisys, j, X_id, j, j); % 获取X_id第(j,j)个元素在决策向量中的索引 c(Xj) = 1; % 设置该位置系数为1 end options = [1e-2, 100, 1e5, 10, 0]; % 优化选项 [copt, xopt] = mincx(lmisys, c, options);

2. 广义特征值问题 (gevp)：最小化标量 λ，使得满足 LMI1(x) < 0 且 LMI2(x) < λ * LMI3(x)。这常用于求解 H∞ 范数（最小化 γ）。
   • 调用格式：[lopt, xopt] = gevp(lmisys, nlfc, options)。nlfc是第一个LMI约束的个数（即LMI1）。

实操心得三：求解器选项与数值稳定性默认的求解器选项可能不适用于所有问题。feasp,mincx,gevp都可以通过options向量调整参数，如初始步长、最大迭代次数、精度等。

options = [1e-2, 200, 1e8, 10, 0]; % 参数含义：[精度, 最大迭代次数, 可行性半径, 迭代显示频率, 目标值] [tmin, xfeas] = feasp(lmisys, options);

如果求解失败（tmin > 0），可以尝试：

• 放宽可行性半径（增大options(3)）。
• 检查LMI问题的标度。如果矩阵A,B的元素数量级差异巨大（如1e-9和1e3），会导致数值问题。尽量对系统模型进行归一化处理。
• 确保问题本身是可行的。有时需要调整性能指标 γ 或初始猜测。

5. 调试技巧与常见问题排查

即使思路正确，构建LMI代码时也极易出错。以下是我在实践中总结的排查清单。

5.1 问题排查速查表

5.2 必备调试工具：lmiedit

lmiedit(lmisys)是LMI工具箱中最强大的调试工具，没有之一。它会打开一个GUI窗口，以矩阵块的形式展示你定义的所有LMI。你可以清晰地看到：

• 每个LMI的维度。
• 每个块是由哪些项组成的。
• 每个项是常数、变量还是变量乘积。
• 变量的位置和系数。

当你的求解结果不符合预期时，第一件事就是打开lmiedit，逐项核对是否与你纸上推导的公式一致。这能解决90%的构建错误。

5.3 数值验证：用解回代验证LMI

求解得到xfeas后，不要仅仅相信tmin < 0。应该将解出的矩阵变量代回原始LMI，计算其最大特征值，确保其为负。

% 假设有两个LMI，解已存储在 xfeas 中 X = dec2mat(lmisys, xfeas, X_id); Y = dec2mat(lmisys, xfeas, Y_id); % 验证第一个LMI: A*X + X*A' + B*Y + Y'*B' < 0 LMI1_left = A*X + X*A' + B*Y + Y'*B'; eig1 = eig(LMI1_left); max_eig1 = max(real(eig1)); fprintf('第一个LMI左端最大特征值: %.4e\\n', max_eig1); % 验证第二个LMI: X > 0 (即最小特征值>0) min_eig_X = min(real(eig(X))); fprintf('矩阵X的最小特征值: %.4e\\n', min_eig_X); if max_eig1 < 0 && min_eig_X > 0 disp('所有LMI约束均严格满足。'); else disp('警告：解不严格满足LMI约束！'); end

这种验证能给你十足的信心，确认整个建模和求解流程是正确的。

6. 从理论到实践：扩展应用场景与思路

LMI的应用远不止于控制器设计。一旦掌握了lmivar和lmiterm的思维，你可以将其应用到众多领域。

1. 观测器（滤波器）设计：对于系统dx/dt = A*x + B*u, y = C*x，设计全维状态观测器d\hat{x}/dt = A*\hat{x} + B*u + L*(y-C*\hat{x})。误差动态系统为de/dt = (A-L*C)*e。通过类似的变量替换W = P*L，可以转化为关于P和W的LMIA*P + P*A' - W*C - C'*W' < 0, P>0，进而解得L = P^{-1}*W。

2. 饱和控制系统分析与设计：考虑执行器饱和sat(u)。利用扇形条件或多面体描述，可以将饱和非线性转化为一组线性微分包含，进而用LMI来估计吸引域或设计抗饱和补偿器。

3. 时滞系统稳定性分析：对于含有时滞的系统，利用Lyapunov-Krasovskii泛函或Lyapunov-Razumikhin函数，可以将稳定性条件转化为包含时滞项的LMI，通常需要引入一些松弛矩阵变量来降低保守性。

4. 组合优化中的半定规划（SDP）：图论中的最大割问题、传感器网络定位等，可以建模为半定规划，这正是LMI工具箱 (mincx) 可以求解的问题类型。

个人体会与最后建议LMI工具箱的学习曲线前期比较陡峭，核心难点在于从连续的数学等式/不等式到离散的、声明式的代码指令的思维转换。我的经验是：永远从最简单的、有解析解的例子开始。例如，对于稳定的矩阵A，Lyapunov方程A'*P+P*A=-Q的解P一定是正定的。你可以先用lyap函数求解，再用LMI工具箱去求解A'*P+P*A < 0, P>0，对比结果是否一致。这种“双轨验证”是建立代码信心的最佳方式。

另外，不要畏惧lmiedit。在初期，每写几行lmiterm就运行一次getlmis和lmiedit查看中间结果，远比一次性写几十行代码然后面对一堆错误要高效得多。当你熟练之后，构建一个中等复杂度的LMI问题就像搭积木一样自然，这种将复杂系统问题转化为可计算模型的能力，会让你在研究和工程中受益匪浅。