本系列文章为线性系统控制的学习笔记。我们将围绕状态空间方程展开,探讨如何分析和计算动态系统。为了聚焦基础推导与数值计算,本系列主要讨论线性时不变系统和离散时间系统。
线性系统控制 1 ———— 状态空间模型
在分析动态系统时,我们通常会面对一个核心问题:如何用数学语言准确描述物理世界的行为?无论是机械系统中的振动控制,还是电路系统中的信号处理,我们最开始接触的工具往往是微分方程。然而,随着系统复杂度的提升,状态空间表示法逐渐成为了现代控制理论的核心。
1 经典视角:高阶微分方程
微分方程是我们描述物理系统动态特性的直观工具。它通过系统变量及其各阶导数之间的关系,揭示了系统是如何随时间演化的。
以一个经典的机械系统——质量-弹簧-阻尼系统为例。其示意图如图1所示。
图1 质量-弹簧-阻尼系统示意图
根据牛顿第二定律,我们可以写出如下的二阶常微分方程:
m y ¨ ( t ) + c y ˙ ( t ) + k y ( t ) = u ( t ) m\ddot{y}(t) + c\dot{y}(t) + ky(t) = u(t)
m y ¨ ( t ) + c y ˙ ( t ) + k y ( t ) = u ( t )
其中:
y ( t ) y(t) y ( t ) 是系统的位移输出
u ( t ) u(t) u ( t ) 是外部施加的力(输入)
m , c , k m, c, k m , c , k 分别代表质量、阻尼系数和弹簧刚度
这个系统是一个经典的单输入单输出系统模型,通过微分方程仅仅描述了输入和输出之间的关系,把系统内部作为一个黑盒。
2 现代视角:状态空间模型
为了打破黑盒,现代控制理论引入了状态变量的概念。状态变量是能够完全确定系统未来行为的最小变量集合。状态空间方程将一个 n n n 阶的高阶微分方程,降维转换成了 n n n 个一阶微分方程的集合,并用矩阵的形式统一表达。它的标准形式包含两个方程:
状态方程: 描述系统内部状态的演化。
x ˙ ( t ) = A x ( t ) + B u ( t ) \dot{\mathbf{x}}(t) = \mathbf{A}\mathbf{x}(t) + \mathbf{B}\mathbf{u}(t)
x ˙ ( t ) = Ax ( t ) + Bu ( t )
输出方程: 描述内部状态与可测量输出之间的关系。
y ( t ) = C x ( t ) + D u ( t ) {\mathbf{y}}(t) = \mathbf{C}\mathbf{x}(t) + \mathbf{D}\mathbf{u}(t)
y ( t ) = Cx ( t ) + Du ( t )
其中:
x ( t ) n × 1 \mathbf{x}(t)_{n \times 1} x ( t ) n × 1 : 状态向量——包含了系统在时刻 t t t 的所有内部状态(比如系统内部所有质量块的位置和速度)。如果系统有 n n n 个状态,它就是一个 n × 1 n \times 1 n × 1 的列向量。
x ˙ ( t ) n × 1 \dot{\mathbf{x}}(t)_{n \times 1} x ˙ ( t ) n × 1 : 状态向量的导数——状态向量随时间的变化率。
u ( t ) m × 1 \mathbf{u}(t)_{m \times 1} u ( t ) m × 1 : 输入向量 / 控制向量——系统从外界接收到的所有激励或控制信号(比如施加在不同位置的力或电压)。如果系统有 m m m 个输入,它就是一个 m × 1 m \times 1 m × 1 的列向量。
y ( t ) p × 1 \mathbf{y}(t)_{p \times 1} y ( t ) p × 1 : 输出向量——我们实际能够测量到或者关心的系统输出。如果系统有 p p p 个输出,它就是一个 p × 1 p \times 1 p × 1 的列向量。
A n × n \mathbf{A}_{n \times n} A n × n : 系统矩阵 / 状态矩阵——它是整个系统最重要的部分,包含了系统内部所有的物理参数(如质量、弹簧刚度、阻尼等)以及它们之间的耦合关系。它决定了系统的自然响应和稳定性。计算矩阵 A \mathbf{A} A 的特征值,就能直接知道系统的极点。
B n × m \mathbf{B}_{n \times m} B n × m : 输入矩阵 / 控制矩阵——它描述了外部输入信号 u ( t ) \mathbf{u}(t) u ( t ) 是如何进入系统内部,并影响各个状态变量的。它决定了系统的可控性。
C p × n \mathbf{C}_{p \times n} C p × n : 输出矩阵 / 观测矩阵——它描述了系统内部的各个状态变量,是如何组合成我们最终在外部传感器上测量到的输出 y ( t ) \mathbf{y}(t) y ( t ) 的。它决定了系统的可观测性。
D p × m \mathbf{D}_{p \times m} D p × m : 直接前馈矩阵 / 传递矩阵——它表示输入信号是否有一条“捷径”可以直接绕过系统内部状态,瞬间影响到输出。在绝大多数实际的机械和物理系统中,输入不可能没有延迟地瞬间改变输出,所以矩阵 D \mathbf{D} D 通常是一个全零矩阵。
3 核心推导:将微分方程转换为状态方程
通过定义新的变量,把高阶导数降为一阶导数。我们依然以上面的二阶系统为例:
m y ¨ + c y ˙ + k y = u m\ddot{y} + c\dot{y} + ky = u
m y ¨ + c y ˙ + k y = u
第一步:定义状态变量
我们将系统的位移和速度定义为状态变量:x 1 = y x_1 = y x 1 = y (位移)x 2 = y ˙ x_2 = \dot{y} x 2 = y ˙ (速度)
状态向量为
x ( t ) = [ x 1 x 2 ] = [ y y ˙ ] \mathbf{x}(t) = \begin{bmatrix} {x}_1 \\ {x}_2 \end{bmatrix} = \begin{bmatrix} {y} \\ \dot{y} \end{bmatrix}
x ( t ) = [ x 1 x 2 ] = [ y y ˙ ]
状态向量的导数为
x ˙ ( t ) = [ x ˙ 1 x ˙ 2 ] = [ y ˙ y ¨ ] \dot{\mathbf{x}}(t) = \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \end{bmatrix} = \begin{bmatrix} \dot{y} \\ \ddot{y} \end{bmatrix}
x ˙ ( t ) = [ x ˙ 1 x ˙ 2 ] = [ y ˙ y ¨ ]
第二步:建立一阶方程组
对于状态方程,我们需要知道状态向量的导数与状态和输入向量之间的关系,其中
x ˙ 1 = y ˙ = x 2 \dot{x}_1 = \dot{y} = x_2
x ˙ 1 = y ˙ = x 2
x ˙ 2 = y ¨ = − k m y − c m y ˙ + 1 m u = − k m x 1 − c m x 2 + 1 m u \dot{x}_2 = \ddot{y} = -\frac{k}{m}y - \frac{c}{m}\dot{y} + \frac{1}{m}u = -\frac{k}{m}x_1 - \frac{c}{m}x_2 + \frac{1}{m}u
x ˙ 2 = y ¨ = − m k y − m c y ˙ + m 1 u = − m k x 1 − m c x 2 + m 1 u
第三步:写成矩阵形式
将上述方程组用矩阵表示,就得到了标准的状态空间方程:
[ x ˙ 1 x ˙ 2 ] = [ 0 1 − k m − c m ] [ x 1 x 2 ] + [ 0 1 m ] u \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \end{bmatrix} = \begin{bmatrix} 0 & 1 \\ -\frac{k}{m} & -\frac{c}{m} \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} + \begin{bmatrix} 0 \\ \frac{1}{m} \end{bmatrix} u
[ x ˙ 1 x ˙ 2 ] = [ 0 − m k 1 − m c ] [ x 1 x 2 ] + [ 0 m 1 ] u
因为我们的输出是位移 y y y ,也就是我们需要测量的值,所以输出方程为:
y = [ 1 0 ] [ x 1 x 2 ] + [ 0 ] u y = \begin{bmatrix} 1 & 0 \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} + [0] u
y = [ 1 0 ] [ x 1 x 2 ] + [ 0 ] u
4 为什么需要使用状态空间方程
这是本篇文章的核心。既然微分方程已经能描述系统,为什么控制工程师和研究人员如此偏爱状态方程?
4.1 方程的求解与计算
经典的传递函数或微分方程处理多输入多输出系统时,由于变量间的耦合,表达式会变得极其冗长和复杂。但在状态空间中,增加输入输出仅仅意味着矩阵 B B B 和 C C C 的维度变大,方程的标准形式 x ˙ = A x + B u \dot{x} = Ax + Bu x ˙ = A x + B u 依然保持简洁不变。并且状态空间方程将复杂的微积分运算,巧妙地转换成了适合计算机的线性代数运算。
4.2 系统的内部结构
微分方程只关心输入和输出的宏观表现。状态方程则让我们能够“看透”系统内部。它不仅能帮助我们判断系统是否稳定,还能让我们分析系统的可控性(我们能否通过输入把系统推到任意状态)和可观测性(我们能否仅通过测量输出反推出所有内部状态)。
4.3 微控制器的适用
现代控制系统几乎都是由微控制器运行的。代码是按离散的采样周期运行的,没有真正的“连续时间”。强行把高阶微分方程写成代码,不仅需要复杂的差分近似,还会引入巨大的截断误差。连续状态方程可以利用状态转移矩阵转化为离散时间状态方程,便于我们在底层硬件上编写控制算法。
4.4 现代控制策略的基石
如果我们想在未来应用更高级的控制理论,比如最优控制中的线性二次型调节器、卡尔曼滤波或者状态反馈控制,这些算法全部都是建立在状态空间模型之上的。
5 状态空间方程优势展现
5.1 内部状态
回到之前推导出的标准状态空间方程。
5.1.1 状态方程
[ x ˙ 1 x ˙ 2 ] = [ 0 1 − k m − c m ] [ x 1 x 2 ] + [ 0 1 m ] u \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \end{bmatrix} = \begin{bmatrix} 0 & 1 \\ -\frac{k}{m} & -\frac{c}{m} \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} + \begin{bmatrix} 0 \\ \frac{1}{m} \end{bmatrix} u
[ x ˙ 1 x ˙ 2 ] = [ 0 − m k 1 − m c ] [ x 1 x 2 ] + [ 0 m 1 ] u
其中,状态变量 x 1 x_1 x 1 是位移,x 2 x_2 x 2 是速度。如果我们把这个矩阵乘法按行拆解开来,就能清晰地看到系统内部是如何演化的。
第一行:运动学约束
x ˙ 1 = x 2 \dot{x}_1 = x_2
x ˙ 1 = x 2
表示位移的变化率,等于当前的速度。作为系统矩阵 A \mathbf{A} A 的第一行,将位置的未来演化与速度的当前状态联系在一起。同时,输入矩阵 B \mathbf{B} B 对应的项是 0(即 + 0 ⋅ u +0 \cdot u + 0 ⋅ u ),这反映了物理系统的基本特性:外力不能直接改变位置,只能通过改变速度(即产生加速度)间接影响位置。
第二行:动力学方程
x ˙ 2 = − k m x 1 − c m x 2 + 1 m u \dot{x}_2 = -\frac{k}{m}x_1 - \frac{c }{m}x_2 + \frac{1}{m}u
x ˙ 2 = − m k x 1 − m c x 2 + m 1 u
我们得到牛顿第二定律的表达形式 a = ∑ F / m a = \sum F / m a = ∑ F / m 。x 2 x_2 x 2 是速度,那么 x ˙ 2 \dot{x}_2 x ˙ 2 就是加速度。这个方程告诉我们,下一瞬间的加速度是由系统内部与外部的多种作用力共同决定的:
位置的反馈项 (− k m x 1 -\frac{k}{m}x_1 − m k x 1 ): 当系统偏离平衡位置时,弹簧产生与位移成正比的恢复力,将系统拉回平衡点。该项体现了系统的刚度特性。
速度的阻尼项 (− c m x 2 -\frac{c}{m}x_2 − m c x 2 ): 阻尼力与速度成正比,并始终阻碍运动,从而不断耗散系统能量。该项反映了系统的能量耗散机制。
外部输入项 (+ 1 m u +\frac{1}{m}u + m 1 u ): 表示外界施加的控制力,是系统能量的输入通道,也是实现控制的关键。
内部演化全貌:输入信号 u u u 通过输入通道作用于系统,直接影响加速度 x ˙ 2 \dot{x}_2 x ˙ 2 。加速度作为速度的一阶导数,在时间积分作用下不断累积,从而改变下一时刻的速度 x 2 x_2 x 2 。更新后的速度 x 2 x_2 x 2 根据第一行的运动学关系,进一步驱动位移 x 1 x_1 x 1 的变化。新的位移与速度状态又共同作用于弹簧和阻尼器,产生新的内部力,并重新决定下一时刻的加速度 x ˙ 2 \dot{x}_2 x ˙ 2 。
由此形成一个连续的动态闭环:输入 → 加速度 → 速度 → 位移 → 内部力 → 加速度。
5.1.2 输出方程
y = [ 1 0 ] [ x 1 x 2 ] = x 1 \mathbf{y} = \begin{bmatrix} 1 & 0 \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} = x_1
y = [ 1 0 ] [ x 1 x 2 ] = x 1
矩阵 C = [ 1 0 ] \mathbf{C} = \begin{bmatrix} 1 & 0 \end{bmatrix} C = [ 1 0 ] 是输出矩阵。我们可以把它纯粹地理解为传感器映射。它决定了内部状态 x \mathbf{x} x 中,有哪些变量能够被我们看到,以及如何被看到。我们之前算出的输出方程 y = [ 1 0 ] [ x 1 x 2 ] = x 1 \mathbf{y} = \begin{bmatrix} 1 & 0 \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} = x_1 y = [ 1 0 ] [ x 1 x 2 ] = x 1 。我们的输出就是纯粹的位移,也就是我们只使用了位移传感器。
如果我们只买了一个测速仪(比如多普勒雷达),那么矩阵 C \mathbf{C} C 就变成 [ 0 1 ] \begin{bmatrix} 0 & 1 \end{bmatrix} [ 0 1 ] 。y = [ 0 1 ] [ x 1 x 2 ] = x 2 \mathbf{y} = \begin{bmatrix} 0 & 1 \end{bmatrix} \begin{bmatrix} x_1 \\ x_2 \end{bmatrix} = x_2 y = [ 0 1 ] [ x 1 x 2 ] = x 2 。我们的输出变成了速度。
如果我们既有位移传感器,又有测速仪,矩阵 C \mathbf{C} C 变成一个 2 × 2 2 \times 2 2 × 2 的单位阵 [ 1 0 0 1 ] \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix} [ 1 0 0 1 ] 。我们的输出 y \mathbf{y} y 就变成了一个包含位移和速度的列向量。
矩阵 C \mathbf{C} C 将物理动态与测量手段彻底解耦。不论底层物理怎么运转,我们都可以通过重新定义 C \mathbf{C} C 矩阵,极其自由地改变系统的输出。
矩阵 D \mathbf{D} D 被称为直接前馈矩阵。它代表了外部输入 u \mathbf{u} u 是否能绕过所有的内部物理状态,瞬间、直接地反映在输出 y \mathbf{y} y 上。
y = C x + D u \mathbf{y} = \mathbf{C}\mathbf{x} + \mathbf{\color{red}{D}\mathbf{u}}
y = Cx + D u
为什么大多数物理系统中 D = 0 \mathbf{D} = \mathbf{0} D = 0 ?
想一想我们的质量-弹簧-阻尼系统。当我们瞬间推了一下质量块(输入 u \mathbf{u} u ),物体的加速度确实瞬间改变了,但它的位置和速度(状态 x \mathbf{x} x )不可能瞬间发生跃变,必须经过时间的积分才能慢慢积累。既然状态不能突变,那么观测到的位移和速度也不可能突变。因此,在绝大多数包含惯性(如质量、电感、电容)的物理系统中,输入无法瞬间变成输出,所以 D \mathbf{D} D 矩阵通常是一个全零矩阵。
什么时候 D \mathbf{D} D 不为零?
如果我们的系统是一个纯电阻网络,输入一个电压 u u u ,根据欧姆定律,输出电流 y y y 会瞬间同比例放大或缩小,没有任何积分延迟。此时 D \mathbf{D} D 就不为零。或者,如果我们强行定义系统的输出 y \mathbf{y} y 不是位移,而是加速度,那么因为推力 u \mathbf{u} u 会瞬间产生加速度,这时候公式里就会出现非零的 D \mathbf{D} D 项。
5.2 多输入多输出系统
下面我们使用双质量-弹簧-阻尼系统来展示使用状态空间方程处理多输入多输出系统,其系统示意图如图2所示。
图2 双质量-弹簧-阻尼系统示意图
根据牛顿第二定律对两个质量块分别进行受力分析,我们可以得到如下的运动学方程:
方程 1
m 1 y ¨ 1 + ( c 1 + c 2 ) y ˙ 1 + ( k 1 + k 2 ) y 1 − c 2 y ˙ 2 − k 2 y 2 = u 1 m_1\ddot{y}_1 + (c_1 + c_2)\dot{y}_1 + (k_1 + k_2)y_1 \mathbf{\color{red}{- c_2\dot{y}_2 - k_2y_2}} = u_1
m 1 y ¨ 1 + ( c 1 + c 2 ) y ˙ 1 + ( k 1 + k 2 ) y 1 − c 2 y ˙ 2 − k 2 y 2 = u 1
方程 2
m 2 y ¨ 2 + c 2 y ˙ 2 + k 2 y 2 − c 2 y ˙ 1 − k 2 y 1 = u 2 m_2\ddot{y}_2 + c_2\dot{y}_2 + k_2y_2 \mathbf{\color{red}{- c_2\dot{y}_1 - k_2y_1}} = u_2
m 2 y ¨ 2 + c 2 y ˙ 2 + k 2 y 2 − c 2 y ˙ 1 − k 2 y 1 = u 2
注意上面方程中标记出来的部分。在方程 1 中,明明是在算 y 1 y_1 y 1 的动态,却混入了 y 2 y_2 y 2 的速度 y ˙ 2 \dot{y}_2 y ˙ 2 和位移 y 2 y_2 y 2 ;方程 2 亦是如此。这意味着,我们无法单独求解任何一个方程。
对比之下,如果我们使用状态空间方程,我们只需要定义四个状态变量:
x 1 = y 1 x_1 = y_1 x 1 = y 1 (位置1), x 2 = y ˙ 1 x_2 = \dot{y}_1 x 2 = y ˙ 1 (速度1), x 3 = y 2 x_3 = y_2 x 3 = y 2 (位置2), x 4 = y ˙ 2 x_4 = \dot{y}_2 x 4 = y ˙ 2 (速度2)
我们就能把一个二阶耦合方程转化为矩阵运算。
[ x ˙ 1 x ˙ 2 x ˙ 3 x ˙ 4 ] = [ 0 1 0 0 − k 1 + k 2 m 1 − c 1 + c 2 m 1 k 2 m 1 c 2 m 1 0 0 0 1 k 2 m 2 c 2 m 2 − k 2 m 2 − c 2 m 2 ] [ x 1 x 2 x 3 x 4 ] + [ 0 0 1 m 1 0 0 0 0 1 m 2 ] [ u 1 u 2 ] \begin{bmatrix} \dot{x}_1 \\ \dot{x}_2 \\ \dot{x}_3 \\ \dot{x}_4 \end{bmatrix} =
\begin{bmatrix}
0 & 1 & 0 & 0 \\
-\frac{k_1+k_2}{m_1} & -\frac{c_1+c_2}{m_1} & \frac{k_2}{m_1} & \frac{c_2}{m_1} \\
0 & 0 & 0 & 1 \\
\frac{k_2}{m_2} & \frac{c_2}{m_2} & -\frac{k_2}{m_2} & -\frac{c_2}{m_2}
\end{bmatrix}
\begin{bmatrix} x_1 \\ x_2 \\ x_3 \\ x_4 \end{bmatrix} +
\begin{bmatrix} 0 & 0 \\ \frac{1}{m_1} & 0 \\ 0 & 0 \\ 0 & \frac{1}{m_2} \end{bmatrix}
\begin{bmatrix} u_1 \\ u_2 \end{bmatrix} x ˙ 1 x ˙ 2 x ˙ 3 x ˙ 4 = 0 − m 1 k 1 + k 2 0 m 2 k 2 1 − m 1 c 1 + c 2 0 m 2 c 2 0 m 1 k 2 0 − m 2 k 2 0 m 1 c 2 1 − m 2 c 2 x 1 x 2 x 3 x 4 + 0 m 1 1 0 0 0 0 0 m 2 1 [ u 1 u 2 ]
此时方程的形式仍然为 x ˙ = A x + B u \dot{x} = Ax + Bu x ˙ = A x + B u ,系统矩阵 A \mathbf{A} A 的维度从单质量系统的 2 × 2 2 \times 2 2 × 2 扩展为了 4 × 4 4 \times 4 4 × 4 。
6 总结
综合上述推导与分析,我们探讨了从高阶常微分方程向状态空间方程过渡的原因,展示如何将一个单输入单输出系统或一个多输入多输出系统的微分方程转换为一个状态空间方程,并尝试通过状态空间方程分析一个系统的内部演化。我们将在之后的文章继续学习状态空间方程的各项特性以及其与传递函数之间的关系。
参考文献
本系列文章主要参考新西兰坎特伯雷大学ENME403-Linear Systems Control and System Identification课程资料。