★ FEATURED ARTICLE 简介这是一份现代控制理论课程第1章“控制系统的状态空间表达式”的配套讲义面向自动化、控制类专业本科生及考研复习者。内容从“状态”“状态变量”“状态空间”等核心概念切入阐明一组最少变量在给定初始状态与输入后如何完整决定系统未来行为进而引出状态方程与输出方程的矩阵化表示。讲义通过R-L-C电路实例完整演示从电路微分方程到状态空间表达式的整理过程让读者直观理解状态矩阵、输入矩阵、输出矩阵和直接传输矩阵各自的作用掌握现代控制理论中最基础的建模方法。同时附带课程要求、主要参考书与章节框架便于把握学习主线、安排复习节奏。讲义还说明了状态变量的物理意义和选取原则例题步骤详细适合课后独立推导。资源共1个PDF文件大小约为2.79MB排版清晰适合打印或配合在线课程使用。已有199人学习可帮助初学者补强线性代数在控制建模中的实际应用。1. 状态空间表达式一套表达式就能描述整个动态系统某工程师调试一套双电机同步驱动系统用传递函数矩阵建模写了整整两页到设计状态反馈时却发现模型根本没法继续推。换用状态空间表达式四个矩阵 A、B、C、D 就把系统完整描述出来仿真、观测器、最优控制全都建立在同一套表达上。这就是为什么控制理论课程把“状态空间表达式”放在第1章它是后续所有设计方法的入口。这篇笔记针对正在跟这门课、或正在做控制算法落地的工程师目标是让你看完就能从微分方程写出状态空间表达式知道参数怎么选、坑在哪里并能用 Python 快速验证自己写的模型对不对。2. 状态变量的选取从微分方程到矩阵微分方程的完整流程2.1 选状态变量的三条规则阶数、物理量、线性无关状态变量的个数是一个硬指标系统独立储能元件的个数也就是微分方程的阶数。弹簧-质量系统有两个储能元件所以需要两个状态变量RLC 电路有电感和电容两个储能元件同样需要两个状态变量。如果数错了后面 A 矩阵的维度一定跟着错仿真器大概率直接抛异常。第二条规则是优先选物理上可测量的量例如位移、速度、电压、电流。原因很直接状态反馈最终要拿传感器测到的量去做控制状态变量选成不可测的量后面还得额外配一个观测器把状态估出来。观测器当然能解决但那是多出来的设计负担。第三条规则是状态变量之间必须线性无关。这里说的不是物理上“看起来不同”而是你选的 n 个变量能张成一个 n 维空间。例如 x1 是位移、x2 是速度这是标准选法如果额外把 x3 也选成位移的倍数那三个变量就只有两个独立维度构造出的 A 矩阵会奇异状态转移过程也存在冗余。实际使用时三条规则会有冲突。优先保证第一条和第三条第二条在条件允许时尽量满足。如果实在找不到物理可测的状态变量就把输出及其各阶导数组成一组相变量这是后面要展开的通用做法。2.2 从弹簧-质量-阻尼系统推导状态空间表达式的完整算例拿一个典型二阶系统来做例子。运动微分方程是m * y(t) c * y(t) k * y(t) f(t)定义状态 x1 y位移x2 y速度。直接对时间求导x1 x2 x2 f/m - c/m * x2 - k/m * x1写成矩阵形式A、B、C、D 就齐了[x1] [ 0 1 ] [x1] [ 0 ] [x2] [-k/m -c/m] [x2] [1/m] [f] y [1 0] [x1] [x2]用代码把矩阵填出来并打印极点是验证模型的第一步import numpy as np m 2.0 # 质量, kg c 0.5 # 阻尼系数, N·s/m k 20.0 # 弹簧刚度, N/m A np.array([[0.0, 1.0], [-k/m, -c/m]]) B np.array([[0.0], [1.0/m]]) C np.array([[1.0, 0.0]]) D np.array([[0.0]]) print(A , A, sep\n) print(极点 , np.linalg.eigvals(A))这段代码里的 m、c、k 就是物理参数改它们能直接观察对极点位置的影响。A(0,1)1 表示 x1 的导数等于 x2这是状态定义决定的不能随便改A(1,0)−k/m 来自弹簧恢复力A(1,1)−c/m 来自阻尼力。如果系统单位不统一例如质量用 kg 而力用 N那么 k/m 的数值就会差好几个数量级这是后文要重点说的量纲问题。2.3 状态变量选取不唯一换一组状态变量极点不变常见的错误是以为状态变量选法是唯一的所以看到别人写的 A 矩阵和自己不一样就怀疑算错了。实际上只要做非奇异线性变换系统的外部特性完全不变。对上面的例子取一组新状态 z1 x1z2 −2*x1 − x2也就是做线性变换 z T xT np.array([[1.0, 0.0], [-2.0, -1.0]]) A_new T A np.linalg.inv(T) B_new T B C_new C np.linalg.inv(T) print(新A极点 , np.linalg.eigvals(A_new)) print(原A极点 , np.linalg.eigvals(A))A_new 是相似变换后的系统矩阵特征值不变说明极点由系统本身决定跟状态变量选择无关。实际工程中把状态变量从物理量换成正则坐标或误差量往往是控制设计的一部分。反馈矩阵 K 的具体数值也会跟着变但闭环极点位置在设计目标下保持不变。2.4 多输入多输出时矩阵维度怎么填不踩空前面的例子是单输入单输出。如果系统有两个输入 u1、u2三个测量输出 y1、y2、y3那 B 矩阵就是 n×2C 矩阵是 3×nD 矩阵是 3×2。这个维度关系其实比矩阵元素值更容易错。最实用的办法是先把端口标在系统方框图上输入箭头数量决定 B 的列数输出箭头数量决定 C 的行数。然后对每个输出写一个 y C x D u 的等式逐个核对。写代码时加断言assert B.shape (A.shape[0], u_dim), B矩阵维度应为 (n, p) assert C.shape (y_dim, A.shape[0]), C矩阵维度应为 (q, n) assert D.shape (y_dim, u_dim), D矩阵维度应为 (q, p)u_dim 是输入通道数y_dim 是输出通道数。不要只看矩阵元素值对不对维度错了仿真器会直接拒绝运行反而是最容易定位的一类错误。3. 从传递函数到状态空间能控标准型、能观标准型与最小实现很多系统建模的起点是频域测量得到的传递函数不是微分方程。这时候要用一个公式化的方法把传递函数转成状态空间表达式。最常见的是能控标准型也叫相变量标准型和能观标准型。3.1 能控标准型把传递函数分母系数直接填进友矩阵对于 G(s) (b1·s b0) / (s² a1·s a0)能控标准型的关系是A [[0, 1], [−a0, −a1]] B [[0], [1]] C [[b0, b1]] D [[0]]为什么这么填因为取状态 x1 yx2 y之后x2 的导数可以写成 −a0·x1 − a1·x2 u所以 A 的最后一行是 [−a0, −a1]而 B 只在最后一行为 1。这就是“友矩阵”的结构很多教材也直接叫它 companion form。适用场景做极点配置状态反馈时这种形式最受欢迎。因为 B 向量非常“干净”设计反馈增益 K [k1, k2] 后闭环系统 A − B·K 的最后一行可以直接写出期望特征多项式的系数。你不需要复杂的坐标变换直接把期望极点对应的多项式系数填进去就能反解出 K。3.2 能观标准型适合观测器设计与能控标准型对偶能观标准型的关系是A [[0, −a0], [1, −a1]] B [[b0], [b1]] C [[0, 1]] D [[0]]它和能控标准型互为对偶Ao AcᵀBo CcᵀCo Bcᵀ。这个关系不是巧合它来自线性系统理论里的对偶原理。如果你正在设计状态观测器用能观标准型可以直接把观测器增益设计问题转换成一个对偶的状态反馈问题。本质上是同一套计算只是换了视角。工程上的选择习惯是做控制器设计用能控标准型做观测器设计用能观标准型。两种标准型只差一个转置但直接按目标场景选型能省掉大量中间推导。3.3 分子分母同阶时D矩阵不能漏初学者最容易漏掉的情况G(s) (b2·s² b1·s b0) / (s² a1·s a0)分子和分母同阶。此时不能直接把 b2 作为 C 的分量丢进去要先做一个长除法G(s) b2 [(b1 − a1·b2)·s (b0 − a0·b2)] / (s² a1·s a0)然后对后面这个真分式部分构造能控标准型直通项 b2 填到 D 矩阵里。工程里做 PID 控制器的状态空间实现时经常会遇到这种情况因为 PID 的传递函数分子阶数至少是分母阶数减一微分项一加就分子分母同阶了。D 矩阵在这个场景下代表“输入直接出现在输出端”是物理上存在的直通路径不是可以随便忽略的项。一个小习惯拿到传递函数后先比较分子分母阶次如果同阶先做除法再构造状态空间表达式。这样 D 矩阵才不会被遗漏C 矩阵里的系数也才是对的。3.4 最小实现先约分再构造否则白多一阶如果传递函数分子分母存在公因子比如 G(s) (s1) / (s² 3s 2)约分之后其实是 1/(s2)。直接按系数构造的表达式阶数会是 2但系统的真正极点只有 −2另一个状态既不可控也不可观是一个“隐藏”在表达式里的冗余结构。用这种非最小实现去设计控制器会出现算出来的反馈增益矩阵实际没有效果的情况。判断最小实现可以计算可控性矩阵 Qc [B, A·B, ...] 和可观性矩阵 Qo [C; C·A; ...] 的秩import numpy as np from scipy.signal import tf2ss num [1, 1] # s 1 den [1, 3, 2] # s^2 3s 2 A, B, C, D tf2ss(num, den) n A.shape[0] Qc np.hstack([B, A B]) Qo np.vstack([C, C A]) print(可控性秩 , np.linalg.matrix_rank(Qc)) print(可观性秩 , np.linalg.matrix_rank(Qo)) print(系统阶数 , n)直接按未约分传递函数构造时跑出来阶数是 2可观性秩往往只有 1说明存在不可观模态约分后 num[1], den[1,2] 时阶数和两个秩都是 1。这才是最小实现。我的建模习惯是先做符号约分再转状态空间最后跑一遍秩检查三步做完才敢把这个模型拿去设计控制器。3.5 反推传递函数验证 A、B、C、D 是否自洽模型构造完可以用公式 G(s) C(sI − A)⁻¹B D 反推传递函数和原始传递函数对比。这一步能同时验证矩阵元素是否填错、D 矩阵是否漏掉、以及状态变量选择是否引入额外模态。from scipy.signal import ss2tf num_back, den_back ss2tf(A, B, C, D) print(反推分子 , num_back[0]) print(反推分母 , den_back)如果反推出来的分子分母和原始传递函数不一致问题大概率出在 A 矩阵的符号上。比如把 −a0 填成了 a0传递函数会完全变样。ss2tf 返回的 num_back 是一个二维数组取第一行即可。4. 状态空间表达式落地必踩的坑维度、量纲、初值与离散化状态空间表达式看着只是四个矩阵但毛刺都在细节里。这一类问题我总结下来集中在四块维度、量纲、初值、离散化。4.1 B矩阵列数写错仿真器直接崩现象一个二输入系统B 矩阵写成了 n×1仿真一跑就报矩阵维度错误。原因没有按输入端口数来确定 B 的列数。系统是单输入还是双输入B 矩阵列数是硬约束写错一列都不行。解决写模型时先画方框图数输入箭头。另一个办法是在代码里直接加断言前面 2.4 节已经写过。这里再强调一遍维度断言建议放在建模函数的第一行每次装配模型都会自动检查。提示矩阵元素值错了仿真器不一定报错但维度错了它一定会报错。把断言写进建模流程等于让编译器帮你检查端口定义。4.2 量纲不统一造成假发散现象理论上稳定的系统仿真输出却快速发散或者响应时间比实测快了几个数量级。原因A 矩阵里填的是 k/m 和 c/m 的数值单位制不统一会导致这些比值偏几个数量级。比如质量用克、力用牛顿算出来的 k/m 数值差了 1000 倍极点位置就会完全变样。解决所有物理量统一用 SI 制至少在构造矩阵前把所有单位换算好。还有一个排查技巧把 A 矩阵的每个元素逐一和量纲公式比对。A(0,1) 的量纲是 1/sA(1,0) 的量纲是 1/s²A(1,1) 的量纲也是 1/s。如果发现 A(1,0) 的量纲不是 1/s²那大概率是单位换算没到位。我从这里踩过的坑不算少现在每次建模都会把 A 矩阵元素连同物理参数一起打印出来人工扫一眼量纲。4.3 初值默认是零实际系统有储能就对不上现象仿真阶跃响应形状对但起始段跟实测曲线差一大截。原因仿真工具默认状态初值为零而实际系统开机时电容可能已经有电压机械系统可能已经有初速或初始位移。状态空间表达式本身是带着初始条件 x(t0) 的但仿真函数不传初值就默认全零。解决用 initial_response 或直接在积分器里设置 x0import control as ct sys ct.ss(A, B, C, D) T np.linspace(0, 10, 1000) t, y ct.initial_response(sys, T, X0[0.05, 0.0])这就是为什么我调试状态反馈控制器时会把阶跃响应和初值响应各跑一次。两个响应不一样才说明模型包含了完整的动态特性。如果只测阶跃响应初值问题很容易被漏掉。实测对比时初值响应常常比阶跃响应更能暴露模型误差。4.4 公因子未约分非最小实现带来隐藏模态现象tf2ss 得到的系统是 2 阶但实际物理系统只有 1 个极点或者状态反馈调参时无论怎么改 K某个模态始终不受控制。原因传递函数分子分母存在公因子直接构造的表达式包含不可控或不可观模态。解决第一步用符号计算对传递函数约分第二步做可控性和可观性秩检查。如果发现秩亏不要硬用这个非最小实现做设计回去先把传递函数约分干净。状态空间下也有自动消去冗余模态的算法但工程上不要依赖这个建模阶段把模型整干净才是根治。约分这个动作看起来很基础但遇到带延迟环节或近似零极点对消的模型时最容易被忽略。4.5 离散化时采样周期太大连续稳定离散发散现象连续系统极点都在左半平面但用零阶保持器离散后离散极点在单位圆外仿真发散。原因采样周期 Ts 大于系统动态的最快时间常数。连续域离虚轴越远的极点离散化后对 Ts 越敏感。Ts 太大时零阶保持器的延迟效应会超过系统的稳定裕度。解决先看连续极点的实部选 Ts 满足 Ts 1/(10·|Re(λmax)|)也就是让采样频率至少是系统带宽的 10 到 20 倍。然后用 c2d 实测一下离散极点from scipy.signal import cont2discrete Ts 0.05 # 采样周期 50ms Ad, Bd, _, _, _ cont2discrete((A, B, C, D), Ts, methodzoh) print(离散极点 , np.linalg.eigvals(Ad))如果离散极点模值大于 1先把 Ts 调小一个数量级再试。这套经验在做嵌入式控制器参数整定时特别常见采样周期不是靠拍脑袋定的而是从连续极点出发反推的。5. 用仿真把状态空间表达式“焊死”响应、极点与李雅普诺夫校验模型写完之后的第一件事不是直接上控制器设计而是先验证你手里的状态空间表达式能真实复现系统的动态。我现在固定用三步响应曲线、极点、李雅普诺夫方程。第一步是阶跃响应。用 python-control 先把连续系统建出来然后跑一条阶跃响应曲线跟实验数据或原始微分方程的数值解做对比。如果状态空间表达式是从微分方程推出来的两者应当完全重合。import numpy as np import control as ct import matplotlib.pyplot as plt A np.array([[0, 1], [-10, -1]]) B np.array([[0], [1]]) C np.array([[1, 0]]) D np.array([[0]]) sys ct.ss(A, B, C, D) t, y ct.step_response(sys) plt.plot(t, y) plt.xlabel(time (s)) plt.ylabel(output) plt.savefig(step_verify.png)第二步是打印极点看是否落在左半平面。这个检查只需要一行 np.linalg.eigvals(A)但每次建模都要做。你会发现不同状态变量组合得到的不同状态空间表达式极点完全一致这是模型正确的必要条件也是验证“状态变量选取不唯一”这句话最快的手段。第三步是我个人的习惯用李雅普诺夫方程做稳定判据而不只看极点。原因是极点判据在线性时不变系统里可靠但遇到参数不确定或非线性系统线性化模型时直接求解李雅普诺夫方程能顺带给出一个能量函数的轮廓这对后续做鲁棒性分析很有用。from scipy.linalg import solve_continuous_lyapunov Q np.eye(2) P solve_continuous_lyapunov(A.T, -Q) print(P特征值 , np.linalg.eigvalsh(P))如果 P 正定系统在李雅普诺夫意义下稳定。注意这里传的是 A.T因为该函数解的是 AᵀP PA −Q 的形式。不同 scipy 版本对参数的解释略有差异建议在本地确认矩阵方程后再用到正式项目里。我自己的习惯是每套状态空间表达式都留一个仿真脚本里面固定放这三段阶跃响应、极点打印、李雅普诺夫检查。模型改动任何参数顺手跑一遍有问题当场暴露。这个习惯救过我好几次因为状态空间表达式的矩阵元素太多肉眼检查总有漏掉的时候而仿真会诚实地告诉你哪里错了。希望帮到你。本文还有配套的精品资源点击获取
阅读完成 · 觉得有帮助？
