第2章 控制系统的时域和频域描述,2.1 状态方程与时域描述 2.2 传递函数与频域描述,2.1 状态方程与时域描述,2.1.1 控制系统的状态空间描述 连续动态系统状态空间的一般形式可以写成,(2.1),(2.2),其中,F(X,t)表示系统所有的可变系数和非线性项。一般而言,系统的输入、输出和状态变量具有不同的维数。为了得到系统的完整描述,定义 G(t)=BU(t) (2.3) B、C、D矩阵一般为非方阵,式（2.1）可以写成,(2.4),式（2.2）可以写成,(2.5),例2.1 将下面的二阶系统表示成标准的状态方程形式。,解：写成标准的状态空间形式,2.1.2 状态方程的创建 假设n维线性微分方程为,(2.6),定义,对于第j个状态变量,(2.7),其中,j定义为,(2.8),对于式（2.7）和（2.8）,可以得到式（2.6）描述的SISO系统的矩阵表示,其中,(2.11),(2.12),例2.2 将下面的三阶线性系统表示成标准的状态空间形式。,解：按照前面介绍的方法,定义,其中,因此,例2.3 将下面的系统转换成标准状态空间形式。(注意到该系统方程右边没有输入的导数项,因此得到的系统矩阵的特征值与给定三阶方程解的特征方程的特征根相同。) 解：定义系统状态,写成矩阵形式,由 可得到,下面我们来验证第二个问题。系统方程的特征根可以写成 3+a12+a2+a3=0 状态矩阵的特征值为,沿矩阵的第一行展开,这里的det(A-I)称为特征方程。可以看出,状态矩阵的特征值与特征方程的根相同。 方程（2.4）的每一个方程中只含有一个导数项,称之为标准形式。然而对于一般的线性系统,在一个方程中可能会包含多个导数项。例如,下面的二阶系统,写成矩阵形式,更为一般的形式是,两边同乘以E-1,有,(2.13),(2.14),(E-1A)称为系统矩阵。,对于方程中含有代数方程的情况,可以通过一系列的代数运算来降低系统的维数。基本步骤包括： （1）重新将方程排序,使得前n1个方程包含导数项,后n2个方程仅包含代数项。 （2）使用矩阵重写原始方程（G=BU）,其中,Xd为包含导数项的状态向量,Xa为没有导数项 的状态向量。,（3）将上一步的矩阵展开，有,（4）得到Xa关于Xd的解,将它代入到微分方程中 ,可得,（5）将新的系统写成标准形式,(2.15),例2.4 将下面的系统表示成标准的状态方程形式。,解：将系统方程写成矩阵形式,上式具有以下形式,方程的解为,其中,2.1.3 非线性系统的线性化 实际上我们接触到的系统都是非线性系统,然而在某个参考状态的某个有限范围内可以采用近似线性化的分析方法。下面将介绍如何对一般的非线性状态方程进行线性化。 假设非线性系统的一般形式为,(2.16),F(X,U,t)包括系统所有的非线性项。系统的状态变 量和输入可以表示成,X0为非线性系统的参考点处的状态和输入。将式（2.17）代入式（2.16）,得到,(2.17),(2.18),将上式在参考点附近进行一阶Taylor展开,对于其中第i个等式,其一阶近似为,(2.19),写成矩阵形式,(2.20),其中,Jx(X0,U0)和Ju(X0,U0)是系统在参考点处的Jacobian矩阵，即,(2.21),系统在参考点附近同样有,(2.22),将式（2.20）和式（2.22）代入式（2.18）,得到,(2.23),式（2.23）表示的是原非线性系统的线性化模型系统。新系统的状态矩阵为A+Jx(X0,U0),同时新的输入矩阵为B+Ju(X0,U0)。,例2.5 计算下列矩阵表示的非线性系统在平衡点( )处的线性化模型。,解：首先计算系统的参考状态。在平衡点处,从而,得到两个平衡点的状态为,计算系统在平衡点处的Jacobian矩阵,在参考状态XT0=0 0处,线性化系统为,而在参考状态XT0=15 -5处,线性化系统为,2.1.4 线性系统的解析解 下面将讨论如何计算一个状态方程描述的线性系统的时域解。在讨论过程中,读者可以看到,矩阵的指数函数在系统解的计算中发挥了重要的作用。 1.无输入的情况 首先考虑一个没有独立输入变量并且只有一个状态变量的最简单情况。系统描述,(2.24),假定解的形式为,(2.25),将其代入原方程,方程的最终解为,(2.26),对于矩阵的情况,即系统含有一个以上的状态变量,(2.27),与标量情况类似,假设解的形式为,(2.28),代入式（2.27），得,由初始条件X(0)=B=X0得,(2.29),以上结果也可以表示成不同的形式。例如,可以以 t=t0为初始条件重新计算,则式（2.29）变为,(2.30),如果采用t=t-t0,有,(2.31),综上所述,无输入线性状态方程的解为,(2.32),其中,(t)称为状态转移矩阵,它是下面方程的唯一解,(2.33),为了验证式（2.23）确实是原系统方程的解,考虑,(2.34),2. 有输入的情况 首先考虑标量（只有一个状态变量）的情况。系统的状态方程为,(2.35),方程两边同时乘以积分因子 ,得到,将方程在t0与t的区间内进行积分,有,最后,两边乘以eat,重新整理方程得,(2.36),如果t0=0,则,(2.37),验证上述解的正确性的方法是将其代入微分方程中去。例如,将式（2.37）代入方程（2.35）中,并根据Leibnitz法则,得到,(2.38),对于矩阵情况,系统描述为,(2.39),与标量的情况类似,方程两边同时乘以积分因子e-At,得到,在t0与t时间区间内进行积分,从而有,最后,两边同时乘以eAt,得到,(2.40),如果t0=0,则,(2.41),同样,为了验证上述解的正确性,将它代入原系统方程,上述解也可以用状态转移矩阵来表示,例如,方程（2.41）采用转移矩阵表示为,2.1.5 线性系统的离散化 正如上面讨论的那样,在对实际系统进行分析和仿真之前,往往需要首先采用计算机计算出系统的解的情况。在高维动态系统的计算机仿真中一般采用两种方法,包括线性系统解析解的离散化和适用于任何系统的数值积分技术。 系统离散化的目标是将线性系统的连续状态方程描述转化成离散形式。假设系统状态方程为,(2.43),将其转化成离散状态方程形式 X(k+1)T=G(T)X(kT)+H(T)U(kT) (2.44) 其中,T是采样周期,G和H是常值矩阵。为了方便,方程（2.44）经常写成下面的形式 Xk+1=GXk+Huk (2.45) 上述方程表示了系统状态随离散时间的迭代关系,如果G和H矩阵都已知,就很容易通过计算机迭代计算系统在各个时刻的状态值。因此,下面的目标就是如何计算离散状态矩阵的值。,为此,首先假设输入U(t)为分段常值的函数,即U(t) Uk。为了推导系统的离散表示,引出前面讲述的连续线性系统的解析解，即,定义t=(k+1)T,并且t0=kT,有,引进新的变量,使得=+kT或者=-kT,从而 d=d (k+1)T-kT=T- U(+kT)=U(kT)=Uk,可以得到,由于,得到最后的解,(2.46),将它与式（2.45）进行比较,最后得到G和H矩阵的计算公式,如果将G按照指数公式展开,有,则可以得到H的展开形式,2.1.6 状态方程的数值积分 在对系统状态方程进行仿真中除了使用上面讲述的矩阵指数函数方法外,还可以采用直接对状态方程进行积分的方法。采用直接积分方法的优点在于它可以很容易地处理时变和相对复杂的非线性系统。 1Euler方法（标量情况） 考虑一般的一阶微分方程,(2.51),在Euler方法中,我们假定方程（2.51）的右边在某个时间段t=tk+1-tk保持常值。这样,方程（2.51）的积分方程可以写成,或者,(2.52),2改进的Euler方法（标量情况） 对积分方程更好的近似方法是假定方程在t时间段内是线性（而不是常值）的,这样得到下面的迭代关系,(2.53),式（2.53）的问题在于计算xk+1时需要计算fk+1,一种可能的解决方法是首先得到xk+1和fk+1的预测值,然后根据情况对xk+1进行修正,以改进第一次得到的预测值。这种方法也称为预测-修正方法。计算步骤为先进行预测计算,(2.54),再进行修正计算,其中,xk+1代表tk+1时刻第一次的预测值, tk+1代表修 正后的最优值。,改进的Euler方法虽然算法简单,却展示了预测-修正数值积分方法的主要思想。我们还可以在计算过程中自适应地调整t的大小,来同时满足计算精度和计算速度的要求。在最新的算法中,只要输入数值积分的容许误差,积分算法将会自动调整计算的步长,以满足计算精度的要求。MATLAB/Simulink中采用的变步长ODE求解算法普遍采用的就是这类算法。 以上关于标量情况的讨论可以很容易地推广到矩阵形式。假设一维微分方程的一般形式是,(2.55),Euler方法（矩阵情况）,其中,改进的Euler方法（矩阵情况） 预测计算,(2.56),(2.57),修正计算,(2.58),其中,(2.59),(2.60),2.1.7 实例 最后通过一个具体的例子来演示前述内容。 例2.6 二阶线性系统的一般形式为,将它写成状态方程形式,有,其中 0=b0 1=b1-a1b0 2=b2-a1b1+a21b0-a2b0,x1(0)和x2(0)为系统的初始条件。 作为一个特殊的例子,图2.1显示的是一个简单的RLC电路。 电路的电压满足平衡方程 ea=eL+eR+eC 其中,图2.1 简单的RLC电路图,i(t)代表电路中的电流。将这些关系代入平衡方程，得,将方程两边进行微分,将上式与二阶系统的一般形式进行比较,可以看出,这样,可以得到标准的状态方程形式,其中,为了得到系统的解析解,假定R=100,L=0.1H,C=0.001F。ea(t)是10V的阶跃信号。下面我们采用三种不同的方法进行时域仿真： 1）用状态转移矩阵表示的系统解析解 对于阶跃输入,LTI系统的解析解可以写成,上式可以写成简化形式,最后一步,我们利用了A-1和eAt可交换的事实。 为了得到最后的表达式,我们使用Sylvester定理来计算状态转移矩阵eAt,首先计算状态矩阵的特征值：,或者,2）LTI系统的离散化 LTI状态方程的离散化方法包含在MATLAB中的控制工具箱（ControlToolbox）中。这个工具箱还包括一整套用于线性系统分析的函数。MATLAB是通过定义标准状态方程中的A、B、C和D矩阵的方法来定义一个LTI状态空间对象： sys=ss(A,B,C,D) 上式将创建一个名为sys的LTI系统对象,接下来可以使用其它函数对其进行各种操作。 例如,对阶跃信号的响应可以通过下式进行计算: Y,T,X=step(sys),或 Y,T,X=step(sys,T) 其中,Y和X分别包含输出和状态向量的时域仿真数据,它是一个三维数组,其第3维对应于输入的维数。 对于第一种用法,T向量中的采样时间数与采样周期是由MATLAB自动确定的。而在第二种用法中,它们是由用户自己定义的。除了使用step函数,MATLAB中具有类似功能的函数还有impulse和lsim,它们分别完成系统的单位阶跃响应和一般输入信号的仿真。对于本例中的RLC电路,调用的格式为,sys=ss(A,10*B,C,D); Y,T,X=step(sys,T); 3）状态方程的数值积分 本例中演示的最后一种方法是通过MATLAB的ode23函数完成的数值积分方法。该函数运用了前面讲述的自适应步长控制算法,采用不同阶的RungeKutta（RK）积分算法来进行误差估计。具体来说,ode23采用二阶和三阶RK算法来对微分方程进行积分。