应用PRO/II进行塔模拟计算时热力学系统的选择应用PRO/II进行塔模拟计算时热力学系统的选择
部分计算模块的讨论
王崇智
石油大学（北京）化学工程研究所（102200）
摘要: 本文总结了应用PRO/II进行炼油、石油化工装置模拟过程中塔模拟计算时热力学系统选择的一些经验，并讨论了塔模拟计算时的一些基本要素和基本技巧。
关键词：PRO/II，塔，模拟，热力学方法，技巧
一、概述
PRO/II软件是美国SIMSCI公司推出的微机版本石油化工工艺流程模拟软件，该软件具备有丰富的物性数据库和热力学方程供用户描述不同状态下的流体热力学过程，对多种炼油、化工工艺过程具有广泛的适应性。该软件不仅可以作为新设计炼油、化工工艺装置的工艺流程模拟软件，同时作为装置标定计算、设备核算的软件。
PRO/II软件在我国的应用十分广泛，其中DOS系统的V3.3、V4.02版本和WINDOWS 操作系统的V4.13 WITH PROVISION V2.0以上版本是比较常用的。PRO/II软件是很多炼油、化工等设计院进行工艺设计的首选工艺模拟软件之一，同时也是炼油、化工等生产单位进行装置标定计算、设备核算的首选工艺模拟软件之一。在实际工作中，有很多时候会遇到解决装置“瓶径”的问题，而塔设备往往是需要进行标定或核算的重要设备之一，这时应用PRO/II软件提供的精馏、吸收、萃取等单元操作过程的严格计算方法进行单塔模拟计算或全流程模拟计算是非常方便的。
本文通过对PRO/II软件提供的热力学方法的初步分析，推荐了PRO/II软件在常见石油化工装置塔模拟计算中的热力学方法。
二、PRO/II热力学方法的初步分析
PRO/II软件提供多种用于流体的气液平衡常数、液液平衡常数、焓、熵、密度和其他传递性能参数等热力学计算方法，由于每种热力学方法有一定的适用范围，在应用PRO/II解决具体问题时，选择合适的热力学方法是能否正确模拟工艺过程的关键。以下分类讨论PRO/II提供的主要的热力学方法。
2.1、普遍化方法
普遍化方法主要包括用于烃类物系计算的SRK方程、PR方程、BWRS方程、GS方程、IGS方程、BK10方程等，各方程的适用范围如下：
热力学方程
适用领域
SRK
气体、炼油过程的烃类物系
SRKKD
炼油过程的烃水物系，尤其高温、高压的气液液过程
SRKM
烃/醇等极性/非极性物系
SRKH
酮/水等极性和高压物系
PR
气体、炼油过程的烃类物系
PRM
烃/醇等极性/非极性物系的气液液过程
PRH
酮/水等极性和高压物系
BWRS
气体、炼油过程的烃类物系
GS
常压以上的炼油物系
IGS
炼油体系的气液液过程
BK10
原油常压、减压蒸馏过程
2.2、液相活度系数方法
液相活度系数方法主要包括用化工、石油化工物系气液、液液、气液液平衡及相关物性参数计算的NRTL（Non-Random Two Liquid）方程、UNIQUAC方程、WILSON 方程、UNIFAC方程、VANLAAR方程、FLORY方程、MARGULES方程等，各方程的适用范围如下：
热力学方程
适用领域
NRTL
有液相活度系数可以利用的化工、石油化工极性物系
UNIQUAC
没有提供气液、液液平衡数据的化工、石油化工极性物系
WILSON
极性物系的气液过程
UNIFAC
任何已知组分结构的物系
VANLAAR
化工、石油化工极性物系的气液、液液过程
FLORY
化工、石油化工极性物系的气液、液液过程
MARGULES
化工、石油化工极性物系的气液、液液过程
2.3、专用数据包方法
PRO/II专用数据包用于计算指定物系的气液、液液平衡及相关物性参数，主要包括GLYCOL数据包、SOUR WATER数据包、ALCOHOL数据包、AMINE数据包等，各专用数据包的适用范围如下：
专用数据包
适用领域
GLYCOL
含有水、乙二醇、三甘醇和气体组分物系的脱水过程。
温度26—204℃；压力≤13.6MPa
SOUR WATER
含有H2S、NH3、CO2、H2O物系的计算。温度20—150℃；压力≤0.345MPa
（原始关联式）；压力≤19.3MPa（修正Van Der Waals方程计算气相逸度）
ALCOHOL
含有醇、H2O和其他极性物系气液、液液、气液液过程。
温度50—110℃；压力≤10.2MPa
AMINE
含有MEA、DEA、DGA、DIPA、MDEA物系的相关计算
三、PRO/II在石油化工装置塔模拟中的热力学方法
石油化工装置种类繁多，以下将分类介绍PRO/II软件在部分装置塔模拟计算推荐采用的平衡常数的热力学计算方法和相应的数据包。
3.1、炼油装置
常见炼油装置塔模拟计算推荐采用的平衡常数的热力学计算方法和相应的数据包如下。
装置名称
塔名称
热力学方法
数据包
备注
常减压装置
初馏塔
GS
SIMSCI
常压塔
GS
SIMSCI
常压汽提塔GS
SIMSCI
减压塔
BK10 SIMSCI
催化裂化装置分馏塔
BK10 SIMSCI
柴油汽提塔BK10 SIMSCI
吸收塔SRK、PR SIMSCI
解吸塔SRK、PR SIMSCI
再吸收塔SRK、PR SIMSCI
稳定塔SRK、PR SIMSCI
催化重整装置脱轻塔SRK、PR
SIMSCI
脱重塔SRK、PR SIMSCI
抽提塔
NRTL ALCOHOL
非芳水洗塔NRTL ALCOHOL
汽提塔
NRTL ALCOHOL
溶剂回收塔NRTL ALCOHOL
苯塔
BK10 SIMSCI
甲苯塔
BK10 SIMSCI
二甲苯塔
BK10 SIMSCI
延迟焦化装置分馏塔
BK10 SIMSCI
柴油吸收塔BK10、GS SIMSCI
加氢裂化装置分馏塔
BK10、GS SIMSCI
脱乙烷塔SRK、PR SIMSCI
脱丁烷塔
BK10、GS SIMSCI
石脑油汽提塔BK10、GS SIMSCI
加氢装置
汽提塔
SRK
SIMSCI
稳定塔
SRK
SIMSCI
MTBE装置
催化蒸馏塔NRTL SIMSCI
催化精馏流程
碳四分离塔NRTL SIMSCI
两器四塔流程
甲醇萃取塔NRTL ALCOHOL
甲醇回收塔
NRTL
ALCOHOL
酸水汽提装置
酸水汽提塔
SOUR
SIMSCI
脱硫装置
干气脱硫塔
AMINE&SRK
SIMSCI
液化气脱硫塔
AMINE&SRK
SIMSCI
溶剂再生塔
AMINE&SRK
3.2、石油化工装置
常见石油化工装置塔模拟计算推荐采用的平衡常数的热力学计算方法和相应的数据包如下。
装置名称
塔名称
热力学方法
数据包
备注
乙烯装置
油洗塔
SRK、PR
SIMSCI
水洗塔
SRK、PR
SIMSCI
SRK、PR SIMSCI
脱乙烷塔SRK、PR SIMSCI
加入Kij
乙烯塔SRK、PR SIMSCI
加入Kij
脱丙烷塔SRK、PR SIMSCI
加入Kij
丙烯塔SRK、PR SIMSCI
加入Kij
脱丁烷塔SRK、PR SIMSCI
环氧乙烷装置吸收塔
NRTL SIMSCI
解吸塔
NRTL SIMSCI
回收塔
NRTL SIMSCI
精制塔
NRTL
乙二醇装置脱水塔NRTL ALCOHOL
乙二醇塔NRTL ALCOHOL
二甘醇塔NRTL ALCOHOL
三甘醇塔NRTL ALCOHOL
乙苯装置吸收塔SRK SIMSCI
稳定塔
SRK SIMSCI
循环苯塔BK10、SRK SIMSCI
脱甲苯塔BK10、SRK SIMSCI
乙苯精馏塔BK10、SRK SIMSCI
多乙苯塔BK10、SRK SIMSCI
苯乙烯装置
苯-甲苯塔
BK10、SRK
SIMSCI
乙苯回收塔
BK10、SRK
SIMSCI
苯乙烯精馏塔
BK10、SRK
SIMSCI
残液精制塔
BK10、SRK
SIMSCI
醋酸乙烯装置
醋酸吸收塔
NRTL、SRKM
SIMSCI
洗涤塔
NRTL、SRKM
SIMSCI
二氧化碳吸收塔 NRTL、SRKM SIMSCI
二氧化碳再生塔 NRTL、SRKM SIMSCI
初馏塔
NRTL、SRKM
SIMSCI
脱水塔
NRTL、SRKM
SIMSCI
脱轻塔
NRTL、SRKM
SIMSCI
脱重塔
NRTL、SRKM
SIMSCI
醋酸装置
醋酸吸收塔
NRTL、SRKM
SIMSCI
洗涤塔
NRTL、SRKM
SIMSCI
二氧化碳吸收塔 NRTL、SRKM
SIMSCI
二氧化碳再生塔 NRTL、SRKM
SIMSCI
四、PRO/II中的重要计算模块
1、RECYCLE模块的应用
很多装置的模拟需要采用RECYCLE进行工艺计算，如：带外部冷回流的分馏过程、萃取蒸馏单元、吸收—解吸过程、洗涤—再生过程等。
在具有循环物流的计算过程中，正确的确定循环物流的初置、参考物流的温度或流量和循环加速对于模拟计算是必要的。
2、CALCULATOR模块的应用
对于有多股物流的工艺过程，PRO/II提供的CALCULATOR模块能够方便的计算各物流的组合工况或物流性质，灵活运用CALCULATOR模块对于完成全流程模拟计算的是一个有效的手段。
3、OPTIMIZER模块的应用
甘醇脱水PR 酸性水汽提SOUR PR
天然气深冷分离PR, PRSV
空气系统PR, PRSV
常压塔PR, PR Options, GS
减压塔PR, PR Option, GS (<10mmHg), Braun K10
乙烯精馏塔Lee Kesler Plocker
高含氢系统PR, ZJ, GS
油藏系统PR, PR Options
水蒸气系统Steam Package, CS, GS
烃类系统，但重点考虑水在烃系统的溶解度Kabadi Danner 抑制水化物PR
极性化工系统活度模型，PRSV
烷基化装置PRSV, NRTL
1、齐次状态方程解 【例】已知线性定常系统的齐次状态方程为 x x ?? ????--=3210 & 试求该状态方程的解。 解 这里我们应用拉氏变换法求系统的状态转移矩阵。 首先计算矩阵 ? ? ????+-=??????---??????=-3213210 00)(s s s s A sI 其次，计算1 ()sI A --及状态转移矩阵()t Φ ??????--++-= --s s s s A sI 2133 211 )(1 ???? ? ???? ??? ++++-+++++=)2)(1() 2)(1(2)2)(1(1)2)(1(3s s s s s s s s s s ?? ?? ???? ??++ +-++ +-+- ++-+=2211221221112 112 s s s s s s s s ?? ? ???+-+---=-==Φ----------t t t t t t t t At e e e e e e e e A sI L e t 222211 2222])[()( 则齐次状态方程的解为 )0(2222)(2222x e e e e e e e e t x t t t t t t t t ??? ???+-+---=-------- 或者 ??? ? ??????? ?+-+---=???? ??--------)0()0(2222)()(21222221x x e e e e e e e e t x t x t t t t t t t t 2、线性变换例题 【例9.15】已知系统具有如下形式 u y y y y 66116')2()3(=+++
第2章 控制系统的状态方程求解 要点： ① 线性定常状态方程的解 ② 状态转移矩阵的求法 ③ 离散系统状态方程的解 难点： ① 状态转移矩阵的求法 ② 非齐次状态方程的解 一 线性定常系统状态方程的解 1 齐次状态方程的解 考虑n 阶线性定常齐次方程 ? ? ?==0)0()()(x x t Ax t x & （2-1） 的解。 先复习标量微分方程的解。设标量微分方程为 ? ??==0)0(x x ax x & （2-2） 对式（2-2）取拉氏变换得 )()(0s aX X s sX =- 移项 0)()(x s X a s =- 则 a s x s X -= )(
取拉氏反变换，得 00 0!)()(x k at x e t x k k at ∑∞ === 标量微分方程可以认为是矩阵微分方程当n=1时的特征，因此矩阵微分方程的解与标量微分方程应具有形式的不变性，由此得如下定理： 定理2-1 n 阶线性定常齐次状态方程（2-1）的解为 00 0!)()(x k At x e t x k k At ∑∞ === （2-3） 式中，∑∞ ==0 !)(k k At k At e 推论2-1 n 阶线性定常齐次状态方程 ???==00 )()()(x t x t Ax t x & （2-4） 的解为 0)(0 )(x e t x t t A -= （2-5） 齐次状态方程解的物理意义是)(0 t t A e -将系统从初始时刻0t 的初始 状态0x 转移到t 时刻的状态)(t x 。故)(0 t t A e -又称为定常系统的状态转移 矩阵。 (状态转移矩阵有四种求法：即定义（矩阵指数定义）法、拉氏反变换法、特征向量法和凯来-哈密顿（Cayly-Hamilton ）法) 从上面得到两个等式 ∑∞ ==0 !)(k k At k At e ])[(11---=A sI L e At 其中，第一式为矩阵指数定义式，第二式可为At e 的频域求法或拉氏反变换法
第三章 系统的分析——状态方程的解 §3-1线性连续定常齐次方程求解 一、齐次方程和状态转移矩阵的定义 １、齐次方程 状态方程的齐次方程部分反映系统自由运动的状况（即没有输入作用的状况），设系统的状态方程的齐次部分为： )()(t Ax t x =& 线性定常连续系统： Ax x =& 初始条件：00x x t == ２、状态转移矩阵的定义 齐次状态方程Ax x =&有两种常见解法：（1）幂级数法；（2）拉氏变换法。其解为 )0()(x e t x At ?=。其中At e 称为状态转移矩阵（或矩阵指数函数、矩阵指数），记为： At e t =)(φ。 若初始条件为)(0t x ，则状态转移矩阵记为：) (0 0)(t t A e t t -=-Φ 对于线性时变系统，状态转移矩阵写为),(0t t φ，它是时刻t ，t 0的函数。但它一般不能写成指数形式。 （１）幂级数法——直接求解 设Ax x =&的解是t 的向量幂级数 Λ ΛΛΛ+++++=k k t b t b t b b t x 2210)( 式中ΛΛ，，， ，，k b b b b 210都是n 维向量，是待定系数。则当0=t 时， 000b x x t === 为了求其余各系数，将)(t x 求导，并代入)()(t Ax t x =&，得： Λ ΛΛΛ&+++++=-1232132)(k k t kb t b t b b t x )(2210ΛΛΛΛ+++++=k k t b t b t b b A
上式对于所有的t 都成立，故而有： ????? ??????======00 3 230 21201!1!31312121b A k b b A Ab b b A Ab b Ab b K K M 且有：00x b = 故以上系数完全确定，所以有： Λ ΛΛΛ+++++=k k t b t b t b b t x 2210)( ΛΛ++++ +=k k t b A k t b A t Ab b 020200! 1 !21 )0()! 1!21(22x t A k t A At I k k ΛΛ+++++= 定义（矩阵指数或矩阵函数）： ∑∞==+++++=022! 1!1!21K k k k k At t A k t A k t A At I e ΛΛ 则 )0()(x e t x At ?=。 （２）拉氏变换解法 将Ax x =&两端取拉氏变换，有 )()0()(s AX X s sX =- )0()()(X s X A sI =- )0()()(1X A sI s X ?-=- 拉氏反变换，有 )0(])[()(1 1x A sI L t x ?-=--
Chapter2状态方程的解 我们要解决的问题是：在系统初始时刻0t t =时，初始状态为00)(x t x =的条件下，对该系统施加控制)(t u ，求出系统状态)(t x 的变化，即求解非齐次方程 （0)(≠t u ）初值问题的解： 00 0)()()()()()(t t x t x t u t B t x t A t x ≥=+=& 或者在系统不加控制)(t u ，（0)(=t u 称为自由系统）的条件下，求出初值)(0t x 对系统状态)(t x 的影响，即求解齐次方程初值问题的解： 00 0)(),()()(t t x t x t x t A t x ≥==& ????离散连续线性定常????离散连续线性时变?? ?? ? ??????数值解解析解非齐次数值解解析解齐次 2.1 线性定常系统状态方程的解 2.1.1 n 阶、线性、定常（无关与时间t A ）连续系统齐次状态方程的解 我们知道：常系数线性微分方程（标量方程）)()(t ax t x =&，0)0(x x =，0≥t 其解为 00 0!)(x k t a x e t x k k k at ∑∞ === 对齐次状态方程（矩阵方程） )()(t Ax t x =&，0)0(x x =，0≥t 很自然，仿照常系数线性微分方程，可得到n 阶线性、定常、连续系统齐次（0)(=t u ）状态方程的解 000! )(x k t A x e t x k k k At ∑ ∞ === 定义矩阵指数：k k k k k At t A k t A At I k t A e ! 1 21!220 ++++=≡∑ ∞ =Λ，它仍是一个矩阵。 若初始时间为0t ，则状态方程的解为 00 00) (!)()(0x k t t A x e t x k k k t t A ∑∞ =--== ∑ ∞ =--=0 0) (! )(0k k k t t A k t t A e 称为定常（连续）系统的状态转移矩阵。 )(0t t A e -物理意义：将系统从初始状态)(0t x 转移到（时刻t 的）状态)(t x 。 2.1.2 矩阵指数At e 的性质
求解系统的状态方程 一、实验设备 PC计算机，MATLAB软件，控制理论实验台 二、实验目的 (1)掌握状态转移矩阵的概念。学会用MATLAB求解状态转移矩阵 (2)学习系统齐次、非齐次状态方程求解的方法，计算矩阵指数，求状态响应； (3)通过编程、上机调试，掌握求解系统状态方程的方法，学会绘制输出响应和状态响应曲线； (4)掌握利用MATLAB导出连续状态空间模型的离散化模型的方法。 三、实验原理及相关基础 (1)参考教材P99~101“3.8利用MATLAB求解系统的状态方程” (2)MATLAB现代控制理论仿真实验基础 （3）控制理论实验台使用指导 四、实验内容 (1)求下列系统矩阵A对应的状态转移矩阵 （a）
（b） 代码： syms lambda A=[lambda 0 0;0 lambda 0;0 0 lambda];syms t;f=expm(A*t) （c） 代码： syms t;syms lambda;A=[lambda 0 0 0;0 lambda 1 0;0 0 lambda 1;0 0 0 lambda];f=expm(A*t) (2) 已知系统
a) 用MATLAB求状态方程的解析解。选择时间向量t，绘制系统的状态响应曲线。观察并记录这些曲线。 （1） 代码： A=[0 1; -2 -3]; B=[3;0]; C=[1 1]; D=[0]; u=1; syms t; f=expm(A*t);%状态转移矩阵 x0=0; s1=f*B*u; s2=int(s1,t,0,t)%状态方程解析解 状态曲线： （2）A=[0 1;-2 -3]; syms t; f=expm(A*t); X0=[1;0]; t=[0:0.5:10]; for i=1:length(t); g(i)=double(subs(f(1),t(i))); end plot(t,g)
第2章 状态空间表达式的解 第1节 线性定常齐次状态方程的解 线性定常齐次状态方程 0(0)x Ax x x ==& 的解为 0()At x t e x = (0)t > 式中，2 2 () 2!!k At k t At e I At A k ∞ ? ==+++=∑ L 证明： 用拉普拉斯变换法。 对 x A x =& 作拉氏变换，得 0()()sX s x AX s -=
1 0()()X s sI A x -=- 11 0()[()]x t L sI A x --=- 因为 2 23111()()sI A I A A I s s s -+++=L 故 1 223111()sI A I A A s s s --=+++L 12023111()[]x t L I A A x s s s -=+++L 2201()2! I At A t x =+++L 0At e x = 顺便可知 ])[(1 1---=A sI L e At 第2节 矩阵指数函数At e 1、At e 的定义和性质
（1）定义 2 2 () 2!!k At k t At e I At A k ∞ ==+++=∑ L 式中 A —线性定常系统系统矩阵，n n ?阶； At e —矩阵指数函数，n n ?阶时变矩阵。 若A 中各元素均小于某定值，At e 必收敛；若A 为实矩阵，At e 绝 对收敛。 （2）基本性质： ◆组合性质： ) (2121t t A At At e e e += 其中21,t t 为相衔接的两时间段。 推论1：I e e e e A t t A t A At ===--0 ) () ( 推论2：) (1 ][t A At e e --=
第三章控制系统状态方程求解 3－1 线性连续定常齐次方程求解 所谓齐次方程解，也就是系统的自由解，是系统在没有控制输入的情况下，由系统的初始状态引起的自由运动，其状态方程为： ………………………………………………………(3 －1) 上式中，X是n×1维的状态向量，A是n×n的常数矩阵。 我们知道，标量定常微分方程的解为： ………………（3 －2） 与（3－2）式类似，我们假设（3－1）的解X(t)为时间t的幂级数形式，即： ………………………………(3 －3) 其中为与X（t）同维的矢量。 将（3－3）两边对t求导，并代入（3－1）式，得：
上式对任意时间t都应该成立，所以变量t的各阶幂的系数都应该相等，即： 即: ……………………………………………（3－4） 将系统初始条件代入（3－3），可得。代入（3－4）式可得： (3) 5） 代入（3－3）式可得（3－1）式的解为：
(3) 6) 我们记： (3) 7） 其中为一矩阵指数函数，它是一个n×n的方阵。所以（3－6）变为： (3) 8） 当（3－1）式给定的是时刻的状态值时，不难证明： (3) 9） 从（3－9）可看出，形式上是一个矩阵指数函数，且也是一个各元素随时间t变化的n×n矩阵。但本质上，它的作用是将时刻的系统状态矢量转移到t时刻的状态矢量，也就是说它起到了系统状态转移的作用，所以我们称之为状态转移矩阵(The State Transition Matrix)，并记： (3) 10） 所以：
【例3－1】已知，求解：根据（3－7）式， 3－2 的性质及其求法 性质1： 【证】根据的定义式（3－7）， 【证毕】 性质2：① ②
2-3 由控制系统的方块图求系统状态空间表达式 系统方块图是经典控制中常用的一种用来表示控制系统中各环节、各信号相互关系的图形化的模型，具有形象、直观的优点，常为人们采用。 要将系统方块图模型转化为状态空间表达式，一般可以由下列三个步骤组成： 第一步：在系统方块图的基础上，将各环节通过等效变换分解，使得整个系统只有标准积分器（1/s ）、比例器（k ）及其综合器（加法器）组成，这三种基本器件通过串联、并联和反馈三种形式组成整个控制系统。 第二步：将上述调整过的方块图中的每个标准积分器（1/s ）的输出作为一个独立的状态变量i x ，积分器的输入端就是状态变量的一阶导数 dt dx i 。 第三步：根据调整过的方块图中各信号的关系，可以写出每个状态变量的一阶微分方程，从而写出系统的状态方程。根据需要指定输出变量，即可以从方块图写出系统的输出方程。 例2-5 某控制系统的方块图如图2-6所示，试求出其状态空间表达式。 解： 该系统主要有一个一阶惯性环节和一个积分器组成。对于一阶惯性环节，我们可以通过等效变换，转化为一个前向通道为一标准积分器的反馈系统。 图2-6所示方块图经等效变换后如下图所示。我们取每个积分器的输出端 信号为状态变量1x 和2x ，积分器的输入端即1x 和2x 。 图2-6 系统方块图
从图可得系统状态方程： ()??? ??? ?+--=-+-==u T K x T x T K K x K u T K x T x x T K x 11211131131121222 2111 取y 为系统输出，输出方程为：1x y = 写成矢量形式，我们得到系统的状态空间表达式： []?????????? ?=???? ??????+? ???????=x y u T K x K K T K x 010********
由传递函数转换成状态空间模型——方法多!!! SISO 线性定常系统 高阶微分方程化为状态空间表达式 SISO ()()()()()()m n u b u b u b y a y a y a y m m m n n n n ≥+++=++++--- 1102211 )(2 211110n n n n m m m a s a s a s b s b s b s G +++++++=--- 假设1+=m n 外部描述 ←—实现问题：有了部结构—→模拟系统 部描述 SISO ? ??+=+=du cx y bu Ax x 实现问题解决有多种方法，方法不同时结果不同。 一、 直接分解法 因为 1 0111 11()()()()()()()() 1m m m m n n n n Y s Z s Z s Y s U s Z s U s Z s b s b s b s b s a s a s a ----?=? =?++++++++ ???++++=++++=----) ()()() ()()(11 11110s Z a s a s a s s U s Z b s b s b s b s Y n n n n m m m m 对上式取拉氏反变换，则 ? ??++++=++++=----z a z a z a z u z b z b z b z b y n n n n m m m m 1) 1(1)(1)1(1)(0 按下列规律选择状态变量，即设)1(21,,,-===n n z x z x z x ，于是有
?????? ?+----===-u x a x a x a x x x x x n n n n 12113 221 写成矩阵形式 式中，1-n I 为1-n 阶单位矩阵，把这种标准型中的A 系数阵称之为友阵。只要系统状态方程的系数阵A 和输入阵b 具有上式的形式，c 阵的形式可以任意，则称之为能控标准型。 则输出方程 121110x b x b x b x b y m m n n ++++=-- 写成矩阵形式 ??????? ? ????????=--n n m m x x x x b b b b y 12101 1][ 分析c b A ,,阵的构成与传递函数系数的关系。 在需要对实际系统进行数学模型转换时，不必进行计算就可以方便地写出状态空间模型的A 、b 、c 矩阵的所有元素。 例：已知SISO 系统的传递函数如下，试求系统的能控标准型状态空间模型。 4 2383)()(2 3++++=s s s s s U s Y 解：直接得到系统进行能控标准型的转换，即
2.5 控制系统的状态空间表达式 2.5 控制系统的状态空间表达式 随着科学技术的发展，被控制的对象越来越复杂，对自动控制的要求也越来越高。面对时变系统，多输入多输出系统、非线性系统等被控量和对控制系统高精度、高性能的严格要求，传统的控制理论已不能适用。同时，计算机技术的发展也要求控制系统地分析，设计中采用计算机技术并在控制系统的组成中使用计算机。因此，适用这些要求的控制系统的另一种数学描述方法----状态空间就应运而生。 2.5.1 状态变量 在对系统动态特性描述中，足以表征系统全部运动状态的最少一组变量，称之为状态变量。只要确定了这组变量在t=时刻的值以及时的输入函数，则系统在任何时刻的运动 状态就会全部确定。状态变量互相间是独立的，但对同一个系统，状态变量的选取并不是唯一的。一个用n 阶微分方程描述的系统，有n个独立变量，这n个独立变量就是该系统的状态变量。 若用表示这n个状态变量，则可以把这n个状态变量看作是向量x(t) 的分量。我们称x(t)为状态变量，它是一个n维向量，记为 分别以状态变量作为坐标而构成的n维空间，称为状态空间。系统在t时刻的状态，就是状态空间的一点。系统在时刻的状态称为初始点，随着时间的变化， x(t)从初始点出发在状态空间描述出一条轨迹，称为状态轨迹。状态魁及表征了系统状态的变化过程。 2.5.2 状态空间表达式 1. 状态方程 由系统的状态变量和输入函数构成的一阶微分方程组，称为系统的状态方程。 对于线性系统，可以写成如下形式
(2.59) 记为 (2.60) 式中x(t)是n维列向量 u(t)是r维输入向量 A是n*n维矩阵，称为系数矩阵 B是n*r矩阵，称为输入矩阵或控制矩阵
例9.20试求下列状态方程在)(1)(t t u =作用下的解。 101111x x u ????=+????????&，1(0)0x ??=???? 解 (1) 积分法 已知非齐次状态方程解的形式为 ?-Φ+Φ=t d Bu t x t t x 0 )()()0()()(τττ ()1 11112101001()1111(1) 1t t t s e s t sI A s te e s s -----?? ?? ??-?? ??-????Φ=-===??????????--?????? ??????--?? L L L 在)(1)(t t u =作用下，为了简化计算，令ττ-=t ' ，有 ?Φ-+Φ=0 )()0()()(t Bd x t t x ’‘ττ ?Φ+Φ=t Bd x t 0 '')()0()(ττ ?Φ+Φ=t Bd x t 0 )()0()(ττ 则： 00101()1t t t t e e Bd d e e te τ τ τττττ????-?? Φ==????????? ??? ?? 所以， 10121()02t t t t t t t e e e x t te e te te ?????? --??=+=???????????????? (2) 拉氏变换法 ()()11 11()(0)()x t sI A x sI A BU s ----????=-+-???? L L ()1111 2 2 110111(0)1 110(1)1(1)t t e s s sI A x te s s s ----???? ????????--??-===? ??????????????? ??????---???? L L L ()11112211 01(1)111()1101(1)1(1)t t s s e s sI A BU s s te s s s ----?????? ??????-??-??-??????-===? ????????? ?????????? ??? ?--??-?????? L L L ()() 1 1 1 1 21()(0)()2t t e x t sI A x sI A BU s te ----?? -????=-+-=???? ? ? ?? L L
状态空间方程 以传递函数为基础的经典控制理论的数学模型适应当时手工计算的局限，着眼于系统的外部联系，重点为单输入-单输出的线性定常系统。伴随计算机的发展，以状态空间理论为基础的现代控制理论的数学模型采用状态空间方程，以时域分析为主，着眼于系统的状态及其内部联系，研究的机电控制系统扩展为多输入-多输出的时变系统。 所谓状态方程是由系统状态变量构成的一阶微分方程组;状态变量是足以完全表征系统运动状态的最小个数的一组变量。状态变量相互独立但不唯一。 状态空间方程可表示成 u x x B A += （状态方程） （2.63） u x y D C += （输出方程） (2.64) 式中， []τ=n x 2x 1x x n 维状态矢量； ??????????=nn 2n 1n n 22221n 11211a a a a a a a a a A n ×n 维系统状态系数矩阵； []τ=r 21u u u u r 维控制矢量； ??????????=nr 2n 1n r 22221r 11211b b b b b b b b b B n ×r 维系统控制系数矩阵； []τ=m 21y y y y m 维输出矢量； ??????????=mn 2m 1 m n 22221n 11211c c c c c c c c c C m ×n 维输出状态系数矩阵； ??????????=mr 2m 1 m r 22221r 11211d d d d d d d d d D m ×r 维输出控制系数矩阵； 系统信号传递方块图如图2-46所示。 D A B C u y + + + + x x ?
第二章状态空间分析法 2－1 状态、状态变量、状态空间、状态方程、动态方程 任何一个系统在特定时刻都有一个特定的状态，每个状态都可以用最小的一组(一个或多个)独立的状态变量来描述。 设系统有n个状态变量x1，x2，…，x n，它们都是时间t的函数，控制系统的每一个状态都可以在一个由x1，x2，…，x n为轴的n维状态空间上的一点来表示，用向量形式表示就是: X = (x 1，x 2 ，…，x n )T X称作系统的状态向量。 设系统的控制输入为：u1，u2，...，u r，它们也是时间t的函数。记： U = (u 1，u 2 ，...，u r )T 那么表示系统状态变量X(t)随系统输入U(t)以及时间t变化的规律的方程就是控制系统的状态方程，如式（2－1）所示。 ………………………………………………………………（2－1） 其中F = (f1，f2，...，f n)T是一个函数矢量。
设系统的输出变量为y1，y2，...，y m，则Y = (y1，y2，...，y m)T 称为系统的输出向量。表示输出变量Y(t)与系统状态变量X(t)、系统输入U(t)以及时间t的关系的方程就称作系统的输出方程，如式2－2所示。 …………………………………………………………. (2－2) 其中G = (g1，g2，...，g m)T是一个函数矢量。 在现代控制理论中，用系统的状态方程和输出方程来描述系统的动态行为，状态方程和输出方程合起来称作系统的状态空间表达式或称动态方程。 根据函数向量F和G的不同情况，一般控制系统可以分为如下四种：?线性定常(时不变)系统（LTI-Linear Time Invariant）； ?线性不定常(时变)系统； ?非线性定常系统； ?非线性时变系统。 在本课程中，我们主要考虑线性定常系统(LTI)。这时，系统的动态方程可以表示如下： …………….（2－3）
第三章 线性系统的运动分析 §3-1线性连续定常齐次方程求解 一、齐次方程和状态转移矩阵的定义 １、齐次方程 状态方程的齐次方程部分反映系统自由运动的状况（即没有输入作用的状况），设系统的状态方程的齐次部分为：)()(t Ax t x = 线性定常连续系统：Ax x = ２、状态转移矩阵的定义 齐次状态方程Ax x = 有两种常见解法：（1）幂级数法；（2）拉氏变换法。其解为)0()(x e t x At ?=。 其中At e 称为状态转移矩阵（或矩阵指数函数、矩阵指数），记为：At e t =)(φ。 若初始条件为)(0t x ，则状态转移矩阵记为：)(00 )(t t A e t t -=-Φ 对于线性时变系统，状态转移矩阵写为),(0t t φ，它是时刻t ，t 0的函数。但它一般不能写成指数形式。 （１）幂级数法 设Ax x = 的解是t 的向量幂级数 +++++=k k t b t b t b b t x 2210)( 式中 ，，， ，，k b b b b 210都是n 维向量，则 +++++=-1232132)(k k t kb t b t b b t x )(2210 +++++=k k t b t b t b b A 故而有： ????? ?? ????== ====003 230 2 12 01!1! 3131 2 121b A k b b A Ab b b A Ab b Ab b K K
且有0)0(b x =。 故 +++++=k k t b t b t b b t x 2210)( ++ +++=k k t b A k t b A t Ab b 02 02 00! 1! 21 )0()! 1!21(22 x t A k t A At I k k ++ ++ += 定义：∑ ∞ == ++ +++=0 2 2! 1! 1!21K k k k k At t A k t A k t A At I e 则)0()(x e t x At ?=。 （２）拉氏变换解法 将Ax x = 两端取拉氏变换，有 )()0()(s Ax x s sx =- )0()()(x s x A sI =- )0()()(1x A sI s x ?-=- 拉氏反变换，有 )0(])[()(11x A sI L t x ?-=-- 则 ])[()(11---==A sI L e t At φ 【例3.1.1】 已知系统的状态方程为x x ?? ? ???=00 10 ，初始条件为)0(x ，试求状态转移矩阵和状态方程的解。 解：（１）求状态转移矩阵 ++ ++ +==k k At t A k t A At I e t ! 1! 21)(2 2φ 此题中： ???? ??=00 10A ， ?? ? ???====00 0032n A A A 所以
第2章 “控制系统的状态空间描述”习题解答 系统的结构如图所示。以图中所标记的1x 、2x 、3x 作为状态变量，推导其状态空间表达式。其中，u 、y 分别为系统的输入、输出，1α、2α、3α均为标量。 3 x 2 x 图系统结构图 解 图给出了由积分器、放大器及加法器所描述的系统结构图，且图中每个积分器的输出即为状态变量，这种图形称为系统状态变量图。状态变量图即描述了系统状态变量之间的关系，又说明了状态变量的物理意义。由状态变量图可直接求得系统的状态空间表达式。 着眼于求和点①、②、③，则有 ①：2111x x x +=α& ②： 3222x x x +=α&③：u x x +=333α& 输出y 为1y x du =+，得 1112223331000100 1x a x x a x u x a x ?? ?????? ????????=+???????????????????????? &&& []123100x y x du x ?? ??=+?? ???? 已知系统的微分方程 (1) u y y y y 354=+++&&&&&& ；(2) u u y y -=+&&&&&&32； (3) u u y y y y 75532+=+++&&&&&&&&& 。试列写出它们的状态空间表达式。 (1) 解 选择状态变量1y x =，2y x =&，3y x =&&，则有：
1223 31231 543x x x x x x x x u y x =??=?? =---+??=?&&& 状态空间表达式为：[]112233123010000105413100x x x x u x x x y x x ????????????????=+????????????????---???????? ????=?????? &&& (2) 解 采用拉氏变换法求取状态空间表达式。对微分方程(2)在零初试条件 下取拉氏变换得： 3222332()3()()() 11()12 23()232 s Y s sY s s U s U s s Y s s U s s s s s +=---==++ 由公式、可直接求得系统状态空间表达式为 1122330100001031002x x x x u x x ?? ????????????????=+? ?????????????????????-?? ?? &&& 123110 2 2x y x x ?????? =- ?????????? (3) 解 采用拉氏变换法求取状态空间表达式。对微分方程(3)在零初试条件 下取拉氏变换得： 323()2()3()5()5()7()s Y s s Y s sY s Y s s U s U s +++=+
实验报告 实验名称利用 MATLAB 求解系统的状态方程 系统的能控性、能观测性分析 系专业班 学号授课老师 预定时间2014-5-28 实验时间实验台号14 一、目的要求 掌握状态转移矩阵的概念。学会用 MATLAB求解状态转移矩阵。 掌握求解系统状态方程的方法，学会绘制状态响应曲线； 掌握线性系统状态方程解的结构。学会用 MATLAB 求解线性定常系统的状态响应和输出响应，并绘制相应曲线。 掌握能控性和能观测性的概念。学会用 MATLAB 判断能控性和能观测性。 掌握系统的结构分解。学会用 MATLAB 进行结构分解。 掌握最小实现的概念。学会用 MATLAB 求最小实现。 二、原理简述 线性定常连续系统的状态转移矩阵为。 函数 step( ) 可直接求取线性连续系统的单位阶跃响应。 函数 impulse( ) 可直接求取线性系统的单位脉冲响应。 函数 lsim( ) 可直接求取线性系统在任意输入信号作用下的响应。 函数 initial( ) 可求解系统的零输入响应。 n 阶线性定常连续或离散系统状态完全能控的充分必要条件是：能控性矩阵
的秩为 n。 线性定常连续或离散系统输出能控的充分必要条件是：矩阵 的秩为m。 n 阶线性定常连续或离散系统状态完全能观测的充分必要条件是：能观测性矩阵 的秩为 n。 三、仪器设备 PC 计算机，MATLAB 软件 四、容步骤 题2.1 A=[0 1;-2 -3];B=[3;0];C=[1 1];D=0; G=ss(A,B,C,D); t=0.5; p=expm(A*t) u1=0;x10=[1;-1]; [y1o,t,x1o]=initial(G,x10,t) t2=0:0.5:10;x20=[0;0];u2=ones(size(t2)); [y2,t2,x2]=lsim(G,u2,t2); plot(t2,x2,':',t2,y2,'-') t3=0:0.5:10;u3=1+exp(-t3).*cos(3*t3);x30=[0;0];
第2章 状态空间表达式的解 第1节 线性定常齐次状态方程的解 线性定常齐次状态方程 0(0)x Ax x x == 的解为 0()At x t e x = (0) t > 式中，2 2 () 2!!k At k t At e I At A k ∞ ? ==+++=∑ 证明： 用拉普拉斯变换法。 对 x A x = 作拉氏变换，得
0()()s X s x AX s -= 1 0()()X s sI A x -=- 1 1 0()[()]x t L sI A x --=- 因为 2 23111()()sI A I A A I s s s -+++= 故 1 223111()sI A I A A s s s --=+++ 12023111()[]x t L I A A x s s s -=+++ 2201()2! I At A t x =+++ 0At e x = 顺便可知 ])[(1 1---=A sI L e At 第2节 矩阵指数函数At e
1、At e 的定义和性质 （1）定义 2 2 () 2!!k At k t At e I At A k ∞ ==+++=∑ 式中 A —线性定常系统系统矩阵，n n ?阶； At e —矩阵指数函数，n n ?阶时变矩阵。 若A 中各元素均小于某定值，At e 必收敛；若A 为实矩阵，At e 绝对收敛。 （2）基本性质： ◆组合性质： ) (212 1t t A At At e e e += 其中21,t t 为相衔接的两时间段。 推论1：I e e e e A t t A t A At ===--0 )() (
·258· 第9章 线性系统的状态空间分析与综合 重点与难点 一、基本概念 1．线性系统的状态空间描述 （1）状态空间概念 状态 反映系统运动状况，并可用以确定系统未来行为的信息集合。 状态变量 确定系统状态的一组独立（数目最少）变量，它对于确定系统的运动状态是必需的，也是充分的。 状态向量 以状态变量为元素构成的向量。 状态空间 以状态变量为坐标所张成的空间。系统某时刻的状态可用状态空间上的点来表示。 状态方程 状态变量的一阶导数与状态变量、输入变量之间的数学关系，一般是关于系统的一阶微分（或差分）方程组。 输出方程 输出变量与状态变量、输入变量之间的数学关系。 状态方程与输出方程合称为状态空间描述或状态空间表达式。线性定常系统状态空间表达式一般用矩阵形式表示： ???+=+=Du Cx y Bu Ax x (9.1) （2）状态空间表达式的建立。系统状态空间表达式可以由系统微分方程、结构图、传递函数等其他形式的数学模型导出。 （3）状态空间表达式的线性变换及规范化。描述某一系统的状态变量个数（维数）是确定的，但状态变量的选择并不唯一。某一状态向量经任意满秩线性变换后，仍可作为状态向量来描述系统。状态变量选择不同，状态空间表达式形式也不一样。利用线性变换的目的在于使系统矩阵A 规范化，以便于揭示系统特性，利于分析计算。满秩线性变换不改变系统的固有特性。 根据矩阵A 的特征根及相应的独立特征向量情况，可将矩阵A 化为三种规范形式：对角形、约当形和模式矩阵。 （4）线性定常系统状态方程解。状态转移矩阵)(t φ（即矩阵指数At e ）及其性质：
·259· (9.8) i ． I =)0(φ ii ．A t t A t )()()(φφφ == iii. )()()()()(122121t t t t t t φφφφφ±=±=+ iv. )()(1t t -=-φφ v. )()]([kt t k φφ= vi. )( ])exp[()exp()exp(BA AB t B A Bt At =+= vii. )( )exp()exp(11非奇异P P At P APt P --= 求状态转移矩阵)(t φ的常用方法： 拉氏变换法 =)(t φL －1])[(1--A sI （9.2） 级数展开法 ++++ +=k k At t A k t A At I e ! 12122 （9.3） 齐次状态方程求解 )0()()(x t t x φ= （9.4） 非齐次状态方程式（9.1）求解 ?-+=t Bu t x t t x 0d )()()0()()(τττφφ （9.5） （5）传递函数矩阵及其实现 传递函数矩阵)(s G ：输出向量拉氏变换式与输入向量拉氏变换式之间的传递关系 D B A sI C s G +-=-1)()( (9.6) 传递函数矩阵的实现：已知传递函数矩阵)(s G ，找一个系统},,,{D C B A 使式（9.6）成立，则将系统},,,{D C B A 称为)(s G 的一个实现。当系统阶数等于传递函数矩阵阶数时，称该系统为)(s G 的最小实现。 传递函数矩阵的实现并不唯一。实现的常用标准形式有可控标准形实现、可观测标准形实现、对角形实现和约当形实现等。 （6）线性定常连续系统的离散化及其求解 对式（9.1）表示的线性定常数连续系统进行离散化，导出的系统离散状态空间描述为 ???+=+=+ )()()( )()()()()1(k D k Cx k y k u T G k x T k x φ