反应器的优化和稳定性

反应器的优化和稳定性1第1页，课件共68页，创作于2023年2月第五章反应器的优化和稳定性A.优化目的：过程利益最大或总费用最小问题分类：静态优化——求稳态操作目标函数的极值动态优化——求非稳态操作目标函数的极值B.稳定性新平衡（或自动恢复）——稳定条件变化——不能建立平衡——不稳定C.控制性：用什么控制系统，多大程度上将装置控制在规定操作条件下——属自动控制问题平衡2第2页，课件共68页，创作于2023年2月5.1目标函数和约束条件5.1.1目标函数对反应器：主要涉及产量问题和效率问题。产量问题：单位时间、单位容积反应器产量最大——“单产”或转化为完成一定量产品反应时间最短——“速度”此类问题对设备费用较大时很重要。效率问题：一定量的原料，获取最多产品，原料费用高时很重要。变量包括三种状态变量如器内温度、浓度、转化率等操作变量流量、原料、冷却水温坐标变量时间、空间目标函数是包含上述变量（部分或全部）的函数。可以转化，如T亦可作操作变量3第3页，课件共68页，创作于2023年2月5.1.2约束条件反应器所特有的及反应过程存在的特定关系式（与上述三种变量间）称为约束条件。分为两种：等式约束：过程方程式（物料平衡、速度式、k式等）等不等式约束：温度极限等、不等式。5.2优化方法优化——根据约束条件求目标函数之极值。

方法有多种：视目标函数与约束条件的具体形式而定。5.2.1微分法操作变量x1

,x2

,x3

,

…,xm

状态变量y1

,y2

,y3

,…,yn目标函数P=P（x1

,x2

,x3

,

…,xm,y1

,y2

,y3

,…,yn）等式约束Fi=Fi（x1

,x2

,x3

,

…,xm,y1

,y2

,y3

,…,yn）=0

4第4页，课件共68页，创作于2023年2月将等式约束写成yi=i（x1

,x2

,x3

,

…

,xm）（写成操作变量的函数）

于是目标函P=P（x1

,x2

,x3

,

…

,xm）（消去yi）若P

和P/xk是连续的，则在满足的xh点存在极值。极值为极大还是极小，由二阶导数符号判别。如:对函数P=P（x1，x2）其极值存在于从而可解得（x10，x20）点若则5第5页，课件共68页，创作于2023年2月下面举例说明应用例1.讨论2级可逆反应的最佳温度。（该问题属于产量问题）(5.11)(5.12)可见，在CA0一定条件下，t由fA及T给出6第6页，课件共68页，创作于2023年2月间歇反应器，t

为操作时间生产速度对连续管式反应器，θ为反应时间生产速度对等温操作，约束条件已写成操作变量的显函数（5.12），代入（5.11）式，然后按fA一定，求T偏导。WmA——为A的一次装料量；V——反应器容积；t，t′——分别为反应时间（操作时间）和装卸料时间。生产量问题化为生产速度问题，即在给定fA条件下使t（或θ）为最小，即相当于生产速度为最大。因而t[(5.11)式]就成为在约束条件(5.12)式下求极值的问题不等式约束0<

fA≤fAe(平衡时转化率)目标函数7第7页，课件共68页，创作于2023年2月可解得TOPT，因很繁，仅写出（5.13）这只是个关系式，尚未给出Topt，须要与原函数联立求解，其解也很麻烦。变通的方法是按照（5.11）式，fA一定时，变T，求t，作图讨论：（5.13）式分母>0

t

tmin

ToptT因因而只有E2＞E1，才存在tmin说明：可逆反应条件0<

fA≤fAe即有rA

≥rC8第8页，课件共68页，创作于2023年2月亦即(开始无其它物料C,D)或即有另外，需说明一点，对∂t/∂T=0所得繁琐方程，可由计算机直接试算迭代求出Topt，及相当的tmin

值对非等温操作问题：随反应进行改变温度，以保持最大反应速度。目标函数对T求偏导9第9页，课件共68页，创作于2023年2月令其为“0”，得（5.15）由（5.15）式，根据CA、CB、CC、CD可求得Topt

在平衡状态时（5.16）Te——平衡温度。由（5.15）、（5.16）式（5.17）（5.16）式有两种用法①已知温度时→平衡浓度②已知浓度时→平衡温度现属于后者k1、k2代入上式10第10页，课件共68页，创作于2023年2月微分法要求从目标函数中消去状态变量yi当yi不能写成操作变量xi的显函数时，不能用微分法，此时，可用拉格朗日乘数法设新的目标函数为：5.2.2拉格朗日乘数法为约束条件根据加上约束条件共2n+m个方程式中，Te可按各时刻CA、CB、CC、CD代入求出，以求Topt

各浓度由推求，可见温度随时间变化，实现困难。其中：11第11页，课件共68页，创作于2023年2月未知量(2n+m)个xh

m个yi

n个λin个2n+m个理论上可求解，但实际上相当麻烦，只有简单问题可得解析解。例2.连续流搅拌槽中的连串反应讨论中间产物R的效率。解：效率等式约束条件由A物料平衡或（5.22）12第12页，课件共68页，创作于2023年2月由R物料平衡或（5.23）新目标函数操作变量：T（含于k1，k2，k3中）及平均停留时间θ状态变量：浓度CA、CR13第13页，课件共68页，创作于2023年2月问题分为：（1）温度T

一定，求最佳θ加上两个约束方程，共5个方程，求解

θ，CA，CR，λ1，λ25个变量，得14第14页，课件共68页，创作于2023年2月

（2）θ一定，求最优温度Topt方法同上（略）此例简单，用微分法可直接求解。例如：应用两个约束方程，可直接得出仅含操作变量的显式，目标函数（5.28）式，对该式直接微分即可求得θopt和Topt若目标函数写成操作变量的函数P=P(x1，x2，，xm)，那么P与xh构成m+1维空间。形象些说，当m=2时，构成三维空间，若以图形表示函数，就为曲面，曲面的最高点（或最低点）即极值点所谓“登山”即从曲面上任意点（初始解）出发向峰点逼近的过程。5.2.3登山法15第15页，课件共68页，创作于2023年2月例如：求y最小值(y显见不能为负值)。设从起点（0,－1）出发，代入函数式得y＝14.4下一步走到（1,1）点，y＝7.11向极值逼近！一步步试算总可找到极值点，但这不科学。应合理确定搜索矢量和步长，才能系统地逐步逼近极值。（1）搜索矢量与步长搜索矢量表示从出发点前进的矢量，步长则表示求前进步幅大小。若以矢量表示出发点、到达点，则与搜索矢量的关系可表示为（5.35）取单位步长，p=1,上例表示为登山方向矢量16第16页，课件共68页，创作于2023年2月即搜索矢量起点、终点Dm必须满足（5.35）式，或由该式协调，然而，都是未知的，需要有一种方法解决这一问题。（2）快速登山法或并且xh

在方向移动距离为

ε

为任意常数，相当于步长符号“”，求极大用“＋”，求极小用“－”。17第17页，课件共68页，创作于2023年2月例3讨论日产10吨化工厂生产费用最小值图e.5.3.a流程图取反应流体压力为x1（kg/cm2）循环比x2（—）各部分费用：反应器分离机循环泵混合器和压缩机目标函数解：以则目标函数可简化为18第18页，课件共68页，创作于2023年2月出发点任选，设x1

=800，x2

=6，P0

=3133333在出发点从而得新的试探点（第一步终点）(x1)1=883，(x2)1=5.78，P1=3109460一般经验，ε使xi变化约10％19第19页，课件共68页，创作于2023年2月以后逐步计算，直到P/xi=0为止，即为极值点。（1000，4）点Pmin=3×106当然，本例变量少，可用微分法直接求解可解得注意：如果用解析法求不出P/xi值，可用xi＋数值法估算（类似定义）得搜索矢量20第20页，课件共68页，创作于2023年2月步长很重要，P曲面复杂时要选小些，曲面简单时可大些；当存在两个以上极值时，一旦到达其一就不能前进，为此应选择几个出发点；到达是极值点还是鞍点，落入鞍点就无法解脱；变步长方法：一次成功，下次取3倍步长；一次失败，下次取1/2步长。为什么搜索失败？如从1点出发，因步长太大，超越极值点2点P值小于1点12起止开始时，步长小些为好，视P增加情况而定。0－1步长太大0－1P增长不多要判别平坦或越过15432021第21页，课件共68页，创作于2023年2月5.2.4线性规划（LinearProgramming）当P为线性函数，且约束（等式、不等式）也为线性时的优化问题就成为线性规划问题。

“实际可行面积”或“可行域”—由线性约束条件所围成的平面面积，优化即围绕此面积周边进行计算,使函数P最大。例4.计算函数最大值约束x≥0,y≥0，y≤x+2,y≤4－x

解：约束条件围成右图面积由K=y+0.5x写出等K值方程y=－0.5x＋K

K＝K0＝0时，y=－0.5x

K＝K1＝1时，y=－0.5x＋1

K＝K2

＝2时，y=－0.5x＋2……12341234xyy=x+2(1,3)K3,5K3K2K1K0线性规划图解“可行域”22第22页，课件共68页，创作于2023年2月不同K值构成一个平行线族，线愈高，则K值愈大，但K必须位于“可行面积”（可行域）内最后得K3.5

是最大值此时y=3x=1由此可见，因K等值线一定要在可行面积内的才有用，故平行移动K线，最后必切（交）于可行面积内的某个顶点（或边），这便是所求极值。第一定理：线性规划目标函数的最佳值必定位于由约束条件规定的凸多边形的一个顶点上。以上是按2个变量讨论的。当变量为N个时，最优值位于N维空间多边形（凸面体）的一个顶点上，不过这种情况需用矩阵代数求解。N维问题：目标函数（5.41）不等式约束（5.42）23第23页，课件共68页，创作于2023年2月引入松弛变量Zi，把不等式化为等式为了符号上取得一致，xi与Zi统改为fj表示（5.43）xi1≤j≤mZj-mm+1≤j≤m

+l=Nfj=由此可改写（5.43）式，得（5.44）24第24页，课件共68页，创作于2023年2月写为矩阵形式或写成向量表示：目标函数写成：（5.45）（5.46）25第25页，课件共68页，创作于2023年2月由于fm+1，…，fm+l不为“0”必须Cm+1，…，Cm+l为“0”，以保证原条件不变上式中，设则，则目标函数可写成：从而，线性规划问题概括为约束（5.48）未知量f1，f2，…，fm+l共m+l个方程l个（不计P=cTf）可见，有无穷多组解，但其中使P达最大者为最优解（5.47）26第26页，课件共68页，创作于2023年2月Dantzig

程序

（计算方法）（运筹学中，迭代或旋转运算）

因方程只有l

个，变量为m+l

个，故先假定l个变量然后求解。（1）（设f1

=f2

=f3

=…=fm=0）求初始可行解根据（5.45）式，f1

…

fm为“0”，又am+1

…

am+l为1，构成单位向量这一点从（5.44）式显而易见即am+1,

am+2

…am+l构成单位系数矩阵上解代入目标函数P=0因f1

…fm都是零，故称从原点出发，为其一解。27第27页，课件共68页，创作于2023年2月解释：

（i）主成分（基变量）由上可见，从原点出发可理解为选定l维坐标系，使符合初始条件。所以，am+1

,

am+2

,

…,am+l应满足条件即am+i＝1都是单位矢量

am+1,am+2,…,am+l称坐标的主成分（基）（ii）剩余成分（非基变量）其余矢量a1，a2

…am都分解到坐标轴上，每个有l个分量，用as表示（s=1,2,

…,m）28第28页，课件共68页，创作于2023年2月（5.52）式中βm+1,s可理解为as在轴m+1轴上的分量，am＋1

为m＋1轴的单位矢量等等。上式可改写为（5.53）实际上即为剩余部分，用单位向量线性表示。（5.54）而ais即为方程组（5.44）中的各系数。（2）从出发点前进，使P增长方法：变换坐标，保持l维，用剩余成分中的一个（一个非基变量）代替主成分中的某一个（某一个基变量

）。非基向量在基向量上的投影。29第29页，课件共68页，创作于2023年2月举例说明：引入松弛变量f3,f4

后约束方程变为：（a）（b）问题有无穷多解设f1＝f2＝0则f3＝8

f4＝6将方程系数列成表（进行基变换）a1a2a3a4b基变量f3a324108f4a442016由方程（a），（b）初始可行解30第30页，课件共68页，创作于2023年2月（a）-（b）/2得（c）（b）/4得（d）此时基变量就成为a3，a1，方程性质未变，列表a1a2a3a4b基变量f3a3031-1/25f1a111/201/43/2此时：f4

称为换出变量；f1

称为换入变量。取f2=f

4=02次可行解P增大然后，再次进行基变换，以f3为换出变量，f2

为换入变量，得

f1=2/3f2=5/331第31页，课件共68页，创作于2023年2月3次可行解此时函数达其极值对于l维坐标系，其基向量由l个单位向量所组成(1)令非基变量为“0”，即可得一初始可行解(2)由此出发，进行基变换，求得可行解和P值(3)每变换一次得一个P，可看出它的增减但是：如何进行变量变换（基变换）？用哪个非基变量as替换基变量am+i才能保证P值增长？(只有满足这一要求，才有实际意义或变换才有意义)。变换方法由符号运算（5.56）32第32页，课件共68页，创作于2023年2月适当选择θ

值，如

这样am+i项被消去，而θas留下，即将as引入基变量（换入变量），am+i成为非基变量（换出变量）。由（5.56）式看出，消去am+i和引入as是任意的，但必须保证变换后的结果能使所有的fm+i都是正值（题给条件），因此，规定θ必须是正值，且i必须使θ值为最小者。这要从基变量变换后P值的变化来看变换前（实质aik最大者）基变量变换后令称为补偿函数由此得33第33页，课件共68页，创作于2023年2月下面举例说明（见p155~p160）例5求函数最大值解：①引入松弛变量（x4，x5，x6）将不等式约束化为等式②求初始可行解，并判断换入基变量和换出基变量可见，只要Cs－Zs>0

P

便增长。（第1表）34第34页，课件共68页，创作于2023年2月a1a2a3a4a5a6bCa41－1010010a542－3010100a651－3001140CS312000ZS0000000CS－

ZS312000（第1表）第1列（f1）Cs–Zs最大，f1为换入基变量第1列θ值θ4值为最小，确定f4

为换出变量因而基变量成为f1,f5,f6

（a1,a5,a6）单位矢量35第35页，课件共68页，创作于2023年2月③第2可行解（将f1

系数向量a1化为单位向量）判别换入换出变量初始可行解为（第2表）第2列Cs–Zs最大，f2为换入基变量a1a2a3a4a5a6bca11－1010013a506－3－41060a6063－50190CS312000ZS3－303003CS－

ZS042－300f5应为换出变量（θ＞0)(原2行—4×1行）(原3行—5×1行）36第36页，课件共68页，创作于2023年2月④求第3可行解，其它同前第2可行解（第3表）第3列Cs-Zs最大，f3为换入基变量上表L3-L2上表L2/6

+L1原L2/6a1a2a3a4a5a6bca110－1/21/31/6023a201－1/2－2/31/6011a6006－1－1130CS312000ZS31－21/32/307CS－

ZS004－1/3－2/3037第37页，课件共68页，创作于2023年2月⑤求第4可行解仅3表3列Cs-Zs

>0，且该列中仅为正β6.3为正，故换出变量为f6。第3可行解（第4表）上表L3/12

+L1上表L3/12

+L2上表L3/6a1a2a3a4a5a6bca11001/41/121/129/43a2010－3/41/121/125/41a3001－1/6－1/61/61/22CS312000ZS312－1/302/39CS－

ZS0001/30－2/3虽所有松弛变量（x5

,x6

,x7）均从基变量中消去了，但此时并非最优解。因第4列Cs-Zs

仍为正，且β1.4

>0

，故目标函数P仍可能增大。故仍需以f4

为换入基变量，f1

为换出变量，求第5可行解。38第38页，课件共68页，创作于2023年2月上表得，第4可行解⑥第5可行解a1a2a3a4a5a6bca440011/31/390a231001/31/381a32/3010－1/92/922CS312000ZS13/31201/97/912CS－

ZS－4/3000－1/9－7/9由Cs-Zs

判别Cs-Zs

<0，此时优化解即为目标函数的最大值为39第39页，课件共68页，创作于2023年2月另外，根据（5.59）式可由当前次的θ

,(Cs-Zs

)值，预测下次的新目标函数值。例如：如上算得第1次：

P=0θ

4=1Cs-Zs=3预测P’=P+θ(Cs-Zs)=0+1×3＝3与计算值相同第2次：

P=3θ5=1Cs-Zs=4预测P’=3＋1×4＝7与计算值相同第3次：

P=7θ6=1/2Cs-Zs=4预测P’=7+½×4＝9与计算值相同因此时表5中Cs-Zs

<0，故目标函数极大值即为12虽然这种方法能预测下次的极值，但仍须以上述表格的运行结果为依据，才能进行继续预测。当然，直接应用表格运算也可解决问题。第4次：

P=9θ1=9Cs-Zs=1/3预测P’=9+1/3×9＝12与计算值相同40第40页，课件共68页，创作于2023年2月目标函数适当选取x(t)，可使目标函数P达极值。P称为泛函数（简称泛函），x(t)称为变函数。5.2.5变分法泛函的极值求解就是变分问题。如：两点A—B间沿何曲线积分时距离最短？yxA(xA,,yA,)B(xB,,yB,)y=y(x,)亦即相当于选择y=y(x)，使积分（泛函）达极值。即求最小问题41第41页，课件共68页，创作于2023年2月泛函的极值条件是变分为0：即得欧拉（Euler）方程或讲义丢掉用上述方法求变函数xq不明显含有变量t当时，可使用另一种形式的欧拉方程或42第42页，课件共68页，创作于2023年2月例6粒子由A→B，时间最短的轨迹方程y=y(x)解：（无阻力）A点速度为零至中间某点（1）(方向切于曲线)设A→B距离s，中间某微段ds（2）粒子由A→B的时间（I相当于目标函数P）（3）此处x为变量yxABy=y(x,)v被积函数（q）不显含变量x，故用方程亦即（4）整理得（5）dxdyds43第43页，课件共68页，创作于2023年2月由三角公式得（6）仍由亦解得（6）式代入上式得44第44页，课件共68页，创作于2023年2月积分可得为摆线方程，R为圆半径A=Rtgθ=dy/dx推论12个以上变函数的变分问题（x为自变量，y、z为x函数）极值满足（5.62）推论2多个自变量的泛函问题（x，y，z为自变量，U为x，y，z函数）45第45页，课件共68页，创作于2023年2月该泛函达极值的必要条件为：（5.63）有约束条件的泛函达极值问题泛函形式约束条件引入拉格朗日乘子，给出新的泛函形式新的泛函极值的必要条件为：（5.64）46第46页，课件共68页，创作于2023年2月例7.圆形管道横断面旋度平方的面积分为最小（条件），平均流速U0，求速度分布。梯度散度旋度圆柱坐标47第47页，课件共68页，创作于2023年2月旋度平方的面积分的泛函形式即为目标函数由流体不可压缩，流量为Const，约束条件组成新的目标函数：圆形断面，速度分布与无关，且VZ=V，故圆管

V仅为r的函数，故有48第48页，课件共68页，创作于2023年2月积分或者为单变量函数，自变量和约束条件的泛函问题，由（5.64）方程则有49第49页，课件共68页，创作于2023年2月可求得由此说明：①与流体力学中层流条件所得速度分布一致;②实际问题一般难以得分析解，数值解为多;③极值须用二阶变分验算。勒让德判据，假如满足

Euler

方程例6（极小）例7（极小）50第50页，课件共68页，创作于2023年2月5.2.6最大原理xi

—状态变量，T—操作变量，t—坐标变量求使P达极值的f0

函数。状态变量速度式：设x0

为虚拟成分，且（5.67）xi

是T的函数反应工程中计算最佳温度是一个重要问题。目标函数（5.66）51第51页，课件共68页，创作于2023年2月现以数学符号引入另一个函数J（5.68）λi(t)称可变拉氏乘子，与t有关，与T无关，又称为时间乘子。λi(t)边界条件（规定）在反应过程中0～θ

时间内，对应的最佳温度为温度分布。如何使P达极值？下面介绍Katz

方法（1960年）设x0

满足52第52页，课件共68页，创作于2023年2月（5.70）对于温度T的变分xi(0)已知λi(0)一定J－P为Const，J增大，P必增大，反之亦然，问题为求J极值。由变分法δT=0，可得极值必要条件。

自变量t，xi…T均为t之函数，J为泛函。同样，xi的变分（5.68）右端积分Ttt53第53页，课件共68页，创作于2023年2月说明：（5.76）（1）求μi（a）可写成（b）（c）54第54页，课件共68页，创作于2023年2月对（c）Talor级数展开（d）比较（a）,（d）式，并考虑定义式（b）,有或（5.75）（2）求λi（5.75）代入（5.76）（5.77）55第55页，课件共68页，创作于2023年2月上式很复杂，可适当选取λi使其简化如（5.78）即取第j个因子可推出由此，式（5.77）右端首末项对消，有（5.79）因

≠0，欲使δJ=0，必有（5.80）上式为使P达极值的条件。问题的边界条件(共轭函数)56第56页，课件共68页，创作于2023年2月说明：以上为问题最优解的必要条件，但不是充分条件。可能有多个解满足上述条件，只有一个为最优。也可能因约束条件限制使过程不能达最佳状态。这时可用彭策根（Pontryagain）最大原理。由哈密顿函数：H与P的极值（极大、极小）相应。现考察（5.83）问题的边界条件57第57页，课件共68页，创作于2023年2月将共轭函数（5.78）和极值必要条件（5.80）代入（5.83）得积分上式且因f0也是常数，故H在最佳点一定是常数。H函数加上λi

把变分问题转化为微分问题，即dH/dt=0得极值条件，但须加上共轭函数条件。后项为“0”，前项必为“0”，即H=Const在最佳条件下，T可由H/T=0求解。例8：连续流反应器中的连串反应A→B→C，设反应为1级试求使产物B产率最大的温度分布。58第58页，课件共68页，创作于2023年2月解：反应速度方程（1）（2）式中：x1、x2

为反应器中任意点A、B浓度（CA、CB），t=l/vl——该点至入口距离，v平均流速，反应常数为：求B产率最大，即必须x2(θ)为最大，目标函数可写为（3）共轭函数为即（4）59第59页，课件共68页，创作于2023年2月本题中而由（4）式，可得（5）（6）温度选取必须使H最大，即60第60页，课件共68页，创作于2023年2月式中且（7）61第61页，课件共68页，创作于2023年2月（1）（2）（5）（6）初始、边界条件则沿反应器的最佳温度分布为宜用数值求解方便起见设，则整个问题化为解下列方程：62第62页，课件共68页，创作于2023年2月5.2.7动态规划以多段横流液—液萃取过程为例。原液流量

Q,浓度

x0,