企业官网建设流程全解析

简介：sssMOR是一套基于MATLAB的状态空间大规模动力学系统模型降阶工具箱，主要面向计算机、电子信息工程、数学等专业的学生和工程师，解决全阶模型计算开销大、难以实时仿真与后续设计的问题。通过低阶近似替代原系统，该工具箱能够在保留主要动态特性的同时显著降低仿真成本，对实时控制和大规模系统分析尤为重要。压缩包共359个文件，大小约8.85MB，以94个m源码、102个html函数帮助、130个png运行结果示例为主体，辅以少量配置、说明和app安装包，覆盖算法实现、使用文档与效果展示。工具兼容matlab2014/2019a/2024a，随附可直接运行的案例数据，代码采用参数化编程、注释明细，调整参数即可适配不同规模系统的降阶需求，也适合二次开发。目前已有56人浏览学习，可支撑课程设计、期末大作业和毕业设计，不仅助益模型降阶理论理解，也能提升实际建模仿真与程序调试能力。

1. sssMOR 是什么：给状态空间大规模系统做模型降阶的 MATLAB 工具箱

几十万阶的状态空间模型，不是每台机器都跑得动，也不是每次仿真都等得起。有限元结构、热传导网络、电力系统、电路模拟里，这类大规模动力学系统很常见。sssMOR 是一个 MATLAB 工具箱，专门对它们做模型降阶（Model Order Reduction，MOR）：输入是稀疏的广义状态空间模型（A、B、C、D、E 矩阵），输出是阶次低得多但输入输出行为足够接近的近似模型，让仿真、频响分析、控制器设计明显加速。标题里的 .rar 就是它的源码包。下面从理论背景讲到安装、参数、方法选型和验证，把"这个工具箱到底能干什么、怎么把它用起来"两件事说透。适合被大矩阵卡住仿真速度的工程师，也适合刚接触降阶的研究生，对已经跑过 balred 的熟手，重点看后面关于稀疏描述系统和误差界的部分。

2. 状态空间与投影降阶：sssMOR 的理论基础与工具箱定位

2.1 广义状态空间模型：为什么 E 矩阵要单独拿出来

教科书里的状态空间模型写作 dx/dt = Ax + Bu，但工程上从有限元或电路方程直接离散出来的，绝大多数是广义形式：

E dx/dt = Ax + Bu, y = Cx + Du

这里 E 是质量矩阵或电容矩阵，A 是刚度矩阵或导纳矩阵。E 可能是奇异的，也就是说系统里含代数约束，属于 DAE（微分代数方程），不能简单写成 dx/dt = E⁻¹Ax + Bu。强行求逆不仅破坏稀疏性，还会把数值性质弄坏。

sssMOR 所依赖的 sss 工具箱（sparse state space）就是为了不展开 E 而设计的：四个矩阵全部按 MATLAB 稀疏格式存储，构造时直接把 E 作为第五个参数传入。10 万阶的问题，稀疏存储只占非零元，内存是几 MB 到几十 MB；一旦误用 full() 展开成稠密矩阵，同样的模型需要约 80 GB 的 double 数组，直接内存溢出。这个差别决定了为什么不能用 Control System Toolbox 的 balred 凑合。

下表列出描述系统与传统表达方式的差异，后面调参数时容易对位：

% 一个体现 E 矩阵作用的简单描述系统 E = speye(3); E(3,3) = 0; % 第三个方程是代数方程，没有动态 A = sparse([0 1 0; -1 0 0; 0 1 1]); B = sparse([0; 1; 0]); C = sparse([1 0 0]); D = sparse(0); sys = sss(A, B, C, D, E, 'name', 'simple_dae');

这里 E 的 (3,3) 位置是 0，表示第三个状态不积分，只随代数方程变化。sss 构造函数按顺序接收 A、B、C、D、E，最后一个参数是名字，用于后续日志和图表识别。注意所有矩阵都用 sparse 类型传入，这是整个模型降阶流程能跑得动大规模系统的前提之一。

2.2 投影降阶：所有方法都在回答"V 和 W 怎么选"

模型降阶的核心想法是：状态 x 维数高，但能量集中在少数方向上。找一个 n×r 的投影矩阵 V，令 x ≈ Vxr，xr 是 r 维降阶状态，r 通常在 10～100。再配合另一个投影矩阵 W，做 Petrov-Galerkin 投影，得到降阶系统：

Er = WᵀEV, Ar = WᵀAV, Br = WᵀB, Cr = CV, Dr = D

从这个公式看，不同降阶方法的差别本质是"V 和 W 怎么构造"。平衡截断取的是可控 Gramian 和可观 Gramian 的主特征方向；Krylov 类方法取的是传递函数在若干插值点上的矩匹配方向；模态截断取的是主导极点对应的特征向量。sssMOR 把这几类选择都实现成独立函数，用户按需求选，而不是一个黑盒。

实现细节上，Ar 是 r×r 小矩阵，可以稠密存储，但 WᵀAV 的组装必须全程用稀疏运算，不能先把 A、V 展开成稠密矩阵再乘。sssMOR 内部用稀疏矩阵乘法加低秩迭代，避免任何一次 O(n²) 的稠密操作，这是它能处理 10 万阶以上模型的关键。

2.3 工具箱定位：和 MATLAB 优化工具箱、控制工具箱不是一回事

搜 MATLAB 工具箱时很容易看到"matlab优化工具箱""matlab机器人工具箱"，容易误以为 sssMOR 是同类"装好就能调"的通用组件。实际上定位差别很大：优化工具箱解决求极值问题；Control System Toolbox 的 balred 虽然也能降阶，但它面向中等规模稠密模型，内部要解稠密 Lyapunov 方程，到 1 万阶就基本跑不动了。sssMOR 走的是"大规模稀疏 + 低秩近似"路线，10 万阶以上才是它的主战场。

另一个区别是，sssMOR 假设你手里已经有状态空间模型。从有限元软件导出质量刚度矩阵、从电路仿真器导出 MNA 方程，这些前处理不属于它的范围。装完之后，先拿它自带的 benchmark 模型把流程走通，再换自己的数据。

3. 安装 sssMOR 并在 MATLAB 中跑通第一个降阶示例

3.1 安装前置条件：先装 sss 工具箱，再装 sssMOR

sssMOR 依赖 sss 工具箱提供的数据结构和基础运算，两者要一起加进 MATLAB 路径。标题里的 .rar 解压后，目录里通常有 sss 和 sssMOR 两个顶层文件夹，有的版本还带 benchmarks 和 demos。Windows 下解压 .rar 用 WinRAR 或 7-Zip；Linux 下如果没有 unrar，先安装再解压，用 bsdtar 也可以：

# Linux 下解压 .rar 源码包 sudo apt install unrar unrar x sssMOR_matlab_code.rar

% 把两个工具箱目录加进 MATLAB 路径，建议写进 startup.m addpath(genpath('/home/user/tools/sss')); % 先加 sss 核心 addpath(genpath('/home/user/tools/sssMOR')); % 再加 sssMOR savepath; % 保存路径预设

参数说明：genpath 会把子目录递归加进去，benchmarks 和 demos 里的 .mat 数据文件之后才能直接 load；savepath 把当前路径表写到默认的 pathdef.m，下次启动 MATLAB 自动生效。路径里尽量不要带中文和空格，addpath 对这类路径的处理历史上有不少兼容问题。

提示：sssMOR 官方要求 MATLAB 不能太老，R2016b 之后的版本基本都能跑，按 R2023b 安装教程的流程操作即可。装完先验证路径是否生效。

help sss % 能弹出帮助说明，说明 sss 装好了 which -all morBalreal % 能找到函数文件，说明 sssMOR 可用

如果 which 返回 empty，说明路径没加对，回到 addpath 那一步检查目录名；如果 help sss 报"未找到"，说明 sss 核心没装进去，sssMOR 一定起不来。

3.2 最小示例：加载模型、降阶、看阶次

sssMOR 自带的 benchmarks 里有几个公开的大规模模型，包括结构动力学的 build 模型、rail 轨道模型、CD player 模型等，规模从几千到十几万阶不等。用它们验证安装最省事：

% 加载 benchmark 模型（benchmarks 目录下的 .mat 文件） load build.mat; % 变量 A,B,C,E 已就位，注意检查是否有 D whos % 先看变量名，有的模型把 E 命名为 M sys_full = sss(A, B, C, D, E, 'name', 'build'); n = size(sys_full.A, 1); fprintf('原系统阶次: %d\n', n); % 降阶到 20 阶，并计时 opts.order = 20; tic; sys_red = morBalreal(sys_full, opts); toc;

说明：sss 构造时五参数顺序是 A、B、C、D、E，不要记错；模型没有 D 时传sparse(0)占位。load 之后先whos看变量名是必要的，有的 benchmark 把质量矩阵命名为 M 而不是 E，语义上要对应上。跑通时控制台大致会打印类似信息：

原系统阶次: 67472 低秩 ADI 迭代 18 步，残差 3.1e-11 Hankel 奇异值: sigma1 = 1.2e+01, ..., sigma20 = 4.7e-06 降阶完成，耗时 1.28s

这段日志里最值得看的是第 20 个 Hankel 奇异值：如果 sigma20 还在 1e-1 量级，说明 20 阶压得太狠，降阶模型误差会偏大，应该调大 opts.order 再看。

3.3 参数传递习惯：opts 结构体是统一入口

sssMOR 的函数基本都收两个参数：第一个是 sss 对象，第二个是配置结构体 opts。不同方法的字段不完全一样，但有几个是公共的：

字段记不全时，直接help morBalreal看注释，或者打开源码看函数开头对 opts 字段的默认值处理，比翻文档快。这个工具箱对使用者友好的地方在于：所有方法共用同一套 opts 风格，换方法时只需要改方法名和少量字段。

4. 核心方法：morBalreal、morIRKA 与可选参数对照

4.1 morBalreal：有误差界的基准方法

平衡截断（Balanced Truncation，BT）是 sssMOR 里理论最扎实的方法。它先算系统的可控 Gramian 和可观 Gramian，做平衡变换让两者相等且对角化，对角元素是 Hankel 奇异值，按大小排序后截掉尾部。优点是：降阶模型保持稳定性，且有严格的 H∞ 误差界：

||G - Gr||∞ ≤ 2 Σσi（i 从 r+1 到 n）

也就是被截掉的 Hankel 奇异值之和的两倍。对需要可信度的场景（控制器验证、安全评估），这个界非常值钱。实现上，sssMOR 用低秩 ADI 或低秩 Krylov 迭代求解稀疏 Lyapunov 方程，把经典 BT 的 O(n³) 稠密成本降到可接受范围。10 万阶以内、单次降阶几秒到几分钟属于正常。

% morBalreal 常用参数 opts.order = 40; % 目标阶次，先看 Hankel 奇异值分布再定 opts.tol = 1e-10; % 低秩 Lyapunov 求解残差容限 opts.maxit = 500; % ADI/Krylov 迭代上限，防死循环 sys_red = morBalreal(sys_full, opts); % 查看 Hankel 奇异值分布，辅助选 r（仅中小规模可用） hsv = hsvd(full(sys_full)); semilogy(hsv, 'o');

参数说明：order 是最直接的目标阶次。tol 控制 Gramian 低秩近似的精度，设太松（比如 1e-4）会让降阶模型实际误差远偏离理论误差界；maxit 是保护参数，迭代不收敛时会提示达到最大迭代次数，此时先调大 maxit，再考虑换预处理。最后一行 hsvd 只对能 full 展开的中小规模模型适用，10 万阶以上不要这么写，想看奇异值分布就依赖 morBalreal 的 verbose 日志。

4.2 morIRKA：迭代有理 Krylov，面向十万阶以上

IRKA（Iterative Rational Krylov Algorithm）是插值类方法的代表。核心思想是让降阶系统与原系统在若干插值点（shifts）上的传递函数矩匹配，再迭代更新插值点位置，直到收敛。它不对 Gramian 做全局求解，每一步只有稀疏 LU 分解和线性求解，因此能处理比 BT 大得多的系统。代价是：没有先验误差界，稳定性不保证，结果依赖初始插值点。

% morIRKA 的参数设置 opts.r = 40; % 降阶阶次（IRKA 里常用 r 表示） opts.iter_max = 200; % 外迭代最大步数 opts.tol = 1e-8; % 插值点相对变化容限 opts.shifts = -logspace(-4, 4, 40); % 初始插值点：负实轴对数均匀分布 opts.verbose = true; % 打印每步插值点移动情况 sys_red = morIRKA(sys_full, opts);

参数说明：shifts 直接决定方法成败。工程经验是先按 logspace 在负实轴对数均匀铺一层，覆盖系统可能的频谱范围；如果模型是振荡型（结构力学、电路），shifts 要带虚部，可以改成-1e2 + 1e2*1i*linspace(-1, 1, 40)这种复平面分布。verbose 打开后观察插值点是否在快速汇聚，如果 50 步还没稳定，大多是初始值离真实主导极点太远，重新铺 shifts。

4.3 方法选择对照表与参数建议

我的习惯是：先跑一次 morBalreal 拿结果和误差界；模型太大跑不动，再降级到 morIRKA。morKrylov 适合批量扫参时的粗近似，不适合做最终交付。另外 r 不要拍脑袋：先看 Hankel 奇异值或 IRKA 收敛后的插值点分布，选在奇异值掉一个量级以上的位置，比盲目设 50 或 100 效果好得多。还有一个容易忽略的点：原系统 E 奇异（DAE）时，各方法对 E 的处理方式不同，换方法后检查 sys_red.E 是否保持了预期结构。

5. 降阶结果的验证技巧与三个常见坑

5.1 用频响做交叉验证，别只看时域

降阶模型时域波形对得上，不代表整个频带内都对得上。最稳的做法是直接对比原系统和降阶系统的频率响应，算相对误差：

omega = logspace(-2, 3, 300); nw = numel(omega); G1 = zeros(1, nw); G2 = zeros(1, nw); for k = 1:nw s = 1i * omega(k); G1(k) = sys_full.C * ((s*sys_full.E - sys_full.A) \ sys_full.B) + sys_full.D; G2(k) = sys_red.C * ((s*sys_red.E - sys_red.A) \ sys_red.B) + sys_red.D; end relErr = norm(G1 - G2, 2) / norm(G1, 2); fprintf('频率范围 [%.0e, %.0e] 内相对误差: %.3e\n', omega(1), omega(end), relErr);

这段代码对 SISO 系统成立，核心是用稀疏反斜杠\求解 (sE-A)\B，不会把大矩阵展开。误差超过 1% 时，先看是否频段选太宽：高频段能量占比低但相对误差大，实际关心哪个频段就只看哪个频段。MIMO 模型把 B、C 的行列循环展开即可，逻辑不变。

5.2 稳定性检查与误差界核对

降阶模型必须做的验收是极点全部在左半平面。BT 方法理论上保持稳定，但数值实现仍要确认；IRKA 不保证稳定，这一步尤其不能省：

if all(real(eig(full(sys_red.A))) < 0) disp('降阶模型稳定'); else disp('降阶模型不稳定，考虑换方法或减小 r'); end

同时把实测误差与 morBalreal 的理论误差界对比。如果实测误差远大于 2 倍截断 Hankel 奇异值之和，先怀疑 opts.tol 设得太松，低秩 Gramian 解本身不够准；其次检查原系统是不是本来就不稳定，BT 的误差界对不稳定系统不成立。

5.3 三个容易踩的坑

第一个坑是漏传 E 矩阵。从有限元导出的模型，把 M 当 E 传给 sss，或者干脆不传，A、B、C 都能凑齐，MATLAB 不会报错，但降阶结果完全错误。每次 load 之后先 whos 核对变量表。第二个坑是 r 设得过大。降阶后的系统是稠密 r×r，仿真复杂度 O(r³)，设到 200 阶可能比原来的稀疏大系统还慢，降阶目的就没了。第三个坑是直接对稀疏模型调 hsvd、bode 这类 Control System Toolbox 函数，sss 对象虽然兼容部分函数，但内部可能偷偷 full 展开，内存瞬间爆炸。碰到这类函数就绕道，用 sssMOR 自带方法或上面这种手写频响循环替代。

本文还有配套的精品资源，点击获取