发明名称 一种基于ARMAX建模的培养基蒸汽喷射灭菌节能控制方法
摘要 本发明公开了一种基于ARMAX建模的培养基蒸汽喷射灭菌节能控制方法,该方法以生物发酵用培养基蒸汽直接加热混合灭菌控制过程为对象,利用生物制药用培养基蒸汽灭菌过程数据建立外部输入自回归滑动(ARMAX)模型,在该模型基础上推导灭菌过程涉及的过程参数对应的传递函数,建立动态矩阵并实现动态矩阵控制(DMC)。为节约控制过程中蒸汽使用量,在达到精确控制的目标的同时,本发明中对控制器的设计上加入了节能指标,使控制过程中的蒸汽使用效率更高,降低生产过程中的蒸汽使用成本;本发明控制方法简单,控制策略和控制器参数设计基于实际工业运行过程数据,对灭菌过程涉及到的工艺环节和设备分析的依赖性小,可实现性强,易于应用。
申请公布号 CN104898416A 申请公布日期 2015.09.09
申请号 CN201510158900.3 申请日期 2015.04.03
申请人 西安交通大学 发明人 曹晖;张士良;张彦斌;贾立新;司刚全
分类号 G05B13/04(2006.01)I 主分类号 G05B13/04(2006.01)I
代理机构 西安智大知识产权代理事务所 61215 代理人 何会侠
主权项 一种基于ARMAX建模的培养基蒸汽喷射灭菌节能控制方法,其特征在于:包括如下步骤:步骤1:利用培养基蒸汽灭菌过程数据,建立灭菌过程ARMAX模型首先收集工业运行现场培养基蒸汽灭菌过程数据,通过对过程数据进行分析,找到灭菌过程所涉及到的过程参数,并确定该控制系统的输入、输出及扰动量,由此确定灭菌控制系统的模型结构;经分析,系统输入为用于加热培养基的蒸汽的单位流量,在实际运行中表现为喷射器蒸汽入口管路的阀门开度,系统输出为喷射器出口温度,即蒸汽和培养基混合后的温度,系统扰动量有:培养基在喷射器入口处温度,培养基流动速度,喷射器入口处蒸汽温度;采用建模的模型结构为ARMAX模型,结合灭菌过程的参数及个数,其形式如下:A(z)y(t)=B(z)u(t)+D<sub>1</sub>(z)e<sub>1</sub>(t)+D<sub>2</sub>(z)e<sub>2</sub>(t)+D<sub>3</sub>(z)e<sub>3</sub>(t) A(z)=a<sub>1</sub>+a<sub>2</sub>z<sup>‑1</sup>+…a<sub>na</sub>z<sup>‑na+1</sup>B(z)=b<sub>1</sub>+b<sub>2</sub>z<sup>‑1</sup>+…b<sub>nb</sub>z<sup>‑nb+1</sup>D<sub>i</sub>(z)=d<sub>i1</sub>+d<sub>i2</sub>z<sup>‑1</sup>+…d<sub>indi</sub>z<sup>‑ndi+1</sup>其中,y(t)为系统输出,u(t)为系统输入量,A(z)为y(t)对应的延迟因子所构成的多项式,B(z)为控制量u(t)对应的延迟因子所构成的多项式,D<sub>1</sub>(z)为第1个扰动量e<sub>1</sub>(t)对应的延迟因子构成的多项式,D<sub>2</sub>(z)为第2个扰动量e<sub>2</sub>(t)对应的延迟因子构成的多项式,D<sub>3</sub>(z)为第3个扰动量e<sub>3</sub>(t)对应的延迟因子构成的多项式,D<sub>i</sub>(z)为第i个扰动量对应的延迟因子所构成的多项式,a<sub>1</sub>,a<sub>2</sub>,……a<sub>na</sub>,b<sub>1</sub>,b<sub>2</sub>,……b<sub>nb</sub>,d<sub>i1</sub>,d<sub>i2</sub>,……d<sub>indi</sub>为待辨识系数;e<sub>1</sub>(t),e<sub>2</sub>(t),e<sub>3</sub>(t)分别为t时刻培养基在喷射器入口处温度,培养基流动速度和喷射器入口处蒸汽温度;步骤2:离散传递函数的推导和控制矩阵的生成利用工业运行数据辨识ARMAX模型中的参数,获得参数模型后直接推导对应过程参数的离散传递函数:<img file="FDA0000693999050000011.GIF" wi="1442" he="163" />故控制量和三个扰动量e<sub>1</sub>,e<sub>2</sub>,e<sub>3</sub>对应的离散传递函数F<sub>u</sub>,F<sub>1</sub>,F<sub>2</sub>,F<sub>3</sub>依次 为:<img file="FDA0000693999050000021.GIF" wi="267" he="151" /><img file="FDA0000693999050000022.GIF" wi="293" he="139" /><img file="FDA0000693999050000027.GIF" wi="356" he="143" /><img file="FDA0000693999050000028.GIF" wi="328" he="155" />为计算各参数对应的动态矩阵,先求各参数在对应离散传递函数下的阶跃响应,其形式如下:对于控制量u,采样时刻1,2,……,n<sub>u</sub>对应的系统输出分别为:g<sub>1</sub>,g<sub>2</sub>,……,<img file="FDA00006939990500000211.GIF" wi="97" he="74" />对于扰动量e<sub>1</sub>,采样时刻1,2,……,<img file="FDA00006939990500000213.GIF" wi="91" he="69" />对应的系统输出分别为:f1<sub>1</sub>,f1<sub>2</sub>,……,<img file="FDA00006939990500000210.GIF" wi="132" he="88" />对于扰动量e<sub>2</sub>,采样时刻1,2,……,<img file="FDA00006939990500000214.GIF" wi="90" he="69" />对应的系统输出分别为:f2<sub>1</sub>,f2<sub>2</sub>,……,<img file="FDA0000693999050000029.GIF" wi="131" he="89" />对于扰动量e<sub>3</sub>,采样时刻1,2,……,<img file="FDA00006939990500000215.GIF" wi="90" he="69" />对应的系统输出分别为:f3<sub>1</sub>,f3<sub>2</sub>,……,<img file="FDA00006939990500000212.GIF" wi="131" he="89" />由此,控制量u,扰动量e<sub>1</sub>,e<sub>2</sub>,e<sub>3</sub>对应的动态矩阵G<sub>u</sub>,G<sub>1</sub>,G<sub>2</sub>,G<sub>3</sub>分别为:<img file="FDA0000693999050000025.GIF" wi="760" he="398" /><img file="FDA0000693999050000026.GIF" wi="814" he="415" /><img file="FDA0000693999050000031.GIF" wi="815" he="413" /><img file="FDA0000693999050000032.GIF" wi="820" he="414" />其中,m为控制步长,N为预测步长;步骤3:控制量的计算基于步骤2所得的各参数动态矩阵,系统预测输出值y′表达为:<img file="FDA0000693999050000033.GIF" wi="818" he="153" />上式中,u为控制量按时间序列构成的向量,E<sub>i</sub>为扰动量e<sub>i</sub>的增量按时间序列构成的向量,G<sub>i</sub>为扰动量e<sub>i</sub>对应的动态矩阵;f<sub>u</sub>为控制量中影响预测值的只与过去时刻数值有关的部分所构成的列向量,f<sub>u</sub>=[f<sub>u</sub>(t+1),f<sub>u</sub>(t+2)……f<sub>u</sub>(t+N)]<sup>T</sup>,<img file="FDA0000693999050000037.GIF" wi="73" he="81" />为扰动e<sub>i</sub>中影响预测值的只与过去时刻数值有关的部分所构成的列向量,<img file="FDA0000693999050000034.GIF" wi="1084" he="82" />其中:<img file="FDA0000693999050000035.GIF" wi="1108" he="144" /><img file="FDA0000693999050000036.GIF" wi="1135" he="160" />其中f<sub>u</sub>(t+k)为f<sub>u</sub>中的第k个元素,<img file="FDA0000693999050000038.GIF" wi="247" he="84" />为<img file="FDA0000693999050000039.GIF" wi="72" he="81" />中的第k个元素,N为预测步长,Δu(t‑i)为t‑i时刻的控制量增量,Δe<sub>i</sub>(t‑j)为扰动量e<sub>i</sub>在t‑j时刻的增量,g<sub>i</sub>为控制量在i时刻的阶跃响应值,fi<sub>j</sub>为扰动量e<sub>i</sub>在j时刻的阶跃响应值;为简便起见,将系统预测输出值y′记为:y′=G<sub>u</sub>u+f其中<img file="FDA00006939990500000310.GIF" wi="608" he="90" />系统建模和动态矩阵获得后,即可确定目标函数进而计算控制量;常规目标函数J形式为:<img file="FDA0000693999050000041.GIF" wi="1337" he="225" />其中λ为平衡目标函数两部分即预测误差和控制耗能权重的系数,N为预测步长,m为控制步长,Δu(t+j‑1)为t+j‑1时刻控制量的变化量,y′(t+j|t)为在t时刻预测的t+j时刻的系统输出值,w(t+j)为t+j时刻的输出参考轨迹,设置为如下形式:w(t+j)=αw(t+j‑1)+(1‑α)r(t+j)w(t)=y(t)其中r(t+j)为系统输出在t+j时刻的设定值,y(t)为t时刻系统输出量,α为0到1之间的常数;求目标函数的最小值,即获得未来控制步长内最优的控制量;但考虑到蒸汽使用成本,为节能降耗起见,目标函数中应当加入对蒸汽使用量的目标约束;蒸汽使用量Q表示为如下形式:Q=T*fs式中,T为时间间隔,fs为蒸汽流速,而fs正比与阀门开度,则控制步长范围内的各个扫描间隔内的蒸汽使用量Q(t+1),Q(t+2),……Q(t+m)为:Q(t+1)=K*[u(t)+Δu(t+1)]*ΔTQ(t+2)=K*[u(t)+Δu(t+1)+Δu(t+2)]*ΔT···Q(t+m)=K*[u(t)+Δu(t+1)+Δu(t+2)……Δu(t+m)]*ΔT式中,K为常系数,ΔT为控制间隔时间,即扫描时间间隔,u(t)为当前时刻控制量输入的数值,m为控制步长;在控制步长内,蒸汽总耗量E为:<img file="FDA0000693999050000042.GIF" wi="403" he="206" />上式表达为如下形式:E=K×m×u(t)+K×M<sup>T</sup>u其中m为控制步长,M为从蒸汽总消耗量的计算过程中推导出的矩阵,u为未来时刻控制量u(t+1),u(T+2),……u(t+m)按时间序列构成的向量,T为转置 符号:u=[u(t+1),u(t+2),……u(t+m)]<sup>T</sup>M=[m,m‑1,……1]<sup>T</sup>如此,加入了对节能指标的考虑后,目标函数J确定为:<img file="FDA0000693999050000051.GIF" wi="1443" he="220" />求目标函数的最小值,即得到未来时刻控制量输入;将目标函数转换成矩阵运算形式后求导,即可得到:<img file="FDA0000693999050000052.GIF" wi="985" he="142" />其中<img file="FDA0000693999050000053.GIF" wi="426" he="93" />w为输出轨迹按时间顺序构成的序列,T为转置符号,I为单位对角矩阵,上标‑1为矩阵逆运算符号;以上即为动态矩阵控制方法得到的未来时刻控制量输入,按照滚动优化的方式,每次计算u都只需要计算第一行第一列的数值,然后送至系统执行器即可。
地址 710049 陕西省西安市咸宁路28号