1、cpsesystem当要求容器的容积一定,求表面积最小,以当要求容器的容积一定,求表面积最小,以使用料最省。使用料最省。x1x22112min32Sxx x3211223xx xVs.tx10,x20一连续反应器如图所示,进行如下反应2AB已知单位体积的液相反应速率为原料A的单位成本折旧、公用工程和其他费用2220d2.02.0(1)dAAAAAcrccxt 1.41A04.0Cc0.620.4CV根据预测,市场只能提供物料A 600单位/h,产品B的市场需求量FB不超过50 单位/h,产品B的价格为C3=2000 元/单位。试确定物料A的进料速度FA0、初始浓度CA0、反应器体积V和转化率各
2、取多大,才能使得该反应器在单位时间内的经济效益是最好的?3B1A021.40.63BA0A0max()(40.4)zC FC FCC FcFV22AAA0Ad2.0(1)dcrcxt A0A0AA0A0A(1)F cr VF cx A0A0A0.550F cx A0600F目标函数或约束条件中有非线性函数的规目标函数或约束条件中有非线性函数的规划问题划问题非线性规划的最优解可能在其可行域中的任非线性规划的最优解可能在其可行域中的任意一点达到意一点达到不一定是全局最优解不一定是全局最优解背景背景理论计算理论计算相对于计算要求,计算能力仍十分有限相对于计算要求,计算能力仍十分有限背景背景为加快计算
3、速度,必须明确各种方法的特点,为加快计算速度,必须明确各种方法的特点,以针对不同问题选择最合适的方法以针对不同问题选择最合适的方法求解思路:迭代迭代从一个选定的初始点x0出发,按照某种特定的迭代规则产生一个点列xkxk有穷点列:最后一个点为最优解xk无穷点列:其中一个点为最优解基本迭代格式 :第k轮迭代点 :第k+1轮迭代点tk:搜索步长pk:迭代方向1kkkkxxt pnkxR1nkxR对于存在0,使则称x*为R上的局部极小点局部极小点,f(x*)称为局部极局部极小值小值严格局部极小点、严格局部极小值严格局部极小点、严格局部极小值():nf xRE*()()xxf xf x若对于任意x,有则
4、x*为R上的全局极小点全局极小点,f(x*)为全局极小值全局极小值严格全局极小点、严格全局极小值严格全局极小点、严格全局极小值*()()f xf x凸集:对于在集合中的每一对点x1和x2,连接两点所形成的直线段上任一点都在此集合内,则该集合为凸集凸集12(1)01xxx凸函数如果函数 满足则称f(x)为F上的凸函数凸函数若则称f(x)为严格凸函数严格凸函数1212(1)()(1)()fxxf xf x()f xF1212(1)()(1)()fxxf xf xyxox1x2 x1+(1-)x2y=f(x)凸函数12(1)fxx12()(1)()f xf xyxox1x2x1+(1-)x2y=f(
5、x)12(1)fxx12()(1)()f xf x非线性规划如f(x)和gi(x)都为凸函数,则称该规划问题为凸规划凸规划可以证明:f(x)的局部极小值也是全局最小值的局部极小值也是全局最小值理想情况理想情况min()()0if xg x 若f(x)有连续的一阶导数,则f(x)为凸函数 对x1、x2 R,有f(x2)f(x1)+f(x1)T(x2-x1)f(x)为严格凸函数 对x1、x2 R,有f(x2)f(x1)+f(x1)T(x2-x1)Hessian矩阵2()()H xHf x 0,00,00,00,00TTTTTxx HxHxx HxHxx HxHxx HxHxx HxH半正定正定半负
6、定负定部分0,不定H为对称矩阵例 判断下列函数的凹凸性(x R)(a)f(x)=3x2(b)f(x)=2x(c)f(x)=5x2(d)f(x)=2x2-x3解(a)f”=6,故f(x)为(严格)凸函数。(b)f”=0,故f(x)既凸又凹(c)f”=-10,故f(x)为(严格)凹函数(d)f”=6-3x,故f(x)既不为凸也不为凹对于多元函数,如何判断H是否正定?特征值特征值f(x)H特征值特征值严格凸函数正定0凸函数半正定0凹函数半负定0严格凹函数负定0 f(x)0对n维函数必要条件:f(x)在x*处一阶可导充分条件H(x*)正定,则x*为极小值,反之为极大值*()0f x例 求函数的所有稳定
7、点解2242121212112()44.54222f xxxxxx xxx x3121121221124.52244044220fxxxx xxfxxxx 解方程组得T1T2T3(1.941,3.854),(1.053,1.028),(0.6117,1.4929)xxx 试判断所得的稳定点是否为最优解2131.7949.764()(1.941,3.854)9.7644H xf 2211.1942.212()(1.053,1.028)2.2124H xf 230.5194.447()(0.6117,1.4929)4.4474H xf 求得各点的H特征值和稳定点类型如下:稳定点f(x)特征值1.9
8、41,3.8540.985537.030.97局部极小点-1.053,1.028-0.513410.503.50(全局)极小点0.6117,1.49292.83007.0-2.56鞍点一维搜索法一维搜索法多项式近似多项式近似求导数方法求导数方法Fibonacci法法0.618法法二次插值法二次插值法三次插值法三次插值法一阶导数一阶导数二阶导数二阶导数最速下降法最速下降法共轭梯度法共轭梯度法牛顿法牛顿法拟牛顿法拟牛顿法一维搜索法一维搜索法步长tk的选定是由使目标函数值沿搜索方向下降最多为依据的,因此这一工作变成了求解以tk为变量的一元函数,故得名一维搜索法一维搜索法。1kkkkxxt p适用于某
9、些不能求得一阶导数解析解的问题如求最小回流比其中ij:组分i对组分j的相对挥发度 xDi:塔顶产品中i组分的组成 :由Underwood公式确定min11ni jDiiijxR11ni jFiiijxq 用经典的微分方法很难求解斐波那契(斐波那契(Fibonacci)法(分数法)法(分数法)0.618法法无需求导,根据函数值判断搜索方向无需求导,根据函数值判断搜索方向适用于求解已知极值区间的单峰函数适用于求解已知极值区间的单峰函数f(x2)f(x1),去掉x1,b0,此时x*a0,x1f(x)xoa0b0 x*x1,x2 在x*的右侧x1x2f(x2)f(x1),去掉a0,x2,此时x*x2,
10、b0f(x)xoa0b0 x*x1,x2 在x*的左侧x1x2f(x2)f(x1):a.去掉x1,b0,此时x*a0,x1 b.去掉a0,x2,此时x*x2,b0 f(x)xoa0b0 x*x1,x2 在x*的两侧x1x2斐波那契数列数列Fn为斐波那契数列 斐波那契分数110 FF,3,2,12nFFFnnnnnFF1n012345Fn112358n67891011Fn1321345589144计算步骤计算步骤选取初始数据,确定单峰区间a0,b0根据缩短率计算Fn,再确定最小n值计算初值t1和t2,计算f(t1)、f(t2)当区间变为a0,t212()()f tf t101221(2)1111
11、(1);()nnaa bt ttFtbbaF反之区间变为t1,b0111012(2)2111(1);()nnat bb ttFtabaF12()()f tf t确定n个搜索点以后,每次的区间缩短率为1231122,nnnnnnFFFFFFFFn次计算能得到的区间长度比为要使精度够大,即:区间缩短的相对精度相对精度1nF1nF如果至某一步则可令121()2ttab211()21()()2tabtaba可以证明对于斐波那契数列,其奇数项和偶数项都各自收敛于同一极限,该极限值等于510.6180339882以0.618作为固定的区间缩短率代替斐波那契法不同的缩短率,就得到了0.618法实施更为容易1
12、(1)0.382tLL210.618,tLLt为 的对称点。计算步骤计算步骤选取初始区间a0,b01000000(1)()0.382()tabaaba2000000()0.618()tabaaba f(t1)f(t2),取a1=a0,b1=t2 t2=t1 t1=a1+0.382(b1-a1)f(t1)f(t2),取a1=t1,b1=b0 t2=a1+0.618(b1-a1)t1=t2 Lt2t1a0b0La1b1t2t1a1b1t2t1搜索n个点后的 区间长度缩短为或者说,迭代k次以后的区间长度变为即 已知搜索的相对精度,迭代次数满足0.618k110000()()0.618nnbaba00
13、()0.618kba例 用0.618法求设初始区间为a0=-1.0,b0=3.0,要求剩余区间长度不大于0.1解 本例可以通过解析法求得精确解2min()2xxx*0.5,()1.75xx第一次搜索选取两个初始试算点比较得,(1,1)000(1,1)(1,2)000(1,2)0.382()0.528,()1.750780.618()1.472,()2.69478xabaxxabax(1,1)(1,2)()()xx可以得到下一次搜索区间11(1,2)(2,2)(1,1)(2,2)(1,1)1,1.472,0.528()()1.75078abxxxxx 第二次搜索第二次搜索计算试算点(2,1)11
14、10.382()10.382(1.4721)0.055696xaba (2,1)(2,2)()()xx确定下一轮搜索区间2(2,1)(3,1)(2,2)21(3,1)(2,2)0.055696,0.5281.472()()1.75078axxxbbxx 迭代迭代次数次数kak-1bk-1x(k,1)x(k,2)(x(k,1)(x(k,2)剩余区剩余区间长度间长度1-130.5281.4721.750782.694782.4722-11.472-0.0556960.5282.058791.750781.5276963-0.0556961.4720.5280.888421.750781.90087
15、0.9441164-0.0556960.888420.3049560.5281.788041.750780.58346450.3049560.888420.5280.6655361.750781.777400.3605860.3049560.6655360.44270.5281.75321.750780.228101.750021.750060.0325不断用低次(不超过三次)多项式来近似目标函数,并逐步用插值多项式的极小点来逼近最优解最速下降法最速下降法共轭梯度法共轭梯度法牛顿法牛顿法拟牛顿法拟牛顿法割线法割线法对于迭代格式使目标函数f下降最快的方向是点xk的负梯度方向称为f在xk的最速下降
16、方向最速下降方向每一轮都从点xk出发沿最速下降方向进行搜索的方法,称为最速下降法最速下降法1kkkkxxt p()()kkkf xpf x 具体步骤选取初始点x0,给定终止误差求梯度向量。计算 ,若 ,停止迭代,输出xk构造负梯度方向进行搜索。求tk,使得令 ,重新求梯度向量()kf x()kf x0()min()kkkkktf xt pf xtpkp1kkkkxxt p例 求解 令x(0)=(1,2)T,1,2为任意实数假设x(0)不是最优点则令221212min()222,0.1f xxxxxT12(22),(22)fxx T(0)12()(22),(22)f x1(0)1(1)(0)(0
17、)(0)2(0)22(1)()2(1)txxtf xt代入目标函数,求最优步长t0有令22(1)1(0)12(0)21(0)12(0)2()min2(1)2(1)22(1)22(1)2f xtttt1(0)(0)d()0,0.5df xtt得T(1)(0)(0)0.5()(1,1)xxf x因故x(1)为最优解,最优值f(x*)=0T(1)(1)()(0,0),()0f xf x对于目标函数的等值线为圆的问题,最速下降法总能一步得到最优解例 用最速下降法求解 取x(0)=(2,2)T2212min()25,0.2f xxxT12()(2,50)f xxxT(0)(0)()(4,100),()1
18、00.08f xf x令代入目标函数,得t(0)=2.003(0)(1)(0)(0)(0)T(0)(0)()()(20.039968,20.9992)f xxxtf xttT(1)(1)(1.92,0.003),()3.68xf xT(1)(1)()(3.84,0.15),()3.843f xf x继续迭代,得kt(k)x1x202.0052.0002100.0810411.8501.920-0.0033.8433.68720.0710.0710.0713.5470.13130.0650.0680.0000.1360.46310-240.0030.0030.0030.1260.16410-2(
19、)()kf x()()kf x可能在任何类型的稳定点终止。通过分析Hessian矩阵来检验是否是极小点。对f(x)的尺度太灵敏,收敛缓慢,容易产生大量摆动(局部最速下降,相邻搜索方向彼此正交)。f(x)为二阶可导函数,且在每个搜索方向上能准确最小化,则较少次迭代即可收敛计算量略大,收敛速度大幅提高搜索方向结合当前梯度方向和前一次梯度方向,搜索方向共轭解大型线性或非线性方程组最有效的方法之一存储空间小,稳定性高1()()()()()2TTkkkkkf xQ xf xf xxxH x若H正定,则Q(x)的稳定点xk+1为最小点,该点可表示为1()()0kkkQ xf xH x 11()kkkxxH
20、f x解得1kkkkxxt p对照可得由此,可反复利用方程直到满足收敛条件1kt 11()kkkxxHf x例 求解解,取x(0)=(0,0)T有22121122min()22f xxxxx xxT(0)1212T()(142),(122)(1,1)f xxxxx 4222H10.50.50.51H得为检验x(1)是否为最优点,1(1)(0)(0)()00.50.51100.5111.5xxHf x T(1)1212T()(142),(122)(0,0)f xxxxx 即x(1)为最优点几何意义()kf x例 用牛顿法求解 该函数的一阶、二阶导数分别为432min()46164,0.01xxx
21、xx322()4(334)()12(21)xxxxxxx若取初始点x(0)=3,则故()52()24xx(0)(1)(0)(0)()5235.167()24xxxx迭代次迭代次数数kx(k)(x(k)(x(k)|(x(k)|03-52245215.167153.352184.33153.35224.33532.302109.44632.30234.0403.38386.8703.38344.0010.005584.0470.0055|(x(k)|0.01收敛速度很快,对于严格凸函数,只需一步即可得到极小值对纯牛顿法模型,每一步步长为1,但如果初始值与局部极小值不是足够接近,通常不收敛需要计算一
22、阶、二阶导数,计算量大,存储空间大对于一阶导数非单调变化的函数,往往会失败为避免计算二阶导数矩阵H及其逆阵,我们设法构造另一个矩阵,用它来逼近二阶导数矩阵的逆阵,称拟牛顿法拟牛顿法构造近似矩阵当f(x)为二阶函数,H为常数矩阵即对于非二阶函数()kH11()()()kkkkf xf xH xx111()()kkkkxxHf xf x111()()kkkkkxxHf xf xBFGS校正:近似Hessian矩阵:()()()()(1)()()()()()()()()()kkTkkkTkkkkTkkTkkxxHGGHHHGxGHG用有限差分代替一阶、二阶导数()()(1)()2()()()()()
23、/2()2()()/kkkkkkkxhxhhxxxhxxhhh:步长例 用拟牛顿法求解 令取x(0)=3,h=10-32min()xxx()()(1)()2()()()()()/()2()()/kkkkkkkxhxhxxxhxxhh3(1)3 23(3.001)(3.0)/103(3.001)2(3.0)(2.999)/(10)0.0050013(10)32.500500.0000020.499500 x7(0.4995)(0.5000)2.5 10最优解析解比牛顿法迭代速度更快适用于写不出目标函数的明确解析式,或导数的解析式对于多参数问题,计算时间和存储空间大函数二阶可导函数二阶可导初值初值几种计算收敛方法的比较 方法方法预测预测搜索搜索方向方向预测预测搜索搜索距离距离是否是否要求要求导数导数是否能是否能保证总保证总是是局部局部下降下降由基准点由基准点预测新方预测新方向时要算向时要算多少次函多少次函数值数值相对收相对收敛速度敛速度直接迭代法直接迭代法 是是 是是 否否 否否 1后期搜后期搜索慢索慢牛顿拉夫牛顿拉夫森法森法是是 是是 是是 是是 n 快快割线法割线法 是是 是是 否否 是是 n相当相当快快拟牛顿法拟牛顿法 是是 是是 否否 否否 1相当快相当快