机械优化设计》讲义

《机械优化设计》讲义》由会员分享，可在线阅读，更多相关《机械优化设计》讲义（260页珍藏版）》请在金锄头文库上搜索。

1、 机械优化设计讲义绪 言 优化设计是 1960 年代初发展起来的一门新学科，它是以电子计算机为工具，使用最优化理论寻求最优设计方案的一种现代设计方法。 最优化理论是一个重要的数学分支，它所研究的问题是讨论在众多的方案中什么样的方案最优以及如何找出最优方案。 这类问题普遍存在于各个领域中。运筹学（Operations Research）用它研究生产、管理、 商业、 军事、 决策等领域中的问题。 优化设计 （Optimal Design）用它处理工程设计领域中的设计问题。 在机械设计领域，传统的设计过程通常按下面步骤进行： 1、在调查分析的基础上，通过估算、经验类比或者实验来选择初始设计参数。 2

2、、对尺寸、强度、刚度、稳定性等各项设计要求进行计算和检查。 3、如果设计要求得不到全部满足，设计人员将调整修改某些设计参数，然后转第 2 步。如此反复，直到所有的设计要求都得到满足为止。 由此可见，传统的机械设计过程本质上是人工反复试凑的过程。用这种方法找到的设计方案，只是众多可行方案中的一个，一般都有再改进的余地。 使用优化设计方法进行机械设计，即用电子计算机的优化计算取代传统设计的人工试凑，不仅能够实现设计计算的自动化，把设计人员从反复检查、反复修改的繁琐计算中解放出来，而且能够获得人工试凑难以得到的、众多可行方案中最优的方案。 一个机械优化设计问题包括两方面内容： 1、把实际的设计问题化

3、为数学规划问题，即建立数学模型。 建立数学模型时，需要应用专业知识来确定设计的限制条件和追求的目标，以确立各设计变量之间的相互关系。 2、求解这个数学规划问题。 根据数学模型的特点，应用优化设计的理论，选择适当的优化算法，使用计算机求解。 第 1 章 优化设计的数学模型 1.1 一个简单的优化设计问题 例 1.1 试设计一个用钢板焊接而成的密封圆筒形容器（图 1.1）。要求其容积为 2 m3，能承受内部 p = 3MPa 的蒸汽压力。受安装空间限制， 要求其外部直径和高度分别为 1 m d 3 m 和 1 m h 3 m。 (d-2t)l p h d (a) 外部尺寸 (b) 受力图 图 1.

4、1 圆筒形容器示意图 t l d-2t d l t l 正应力 所产生的内力： 2tl 蒸气压力 p 所产生的外力： (d-2t)l p lptdtl)2(2 由此可得该容器的强度条件： 2)2(tptd 如果选用厚度t在 1mm 20mm 之间的 Q235 钢板，那么其许用应力为 160 MPa，强度条件被整理为 3d - 326t 0 。因此该容器的全部设计要求为： 201300010003000100003263102)2(4)2(92thdtdthtd (1.1) 先选 t = 8 mm ，然后根据 3d - 326t 0 ，得 d 837 mm 。 与直径约束 d 1000 mm 相

5、抵触， 于是修改设计。 这次选 t (d-2t)l p h d (a) 外部尺寸 (b) 受力图 图 1.1 圆筒形容器示意图 t l d-2t d l t l = 13 mm ，得 d 1413 mm 。根据直径条件选 d = 1100 mm ，并代入容积条件，获得高度 h = 2233.7 mm 。由于 2233.7 mm 在规定高度范围之内，所以就得到一个能满足全部设计要求的方案： t = 13 mm，d = 1100.0 mm，h = 2233.7 mm 使用这种方法还可以设计出另外的方案。式(1.1)中只有 1 个方程却有 3 个未知量，因此存在无穷多个解，即存在无穷多个满足设计要求

6、的方案。 如果根据容器耗费材料的多少来衡量设计方案的优劣，那么需要把例 1.1 的问题表述为： 201300010003000100003263102)2(4)2(s.t.)2(2),(min922thdtdthtdtthdtdhdtf (1.2) “min”是英文“minimize”的缩写，它后面的函数称为目标函数，意为使目标函数的值最小化。 “s.t.”是英文“subject to”的缩写，意为“受约束于” 。它后面的等式和不等式是使目标函数最小化时的约束条件，分别称为等式约束和不等式约束。 1.2 数学模型的建立 优化设计的数学模型由设计目标和设计约束两部分组成： 设计目标：“min”及

7、其后面的目标函数。 设计约束：“s.t.”及其后面的等式和不等式。 实际设计问题数学模型的抽象。 确切：模型不能失真。 简洁：不要太复杂。 需要注意：适当确定设计变量，合理构造目标函数，正确列出约束条件 1.2.1 适当确定设计变量 机械设计方案参数：几何量（结构尺寸、位置关系等）；物理量（材料的弹性模量、许用应力，零件的工作速度、加速度等）；工程参数（比如齿轮的模数和齿数、轴的挠度和自振频率等等）。 设计变量：其变动会直接或间接地影响目标函数值大小的那些设计参数。 为了降低问题求解的难度，设计者应尽量减少设计变量的数目。 1、例：圆筒形容器：储油罐 钢板厚度t可以不作为设计变量出现在数学模型

8、中（尽管设计说明书需要这个参数）。这时若用 d 和 h 表示容器内壁的尺寸，其变化范围为 1 m3 m，则其数学模型为 30001000300010001024s.t.2),(min922hdhddhdhdf （求解结果： d = h = 1365.6mm）。 2、 如果根据需要和经验事先选定一种材料。 那么与材料相关的弹性模量、许用应力等参数都可取为常数，使数学模型变得简洁。 3、离散变量（齿轮的齿数、模数、弹簧钢丝直径、板材厚度等）先处理成连续变量，降低求解难度。求得结果后，再让它们取结果附近的离散值。在一定程度上损害了数学模型的确切性，但很多时候可近似使用。 20130001000300

9、0100003263102)2(4)2(s.t.)2(2),(min922thdtdthtdtthdtdhdtf例如式(1.2)中的优化结果为： t = 9.2 mm，d =1000.0 mm，h = 2661.3 mm 9.2 mm 不属于国家标准的钢板厚度系列。t 应当取 9.2 mm 附近的系列值。现取 t = 10mm，把式(1.2)中关于 t 的不等式约束替换为等式约束，问题成为 10300010003000100003263102)2(4)2(s.t.)2(2),(min922thdtdthtdtthdtdhdtf 优化结果： t = 10.0 mm，d = 1086.7 mm，h

10、 = 2258.1 mm 。 1.2.2 合理构造目标函数 目标函数：设计者度量设计方案优劣程度的量化指标。 1、许多机械设计问题，都以质量最小（它通常与材料体积或表面积最小等价）作为设计目标。 减小质量 = 节省材料、节省加工费、减小惯性力、降低能耗。 2、对于应力集中现象严重的构件：应力集中系数最小。 对于精密仪器：测量误差最小 对于再现运动轨迹的机构：轨迹误差最小。 总之， 应当从设计对象的用途出发， 以最重要最具代表性的指标作为目标函数。 例：圆筒形容器：饮料罐（300ml）。 设计目标：外观漂亮。d / h = 0.618 时，外观最匀称。因此可以用 d / h 与 0.618 之间

11、的正负偏差最小（尽量符合黄金分割率）作为设计目标。其数学模型为 06001034s.t.618. 0),(min52hdhdhdhdf 结果： d = 60.0mm，h = 106.1mm。 3、多目标优化问题化为单目标优化问题。 把最重要的某个设计目标作为单目标优化问题的设计目标，把其它设计目标处理成约束条件。 例：圆筒形容器：油漆桶（20000ml） 双目标：制造油漆桶的材料最少以降低成本，外观匀称以吸引顾客。 可以把用料最少作为设计目标， 把外观匀称用 d/h 与 0.618 之间的正负偏差不大于某个允许值的约束条件来表示。 假设该允许值为 0.1，这个双目标优化问题就转化成单目标优化问

12、题： 001 . 0618. 01024s.t.2),(min722hdhdhddhdhdf 结果：d = 263.4mm，h = 366.9mm。 4、设计变量不一定在目标函数的表达式里出现。 例如式(1.2)的问题，材料体积也等于容器所占空间的总体积减去其容积： 921024),(hdhdtf (1.7) 因此式(1.2)中的目标函数 tthdtdhdtf)2(2),(2 完全可以用式(1.7)来代替。这时的目标函数中虽然不出现壁厚 t，但它仍然是 t 的函数，是壁厚 t 的隐函数。 5、原目标函数上乘（除）一个常数或加（减）一个常数，只会改变目标函数值的大小，并不影响最优方案本身。 为了

13、使数学模型简洁，可以把原目标函数中的这些常数去掉。因此，例 1.1 的优化设计数学模型也可写为 201300010003000100003263102)2(4)2(s.t.),(min922thdtdthtdhdhdtf 结果仍然是 t = 9.2 mm， d = 1000.0 mm， h = 2661.3 mm 。 1.2.3 正确列出约束条件 约束条件是设计方案必须满足的各种设计限制。 约束条件分两类，一类称为性能约束，一类称为边界约束。 性能约束 零件设计：强度条件、刚度条件、振动稳定性条件、耐热性条件等。 机构设计：装配条件、邻接条件、传动比条件、速度条件、加速度条件等。 许多性能约束

14、实际就是设计规范的计算公式， 它们可以根据力学、机械学、几何学的知识推导出来。 边界约束： 对设计变量取值范围的限制， 它给出设计变量的上边界和下边界。 边界约束根据设计对象的结构需要或经验给出。注意既不能遗漏必要的边界约束，也不能无根据地缩小边界约束的范围。 1、数学模型必须列出全部必要的约束条件，不能遗漏。 例：在锅炉设计的式(1.2)中，如果遗漏了最后三个边界约束，成为 03263102)2(4)2(s.t.)2(2),(min922tdthtdtthdtdhdtf (1.9) 那么它的解是 t = 2.4 mm，d = 260.6 mm，h = 38909.2 mm。它高达近 40米，

15、与其说是一个锅炉不如说是一根旗杆。 2、注意约束条件不能互相矛盾。 例：如果把式(1.2)的钢板厚度的上边界减小为 8 mm，使它变为 81300010003000100003263102)2(4)2(s.t.)2(2),(min922thdtdthtdtthdtdhdtf (1.10) 那么该问题无解。之所以发生这种情况，是因为强度约束和直径约束隐含地要求 t 9.2，但修改后的厚度约束却要求 t 8，两者互相矛盾。所以，不要没有根据地随意缩小设计变量的取值范围。 3、 当某个设计变量可以通过某个等式约束表示成其它设计变量的函数时， 目标函数和约束函数中的这个设计变量便可以用这个函数替换。这

16、样不仅能减少数学模型的设计变量的个数，还能减少约束条件的个数。 例如式(1.3)的问题， 由于储油罐的高 h 可以根据等式约束表示成直径 d 的函数：30001000300010001024s.t.2),(min922hdhddhdhdf 29108dh 所以式(1.3)可以简化为 3000108100030001000s.t.1082)(min2992dddddf (1.11) 4、数学模型中的有些约束条件可能是多余的。多余的约束条件的存在不影响求解的结果。把它们从数学模型中剔除可以简化数学模型，加快求解速度。 例如式(1.11)中有 4 个不等式约束，其中最后两个可改写为300010810

17、0030001000s.t.1082)(min2992dddddf 661083108d 由于 30001081000310866 所以可以把式(1.11)的数学模型进一步简化为 6921081000s.t.1082)(minddddf (1.12) 当式(1.12)中的约束条件被满足时，式(1.11)中的另两个约束条件自然被满足，它们是多余的约束，所以可以删去它们。 问题(1.12)应用微分学的极值理论就可求解。令式(1.12)中的目标函数的一阶导数为零，求其可能的极值点： 3392932920001080108108)(ddddddf (1.13) 因为32000d满足两个约束条件，而且这

18、时目标函数的二阶导数 01016)(39 ddf (1.14) 所以32000d是目标函数 f(d)的极小点。 由此得储油罐的高 h 和直径 d 都为32000mm。 一般的优化设计问题无论怎样简化， 也不可能化成式(1.12)那样简单的形式，更不是只用微分学的极值理论就可轻易求出它的解析解的。要求解一般的复杂优化设计问题，必须借助电子计算机的力量。 本课程的主要内容就是讨论如何编制计算程序来求解优化设计的数学模型。 1.3 数学模型的一般形式 1.3.1 数学模型的一般形式 优化问题的数学模型可以归纳整理成下面的一般形式： qvgnpuhfvun, 2 , 10)(, 2 , 10)(s.t

19、)(minxxRxx Tnnxxxxxx,2121x 式中，x 是 n 维欧氏空间 Rn 中的一个列向量：它称为设计向量。设计变量 x1、x2、xn 是设计向量 x 的 n 个分量。 1、目标函数 f(x) 取最大值的方案为最优方案时，设计目标用 max f(x) 表示。由于 max f(x) = min (- f(x) ) 所以两者可以互相转换，式(1.15)不失其一般性。 2、 由于不等式的两端同乘以 “-1” 可以使不等号 “” 和 “” 互相转换，所以优化问题的不等式约束可以统一为 gv(x) 0 的形式（也可以统一为 “ 0” 的形式）。 3、 hu (x) 和 gv(x) 称为约束

20、函数。 例：油漆桶设计问题，令 x = (x1 , x2)T = (d , h)T，则 001 . 0618. 01024s.t.2),(min722hdhdhddhdhdf 0)(0)(0518. 0)(0718. 0)(0108)(s.t.)5 . 0(),(),(min241312221172211212121xgxgxxgxxgxxhxxxxxfhdfxxxxx 满足全部约束条件的 x 称为优化设计问题的可行解。一个可行解就是一个可行的设计方案。全体可行解构成的集合称为可行域。所谓优化设计，就是要寻找使目标函数 f(x) 取最小值的可行解 x*，这种可行解称为最优解，即最优设计方案。最

21、优解 x*所对应的目标函数值 f(x*)称为最优值。 式(1.4)的最优解 x* = (x1, x2, x3)T = (t, d, h)T = (10, 1086.7, 2258.1)T (mm)，最优值 f(x*) = 94954.9 cm3。 把可行域记为 D： vuvgpuhDxxx, 2 , 1, 0)(;, 2 , 1, 0)( (1.17) 数学模型的更简洁的形式： )(minxxfD (1.18) 它表示在 n 维设计空间的可行域 D 上，使目标函数 f (x)最小化。 1.3.2 优化问题的分类 当数学模型中不包含约束条件，即 p = q = 0 时，称为无约束优化问题，否则称

22、为约束优化问题。 当目标函数和所有的约束函数都是设计变量的线性函数时，称为线性规划问题，否则称为非线性规划问题。 生产计划管理的优化问题：线性规划。工程设计中的优化问题：非线性规划。 优化设计问题的复杂程度：设计变量和约束条件的个数。 小型优化问题：设计变量和约束条件的个数都不超过 10。 大型优化问题：设计变量和约束条件的个数超过 50。 中型优化问题：介于二者之间。 第 2 章优化问题有关的数学概念 2.1 n 维欧氏空间 1、当 n 维实向量空间 Rn中的任意两个列向量 TnTn,2121 的内积运算被定义为 niiiTT1 (2.1) 时，这个 n 维向量空间称为 n 维欧氏空间，本书

23、仍然用 Rn表示。 2、n 维欧氏空间中向量 的长度（模或范数）被定义为 niiT12 (2.2) 3、n 维欧氏空间中两个向量 , 之间的夹角,被定义为 Tarccos, (2.3) 因此 ,cosT (2.4) 2.2 正定二次型 含有 n 个变量 x1、x2、xn的二次多项式 )(2222),(11,1, 1113113211222222211121RijjiijnjijiijnnnnnnnnnnaaaxxaxxaxxaxxaxxaxaxaxaxxxP (2.5) 称为实二次型，简称二次型。它可用一个向量Tnxxx,21x和一个实对称矩阵 A 表示成矩阵乘积的形式： AxxxTnnnnn

24、nnnxxxaaaaaaaaaxxxP2121222211121121,)( (2.7) 如果对任何实向量 x0， 都有 P(x) = xTAx 0， 则称二次型 P(x) = xTAx 是正定的。 此时矩阵 A 称为正定矩阵， 或称矩阵 A 是正定的。 实对称矩阵A为正定的充分必要条件是A的各阶顺序主子式都大于零。即 0, 0, 02122221112112221121111nnnnnnaaaaaaaaaaaaaa (2.8) 如果对任何实向量 x0，都有 P(x) = xTAx 0，且有 x(0)0，使 x(0)TAx(0) = 0，则称二次型 P(x) = xTAx 是半正定的。此时矩阵

25、 A称为半正定矩阵，或称矩阵 A 是半正定的。 实对称矩阵A为半正定的充分必要条件是A的各阶主子式都大于或等于零。即 )1 (021212221212111niiiaaaaaaaaakiiiiiiiiiiiii ii ii ikkkkkk ( 2.9) 如果对任何实向量 x0，都有 P(x) = xTAx 0，0 = 0；二阶主子式为 0，可知)(12xf在整个平面上是半正定的，因此 f1(x)是凸函数； 由)(22xf的一阶顺序主子式为 2 0，二阶顺序主子式为 110，可知)(22xf在整个平面上是正定的，因此 f2(x)不仅是凸函数，而且是严格凸函数。 凸函数具有下面重要的性质： 1、如

26、果函数 f(x)是定义在凸集 D 上的凸函数，实数 a 0，那么函数 a f(x)也是定义在凸集 D 上的凸函数。 2、如果函数 f1(x)、f2(x) 是定义在凸集 D 上的凸函数，那么函数 f1(x) + f2(x) 也是定义在凸集 D 上的凸函数。 由此推论：， 如果函数 f1(x)、f2(x) 、fk(x)是定义在凸集 D 上的凸函数，实数 a 1、a 2、a k 0，那么函数 kiiifa1)(x 也是定义在凸集 D 上的凸函数。 例如， 因为例 2.1 中的211)(xfx和5103)(212122212xxxxxxf x都是整个平面上的凸函数，所以 51033)()(2)(212

27、1222121xxxxxxfffxxx 也是整个平面上的凸函数。 关于凸函数的重要定理： 定理 2.4 如果函数 g(x) 是定义在凸集 S 上的凸函数，那么集合 SgDxxx, 0)( 是凸集。 证明 在集合 D 中任取两点 x (1), x (2)，则对任何数值 0, 1，都有 g(x (1)0 和 (1- )g(x (2)0。由于 g(x)为凸函数，所以 0)1 ()1 ()2() 1 ()2() 1 (xxxxggg 即点)2()1 ()1 (xxx满足集合 D 的条件 g(x)0，在集合 D 之内。根据凸集定义，D 是凸集。 利用这个定理，可以通过对优化问题中约束函数是否为凸函数的判

28、别，来判断其可行域是否为凸集。 2.5.3 凸规划 对于优化问题 )(minxxfD (2.34) 如果 D 是凸集，f(x)是 D 上的凸函数，则称该问题为凸规划。 对于优化问题 qvgnpuhfvun, 2 , 10)(, 2 , 10)(s.t)(minxxRxx (2.36) 当 f(x)和全部 gv(x)都是可行域 D 上的凸函数，全部 hu(x)都是可行域 D 上的线性函数时，该问题是凸规划。 通过对优化问题是否为凸规划的判别，可以知道 f(x)的局部极小点是否是最优点。f(x)的局部极小点是这样来定义的： 如果在 x*D 的充分小的邻域内有 f(x) - f(x*)0， 则称 x

29、*为 f(x)的局部极小点。 如果该不等式去掉等号， 则称 x*为 f(x)的严格局部极小点。 目标函数 f(x)可能有若干个相同或不同的局部极小点，但是只有其中的一个(或一些)取最小值的才是最优点。 凸规划问题具有下面两个定理所示的重要性质。 定理 2.5 如果凸规划存在局部极小点，则局部极小点就是最优点。 证明 用反证法。 设 x*D 是 f(x)的局部极小点。 如果 x*不是 f(x)的最优点，即存在另外的点 x (k) D，使 f(x (k) f(x*)，则对任意0, 1，有 x = x* + (1- ) x (k) D。由于 f(x)是凸函数，所以 0 f (x) x (a) 图 2

30、.10 多个最优点与唯一最优点 x(1) x (2) x (3) f (x) (b) 0 x x* 0 f (x) x x(3) x(1) x(2) 图 2.9 局部极小点与最优点 f(x) = f x*+(1- ) x (k) f(x*) +(1- ) f(x (k) f(x*) + (1- ) f(x*) = f(x*) 令 1，则 xx*，即在 x*的邻域内存在点 xD，使 f(x) f(x*)。这与 x*是 f(x)的局部极小点矛盾。因此 x*是 f(x)的最优点。 定理 2.6 如果凸规划的目标函数 f(x) 是严格凸函数， 又存在局部极小点，那么它的局部极小点是唯一的。 证明 用反

31、证法。设存在不同的两点 x (1), x (2)D，它们取相同的极小值 f(x (1) = f(x (2) = C，则对任意 0 1，有 x = x (1) + (1- ) x (2) D。由于 f(x)是严格凸函数，所以 f(x) = f x (1)+(1- ) x (2) f(x (1) +(1- ) f(x (2) = C + (1- )C = C 这与 C 是 f(x)在 D 上的极小值矛盾。因此 f(x)的局部极小点是唯一的。 凸规划的局部极小点 = 最优点。（目标函数是严格凸函数时，最优点只有一个） 优化问题的求解就是寻找局部极小点，并且从中找出最优点。 第 3 章 优化问题的极值

32、条件 3.1 无约束优化问题的极值条件 无约束优化问题 nfRxx)(min (3.1) 3.1.1 一阶必要条件 定理 3.1 设函数 f(x) 在 x*点处连续可微，若 x*是 f(x)的局部极小点，则必有 x*点处的梯度0x)(f。 证明 把 f(x)在 x*点处作一阶泰勒展开，得 xxxxxxxxxxxxxxofffofffTT)()()()()()()()( (3.2) 因为 x*是极小点， 所以在 x*的充分小的邻域内有 f(x) - f(x*)0， 因此 0)()()()(xxxxxxxxxxxxxofffT (3.3) 令 xxxxd，于是 0)(xxxxdxofT (3.4)

33、 因为 d 是一个单位向量，所以当 xx*时 0)(dxTf (3.5) d 为任意方向的单位向量时该不等式都应成立，设 )()(xxdff (3.6) 并把它代入式(3.5)，得 0)(0)()()()()()(2xxxxxxxfffffffT (3.7) 所以 0x)(f。 0 f (x) x x* x (1) x (2) x (3) 图 3.1 驻点、鞍点、极小点和最优点 满足一阶必要条件的点称为 f(x) 的驻点。 驻点分为三种类型：局部极小点、局部极大点和鞍点。 鞍点是沿某些方向为极小，沿另一些方向为极大的驻点。 3.1.2 一阶充要条件 定理 3.2 设 f(x)是定义在 Rn上的

34、可微凸函数，x*Rn，则 x*为问题(3.1)的最优点的充要条件是0x)(f。 证明 必要性是显然的。如果 x*是最优点，自然是局部极小点，由定理 3.1 知0x)(f。 现在证明充分性。设0x)(f，则对任意的 xRn，有 0)()(xxxTf (3.8) 由于 f(x)是可微凸函数，根据定理 2.2，必有 )()()()()(xxxxxxffffT (3.9) 即 x*是最优点。 3.1.3 二阶必要条件 定理 3.3 设函数 f(x)在 x*点处二次连续可微， 若 x*是 f(x)的局部极小点，则必有0x)(f，并且海赛矩阵)(2xf正定或半正定。 证明 把 f(x)在 x*点处作二阶泰

35、勒展开，并考虑到0x)(f，得 2222)()(21)()()()(21)()()()(xxxxxxxxxxxxxxxxxxxxxofffoffffTTT (3.10) 因为 x*是极小点，所以在 x*的充分小的邻域内有 f(x)-f(x*)0，因此 0)()(21)()(22222xxxxxxxxxxxxxxxofffT (3.11) 令 xxxxd，于是 0)(21222xxxxdxdofT (3.12) d 是任意方向的单位向量，当 xx*时 0)(2dxdfT (3.13) 因此)(2xf必是正定或半正定的。 3.1.4 充分条件 定理 3.4 设函数 f(x)在 x*点处二次连续可微

36、，若0x)(f，且海赛矩阵)(2xf正定，则 x*是 f(x)的严格局部极小点。有 证明 把 f(x)在 x*点处作二阶泰勒展开，并考虑到0x)(f，得 2222)()(21)()()()(21)()()()(xxxxxxxxxxxxxxxxxxxxxofffoffffTTT (3.14) 用反证法，假设 x*不是严格极小点，那么 f(x)f(x*)0，因此 0)()(21)()(22222xxxxxxxxxxxxxxxofffT (3.15) 令 xxxxd，于是 0)(21222xxxxdxdofT (3.16) d 是单位向量，当 xx*时 0)(2dxdfT (3.17) 这个结果与)

37、(2xf是正定矩阵矛盾，因此 x*必是严格局部极小点。 总之，对于驻点 x*，如果)(2xf正定， x*是 f(x)的严格局部极小点； 如果)(2xf既不是正定的又不是半正定的， x*不是 f(x)的局部极小点；如果 f(x)在 Rn上是凸函数，x*是最优点；如果 f(x)在 Rn上是严格凸函数，x*是唯一最优点。 例 3.1 利用极值条件对问题 12232313131)(minxxxxfx 求解。 解 由一阶必要条件 0xxx222212121)()()(xxxxfxff 解得驻点 21012101)4()3()2()1(xxxx 由 f(x)的海赛矩阵 22002)(212xxf x 得各

38、驻点的海赛矩阵： 2002200220022002)4(2)3(2)2(2)1(2xxxxffff 由于)()()4(2)1(2xxff，不定，所以 x (1), x (4)不是极小点；)()3(2xf负定，所以 x (3)不是局部极小点而是局部极大点；)()2(2xf正定，所以点 x (2)是严格局部极小点。因此最优点 x* = x (2) = (1, 2)T。 3.2 等式约束优化问题的极值条件 对于等式约束优化问题 npuhfun, 2 , 10)(s.t)(minxRxx (3.18) 为了讨论方便，先引入向量约束函数 Tuhhh)(,),(),()(21xxxxh 把式(3.18)表

39、示成 0xhRxx)(s.t)(minnf (3.19) 由于存在 p 个等式约束方程，因此等式约束优化问题的 n 个变量并不是完全独立的。 不失一般性， 可以把前 p 个变量看作非独立变量，把剩下的(n - p)个变量看作独立变量。于是设计向量 x 可分解为两部分： DFxxx (3.20) 其中 xF是非独立变量组成的向量，xD是独立变量组成的向量： TnppDTpFxxxxxx),(,),(2121xx 因为 xF可看作是由等式约束方程组所决定的 xD的隐函数，所以目标函数 f(x)可看作是仅以 xD为自变量的函数： ),()(DDFDfFxxxx (3.21) 因此等式约束问题(3.1

40、8)与下面无约束优化问题是等价的： pnDDFRxx)(min (3.22) 为了研究问题(3.22)的极值条件（也就是问题(3.18)的极值条件），引入向量约束函数 h(x)的一阶偏导数矩阵： nppppppnppxhxhxhxhxhxhxhxh)()()()()()()()()(11111111xxxxxxxxxxh (3.23) 并把它分解成与 xF和 xD对应的两部分： DFxxhxxhxxh)(,)()( (3.24) 其中 npppnpDppppFxhxhxhxhxhxhxhxh)()()()()(,)()()()()(11111111xxxxxxhxxxxxxh (3.25) 现

41、在用微分法求 F(xD)的梯度。F(xD)的微分为 DTFTDDFDdfdfdfdFDFxxxxxxxxxx)()(),()( (3.26) 其中 TnppDTpFTnppTpdxdxdxddxdxdxdxfxfxffxfxfxffDF,)(,)(,)()()(,)(,)()(21212121xxxxxxxxxxxx (3.27) 为了把式(3.26)中的 dxF用 dxD来表示，把约束方程的两端微分： 0)()()()()(0)()()()()(1111111111111nnpppppppppnnppppdxxhdxxhdxxhdxxhdhdxxhdxxhdxxhdxxhdhxxxxxxxx

42、xx(3.28) 并按式(3.24)分解，得 0xxxhxxxhDDFFdd)()( (3.29) 设 Fxxh)(非奇异，则有 DDFFddxxxhxxhx)()(1 (3.30) 把它代入式(3.26) DTFTDDFDdfdfdfdFDFxxxxxxxxxx)()(),()( 得 DDFTTDdffdFFDxxxhxxhxxxxx)()()()()(1 (3.31) 因此 F(xD)的梯度为 TDFTTDffFFDxxhxxhxxxxx)()()()()(1 (3.32) 设 x*为无约束优化问题(3.22)的极小点，则有0x)(DF，即有 0xxhxxhxxxxDFTTffFD)()(

43、)()(1 (3.33) 引入一个 p 维列向量 Tp),(21，令 1)()(FTTfFxxhxx (3.34) 式(3.33)成为 0xxhxxDTTfD)()( (3.35) 把式(3.34)两边右乘 Fxxh)(并整理后，得 0xxhxxFTTfF)()( (3.36) 因为目标函数f(x)在x*处的梯度也由与xF和xD对应的两部分所组成： TTTfffDF)(,)()(xxxxx (3.37) 所以式(3.35)和式(3.36)可合并为 0xxhx)()(TTf (3.38) 把它的两端转置后为 0xxhxTf)()( (3.39) 即 0xx0xxxxpuuupphfhhhf121

44、21)()()()()()( (3.40) 这就是 x*为等式约束问题的极小点的一阶必要条件。由于函数 puuuhfL1)()(),(xxx (3.41) 在点（x*, *）处的梯度等于零时，可以写为 0xhxxxxxx)()()(),(),(),(1puuuhfLLL (3.42) 正好是式(3.40)和 x*所满足的等式约束条件，所以 L(x, )是研究约束优化问题的重要函数， 称为拉格朗日(Lagrange)函数。 = (1,1,p)T称为拉格朗日乘子向量。 3.2.1 一阶必要条件 定理 3.5 设函数 f(x)和 hu(x), (u =1, 2, p)在可行点 x*处连续可微，)(x

45、uh, (u =1, 2, p)线性无关，若 x*是等式约束优化问题(3.18)的局部极小点，则必存在一个实向量Tp,21，使得下式成立： 0xxpuuuhf1)()( (3.43) 定理 3.5 表明：目标函数在 x*处的梯度必为约束函数梯度的线性组合。 式(3.43)和原约束方程构成含有n+p个未知量及n+p个方程的方程组，求解这个方程组可以获得满足一阶必要条件的点 x*。这就是著 名的求解等式约束问题的拉格朗日乘子法。对凸规划问题，满足一阶必要条件的点 x*就是最优点。 3.2.2 二阶必要条件 定理 3.6 设函数 f(x)和 hu(x), (u =1, 2, p)在可行点 x*二次处

46、连续可微， 若 x*是问题(3.18)的局部极小点， 则必存在一个实向量 *Rp，使得一阶必要条件成立；并且对于满足条件 puhuT, 2 , 1, 0)(xd (3.44) 的任意向量 dRn，均有 0,2dxdxLT (3.45) h1(x)=0 h2(x)=0 D )( xf)(2xh)(1xh x* 图 3.2 等式约束一阶必要条件的几何意义 成立。 式(3.45)中的xx,2L是拉格朗日函数 L(x, )在 x*处关于 x 的海赛矩阵： puuuhfL1222)()(,xxxx (3.46) 3.2.3 充分条件 定理 3.6 设函数 f(x)和 hu(x), (u =1, 2, p

47、)在可行点 x*二次处连续可微，若存在一个实向量 * Rp，使得其一阶必要条件成立；并且对于满足条件 puhuT, 2 , 1, 0)(xd (3.47) 的向量 dRn，如果 d = 0，或者虽然 d 0 但是 0),(2dxdxLT (3.48) 成立，则 x*是问题(3.18)的严格局部极小点。 值得注意的是，约束优化问题的二阶必要条件和充分条件比无约束优化问题的条件宽松： 它只要求拉格朗日函数 L(x, ) 关于 x 的海赛矩阵在约束超曲面的切平面上（半）正定，并不要求它在整个欧氏空间中（半）正定。如果它在整个欧氏空间中是（半）正定的，那么在其子空间的约束超曲面的切平面上必然是（半）正

48、定的。 例 3.3 求解问题： 01s.t.)(min2121xxxxfx 解 计算目标函数和约束函数的梯度与海赛矩阵： 0xxxx)(0110)(11)()(2212hfhxxf 由一阶必要条件 0xx12)()(xxhf 和等式约束方程，解得拉格朗日函数的驻点：2121xx。 求满足式0)(xdhT的非零向量 d =(d1, d2)T。由 011,2121dddd 得 d2 = -d1，所以 d = (d1, - d1)T，d1可取非零的任意实数。因为 020110,)()(21111122dddddhfTdxxd 所以T21,21*x是严格局部极小点。 本问题只有这一个局部极小点， 因此

49、它就是最优点。 从这个例子可以看出，拉格朗日函数关于 x 的海赛矩阵在整个欧式空间中虽然是不定的，在约束超曲面的切平面上却可能是正定的。 3.3 不等式约束优化问题的极值条件 对于不等式约束优化问题 qvgfvn, 2 , 10)(s.t)(minxRxx (3.49) 可以引入 q 个松弛变量 z1, z2 , , zq，把不等式约束 gv(x)0 化为等式约束： qvzgvv, 2 , 10)(2x (3.50) 从而把原问题化为等式约束优化问题 qvzgfvvn, 2 , 10)(s.t)(min2xRxx (3.51) 因此可以用处理等式约束问题的方法，构造广义拉格朗日函数： qvvv

50、vzgfL12)()(),(xxzx (3.52) 其中 = (1, 2 , , q ) T 为广义拉格朗日乘子向量，z = ( z1, z2 , , zq ) T为松弛向量。 由广义拉格朗日函数(3.52)，即可得出点 (x*, *, z*) T为驻点的一阶必要条件： qvvvgfL1)()(),(0xxzxx (3.53) ), 2 , 1(0)()(),(2qvzgLvvxzx (3.54) ), 2 , 1(02),(qvzLvvzxz (3.55) 由式(3.54)知，gv(x*)与vz同时为 0 或同时不为 0，所以式(3.55)与 ), 2 , 1(0)(qvgvvx (3.57

51、) 完全等价，式(3.55)的条件可替换成式(3.57)。而式(3.54)写成 ), 2 , 1(0)()(2qvzgvvx (3.58) 后，就是原问题的约束条件，所以问题(3.49)的一阶必要条件中并不出现松弛向量 z*。 3.3.1 一阶必要条件(K-T 条件) 定理 3.8 设函数 f(x)和 gv(x), (v =1, 2, q)在可行点 x*处连续可微，)(xvg, (vV)线性无关，若 x*是不等式约束优化问题(3.49)的局部极小点， 则必存在一个实向量Tq,21， 使得下面各式成立： qvvvgf1)()(0xx (3.59) ), 2 , 1(0)(qvgvvx (3.60

52、) ), 2 , 1(0qvv (3.61) 这个条件称为库恩-塔克(Kuhn-Tucker)条件，简称为 K-T 条件，满足 K-T 条件的可行点称为 K-T 点。 K-T 条件中的前两个条件已在前面导出，下面用反证法证明第三个条件。 证明 设 dRn是 x*处的一个可行方向，即设 VvgvT0)(xd (3.62) 当Vv 时（V 表示在 x*处起作用约束的下标集合） ，0)(xvg，所以由式 ), 2 , 1(0)(qvgvvx，有 Vvv0 (3.63) 因此 K-T 条件的第一、二两个条件还可以写成下面的形式： Vvvvgf0xx)()( (3.64) 把式(3.64)两端左乘 dT

53、，得 VvvTvTVvvTvTgfgf)()(0)()(xdxdxdxd (3.65) 如果0v，即0v，那么由式(3.62) VvgvT0)(xd，有 VvvTvTgf0)()(xdxd (3.66) 式(3.62)和式(3.66)说明 d 既是 x*处的一个可行方向， 又是目标函数 f(x*) 在 x*处的下降方向。如果沿 d 方向移动微小距离，所得新的可行点的目标函数值将比 f(x*)小，这与 x*是局部极小点矛盾，所以必有 Vvv,0 (3.67) 把式(3.67)与式(3.63)综合起来，就得 K-T 条件的第三个条件。 从几何意义上看，在不等式约束问题的最优点处，目标函数的梯度向量

54、必须处在由起作用的约束曲面的负法向量所张成的凸多面锥中（图 3.3 a）。 使用式(3.64)和式(3.67)，可以把 K-T 条件写成 VvgfvVvvv0)()(0xx (3.68) g1(x)=0 g2(x)=0 x* )(1xg)(2xg)(xfg3(x)=0 )()1(3xg)()1(xf)()1(1xgx(1) 图 3.4 K-T 点与非 K-T 点 )()(3kgx)()(1kg x)()(2kgx)()(kf xx(k) )(3xg)(1xg)(2xg)(xfx* (a) x*是极小点 (b) x(k)不是极小点 图 3.3 K-T 条件的几何意义 的形式。即 目标函数在 x*

55、处的梯度为起作用的各约束函数负梯度的非负线性组合。 如果给定一个可行点 x(k)， 验证它是否为 K-T 点， 只需把该点代入方程(3.68)，看有无满足该方程的v存在就行了。如果欲求问题的K-T点， 就需要在可行域中求解由(3.59)和式(3.60)组合而成的方程组，方程组正好含有 n+q 个未知量和 n+q 个方程。 如果 f(x)和 gv(x) (v=1,2,q)都是凸函数（凸规划问题），那么满足 K-T 条件的点 x*就是最优点。对于非凸规划问题，局部极小点x*还必须满足下面的二阶必要条件。 3.3.2 二阶必要条件 定理 3.9 设函数 f(x)和 gv(x), (v =1, 2,

56、q)在可行点 x*处二次可微，若 x*是不等式约束优化问题(3.49)的局部极小点，则必存在一个实向量Tq,21，使得 K-T 条件成立；并且对于满足条件 )0,(, 0)()0,(, 0)(vvTvvTVvgVvgxdxd (3.69) 的任意向量 dRn，均有 0,2dxdxLT (3.70) 成立。式中 VvvvgfL)()(,222xxxx (3.71) 3.3.3 充分条件 定理 3.10 设函数 f(x)和 gv(x), (v =1, 2, q)在可行点 x*处二次可微，若存在一个实向量Tq,21，使得 K-T 条件成立；并且对于满足条件 )0,(, 0)()0,(, 0)(vvT

57、vvTVvgVvgxdxd (3.72) 的向量 dRn，如果 d = 0，或者虽然 d 0 但是 0),(2dxdxLT (3.73) 成立，则 x*是问题(3.49)的严格局部极小点。 例 3.4 验证点 x(1) = (0, 0)T 和 x(2) = (1, 1)T 是否为问题 0)(0)(s.t)2()(min21222112221xxgxxgxxfxxx 的 K-T 点。 解 目标函数和约束函数的梯度分别为： 11)(21)(2)2(2)(22121xxxgxgxxf 先验证点 x(1) = (0, 0)T。该点是可行点，起作用的约束函数下标集V = 1, 2。如果该点满足 K-T

58、条件，则有 0011010421 即 004221 解得 1 = 4，2 = 0。由于 1 0，所以 x(1)不是 K-T 点。 再验证点 x(2) = (1, 1)T。该点也是可行点，起作用的约束函数下标集 V= 1, 2。如果该点满足 K-T条件，则有 0011212221 即 022022121 解得 1= 0，2 = 2。由于两个拉格朗日乘子都非负，所以 x(2)是 K-T 点（图 3.5） 。 例 3.5 求解问题 x1 x2 x(1) x(2) 图 3.5 例 3.4 的两个可行点 0)(02)(s.t) 1()(min22211221xgxxgxxfxxx 解 目标函数和约束函数

59、的梯度分别为： 10)(11)(1) 1(2)(211xxxggxf 根据 K-T 条件 2 , 102 , 10)()()(21vvggfvvvvvvx0xx 有 0) 1(211x (1) 0121 (2) 0)2(211 xx (3) 022x (4) 0,21 (5) (1)至(4)是一个非线性方程组，问题归结为求这个方程组满足(5)和原约束不等式的解。一般说来，求解非线性方程组比较复杂，但本问题比较简单。 在(4)中，若 2 = 0，则由(2)得到 1 = 1 0，不满足(5)。因此只能有 x2 = 0，这时(3)成为 0)2(11x (6) 在(6)中，若 x1 2 = 0，则由(

60、1)得到 1 = 2 0，可见它满足 K-T 条件的(3.79)式，因此 x*是 K-T 点。 例 3.7 检验 x* = (1, 5)T 是否为约束优化问题 026)(01)(06)(s.t)(min2221211211221xxgxgxxhxxfxxxx 的 K-T 点。 解 把该点代入全部约束条件知道它是可行点，而且全部约束都起作用，V= 1, 2。目标函数和起作用的约束函数的梯度分别为： 21211122)(01)(11)(12)(xxgghxfxxxx 代入 x*处的 K-T 条件： 0101022010201111221211211 这是 3 个未知数 2 个方程的方程组，只要存在

61、一组满足 K-T 条件的 1, 1, 2，x*就是K-T 点。设 1 = 0，得，2 = 1/10 0，1 = 22/10 0 所以 x*是 K-T 点。 这个例子说明，K-T 点所对应的拉格朗日乘子的值不一定是唯一的。在本例中，只要1-11/4, 1，就能满足 K-T 条件。 另外，h1(x)和 g1(x)是线性函数，因而是凸函数；由 2002)(0002)(222xxgf 可知 f(x)和 g2(x)也是凸函数，所以原问题是凸规划问题。根据一阶充要条件，x* = (1, 5)T就是其最优点。 例 3.8 求解约束优化问题 0)3(3)(010)3()(s.t)(min2211222111x

62、xgxxhxfxxx 解 目标函数和约束函数的梯度与海赛矩阵分别为： 1) 3(6)(2) 3(2)(01)(11211xgxxhfxxx 0006)(2002)()(12122xx0xghf 根据 K-T 条件 00)3(3)(001)3(62)3(201)()()(1221111112111111xxgxxxghfxxxx 和等式约束条件，得下列方程组 0) 3(6) 3(211111xx (1) 02121x (2) 0) 3( 32211xx (3) 010) 3(2221xx (4) 当 10 时，根据式(3)，有 3)3(221xx (5) 把(5)代入(4)，得 0303222

63、xx (6) 从中解得 x2 = 10/3 或 x2 = - 3。 把 x2 = 10/3 代入(5)， 无解； 把 x2 = - 3 代入(5)， 得 x(1) = (2, -3)T 和 x(2) = (4, -3)T 两点。这两点都是可行点。 把 x(1) = (2, -3)T代入(1),(2)，得1 = 1/38，1 = -3/19。由于10，所以 x(2)是 K-T 点。 当 1 = 0 时，式(1)成为 0)3(2111x (7) 可知10。这时根据(2)，得 x2 = 0。把 x2 = 0 代入(4)，得 10) 3(21x (8) 解之得1031x，因此得TT)0,103(,)0

64、,103()4()3(xx两点。这两点都是可行点。 把 x(3)， 1 = 0 代入(1), 得10211。 可知x(3)是 K-T 点。 把 x(4)，1 = 0 代入(1), 得 10211。可知x(4)也是 K-T 点。 下面用二阶条件判别 x(2), x(3), x(4)三个 K-T 点是否是局部极小点。 1、对 x(2)点，先求二阶条件中的向量 d = (d1, d2)T。由于两个约束都起作用，而且 10，因此由 0616,)(06262,)(2121)2(12121)2(1ddddgddddhTTxdxd 得 d1= d2= 0，即 d = 0。因此 x(2) 满足二阶充分条件，它

65、是严格局部极小点。 2、对 x(3)点，只有等式约束是起作用的约束，根据 01020102,)(121)3(1dddhTxd 得 d1= 0，即 d= (0, d2)T，d2为任意实数。这时有 0100101001010)()(),(2222)3(121)3(211)3(2dddhfLTTdxxddxdx 它不满足二阶必要条件，因此点 x(3)不是局部极小点。实际上，它是一个极大点。 3、对于 x(4)点，只有等式约束是起作用的约束。在该点处 01010100101)()(),(22212121)4(121)4(211)4(2ddddddhfLTTdxxddxdx 即对任何实向量 d0，拉格朗

66、日函数在 x(4)处关于 x 的海赛矩阵都是正定的，因此该点满足二阶充分条件，x(4)是严格局部极小点。 x1 x2 x(1) x(2) 图 3.7 例 3.8 的 4 个可行点 x(4) x(3) 在两个严格局部极小点中，因为4)(103)()2()4(xxff，所以T)0,103()4(x是该问题的最优点。 第 4 章 优化问题解法概述 优化设计问题的求解方法：解析法、图解法和数值解法。 4.1 解析法和图解法 解析法根据优化问题的极值条件，通过目标函数和约束函数的一阶导数和二阶导数的分析计算来求得问题的最优点。 实际的优化设计问题：设计变量多、约束条件多，非线性。 求解复杂、困难、不可能

67、。不实用。 图解法：简便、直观、结果一目了然。但是只能求解两个设计变量的问题， 精确绘出目标函数的等值线和约束函数的约束曲线也是很困难的，所以图解法也不实用。 x1 x2 x(1) x(2) 例 3.8 的附图 x(4) x(3) 数值解法：依靠电子计算机的大量数值计算来求得优化问题的最优解。 4.2 数值解法概述 4.2.1 数值解法的基本思想 先在设计空间中选取一个初始设计点 x(0)。 然后按某种规则对设计方案进行改进，得到一个新的设计点x(1)。 如果新的设计点尚未令人满意，则对新的设计方案进行再改进。如此反复，直到获得满意的设计点为止。 设计方案的改进过程将在设计空间中产生一个点列：

68、 x(k), (k = 0,1,2,) 如果对设计方案进行改进的规则是合理的，那么点列x(k)将收敛于优化问题的极小点 x*，即 g1(x)=0 图 1.11 的最优点 x2 0 x1 x* d(2) x(k) d(1) f (x)=ci 图 1.12 的最优点位置 g1(x)=0 x2 g2(x)=0 f (x)=ci 0 x1 x* (b) xx)(limkk (4.1) 对设计方案进行改进的规则称为算法。 最常见的算法是迭代下降算法。 迭代：从一个点 x(k)求得一个后继点 x(k+1)后，把计算公式中的 k用 k+1 替代，形成计算循环。 下降： 每次计算循环都让目标函数值有所减小，

69、从而保证点列x(k)收敛于局部极小点。 数值计算法的求解过程是一个逐渐逼近极小点 x*的过程，通常把这种逼近称为对 x*的“搜索”。 x1 x2 0 x(2) x(0) x(1) x* 图 4.2 数值解法的搜索过程示意图 常用的非线性规划的迭代公式为 ), 2 , 1 , 0()()()1(kkkkkSxx (4.2) S(k)：第 k+1 次迭代计算时的搜索方向； k：第 k+1 次迭代计算时的步长因子； kS(k)： 是 x(k)与 x(k+1)之间的距离， 即从 x(k)到 x(k+1)的搜索步长。 由 x(k)出发沿方向 S(k)求步长因子k的过程叫做一维搜索， 它是求解最优化问题的

70、基本步骤之一。 4.2.2 算法的收敛性 1、如果某种算法构造出的点列x(k)能够在有限步内得到优化问题的极小点 x*，那么这种算法便称为是收敛的算法。 如果只有当 x(0)充分接近极小点 x*时，算法产生的点列x(k)才收敛于 x*，则称该算法为具有局部收敛性的算法。 0 x2 x3 x1 x(k) x(k+1) S(k) )(kkS图 4.3 从一个设计点搜索到 另一个设计点的几何含义 x2 0 x* x* x(0) x(0图 4.4 多个极小点的情形 如果对于可行域内的任何初始点 x(0)，算法产生的点列x(k)都收敛于 x*，则称该算法为具有全局收敛性的算法。 2、一个好的算法必须以较

71、快的速度收敛到极小点 x*。 算法收敛速度的快慢，用算法的收敛级来衡量。 定义 4.1 设点列x(k)收敛于 x*，如果对于某个实数 p1，有 0lim)()1(pkkkxxxx (4.3) 成立，那么 p 就称为点列x(k)的收敛级，或者称点列x(k)是 p 级收敛的。 收敛级 p 越大，点列x(k)收敛得越快。当收敛级 p 相同时，收敛比 越小，点列x(k)收敛得越快。 3、一个函数在某点附近可以用一个二次函数来近似。因此，一个算法如果对于二次目标函数的收敛速度快， 那么对一般的目标函数也会有较快的收敛速度。 评价一个算法的收敛速度，一般看它用在二次目标函数时的表现如何。 反过来，为了建立

72、一个好的算法，往往先针对二次目标函数建立一个收敛速度快的算法，然后把它推广到一般目标函数的问题上去。 4.2.3 迭代的终止准则 每一种算法，都要给出终止迭代过程的准则，以便迭代进行到一定程度时能够停止计算。常用的终止准则有 3 种： 1、相邻两点的距离充分小，即 )()1(kkxx (4.4) 或者 )()()1(kkkxxx (4.5) 2、相邻两点的目标函数值的下降量充分小，即 )()() 1()(kkffxx (4.6) 或者 )()()()()1()(kkkfffxxx (4.7) 3、在无约束优化问题中，最新点的梯度充分接近于零，即 )()(kf x (4.8) 4.2.4 梯度的

73、数值计算 用数值解法求解优化问题的时候，一般并不使用梯度的解析式来计算某点的梯度值，而是通过对该点邻近的函数值的计算，获得它的近似值。 很小时， 一元函数 f(x) 在点 x(k) 处的一阶导数可以用它的一阶差商来近似。 向前差商公式： hxfhxfhxfhxfdxxdfkkkkhk)()()()()()()()()(0)(lim (4.9) h ：差分步长 向后差商公式： hhxfxfhhxfxfdxxdfkkkkhk)()()()()()()()()(0)(lim(4.10) 中心差商公式： hhxfhxfhhxfxfhxfhxfdxxdfkkkkkkk2)()()()()()(21)()

74、()()()()()()( (4.11) n 元函数 f(x)在点 x(k)处的一阶偏导数也可以用它的一阶向前差商公式 ), 2 , 1()()()()()()(nihfhfxfkikikxexx (4.12) 或一阶向后差商公式 ), 2 , 1()()()()()()(nihhffxfikkikexxx (4.13) 或中心差商公式 ), 2 , 1(2)()()()()()(nihhfhfxfikikikexexx (4.14) 来近似计算。式中 ei是第 i 个分量为 1，其余分量为 0 的 n 维列向量，即 Ti00100e (4.15) 向前差商近似计算 n 元函数 f(x) 在点

75、 x(k)处的梯度的程序框图。图中的 ri 是梯度的第 i 个分量。 框图中的符号“:=”表示把它右边表达式的值赋给它左边的变量。在一般程序设计语言（C, Fortran, Basic）的中，它被直接写成等号“=” 。本书程序框图中的该含义全部使用“:=” 。 图 4.6 梯度的数值计算框图 i := i+1 yes no 开 始 给定 x(k), h, n i := 1 y1 := f(x(k) 结 束 hyyrcxfyhxxxcikikkikiki/)(:)(:12)()(2)()()(xTnkrrrf,:)(21)(xi = n ? 第 5 章 一维搜索 5.1 一维搜索的概念 1、当点

76、 x(k)与搜索方向 S(k)确定后，从 x(k)出发，沿 S(k)方向的射线上任意点 xs的计算公式为 )()(kksSxx (5.1) 这里， 是它的步长因子。 当 取不同数值时，点 xs在射线上移动，其目标函数值 f(xs)也发生变化。 目标函数值最小的那一个点所对应的步长因子称为最佳步长因子。 x(k) S(k) )(kkS图 5.1 一维搜索示意图 xs )(kSx(k+1) )()1( kf x 所谓一维搜索，就是在已知 x(k)与 S(k)的情况下，寻找最佳步长因子k ，使迭代产生的新点 x(k+1)处的目标函数值为最小。 一维搜索就是求以 为变量的一维优化问题 1)()() 1

77、()()(minRSxxkkkff (5.2) 的最优解 k 。而最佳步长因子即为 )(min|)()(1kkkfSxR (5.3) 2、因为 x(k+1)是搜索方向 S(k)上目标函数值最小的点，所以它是向量 S(k)与目标函数等值面的切点，即搜索方向 S(k) 与 x(k+1)处的梯度向量 f(x(k+1)垂直： 0)() 1()(kTkf xS (5.4) 3、如果目标函数 f(x)是一个二次函数，那么它在 x(k)处的二阶泰勒展开式是一个等式，即 x(k) S(k) )(kkS图 5.1 一维搜索示意图 xs )(kSx(k+1) )()1( kf x )()(2)(2)()()()(

78、)(2)()()()()()()1()()(21)()()()(21)()()()()(kkTkkTkkkkTkkTkkkkkffffffffSxSSxxSxSSxxSxx (5.5) 这时，问题(5.2)的最优解 k 可根据它的一阶必要条件求出： )()(2)()()()()(2)()()()1()()()(0)()()()(kkTkkTkkkkTkkkTkkfffffSxSSxSxSSxx (5.6) 例 5.1 试对目标函数 2322212121)(xxxfx 求在 x(k) = (1, 1, 1)T 处沿方向 S(k) = (-1, -2, 0)T 的最佳步长因子 k。 解 f(x)是

79、二次函数，因此最佳步长因子k可直接用式(5.6)求出。由 100010002)(112)(2)(2)(321xxxffxxxfk 得 3242402102222021100010002021021112k 4、对于非二次函数，式(5.6)不成立，所以式(5.6)并没有实用价值。 实际求取最佳步长因子 k 时， 使用一维优化问题的数值解法。 5、由于每次迭代都要进行一维搜索，所以一维搜索是数值计算法最基本的内容，它在程序运行中不仅占用的机时较多，而且对一些算法的计算效果有重要影响。 维搜索一般分两步进行。第一步确定k所在的大致区间，该区间称为搜索区间；第二步把包含k的搜索区间不断缩小，最终找出k

80、 。 5.2 搜索区间的确定 确定搜索区间：在 数轴上，寻找 a, b 两点，保证 k 落在a, b 之间。 假定已经找到一个搜索区间a, b，并且其中只有一个局部极小点k ， 那么目标函数)()()(kkfSx在 a, b 内应呈 “高-低-高” 形状。 a(1) b(3) k 2 f(x(k)+ S(k) 0 图 5.2 只有一个极小点的搜索区间 每向前试探一步，就把 1, 2, 3 三点依次向前移动一步。直到目标函数在 1, 3区间内呈“高-低-高”形状为止。 如果试探过程中所用步长h为负数， 那么最后的 1, 3 是沿 的反方向试探获得的，这时 31，局部极小点 k 位于 3, 1 之

81、中。 这种方法依靠步长 h 的正负来控制试探的方向，既可进又可退，因此称为进退法。 例 5.2 试对目标函数 2322212121)(xxxfx 求在x(k) = (1, 1, 1)T处沿方向S(k) = (-1, -2, 0)T进行一维搜索时的搜索区间a, b。（初始步长h = 0.2） 解 搜索从 x(k)处开始，因此给定0 = 0，h = 0.2。 24321)21 (21)1 ()(021111222)()()()(kkkkfSxSx 按图 5.3 的计算步骤，先计算1, 2处的目标函数值： 1 =0 =0 y1=f (x(k)+1S(k) ) = 2 2 =1+ h = 0.2 y2

82、=f (x(k)+2S(k) ) = 30.22 - 40.2 + 2 = 1.32 因为 y2 y1，所以接着计算 3 =2 + h = 0.4 y3=f (x(k)+3S(k) ) = 30.42 - 40.4 + 2 =0.88 因为 y3 y2，所以让 h =2h = 0.4 1 =2 = 0.2 y1 = y2 = 1.32 2 =3 = 0.4 y2 = y3 = 0.88 并继续向前试探： 3 =2 + h = 0.4+0.4 = 0.8 y3 = f (x(k)+3S(k) ) = 0.72 因为 y3 0 所以得结果： a, b= 1,3 = 0.4, 1.6。 5.3 区间

83、消去法原理 对于搜索区间a, b内的任意两点 a1, b1， a1 y2? |b-a| ? yes 令 xR1，x(k) = 0，S (k) = 1， = x，则 f(x(k) + S(k) = f(x)。先用进退法求搜索区间：取 0 = 0，h = 0.2，搜索 4 次获得搜索区间 a, b = - 1.4, - 0.6。然后按图(5.7)所示的步骤在a, b中找最优点。 表 5.1 为本问题的计算过程及结果。本问题的解析解为 x* = - 1，与数值解结果 x* = - 1.0 相同。 表 5.1 例 5.3 计算过程及结果 次 数 a a1 b1 b | b-a| y1 比较 y2 0

84、-1.4 -1.0944 -0.9056 -0.6 -0.9911 = -0.9911 1 -1.4 -1.2111 -1.0944 -0.9056 0.49 -0.9554 -0.9911 2 -1.2111 -1.0944 -1.0223 -0.9056 0.31 -0.9911 -0.9995 3 -1.0944 -1.0223 -0.9056 0.190.2 结 果 x* = k = 0 .5 (a+b) = -1.0 例 5.4 试用黄金分割法求目标函数 2322212121)(xxxfx 在 x(k) = 1, 1, 1T处，沿方向 S(k) = - 1, - 2, 0T的最佳步长

85、因子 k。 解 例 5.2 已求得了本问题的搜索区间 a, b = 0.4, 1.6，并导出了从 x(k) = (1, 1, 1)T出发，沿方向 S(k) = (- 1, - 2, 0)T上各点的目标函数值的计算公式 243)(2)()(kkfSx 表 5.2 是计算过程及结果（ = 0.05） 。数值解结果为 k = 0.67，与解析解 k = 2/3 相同。 表 5.2 例 5.4 的计算过程及结果 次 数 a a1 b1 b | b-a| y1 比较 y2 0 0.4 0.85840 1.14160 1.6 0.77695 1.34335 1 0.4 0.68329 0.85840 1.

86、14160 0.74 0.66750 0.66750 3 0.57511 0.68329 0.75018 0.85840 0.28 0.66750 0.66750 5 0.64200 0.68329 0.70886 0.75018 0.17 0.66750 0.67201 6 0.64200 0.66754 0.68329 0.70886 0.07 0.66667 0.66750 7 0.64200 0.66754 0.68329 0.04 0 ?”的条件不成立，说明 p 点落在a, b之外。这两种情况都是计算机的舍入误差造成的，只能在a, b区间收缩得非常小、a, c, b 三点很接近时才会

87、发生。这时计算可以结束，区间中的 c 点即为极小点。 yes no 结 束 k := c yes 图 5.10 二次插值法程序框图 no yes no 开 始 给定 a, b, y1 := f(x(k)+ aS(k) y3 := f(x(k)+ bS(k) c := 0.5 (a + b) y2 := f(x(k)+ cS(k) yp := f(x(k)+ pS(k) C1 := (y3 y1)/( b a) C2 := (y2 y1)/( c a) - C1 / ( c b) p := 0.5(a+b-C1/C2) k := p no yes b := p y3:= yp a := c y1

88、:= y2 c := p y2:= yp b :=y3:=c :=y2:=no yes yes ye (p-a)(b-p)0 ? |c-p| ? y2 yp ? p c ? yp c y2 yp a b c p y y p c y2 yp a b c p y y p c y2 yp a b c p y2 y p c y2 yp 例 5.5 试用二次插值法求目标函数 2322212121)(xxxfx 在 x(k) = (1, 1, 1)T处，沿方向 S(k) = (-1, -2, 0)T的最佳步长因子k。 解 根据例5.2求得的搜索区间a, b = 0.4, 1.6和目标函数值的计算公式 24

89、3)(2)()(kkfSx 令 = 0.0001，用图 5.13 的步骤计算，其过程及结果列于表 5.3 中。数值解结果为k = 0.6667，与解析解k = 2/3 相同。 表 5.3 例 5.5 的计算过程及结果 次 数 a c (p) p (c) b | c-p| y2 比较 yp 1 0.4 p: 0.66667 c: 1.00000 1.6 0.3333 1.000000 0.666667 0.4 p: 0.66667 c: 0.66667 1.00000 0.0000 0.666667 0.666667 | c-p| ，故结束计算。 结 果 因 y2yp，故 k = c = 0.6

90、667 例 5.6 试用二次插值法求解一维优化问题 64s.tsin1)(minxxxxf 解 对于任意 x4, 6，有 sinx -0.788494 2 4 p: 4.7543 c: 4.7807 5.0000 0.0264 -0.788494 -0.788786 3 4 c: 4.7543 p: 4.7563 4.7807 0.0020 -0.788786 -0.788789 4.7543 c: 4.7563 p: 4.7566 4.7807 0.0004 -0.7887886 -0.7887887 | c-p| ，故结束计算。 结 果 因 yp f(x(k)的情况，这时的牛顿法是不收敛的

91、。所以牛顿法只具有局部收敛性，不具备全局收敛性；不是一个实用的算法。 修正牛顿法在牛顿法的迭代过程中增加一维搜索的步骤： )()()(1)(2)(kkkffxxS (6.52) ), 2 , 1 , 0()()() 1(kkkkkSxx (6.53) 6.3.2 变尺度法 一、变尺度法的基本思想和迭代公式 基本思想： 用一个矩阵 H(k)来代替修正牛顿法中的1)(2)(kf x， 在迭代过程中， H(k)逐渐逼近1)(2)(kf x。 迭代公式： )()()()(kkkf xHS (6.54) ), 2 , 1 , 0()()() 1(kkkkkSxx (6.55) 矩阵 H(k)在迭代过程中

92、不断变化，称为变尺度矩阵。 二、变尺度矩阵 H(k)的计算公式 1、为了保证搜索方向S(k)是使目标函数值下降的方向，要求 S(k)与)()(kf x之间的夹角小于 90，即 0)()()(kTkf xS (6.56) 把式(6.54)代入(6.56)，得 0)()()()()(kTkTkffxHx (6.57) 当 H(k)是对称矩阵时，H(k)T= H(k)，这时 0)()()()()(kkTkffxHx (6.58) 所以 H(k)应当是一个对称正定矩阵。 2、为了保证 H(k)在迭代过程中逐渐逼近1)(2)(kf x，H(k)必须满足拟牛顿条件。 设在第 k 次迭代后，得到点 x(k+

93、1)。由于目标函数 f(x)在点 x(k+1)的二阶泰勒近似式为 )()(21)()()()()1()1(2)1()1()1()1(kkTkkTkkffffxxxxxxxxxx (6.59) 所以在点 x(k+1)附近，f(x)的梯度为 )()()() 1() 1(2) 1(kkkfffxxxxx (6.60) 令 x = x(k)，则 )()()() 1()() 1(2) 1()(kkkkkfffxxxxx (6.61) 记 )() 1()()()1()(),()(kkkkkkffxxxxxg (6.62) 于是有 )() 1(2)()(kkkfxxg (6.63) )()(1)1(2)(k

94、kkfxgx (6.64) 这样， 在获得 x(k)和 x(k+1)点以及目标函数在这两点处的梯度后， 便可根据式(6.64)估计在 x(k+1)处的海赛矩阵的逆。因此，如果 H(k+1)满足 )()() 1(kkkxgH (6.65) 那么它就能在迭代过程中逼近矩阵1) 1(2)(kf x。该式称为拟牛顿条件。 3、H(k)的递推公式 为了适应迭代计算的需要，可以把 H(k)构造成递推的形式： ), 2, 1, 0(,)()() 1(kkkkHHH (6.66) )(kH：第 k 次迭代时的校正矩阵。H(0) = I。假设校正矩阵取下列形式 TkkkTkkkkvvuuH)( (6.67) 其

95、中 uk、vk是 n 维待定列向量，k、k是待定常数。)(kH是一个形式最简单的待定对称正定矩阵。 把式(6.66)和(6.67)代入(6.65)的拟牛顿条件： )()()()()()()()()()()()()()()(kkkkTkkkkTkkkkkkkkkkkkgHxgvvguugHxgHxgHH (6.68) 满足这个等式的待定向量有多种取法，不妨取 )()()()()(kkkTkkkkkTkkkgHgvvxguu (6.69) 注意到)()(,kTkkTkgvgu都是数量，因此可把待定向量取为 )()()(,kkkkkgHvxu (6.70) 把它们代入(6.69)，得待定常数 )()

96、()()()(1,1kkTkkkTkkgHggx (6.71) 从而获得一个校正矩阵公式： ), 2, 1 , 0(,)()()()()()()()()()()()(kkkTkkTkkkkTkTkkTkkkTkkkkgHgHggHgxxxvvuuH (6.72) 因为初始变尺度矩阵 H(0)和校正矩阵)(kH都是对称正定矩阵， 所以式(6.66)的 H(k+1)也是对称正定矩阵，这就保证每次迭代都能使目标函数值有所下降。 使用式(6.66)和(6.72)计算变尺度矩阵的变尺度法称为 DFP 法，（1959 年由 Davidon 首先提出，1963 年被 Fletcher 和 Powell 加以

97、改进）。 可以证明，对于二次凸函数 cfTTxbAxxx21)(， DFP 法构造出的搜索方向是一组 A 共轭方向，经过 n 次迭代一定能够达到极小点。所以 DFP 法也是一种共轭方向法，它具有二次收敛性。对于一般凸函数，DFP 法是全局收敛的。 三、变尺度法的迭代步骤 (1) 任取初始点 x(0)，给定允许误差 ； (2) 计算)()0(xf，如果)()0(xf，得 x* = x(0)，停止迭代，否则转(3)； (3) 令 k = 0，取 H(0) = I； (4) 计算搜索方向 )()()()(kkkf xHS。 (5) 一维搜索，求得 x(k+1) = x(k) + kS(k)； (6)

98、 如果)()1(kf x，得 x* = x(k+1)，停止迭代；否则转(7)； 图 6.5 变尺度法程序框图 x* := x(k+1) 开 始 给定 x(0), 结 束 no no yes yes k := k+1 )()()1()()()()()(:)(min|:)(:1kkkkkkkkkkffSxxSxxHSRIH:0:)0(k?)()1(kf x)()()1(:kkkHHH?)()0(xfx* := x(0) no yes )(:)(:)1()0()1()0(kkffxxxxk = n-1 ? (7) 若 k = n-1，令 x(0) = x(k+1)，)()()1()0(kffxx，转

99、(3)；否则使用式(6.65) )和(6.71)计算下次迭代需用的变尺度矩阵H(k+1)， 令k = k+1，转(4)。 四、其它变尺度法 满足拟牛顿条件的校正矩阵不是唯一的，选择不同的校正矩阵就得到不同的变尺度法。 比较著名而且实用效果很好的变尺度法， 除 DFP 法外， 还有 BFGS法。BFGS 法是由 Broyden，Fletcher，Goldfarb 和 Shanno 在 1970 年提出的。它的迭代步骤及程序框图与 DFP 法的完全相同，但是使用的校正矩阵为 ), 2, 1, 0(,1)()()()()()()()()()()()()()()()()()(kkTkkTkkTkkkk

100、TkTkkkTkkkTkkgxHgxxgHgxxxgxgHgH (6.73) BFGS 法的变尺度矩阵 H(k+1)也满足拟牛顿条件， 与 DFP 法相比具有更好的计算稳定性。BFGS 法是目前公认的最好的拟牛顿算法。 例 6.5 用 DFP 法求解下面问题 22212)(minxxfx 初始点 x(0) = (1, 1)T ， = 1/10 。 解 本题求解对象与例 6.1 相同。 2124)(xxf x 2)(22)(1)(2)(2)(1)(1)()(22kkkkkkkssxsxs 按程序框图 6.5，其迭代过程如下： 开始时 k = 0， Tf24)()0(x 因为 10/152)2()

101、4()(22)0(xf 所以需进行第一次迭代，这时 24241001)()0()0()0(Tf xHS 185)2()4(21)2(1)4(2)()(22222)0(22)0(1)0(2)0(2)0(1)0(10ssxsxs 41912418511)0(0)0()1(Sxx 因为 101954)(2194)()1()1(xxff 而且 k = 0 n-1 = 2-1 = 1，所以继续迭代，计算变尺度矩阵： 14910242194)()(1295114191)0()1()0()0()1()0(xxgxxxff 3053838863061144161711224181100114100114100

102、11414100114121212211001)0()0()0()0()0()0()0()0()0()0()0()0()1(gHgHggHgxxxHHTTTT 当 k =1 时， 4117421943053838863061)()1()1()1(xHSf 3617)4() 1 (24)4() 1(1241791)()(22222)1(22)1(1)1(2)1(2)1(1)1(11ssxsxs 004117436174191)1(1)1()2(Sxx 这时， 1010)(00)()2()2(xxff 已满足精度要求，因此停止迭代，输出数值解： 00)2(xx 该问题的解析解也是 x*=0, 0T

103、。 本例题验证了对二次凸函数， DFP法经 n 次迭代就可得到极小点的论断。事实上，由 04141617441200424174)1()0(ASST 可知，例题两次迭代的搜索方向是 A 共轭的，DFP 法也是一种共轭方向法。 6.4 坐标轮换法和鲍威尔法 6.4.1 坐标轮换法 坐标轮换法非常简单，其搜索方向就是坐标轴的方向 e1 , e2 , , en，迭代公式为 k = n-1 ? ?)0()1(xxkx* := x(k+1) 图 6.6 坐标轮换法程序框图 开 始 给定 x(0), 结 束 no no yes k := 0 yes k := k+1 x(0) := x(k+1) S(k)

104、 := ek+1 x(k+1) := x(k)+kS(k) )(min|)()(1kkkfSxR ) 1, 2, 1, 0()()()1(1)(nkkkkkkkSxxeS (6.76) 迭代步骤为 (1) 任取初始点 x(0)，给定允许误差 ； (2) 令 k = 0; (3) 以 ek+1为搜索方向作一维搜索，求得 x(k+1) = x(k) + k ek+1； (4) 若 k = n-1，转(5)；否则令 k = k+1，转(3)。 (5) 如果)0()1(xxk，得 x* = x(k+1)，停止迭代；否则令 x(0) = x(k+1)，转(2)； 坐标轮换法的锯齿现象有时很严重，收敛速度

105、慢。 如果目标函数的等值线存在脊线，那么当搜索点 x(k)来到脊线上后，由于坐标轮换法不能提供使目标函数值下降的搜索方向，迭代计算将无法继续进行。 图 6.8 等值线存在脊线的情况 x1 x2 0 x(k) 脊线 图 6.7 坐标轮换法的搜索过程 x1 x2 0 e1 e2 e1 e2 x(b) x(a) 6.4.2 鲍威尔基本算法 如果沿着 e1、e2作一维搜索获得点 x(b)，再沿着 e1、e2作一维搜索获得点 x(a)之后，接着沿 x(a)- x(b)的方向进行搜索，收敛速度就可能加快。 对于 n 维二次凸函数 cfTTxbAxxx21)( (6.77) 存在如下定理： 定理 6.3 对

106、于式(6.77)的二次凸函数，设 S(0), S(1), , S(n-1)关于 A共轭，如果从两个不同的点出发，分别依次沿 S(0), S(1), , S(n-1)作一维搜索， 最后得到不同的点x(a)和点x(b)， 则向量 x(a)- x(b) 与S(0), S(1), , S(n-1)关于 A 共轭，即 ) 1, 2, 1, 0(0)()()(nkbaTkxxAS (6.78) 证明 根据定理的条件和式(6.26) ) 1, 2, 1, 0(0)()()(nkfnTkxS (6.26) 有 ) 1, 2, 1, 0(0)() 1, 2, 1, 0(0)()()()()(nkfnkfbTka

107、TkxSxS (6.79) 两式相减，得 ) 1, 2, 1, 0(0)()()()()(nkffbaTkxxS (6.80) 由 bAxxf，知 )()()()()()(babaffxxAxx (6.81) 所以 ) 1, 2, 1, 0(0)()()(nkbaTkxxAS (6.82) 这个定理表明，从两个不同的点出发，分别沿相同的方向作一维搜索，也可以产生关于 A 的共轭向量组。 如果图 6.7 中的目标函数是二次凸函数，那么 x(a)- x(b) 与 e2是 A共轭的，因此沿 x(a)- x(b) 的方向进行一维搜索，收敛速度就能加快。 鲍威尔法的迭代步骤： (1) 任取初始点x(0)

108、， 给定允许误差 ， 令S(k) = ek+1 , (k = 0, 1, 2, , n-1) ； 图 6.7 坐标轮换法的搜索过程 x1 x2 0 -e1 e2 e1 e2 x(b) x(a) (2) 依次沿 S(k)方向作一维搜索，求得 x(k+1) = x(k) + k S(k) , (k = 0, 1, 2, , n-1)； (3) 令 S(n) = x(n)- x(0)， 沿 S(n)方向作一维搜索， 求得 x(n+1) = x(n) + nS(n)； (4) 如果)0()1(xxn，得 x* = x(n+1)，停止迭代；否则令 x(0) = x(n+1)，S(k) = S(k+1)

109、, (k = 0, 1, 2, , n-1)，转(2)； 对于二次凸函数，鲍威尔法各阶段的始点和终点所确定的向量是A 共轭的，因此至多经过 n 个阶段的迭代就可以求得极小点。 对于非二次目标函数，一般经 n 个阶段的迭代达不到极小点，这时可将所得的点作为初始点继续搜索下去，直到满足精度要求为止。 图 6.9 鲍威尔法的搜索过程 (上标中第一位数字表示阶段数) x1 x2 0 x(1.0) x(1.2) x(1.1) x(2.3) x(2.0)= x(1.3) x(2.1) S(1.1) S(1.0) S(1.2) S(2.0) S(2.1) x(2.2) S(2-2) 6.4.3 鲍威尔改进算

110、法 鲍威尔基本算法产生的 n 个向量有可能线性相关或者近似线性相关，这时张不成 n 维空间，得不到真正的极小点。 在鲍威尔改进算法中，鲍威尔提出了新的方向替换原则：在一个阶段的迭代结束后，对 S(0), S(1), , S(n)等 n+1 个方向进行比较，选出n 个最接近共轭的方向作为下一阶段的搜索方向。 具体做法：计算函数值 f* = f( 2x(n)- x(0) ) 和各相邻点的函数值之差 k = fk - fk+1 = f(x(k)- f(x(k+1), (k = 0, 1, 2, , n-1)，找出其中最大的一个，即找出 m = maxk；如果 f*- f0 0 与 20200*22f

111、ffffffmmnn 同时成立，则构造新的搜索方向组，否则下一阶段仍使用原来的搜索方向。 鲍威尔改进算法虽然不再具有二次收敛性，但是确实消除了搜索方向线性相关的现象，它是一个比较有效的直接算法。 鲍威尔改进算法的迭代步骤为： (1) 任取初始点 x(0)， 给定允许误差 ， 令 f0 = f(x(0)， S(k) = ek+1 , ( k = 0, 1, 2, , n-1) ； (2) 依次沿 S(k)方向作一维搜索，求得 x(k+1) = x(k) + k S(k) , fk+1 = f(x(k+1)， (k = 0, 1, 2, , n-1)； ?)0()( xxnx* := x(n) 图

112、 6.10 鲍威尔改进算法程序框图 开 始 给定 x(0), k = n-1 ? 结 束 no no yes k := 0 yes k := k+1 S(n) := x(n) - x(0) x(0) := x(n) +nS(n) )(min|:)()(1nnnfSxR x(k+1) := x(k) +kS(k) fk+1:= f(x(k+1) )(min|:)()(1kkkfSxRf0 := f(x(0) S(k) := ek+1 f* := f(2x(n) - x(0) ) x(0) := x(n) f := f(x(n) 是否满足 判别条件? no yes (3) 如果)0()(xxn，得

113、 x* = x(n)，停止迭代；否则令)(max1101kknkmmmffff，f* = f(2x(n) - x(0) )； (4) 如果 f*- f0 0，并且20200*22fffffffmmnn，则转(5)；否则搜索方向组不变，令 x(0) = x(n)，f0 = f(x(n)，转(2)； (5) 令 S (n) = x(n)- x(0) ， 沿 S(n)方向作一维搜索， 求得 x(0) = x(n) + nS(n)， f0 = f(x(0)； (6) 令 S(k) = S(k+1) , (k = m, m+1, , n-1)，转(2)。 例 6.6 用鲍威尔改进算法求解下面问题 121

114、222122123)(minxxxxxfx 初始点 x(0) = (2, 2)T ， = 1/100 。 解 迭代过程中的最佳步长因子暂且使用下面的解析式计算： 1113)(23)(21221xxfxxxxf )(2)(1)(2)(1)(2)(1)(1)(2)(2)(1)()(2)()()(111323)()(kkkkkkkkkkkkTkkTkkssssssxxxxffSxSSx 2)(2)(2)(12)(1)(1)(2)(2)(2)(1)(1)(2)(323kkkkkkkkkkssssxxsxxs 按程序框图，其迭代过程如下： 开始时，f0 = f(x(0) = 0，S(0) = e1 =

115、(1, 0)T，S(1) = e2 = (0, 1)T 第 1 阶段 32)(2)(3232)0(2)0(2)0(12)0(1)0(1)0(2)0(2)0(2)0(1)0(10ssssxxsxxs 32)(3232013222)1(1)0(0)0()1(xSxxff 32)(2)(3232)1(2)1(2)1(12)1(1)1(1)1(2)1(2)1(2)1(1)1(11ssssxxsxxs 98)(113410323232)2(2)1(1)1()2(xSxxff 因为 10012321121134)0()2(xx 所以需进行第 2 阶段迭代，这时 98)2(*113211211382)0()

116、2()0()2(xxxxff 329832,320max,max21100ffffm 因为 2022020202002200*243128729643298989162*22*222431289832*098*ffffffffffffffffmmnnm 所以需构造新的搜索方向组。令 11321121134)0()2()2(xxS 21)(2)(3232)2(2)2(2)2(12)2(1)2(1)2(2)2(2)2(2)2(1)2(12ssssxxsxxs 得第 2 阶段的初始点和搜索方向： 1)(111132211134)0(0)2(2)2()0(xSxxff 113210)2()1()1()

117、0(SSSS 第 2 阶段 0)(2)(3232)0(2)0(2)0(12)0(1)0(1)0(2)0(2)0(2)0(1)0(10ssssxxsxxs 11)0(0)0()1(Sxx 0)()(3232)1(2)1(2)1(12)1(1)1(1)1(2)1(2)1(2)1(1)1(11ssssxxsxxs 11)1(1)1()2(Sxx 因为 100101111)0()2(xx 所以得x* = x(2) = (1, 1)T，结束迭代计算。 第 7 章 约束优化方法（线性规划） 当优化问题 qvgnpuhfvun, 2 , 10)(, 2 , 10)(s.t)(minxxRxx (7.1) 的

118、目标函数 f(x) 和所有的约束函数 hu(x)、gv(x) 都是线性函数时，称为线性规划问题。因为线性函数是凸函数，其可行域是凸集，所以线性规划问题是凸规划问题。 线性规划是最优化理论中方法最成熟、应用最广泛的一个分支。 7.1 线性规划的几何含义 例 7.1 用图解法求解线性规划问题 0)(0)(01223)(042)(s.t.)(min241321221121xgxgxxgxxgxxfxxxxx 解 可行域是一个凸多边形 ABCD。 目标函数的等值线是一个平行直线族，这些互相平行的等值线沿着其负梯度方向 Tf) 1, 1 ()(x 移动时，目标函数值是下降的。 显然，在可行域中，经过 A

119、 点的目标函数等值线具有最小值，所以本题的最优点就是 A 点。A 点是直线 g1(x) = 0 与 g2(x) = 0 的交点，其坐标为(2, 3)，因此 x* = (2, 3)T。 对于 n 维的线性规划问题来说， 约束函数和目标函数都是线性函数，其可行域一般是 n 维空间中由超平面围成的凸多面体，目标函x1 x2 g1(x)=0 g2(x)=0 g3(x)=0 g4(x)=0 f(x)= -5 A B C D 图 7.1 线性规划问题的图解法 )(xf (2, 3)T (0, 2)T (4, 0)T (0, 0)T 数的等值面则是一个平行超平面族。如果它存在唯一最优点，那么它也可以在凸多面

120、体的某个顶点上找到。 线性规划问题还可能出现三种特殊情形： 1、有多重解，2、目标函数无下界（无最优解），3、无解（可行域为空集）。 无论是有唯一解还是有多重解，只要线性规划问题存在最优解，就可以在凸多面体的顶点中找到它（或者找到它们当中的一个）。因此，寻找线性规划问题的最优点，没有必要在整个可行域中搜索，只需在凸多面体的有限个顶点中搜索即可。 7.2 线性规划的标准形式 线性规划问题的标准形式称为线性规划的标准型： nnxcxcxcf2211)(minx (7.2) 图 7.2 线性规划的三种特殊情形 x1 x2 0 f(x) = c x1 x2 0 a) 多重解 b) 无最优解 c) 无可

121、行解 x1 x2 0 )(xf)(xf mnmnmmnnnnbxaxaxabxaxaxabxaxaxa22112222212111212111s.t. (7.3) ), 2 , 1(0njxj (7.4) 其中 aij，cj， bi为常数，bi0, ( i =1, 2, , m，j =1, 2, , n)。 一般称 cj 为价格系数，称式(7.3)为约束方程组，称式(7.4)为非负约束。线性规划的标准型也可写为 ), 2 , 1(0), 2 , 1(s.t.)(min11njxmibxaxcfjinjjijnjjjx (7.5) 若记 c = (c1, c2, , cn)T， b = (b1,

122、 b2, , bm)T， mnmmnnnmijaaaaaaaaaa212222111211)(A， 线性规划的标准型可以写成矩阵形式： 0xbAxxcxs.t.)(minTf (7.6) 矩阵 A 称为约束矩阵，b 称为右端向量。 若记 A = p1, p 2, , p n，其中 pj = (a1j, a2 j, , am j)T，线性规划的标准型还可表示成 0xbpxcxjnjjTxf1s.t.)(min (7.7) 当原问题不是标准形式时，总能用下列方法化成标准形式： (1) 若第 i 个等式约束中 bi0，则用(- 1)乘该等式的两端。 (2) 若第 j 个变量 xj没有非负限制（这时称

123、 xj为自由变量），则引入两个非负变量0, jjxx，令 jjjxxx 把它代入目标函数和约束条件中去。 (3) 若第 i 个约束是不等式约束 injjijbxa1， 则增加一个非负松弛变量 xn+i0，把它化为等式约束 iinnjjijbxxa1， 与之类似，若第 i 个约束是不等式约束 injjijbxa1， 则增加一个非负松弛变量 xn+i0，把它化为 iinnjjijbxxa1 这时把目标函数中与 xn+i 对应的价格系数 cn+i 取为 0。 例 7.2 把下面的线性规划问题化为标准型 0,52327s.t.32)(min21321321321321xxxxxxxxxxxxxxf x

124、 解 增加两个非负松弛变量 x4, x50，并令333xxx ，原问题就化为标准型： 0,5002232070s.t.00332)(min543321543321543321543321543321 xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxf x 7.3 线性规划的基本概念和基本定理 设约束矩阵 A = (p1, p 2, , p n)的秩为 m，且 mn，则 A 中必存在非奇异的 m 阶子矩阵 B。为记述方便，不妨设 B 由 A 的前面 m 个列向量组成（可以通过列调换做到这一点），即 B = p1, p 2, , pm 方阵 B 称为线性规划问题的一个基矩阵。 基矩阵中

125、包含的列向量 p1, p 2, , p m称为基向量，与基向量对应的变量 x1, x 2, , x m称为基变量。 其余的列向量 pm+1, p m+2, , p n称为非基向量，与之对应的变量 xm+1, x m+2, , x n称为非基变量，非基向量构成的矩阵 N = pm+1, p m+2, , pn 称为非基矩阵。 在约束方程组 mnmnmmnnnnbxaxaxabxaxaxabxaxaxa22112222212111212111s.t. 中取定一个基矩阵 B 之后，令相应的非基变量全部为 0，约束方程组就成为 p1x1 + p 2x2 + p 3x3 + + pmxm = b 这时方

126、程个数与变量个数都是 m，方程组有唯一解。线性规划的这种解称为与基矩阵 B 对应的基解。因为从约束矩阵 A 中最多可取出mnC个基矩阵，所以线性规划最多有mnC个基解。基解中的非零分量不一定全部非负，满足非负约束条件的那部分基解称为基可行解。显然，基可行解的个数不会超过 )!( !mnmnCmn。 基可行解的重要性质： 定理 7.1 x(k)为标准型线性规划问题的基可行解的充分必要条件是 x(k)为可行域的顶点。 定理 7.2 若标准型线性规划问题存在最优解， 则必有基可行解是最优解。 如果线性规划问题有最优解，那么只需从构成可行域的凸多面体的顶点中寻找即可。 例 7.3 对例 7.1 的线性

127、规划问题 0)(0)(01223)(042)(s.t.)(min241321221121xgxgxxgxxgxxfxxxxx 求其基可行解和最优解。 解 先把问题化为标准型： 0,122342s.t.)(min432142132121xxxxxxxxxxxxfx 然后从约束矩阵 10230121A 中取不同的基矩阵 B(k), ( k = 1, 2, )，求它的全部基解。对 x1 x2 g1(x)=0 g2(x)=0 g3(x)=0 g4(x)=0 f(x)= -5 A B C D 图 7.1 线性规划问题的图解法 )(xf (2, 3)T (0, 2)T (4, 0)T (0, 0)T 23

128、2121)1(ppB， 解方程组 1242321)1(2)1(1xx， 得基解 TTxx0, 0, 3, 20, 0,) 1 (2) 1 (1) 1 (x 对 1001,1202,0212,1301,031143)6(42)5(32)4(41)3(31)2(ppBppBppBppBppB 等基矩阵，用同样方法可求得与之对应的基解： TTTTTTTTTxxxxxxxxxx12, 4, 0, 0, 0, 0,8, 0, 2, 0, 0, 0,0, 8, 6, 00, 024, 0, 0, 4, 0, 0,0, 8, 0, 40, 0,)6(4)6(3)6()5(4)5(2)5()4(3)4(2)4

129、()3(4)3(1)3()2(3)2(1)2(xxxxx其中x(1), x(2), x(5), x(6)满足非负约束，因此是基可行解，它们正好与图 7.1 中的 A, D, B, C 各点相对应。 因为 f(x(1) = min f(x(1), f(x(2), f(x(5), f(x(6) = min-5, -4, -2, 0 ，所以对该问题的标准型，最优解为 x*= x(1) = (2, 3, 0, 0)T；对原问题，最优解为 x*= (2, 3)T。 当 n = 30，m = 15 时，基解的个数可能达到1530C1.55108之多。因此尽管基可行解的个数是有限的，仍然需要寻找一种算法，在

130、基可行解中有选择性地进行搜索，尽快地获得最优解。 1947 年，美国数学家丹茨格（Dantzig）提出了著名的求解线性规划问题的单纯形法（Simplex method）。单纯形法构造的迭代格式，在使目标函数值有所下降的前提下， 从一个基可行解求得另一个基可行解，因而能够较快地搜寻到最优解。 单纯形法经过五十多年的发展，形成许多变种。它使用方便、行之有效，容易掌握，是求解线性规划问题的权威性算法。 7.4 单纯形法 7.4.1 线性规划的规范式0,122342s.t.)(min432142132121xxxxxxxxxxxxfx 对于线性规划的标准形式 0xbAxxcxs.t.)(minTf (

131、7.8) 把 约 束 矩 阵A用m m基 矩 阵B(k)和m (n-m) 非 基 矩 阵N(k)分 解 为314212104110032011210120132Txxxx BNAbcxxBN A = B(k) N(k) 把列向量 x 和 c 也作相应的分解： )()()()(kNkBkNkBcccxxx，3142113142112212101240132121041012011201322412432123212xxxxxxxxxxxxxx 这时约束方程组 Ax = b 可表示为 bxxNB)()()()(kNkBkk 即 bxNxB)()()()(kNkkBk 解得 )()(1)(1)()(

132、)()(kNkkkkBxNBbBx 把它代入目标函数，有 )()(1)()()(1)()()()()()(1)(1)()()()()()()()()()()(kNkkTkBTkNkTkBkNTkNkNkkkTkBkNTkNkBTkBTfxNBccbBcxcxNBbBcxcxcxcx 若令0,)23(12)2(4s.t.)(min432121421321xxxxxxxxxxxxfx )()()2()()1()(1)()2(1)()1(1)()()(1)()()(2)(1)(1)()()()()()(2)(11)()(1)()()()()(,)(,)(,)(kmnkmmkmmknkmkmmnmki

133、jkkTknkmkmkkTkBTkNkNTkmkkkkkTkBkaaaaaaabbbfNBNBccbBbbBc 可把原问题表示成0,)23(12)2(4s.t.)(min432121421321xxxxxxxxxxxxfx njxmixabxxffjnmjjkijkiinmjjkjk, 2, 1, 0, 2, 1,s.t.)(min1)()(1)()(x (7.10) 线性规划的这种形式称为与基变量 x1, x 2, , x m对应的规范式。 如果基变量不是由前 m 个变量组成，上述过程同样可以进行。若用 S(k)表示基变量的下标集，T(k)表示非基变量的下标集，可写出与任意一组基变量对应的规

134、范式：0,)23(12)2(4s.t.)(min432121421321xxxxxxxxxxxxfx )()()()(minkTjjkjkxffx (7.11) )()()(,s.t.)(kTjjkijkiiSixabxk (7.12) njxj, 2, 1, 0 (7.13) 规范式中的)(kj称为变量 xj的判别数。 7.4.2 基可行解的转换0,)23(12)2(4s.t.)(min432121421321xxxxxxxxxxxxfx 在线性代数中，线性规划的式 )()()(,)(kTjjkijkiiSixabxk 是线性方程组的通解，非基变量 xj, ( jT(k)是通解中的自由变量。

135、令非基变量 xj为0，可得到它的一个特解，这就是与基矩阵 B(k)对应的线性规划的基解： Tknkkkxxx)()(2)(1)(,x0,)23(12)2(4s.t.)(min432121421321xxxxxxxxxxxxfx )()()()(, 0,kkkjkjTjSjbx 当0)(kjb时，x(k)是一个基可行解，这时的目标函数值 )()()()()()()(kTjkjkjkkfxffkx (7.14) 假设已经获得一个与基矩阵 B(k)对应的基可行解 x(k)，这时如果让一个非基向量 pr, ( rT(k)进入 B(k)，同时让一个基向量 pl, ( lS(k)离开 B(k)，就可以生成

136、一个新的基矩阵 B(k+1)，从而获得与 B(k+1)对应的另一个基解 x(k+1)。 为了保证x(k+1) 仍然是一个基可行解，并且保证目标函数值 f(x(k+1)较 f(x(k)有所下降，必须对进入 B(k) 的非基向量 pr和离开 B(k)的基向量 pl进行严格挑选。 进入 B(k)的非基向量 pr 称为进基向量，对应的非基变量 xr 称为进基变量；离开 B(k)的基向量 pl 称为离基向量，对应的基变量 xl 称为离基变量。 (1) 选择进基变量 xr 的方法 由式(7.11)： )()()()(minkTjjkjkxffx 得 x(k)和 x(k+1)的目标函数值0,)23(12)2

137、(4s.t.)(min432121421321xxxxxxxxxxxxfx )()()1()()()1()()()()()()(kkTjkjkjkkTjkjkjkkxffxffxx (7.15) 因为在 x(k)中，全部非基变量为零；在 x(k+1)中，除 xr 之外的非基变量为零，所以 0,)23(12)2(4s.t.)(min432121421321xxxxxxxxxxxxfx ) 1()()() 1()()()()(krkrkkkkxffffxx 为了使 f(x(k+1) f(x(k) 必须有 0)(kr 因此，为了使新的基可行解的目标函数值有所下降，xr 的判别数应当是负数。 通常选其

138、判别数为 0min)()()(kjkjkr (7.20) 的 xr为进基变量。 (2) 选择离基变量 xl 的方法 对于新的基可行解 x(k+1)，根据 )()()(,)(kTjjkijkiiSixabxk0,)23(12)2(4s.t.)(min432121421321xxxxxxxxxxxxfx rjTjxTjxkkjkkj, 0, 0)()1()()( 有 )()1()()()1()()()1(,)(kkrkirkiTjkjkijkikiSixabxabxk (7.21) 因为xl, ( lS(k)是 x(k)的离基变量，所以它在变为 x(k+1)的非基变量后，应该取0： )()1()(

139、)1()()()1(0klkrklrkrklrklklbxaxabx (7.22) 由于) 1( krx和)(klb都是非负的， 所以0)(klra。 这就是说， 离基变量的下标 l 必须在满足 )()(, 0kkirSia (7.23) 的那些下标 i 当中进行选择。0,)23(12)2(4s.t.)(min432121421321xxxxxxxxxxxxfx 设0)(klra，得 xr 在成为 x(k+1)的基变量之后所取的值： )()()1(klrklkrabx (7.24) 为了保证 x(k+1)的全部基变量非负，根据式(7.21)，应该有 liSixabxkkrkirkiki, 0)

140、() 1()()() 1( (7.25) 当0)(kira时，式(7.25)是成立的。当0)(kira时，要让式(7.25)成立，必须有 liSiabxkkirkikr,)()()()1( (7.26) 由式(7.24)和(7.26)，得0,)23(12)2(4s.t.)(min432121421321xxxxxxxxxxxxfx liSiababkkirkiklrkl,)()()()()( (7.27) 显然，为使式(7.27)成立，只要取 0min)()()()()()(kirkirkiSiklrklaababk (7.28) 即可。 根据该式选择离基变量 xl 的下标 l， 就能保证新的

141、基解 x(k+1)是一个基可行解。 7.4.3 最优解的判定准则 单纯形法借助判别数)(kj来判断当前的基可行解 x(k)是否是最优解。 定理 7.3 设 x(k)是线性规划的一个基可行解，在其规范式中，若 0min)()(kjTjk，则 x(k)是线性规划的唯一最优解；若 0min)()(kjTjk，则 x(k)是线性规划的多重解中的一个最优解。 定理 7.4 在线性规划的规范式中，若存在判别数0)(kr，且全部0)(kira, ( iS(k)，则线性规划无最优解。 证明 以 xr作为进基变量，求取新的基可行解 x(k+1)时必须满足的条件为 )()(, 0kkirSia 在本定理中，0)(

142、kira，因此必须取0)(kira。由 liSixabxkkrkirkiki, 0)() 1()()() 1( 可知，这时xr取任何正数，都能保证 x(k+1)是基可行解。若取 )1(krx，则有 ) 1()()() 1()(krkrkkxffx 即目标函数值在可行域上无界，线性规划无最优解。 7.4.4 单纯形法的计算公式 7.4.5 单纯形法的迭代步骤 单纯形法的迭代步骤如下： (1) 给出初始基可行解 x(0) ，即给出基变量与非基变量的下标集合 S(0), T(0)，以及系数)0()0(,ijiab, (i S(0), j T(0)，计算判别数)0(j, ( j T(0)； (2) 若

143、 0min)0()0(jr，则 x(0)就是最优解，停止计算。否则选择进基变量 xr，令 k = 0，转(3)； (3) 若)( , 0)()(kkirSia，则线性规划无最优解，停止计算。否则转(4)； (4) 根据0min)()()()()()(kirkirkiSiklrklaababk选择离基变量 xl，令 S(k+1) = S(k)+r-l， T(k+1) = T(k) - r+l； (5) 用式(7.45)和(7.47)计算)1()1(,kjkib。 (6) 若0min) 1() 1(kjkr，则获得了最优解 x* = x(k+1)，结束计算。否则选择进基变量 xr，并用式(7.46

144、)计算) 1( kija，令 k = k+1，转(3)。 作为计算结果，最优解 x* = x(k+1)的各分量为 )1()1()1()1(, 0,kkkikiTiSibx 7.4.6 表格形式的单纯形法 若令 ijSiijSiaSjkkkijkkj, 0, 1, 0)()()()()( 则可把线性规划的规范式(7.11)(7.13)： )()()()(minkTjjkjkxffx njxSixabxjkTjjkijkiik, 2 , 1, 0,s.t.)()()()( 0,)23(12)2(4s.t.0)(min432121421321xxxxxxxxxxxxfx 表示为0,12023402s

145、.t.0)(min43214321432121xxxxxxxxxxxxxxfx njxSibxaxffjkkinjjkijnjjkjk, 2, 1, 0,s.t.)(min)()(1)(1)()(x 引入一个变量 f ，式(7.49)进一步化为下面的等价形式： 0,120234020s.t.min43214321432121xxxxxxxxxxxxxxff njxSibxaxfffjkkinjjkijnjjkjk, 2, 1, 0,s.t.min)()(1)(1)()( (7.50) 把式(7.50)的约束方程组整理为 0,1202304020000s.t.min432143214321432

146、1xxxxxxxxfxxxxfxxxxff )()(1)()(1)(,0kkinjjkijknjjkjSibxaffxf (7.51) 然后把方程组的增广系数矩阵置于一张表中，就得到如表 7.1 那样的单纯形表。 用单纯形法求解线性规划问题，实际上就是求解线性方程组(7.51)。只是在对它的增广系数矩阵做行的初等变换时，要按一定的规则选择自由未知量（非基变量），以保证所得的基解是基可行解 （0)(kjb） ， 并保证每次所得基可行解的目标函数值要比上次的小 （f (k) 0，如果原问题存在可行解xx的话，辅助问题将存在可行解0xxxa, ，其目标函数值00eT，将比最优值aTxe还小，这是不可

147、能的。 (2) 0ax， 而且基变量中已不存在人工变量。 这时可以把xx 作为原问题的初始基可行解， 因为此时 Abx ， 说明xx 是原问题的基可行解。 (3) 0ax，但是基变量中还存在人工变量。这时xx 虽然是原问题的一个可行解， 但还不是基可行解， 因此需要继续进行迭代计算，直到基变量全部变成 非人工变量 x 的分量为止。这时的迭代仍然用主元消去法， 但是由于此时与基变量 xai 对应的判别数和右端向量的分量 bi 都是零，所以主元的选取无须遵守式(7.52)和(7.53)的准则，只要能把人工变量转换成非基变量就行。 当基变量中已不存在人工变量时，原问题的初始基可行解就得到了。 例 7

148、.9 用单纯形表求解问题 0,312s.t.2)(min211212121xxxxxxxxxfx 解 先把问题化为标准型： 0,312s.t.2)(min543215142132121xxxxxxxxxxxxxxxfx 由于标准型的约束矩阵并不包含 3 阶单位矩阵，因此还要引进人工变量 x6和 x7，构造一个辅助问题： 7 , 1, 0312s.t.)(min517421632176jxxxxxxxxxxxxxfjax3176421732162312xxxxxxxxxxxxx 第一阶段用单纯形法解辅助问题。求解之前需要先把辅助目标函数化为规范式的形式。由于此时的基矩阵 B(k)是一个单位矩阵，

149、因此根据式(7.9) )(1)()()()(1)()()()(,)(kkTkBTkNkNkTkBkfNBccbBc 有 )()()()()()(,kTkBTkNkNTkBkfNccbc (7.58) 由上式，得 3001000131100101) 1 () 1 (2010011123000110) 2(07654321jixxxxxxxk21011010211001011/111011) 2 (0) 2/ 1 (120011) 2(017654321jixxxxxxxk2/ 32/ 12/ 112/ 12/ 100/2/ 32/ 12/ 102/ 12/ 101/2/ 12/ 12/ 102/

150、 12/ 110/0110000027654321jixxxxxxxk表 7.7a 第一阶段的单纯形表 11020001101101110113312011)0(4)0(3)0(2)0(1)0(5)0(7)0(6)0(4)0(3)0(2)0(1)0(4)0(3)0(2)0(1)0(5)0(7)0(6)0(5)0(7)0(6)0(ppppcccccccbbbcccfa 这样就可写出辅助问题的初始单纯形表。表 7.7a 给出了它的迭代情况，迭代两次获最优解。这时人工变量 x6和 x7都等于零，并且都是非基变量，因此得原问题的初始基可行解： x(0) = (3/2, 1/2, 0, 0, 3/2)T

151、 在进行第二阶段计算之前， 先把表中与人工变量所对应的列删除，然后把表中的系数代入式(7.58) )()()()()()(,kTkBTkNkNTkBkfNccbc 43214324312321252212121212123xxxxxxxxxx 求原问题 0,312s.t.2)(min543215142132121xxxxxxxxxxxxxxxfx 的目标函数值和判别数： 2/312/ 12/ 10032/302/ 12/ 101/2/ 10) 2/ 1 (2/ 110) 1 (2/50) 2/3(2/ 100054321jixxxxxk110) 1 (10) 1 (200111/101120/

152、400) 2(30154321jixxxxxk110110/310001/211010/620010254321jixxxxxk表 7.7b 第二阶段的单纯形表 23212/12/12/12/12/12/1021252/32/32/1021)0(4)0(3)0(5)0(1)0(2)0(4)0(3)0(4)0(3)0(5)0(1)0(2)0(5)0(1)0(2)0(ppcccccbbbcccf 表 7.7b 给出了原问题标准型的初始单纯形表和它的迭代过程。 最后获得原问题的最优解 x* = (3, 0)T和最优值 f(x*) = -6。 例 7.10 用单纯形表求解问题 0,044422s.t.

153、)(min3213132132121xxxxxxxxxxxxxfx 解 引入松弛变量 x4，把问题化为标准型，再引入人工变量 x5和x6，得辅助问题： 6 , 1, 0044422s.t.)(min6315321432165jxxxxxxxxxxxxxxfjax 321653163215243400444xxxxxxxxxxxx 由式(7.58)写出辅助目标函数的规范形式： 321652434)(xxxxxfax 辅助问题只迭代一次就获得最优解。这时人工变量 x5和 x6虽然等于零，但是仍然是基变量，还需把它们变成非基变量才能进入第二阶段的计算。这时主元的选取与j和i无关，只要它位于人工变量所

154、在的行，原来变量所在的列即可（见表 7.8a 中 k = 1, 2 时的主元） 。当人工变量 x5和 x6都变成非基变量时，就得到原问题标准型的一个初始基可行解： 0100101/4010144120011) 2 (1) 1 (40002) 4(30654321jixxxxxxk010010) 1 (/0012302/1002/ 12/ 112/ 1/00024011654321jixxxxxxk0100101/0212) 5(00/12/ 102/ 1010/2654321jixxxxxxk05/ 35/ 15/ 2001/05/ 25/ 15/ 2100/12/ 102/ 1010/365

155、4321jixxxxxxk 表 7.8a 第一阶段的单纯形表 x(0) = (0, 1, 0, 0)T 在进行第二阶段计算之前，先把表中人工变量的列删除，并由式(7.58) 求原问题的目标函数值和判别数： 1015/25/22/11011001101)0(4)0(1)0(3)0(2)0(4)0(4)0(1)0(3)0(2)0(1)0(3)0(2)0(pccccbbbcccf 4214241101121152xxxxxxx 这时唯一的非基变量的判别数大于零， 所以初始基可行解已经是原问题标准型的唯一最优解（表 7.8b ） ，即原问题的最优解 x*= (0, 1, 0)T 和 最优值 f(x*)

156、 = -1。 表 7.8b 第二阶段的单纯形表 05/2001/05/2100/12/1010/110/100004321jixxxxk 第 8 章 约束优化方法（非线性问题） 对于非线性约束优化问题 qvgnpuhfvun, 2 , 10)(, 2 , 10)(s.t)(minxxRxx (8.1) 其数值解法可以分为两类： 间接解法：把原问题化为一系列无约束优化问题（惩罚函数法、增广乘子法等）。 直接解法：在可行域内逐步改善设计点（随机方向法、复合形法等）。 8.1 惩罚函数法 8.1.1 外罚函数法（外点法） 对于等式约束问题 npuhfun, 2 , 10)(s.t)(minxRxx

157、(8.2) 可以定义一个辅助函数 puuhff12)()(),(xxx (8.3) 构造一个无约束的辅助问题 puuhff12)()(),(minxxx (8.4) 其中 参数 是一个充分大的正数。 辅助问题的极小点x 必使得 hu(x ) 接近于零，否则它的第二项将是很大的正数，x 就不是极小点。hu(x )接近于零说明该点近似满足约束条件，这时辅助函数 f(x, ) 近似等于 原目标函数 f(x)。 因此，通过求解辅助问题(8.4)，能够获得原问题的一个近似解。辅助问题对原问题的近似程度，与参数 的大小有关。 越大，近似程度越高。当 +时，hu(x ) 0（否则 f(x , ) +），辅助

158、问题的极小点x 逼近于原问题的极小点 x*。 对于不等式约束问题 qvgfvn, 2 , 10)(s.t)(minxRxx (8.5) 构造一个无约束的辅助问题 qvvgff12)(, 0max)(),(minxxx (8.6) 其中参数 也是一个充分大的正数。当 x 为可行点时（gv(x) 0） ， 0)(, 0maxxvg 当 x 为不可行点时（gv(x) 0） ， )()(, 0maxxxvvgg 所以辅助问题的极小点将近似于原问题的极小点。当 时，辅助函数的第二项逼近于零，辅助问题的极小点逼近于原问题的极小点。 对于一般的约束优化问题，综合上面等式约束和不等式约束两种情况，无约束的辅助

159、问题为 )()(),(minxxxPff (8.7) 其中 是一个充分大的正数，P(x)为 qvvpuughP1212)(, 0max)()(xxx (8.8) 当 x 为可行点时，P(x) = 0，从而有 f(x, ) = f(x)；当 x 为不可行点时， P(x)是一个很大的正数，它的存在是对 x 不在可行域内的一种惩罚。 通常把 P(x) 称为惩罚项， 把 称为罚因子，把P(x)称为罚函数，把 f(x, )称为增广目标函数。 例 8.1 求解下列非线性约束优化问题 01)(s.t.)(min112221xgxxfxx 解 定义增广目标函数 01,) 1(01,)(, 0max),(121

160、222112221212221xxxxxxxgxxfxx 当011x时，令 02, 022211xxfxxf 得点(0, 0)T，该点不满足不等式011x，不是它的极小点。 当011x时，令 02, 0) 1(2222111xxfxxxf 由此得辅助问题的极小点和极小值： Txx0,121x 1111) 1(), (22212221xxxfx x 是原问题的近似解， 越大，x 越接近原问题的极小点 x*。当 + 时，x x*，因此 1) (lim)(,01limxxxxff。 进行数值计算的时候， 不可能取无穷大，只能取有限的正数。这时，上例原问题的约束函数 01111) (1xg 不满足约束

161、条件，点x 落在可行域之外。因此P(x) 称为外罚函数，称为外罚函数法或外点法。 为了求得足够精度的x ，需要把罚因子 取得足够大才行。一般先取不太大的罚因子，求得一个近似程度不高的极小点，然后将近似极小点作为后续计算的初始点，形成计算循环。 在计算循环中，罚因子逐步增大，最终使近似极小点的近似程度达到给定的精度要求。 可以证明，如果取一个趋向无穷大的严格递增的正数序列k ： kkklim,0321 从某个正数 1 开始，对每个 k 求解辅助问题 )()(),(minxxxPffkk (8.9) 能够得到一个极小点的序列)(kx，那么这个序列将收敛于原约束问题的极小点 x*。 通过求解一系列无

162、约束问题来得到约束问题极小点的方法，称为系列无约束极小化方法，简称 SUMT(Sequential Unconstrained Minimization Technique)方法。因此外罚函数法又称为 SUMT 外点法。 外罚函数法的迭代步骤如下： (1) 给定初始点 x(0)、初始罚因子1、罚因子放大系数 C 1 和允许误差 ，置 k = 1； (2) 以 x(k-1)为无约束问题(8.9)的初始点， 调用无约束优化程序， 求解该辅助问题，获其极小点 x(k)； (3) 如果 k = kP(x(k) 1。 例 8.2 用外点法求解下列非线性约束优化问题 001s.t.) 1(31)(min2

163、1231xxxxfx 解 构造增广目标函数 2221231, 0max1, 0max) 1(31),(xxxxfkkx 取初始点 x(0) = (- 4, - 4)T，初始罚因子 1= 1，罚因子放大系数 C =10，允许误差 = 0.001，求解无约束问题)(minxf，迭代情况如表 8.1。 第 5 次迭代结束时满足迭代终止条件， 因此得极小点 x* x(5) = (1.000, 0.000)T。 表 8.1 外点法迭代情况 k k x1 (k) x2 (k) f (x(k) k =kP(x(k) 1 1 0.2361 - 0.5000 0.1295 0.8336 2 10 0.8322

164、- 0.0500 2.0001 0.3067 3 102 0.9804 - 0.0050 2.5840 0.0410 4 103 0.9980 - 0.0005 2.6582 0.0042 5 104 0.9998 - 0.0001 2.6658 0.0004 8.1.2 内罚函数法（内点法） 外点法把迭代点定义在可行域之外， 当迭代点落在约束边界上时，惩罚项为零，从而“诱逼”迭代点（外点）朝着约束边界移动。 内点法与外点法相反，它把迭代点定义在可行域之内，当迭代点落在约束边界上时，让惩罚项为无穷大，从而“阻截”迭代点越过约束边界。 由于内点法的迭代点始终是内点， 所以只适用于不等式约束问题。

165、对于问题 qvgfvn, 2 , 10)(s.t)(minxRxx (8.10) 内点法构造一个增广目标函数 )()(),(xxxBff (8.11) 把原问题转化为无约束的辅助问题 )()(),(minxxxBff (8.12) 其中 是一个足够小的正数，B(x)为 qvvqvvgBgB11)(ln)(,)(1)(xxxx或 (8.13) 因为当 x 从可行域内部趋近约束边界时， B(x) +， 从而使增广目标函数 f(x, ) +，所以辅助问题(8.12)的极小点不会越过约束边界。由于 B(x)是防止辅助问题的极小点穿越约束边界而设置的障碍，通常把 B(x) 称为障碍项，把 称为障碍因子或

166、罚因子，把B(x)称为障碍函数或内罚函数。 障碍因子 非常小时，辅助函数 f(x, ) 近似等于原目标函数 f(x)。辅助问题的极小点是原问题的一个近似解。辅助问题对原问题的近似程度，与 的大小有关。 越小，近似程度越高。 为了求得足够精度的近似极小点，采取系列无约束极小化方法，取一个严格递减而且趋于零的正数序列k： 0lim,321kkk 从某个正数1开始，对每个k求解辅助问题 )()(),(minxxxBffkk (8.14) 可以证明，由此获得的极小点序列)(kx必将收敛于原约束问题的极小点 x*。 内点法的迭代步骤如下： (1) 给定初始内点 x(0)、 初始障碍因子1、 障碍因子缩小

167、系数 C (0, 1)和允许误差 ，置 k =1； (2) 以 x(k-1)为无约束问题(8.14)的初始点，调用无约束优化程序，求解该辅助问题，获其极小点 x(k)； (3) 如果 k = kB(x(k) ，则得原问题的近似极小点 x* = x(k)，停止计算；否则令 k+1 = Ck，置 k = k+1，转(2)。 内点法的初始点 x(0)必须是内点，最好离约束边界远一些。一般取 1=1，C = 0.1 0.7 1。 例 8.3 用内点法求解例 8.2 的约束优化问题 001s.t.) 1(31)(min21231xxxxfx 解 定义增广目标函数 21231111) 1(31),(xxx

168、xfx 先用解析法求解),(minxf。设x 是辅助问题的极小点，则有 01), (, 0) 1() 1(), (22221211xxfxxxfxx 由此得辅助问题的极小点和极小值： 3211131), (,1xxfxx x 是原问题的近似解， 越小，x 越接近原问题的极小点 x*。当 0 时，x x*，因此 38), (lim)(,01lim00xxxxff。 现在用系列无约束极小化方法求其数值解。其增广目标函数为 21231111) 1(31),(xxxxfkkx 取初始内点 x(0) = (4, 4)T，初始障碍因子1 = 1，障碍因子缩小系数 C = 0.1，允许误差 = 0.001，

169、得表 8.2 的迭代结果。迭代到第 8 次时满足迭代终止条件，因此得极小点 x* x(8) = (1.000, 0.000)T。 表 8.2 内点法迭代情况 k k x1 (k) x2 (k) f (x(k) kB(x(k) 1 1 1.4142 1.0000 5.6904 3.4142 2 10-1 1.1473 0.3162 3.6164 0.9953 3 10-2 1.0488 0.1000 2.9667 0.3049 4 10-3 1.0157 0.0316 2.7615 0.0954 5 10-4 1.0050 0.0100 2.6967 0.0300 6 10-5 1.0016 0

170、.0032 2.6762 0.0095 7 10-6 1.0005 0.0010 2.6697 0.0030 8 10-7 1.0002 0.0003 2.6676 0.0009 8.1.3 混合罚函数法（混合法） 混合罚函数法能够处理既有不等式约束又有等式约束的优化问题。 对于一般的约束优化问题 qvgnpuhfvun, 2 , 10)(, 2 , 10)(s.t)(minxxRxx 混合罚函数法构造的系列无约束极小化问题为 0lim,)()(1)(),(min321kkkkkkBPffxxxx (8.15) 其中 qvvpuugBhP112)(1)(, )()(xxxx (8.16) 混合

171、罚函数法具有内点法的求解特点， 迭代过程在可行域内进行，计算步骤和罚因子的选取都与内点法类似。 (1) 给定初始内点 x(0)、 初始罚因子 1、 罚因子缩小系数 C (0, 1)和允许误差 ，置 k = 1； (2) 以x(k-1)为无约束问题(8.15)的初始点，调用无约束优化程序，求解该辅助问题，获其极小点 x(k)； (3) 如果 k =)()(xxBPkk， 则得原问题的近似极小点x* = x(k)，停止计算；否则令 k+1 = Ck，置 k = k+1，转(2)。 外点法、内点法和混合法的程序框图可以用图 8.1 统一表示，但三种方法使用的计算式不同。 例 8.4 用混合法求解约束

172、问题 01052s.t.3) 1()(min1222221xxxxxfx 解 定义增广目标函数 1)52(3) 1(),(12222221xxxxxfkkkx 取初始内点 x(0) = (0, 0)T，初始罚因子1 = 1，罚因子缩小系数 C = 0.1，允许误差 = 0.002，迭代到第 8 次时满足迭代终止条件，得极小点 x*x(8) = (1.00, 2.50)T。迭代情况如表 8.3 所示。 表 8.3 混合法迭代情况 k k x1 (k) x2 (k) f (x(k) k 1 1 1.7937 3.8333 - 25.5645 8.3710 2 10-1 1.3684 2.8433

173、- 16.4792 1.7628 3 10-2 1.1710 2.6026 - 14.5518 0.4793 4 10-3 1.0794 2.5319 - 13.9997 0.1411 5 10-4 1.0368 2.5100 - 13.8289 0.0429 6 10-5 1.0171 2.5032 - 13.7750 0.0133 7 10-6 1.0079 2.5010 - 13.7579 0.0041 8 10-7 1.0036 2.5003 - 13.7525 0.0013 1、 收敛速度判别数 (0, 1) 和允许误差 ， 置 k = 1； (2) 以 x(k-1) 为无约束问题(

174、8.26)的初始点，调用无约束优化程序求解，获其极小点 x(k)； (3) 如果k，则得原问题的近似极小点 x* = x(k)，停止计算；否则用式(8.26)计算乘子向量)1( k，转(4)。 (4) 如果1kk，则令 k+1 = k ，转(5)；否则令 k+1 = C k ，转(5)。 (5) 置 k = k+1，转(2)。 在没有其它信息的情况下，初始乘子向量)1(可以取为单位向量，初始罚因子1一般可取为 1，罚因子放大系数 C 可取为 10，收敛速度判别数 可取为 0.25。 例 8.6 用增广乘子法求解约束优化问题 012)(s.t.4)(min2122121xxhxxxfxx 解 构

175、造增广拉格朗日函数 22121)(22121)() 12(2) 12(4),(xxxxxxxLkkkkx 设 x(k)是),()(kkLx的极小点，则有 0) 12(20) 12(2242)(2)(1)()(22)(2)(1)()(11kkkkkkkkkkxxxxLxxxxL 由此得 Tkkkkkkkkxx233,3224)()()(2)(1)(x 把 x(k)代入乘子迭代公式(8.25) ), 2, 1(),()()() 1(puhkukkukux 有 123332)24(2)()()()1(kkkkkkkk 得 kkkkk326322)()1( 当 k 4/3 时 )(k收敛，而且 k 越

176、大收敛越快。设 k 取大于 4/3 的某个固定数值，并且让)(k，即对上式取极限： kkk326322 得2。这时 x(k)x*，因此本问题的极小点为 ()( )03324244lim,123322332kTTkkkkkkk xx 现在用系列无约束极小化方法求其数值解。取初始点 x(0) = (0, 0)T，初始乘子1)1(，初始罚因子1 = 1，罚因子放大系数 C = 10，收敛速度判别数 = 0.25，允许误差 = 0.001，求解无约束问题 ),(min)(kkLx，迭代情况如表 8.4。第 5 次迭代结束时满足迭代终止条件，因此得最优解 x* x(5) = (0.000, 1.000)

177、T。从表中可以看到，增广乘子法的罚因子k不必很大，而且收敛较快。 表 8.4 增广乘子法迭代情况 k k )(k x1 (k) x2 (k) f (x(k) 1kk k 1 1 1 -6.0000 4.0000 -44.000 9.000 9.0000 2 10 -8.0000 0.4286 0.7857 2.1480 0.071 0.6429 3 10 -1.5714 -0.0306 1.0153 0.9074 0.071 0.0459 4 10 -2.0306 0.0022 0.9989 1.0066 0.071 0.0033 5 10 -1.9978 -0.0002 1.0001 0.9

178、995 0.071 0.0002 1、 收敛速度判别数 (0, 1)和允许误差 ， 置 k = 1； (2) 以 x(k-1)为无约束问题(8.42)的初始点，调用无约束优化程序求解，获其极小点 x(k)； (3) 如果k，则得原问题的近似极小点 x* = x(k)，停止计算；否则用式(8.43)计算乘子向量) 1() 1(,kk，转(4)。 (4) 如果1kk， 则令 k+1 = k ， 转(5)； 否则令 k+1 = Ck ，转(5)。 (5) 置 k = k+1，转(2)。 增广乘子法的罚因子可取某个有限值，克服了罚函数法使增广目标函数出现病态的缺点。计算经验表明，增广乘子法的适用性比罚

179、函数法好，收敛速度也比罚函数法快，是目前求解约束优化问题最好的算法之一。 例 8.7 用增广乘子法求解约束优化问题 01s.t.2)(min212221xxxxf x 解 构造增广拉格朗日函数 2)(221)(2221)()1 (, 0max212),(kkkkkkxxxxLx 当 0)1 (21)(xxkk 时，令 04, 022211xxLxxL 得 x(k) = 0。由0, 0)(kk，知此点不满足不等式0)1 (21)(1xxkk，因此 x(k) = 0 不是极小点。 当 0)1 (21)(xxkk 时，令 0)1 (40)1 (221)(2221)(11xxxxLxxxxLkkkk

180、得 Tkkkkkkk34,34)(2)()()(x 对于任何0, 0)(kk，该点都满足不等式0)1 (21)(xxkk，因此它是极小点。 把 x(k)代入乘子迭代公式(8.38)，得 kkkkkkkkk344344344344, 0max)()()1( 当 k 0 时)(k收敛，而且k越大收敛越快。设 k = 2，得 5452)()1(kk 令 k，得34)(k。这时 x(k)x*，因此 3/13/223423/4,234)23/4(2lim)(Tkkxx 表 8.5 给出了数值计算的迭代情况。 计算时取初始点 x(0) = (0, 0)T，初始乘子1) 1 (，初始罚因子1= 1，罚因子放

181、大系数 C = 10，收敛速度判别数 = 0.25，允许误差 = 0.002。从表中可以看到，增广乘子法收敛较快，经 4 次迭代就已满足精度要求，其最优解 x* x(4) = (0.67, 0.33)T。 表 8.5 增广乘子法迭代情况 k k )(k x1 (k) x2 (k) f (x(k) 1kk k 1 1 1 0.5714 0.2857 0.4898 0.1429 0.1429 2 1 1.1429 0.6122 0.3061 0.5622 0.5714 0.0816 3 10 1.2245 0.6603 0.3301 0.6539 0.1177 0.0096 4 10 1.3205

182、 0.6659 0.3330 0.6652 0.0011 0.00111 时,)()(kk都收敛，而且 k 越大收敛越快。设k = 2，得 31631)()1(kk )()1(21kk 令 k，得0, 4)()(kk。这时 x(k)x*，因此由 Tkkkkkkk243210,22)()()(x 得 2/51224342210,2222lim)(Tkkxx 表 8.6 是数值计算的迭代情况。取初始点 x(0) = (0, 0)T，初始乘子1, 1) 1 () 1 (，初始罚因子 1 = 1，罚因子放大系数 C =10，收敛速度判别数 = 0.25，允许误差 = 0.002。把它与表 8.3 对比

183、可以看出，增广乘子法优于混合法。增广乘子法经 4 次迭代就得其最优解 x* x(4) = (1.00, 2.50)T， 而混合罚函数法需迭代 8 次才能获得同样精度的结果。 表 8.6 增广乘子法迭代情况 k k )(k )(k x1 (k) x2 (k) f (x(k) 1kk k 1 1 1 1 1.3333 5.5000 -46.639 7.1255 36.3333 2 10 7.0000 0.6667 1.0556 2.3421 -12.509 0.0043 0.1553 3 10 3.8421 0.1111 1.0093 2.5083 -13.816 0.0614 0.0095 4

184、10 4.0083 0.0185 1.0015 2.4996 -13.747 0.1619 0.0015 n）个向量，然后从中选取一个相对好的向量，作为迭代计算所使用的搜索方向。它的具体操作步骤如下： (1) 给出一个步长因子 0。使用随机函数产生 mn 个服从均匀分布的伪随机数 rij (0, 1)，并用下式把它转化为(-1, 1)区间内的伪随机数 R ij 。 ), 2 , 1;, 2 , 1(, 12njmirRijij (8.49) (2) 用下式生成 m 个单位向量： ), 2 , 1(,12112)()(2)(1)(miRRRReeeiniinjijiniiie (8.50) (3

185、) 从 x(0)出发，用下式生成 m 个随机点 x(i)： ), 2 , 1()()0()(miiiexx (8.51) 这 m 个随机点分布在以点 x(0)为中心，以 为半径的超球面上（图 8.3） 。 图 8.3 随机方向法搜索示意图 x(L) x1 x2 x* x(0) x(0) x(0) x(L) 0 x(1) (4) 检查 x(i)是否为可行点。对其中的可行点，比较它们的目标函数值， 并挑选出目标函数值最小的那一个随机点x(L)。 接着比较x(L), x(0)两点的目标函数值，如果 f(x(L) f(x(0)，就可以判定 e(L) 是一个下降可行方向。总之，由计算机生成的 m 个随机

186、点中，如果其中的点 x(L) 满足条件 )()(, 2 , 1;, 2 , 1, 0)()(min)()()()()()(kLiviLffmiqvgffxxxxx (8.52) 那么该点对应的向量 e(L)就是一个下降可行方向，因此可把 e(L) 选取为搜索方向。如果 m 个随机点都不满足条件(8.52)，那么把步长因子 减半， 重新进入步骤(3)， 生成新的随机点。 重复这个过程， 直到获得满足条件(8.52)的点 x(L)，或者 经过反复减半而变得充分小时为止。 当 变得充分小时（ ），说明 x(0) 的可行邻域内已无目标函数值更小的点，因此 x(0)就是一个极小点（图 8.4）。 x(0

187、) 0.5 x* 图 8.4 极小点邻域示意图 8.3.3 后继点的确认 确定搜索方向 e(L)之后，随机方向法以 x(L)作为新的初始点，使用迭代公式 ), 2 , 1 , 0()()()1()()0(kLkkLexxxx (8.53) 计算它的后继点 x(k+1)，然后检查 x(k+1)是否满足下降可行条件： ), 2 , 1(, 0)()()()1()()1(qvgffkvkkxxx (8.54) 如果条件(8.54)能够满足，说明后继点 x(k+1)优于当前点 x(k)，因此继续使用迭代公式(8.53)，搜索这个方向上更好的点。如果条件(8.54)不能满足，说明当前点 x(k)已经是这

188、个方向上较好的点了，因此结束这个方向的搜索，并把当前点 x(k)作为新一轮搜索的初始点 x(0)，开始新一轮计算。 图 8.3 随机方向法搜索示意图 x(L) x1 x2 x* x(0) x(0) x(0) x(L) 0 x(1) x(0) 0.5 x* 图 8.4 极小点邻域示意图 8.3.4 算法的步骤 随机方向法的程序框图如图 8.5 所示。它的主要计算步骤如下： (1) 给定随机方向的个数 m，搜索的步长因子 ，允许误差 。 (2) 使用式(8.48)生成一个随机点 x(0)，检查所得随机点 x(0)是否满足全部约束条件 gv(x(0)0, (v=1,2, q)，若不满足，则重新生成随

189、机点，直到出现一个满足全部约束条件的初始可行点 x(0)为止。 (3) 使用式(8.49)和(8.50)生成 m 个随机单位向量 e(i)，并使用式(8.51)获得 m 个随机点 x(i)，这些随机点位于以 x(0)为球心、以 为半径的超球面上。 (4) 使用式(8.52)检查 m 个随机点 x(i) 中是否存在优于 x(0)的可行点 x(L)。如果不存在，则转(6)；否则令 x(0) = x(L)，k = 0，转(5)。 (5) 使用式 )()()1(Lkkexx 计算后继点 x(k+1)，并使用式(8.54)检查 x(k+1)是否为优于 x(k)的可行点。如果 x(k+1)优于 x(k)，

190、则令 k = k+1，继续使用该式搜索这个方向上更好的点。否则令 x(0) = x(k)，转(3)。 (6) 如果 ，则获得近似极小点 x*x(0)，结束迭代；否则令 = 0.5，转(3)。 由式(8.51)和迭代终止条件 ( )可知，所得结果的误差为 exx)0( 8.4 复合形法 对于不等式约束问题 qvgfvn, 2 , 10)(s.t)(minxRxx (8.55) 复合形法首先在可行域内构造一个不规则的多面体，然后比较这个多面体各顶点处的目标函数值，通过去掉目标函数值最大的顶点、增添一个目标函数值较小的顶点的方法，不断改变多面体的形状和位置， 最终让多面体的一个顶点逼近所解问题的极小

191、点。 这个多面体就是所谓的复合形，复合形的顶点一般取为 n + 1 2n 个。 8.4.1 初始复合形的产生 为了保证初始复合形位于可行域内，必须使它的顶点全部为可行点。设复合形有 m 个顶点，可用下面的算法产生一个初始复合形： 一、生成 m 个随机点，其中至少有一个是可行点。 (1) 令 i =1，给出变量 x 的下限值 j 和上限值 j ，即让 ), 2 , 1(,njxjjj (2) 调用程序设计语言的随机函数，产生 n 个服从均匀分布的伪随机数 rj (0, 1), (j = 1, 2, , n)。 (3) 生成一个随机点 x(i)，它的各个分量为 ), 2 , 1(),()(njrx

192、jjjjij (8.56) (4) 如果 x(i)满足全部约束条件 gv(x(i)0, (v = 1,2, q)，则把它记入可行点的下标集合 I，否则记入不可行点的下标集合 T。 (5) 若 i = m，则转(6)；否则令 i = i+1，转(2)。 (6) 若下标集合 I 为空集，则转(2)；否则至少已经获得一个可行点，此阶段任务完成。 二、把 m 个随机点中的不可行点全部移到可行域中来。 (1) 设下标集合 I 中有 s 个元素，用下式计算 x(i)中所有可行点的中心 x(c)： Iiics)()(1xx (8.57) (2) 反复使用下式把非可行点朝中心 x(c)移动，直到所有的非可行点

193、都成为可行点时为止： Ticici),(5 . 0)()()()(xxxx (8.58) 有了 m 个可行点，就获得了初始复合形（图 8.6）。 8.4.2 复合形法的原理和计算步骤 根据复合形各顶点目标函数值的大小， 可以找出三个特殊的顶点：最坏点 x(H)、次坏点 x(G)和最好点 x(L)。即这三个顶点满足下列条件 miffHimiffmiffiLiGiH, 2 , 1)(min)(;, 2 , 1)(max)(, 2 , 1)(max)()()()()()()(xxxxxx (8.59) 图 8.7 复合形及其搜索方向示意图 x1 x2 0 S x(L) x(H) x(c) x(G)

194、x(T) 图 8.6 非可行点两步移入可行域内生成初始复合形的示意图 x(1) x(2) x(3) x(c) x(3) x(3) 显然，用一个目标函数值更小的可行点 x(T)去替换最坏点或次坏点，就可改变复合形的形状，并使新复合形的中心更靠近所解问题的极小点。因此，复合形法的主要工作就是寻找这种替换点。 在获得一个复合形后， 复合形法总是先计算除去最坏点 x(H)之后的各顶点的点集中心 x(c) mHiiicm, 1)()(11xx (8.60) 一般情况下，最坏点 x(H)与 x(c)的连线方向 S = x(c) x(H) (8.61) 是目标函数值下降的方向（图 8.7） 。如果可行域是凸

195、集，x(c)必为可行点，S 是一个下降可行方向。如果可行域是非凸集，x(c)可能为非可行点（图 8.8） 。为了获得可行方向，当 x(c)为非可行点时需要重新生成复合形。 生成新复合形时所用变量的下限 = (1,2, ,n)T和上限 = (1, 2, , n)T改用式 ), 2 , 1(,max,min)()()()(njxxxxLjcjjLjcjj (8.62) 计算，以求新的 x(c)为可行点。 获得下降可行方向 S 后，复合形法使用下式计算替换点 x(T)： x(c) 图 8.8 x(c)为非可行点的情况 x(L) x(H) x(G) Sxx)()(cT (8.63) 式中 为步长因子。

196、 是可调整的，以保证 x(T)既是一个可行点，其目标函数值又能有所下降。在接近极小点的地方， 取很小的值才能做到这一点，因此当迭代计算获得近似极小点时，复合形将变得非常小，以至各顶点的目标函数值的差别几近于零： miLiffm12)()()()(1xx (8.64) 式中 为允许误差。 复合形法以各顶点与最好点的目标函数值之差的均方根值 不大于 为迭代终止条件。 如果在上述方向上找不到优于 x(H)的替换点， 复合形法就改为寻找次坏点 x(G)的替换点。这时，只需把公式(8.60)和(8.61)中的 x(H)换成x(G)就行了。一旦找到 x(T) ，就可用它替换最坏点或次坏点，构成一个其中心

197、x(c)更靠近极小点的新复合形。重复上述过程，复合形将逐渐向极小点移动，形状变小，最终可取它的最好点 x(L)作为极小点 x*的近似值。 复合形法的计算步骤如下： (1) 给定复合形的顶点数 m，允许误差 。 (2) 按前述方法生成初始复合形，计算复合形各顶点的目标函数值。 (3) 确定复合形的最坏点 x(H)、次坏点 x(G)和最好点 x(L)。如果不等式(8.64)成立，则已获得近似极小点 x* x(L)，结束计算；否则转(4)。 (4) 使用公式(8.60)计算除最坏点 x(H)之外的(m-1)个顶点的中心x(c)。检查 x(c)是否满足全部约束条件 gv(x(0) 0, (v=1,2,

198、 q)，若不满足，则用式(8.62) 重新设定变量的下上限，转(2)生成新复合形；否则转(5)。 (5) 令 S = x(c) x(H)， = 1.3。 (6) 使用下式计算 x(T) Sxx)()(cT (7) 用不等式约束条件检查 x(T)是否为可行点， 若是则转(8)； 否则令 = 0.5，转(6)。 (8) 比较 f(x(T) 和 f(x(L) 的大小，若 f(x(T) f(x(H)，则令 x(H) = x(T)， f(x(H) = f(x(T)，构成新的复合形，转(3)；否则转(9)。 (9) 如果 10-5，说明 S 不是下降可行方向，应该改变搜索方向，因此令 x(H) = x(G

199、)，f(x(H) = f(x(G)，转(5)；否则还可缩短步长，因此令 = 0.5，转(6)。 图 8.9 是按上述步骤进行迭代计算的程序框图。 8.5 简约梯度法及广义简约梯度法 1963 年，Wolfe 把线性规划的单纯形法推广到带线性约束的非线性规划问题，提出了简约梯度法。1969 年，Abadie 等人把它发展成求解一般约束优化问题的广义简约梯度法。 8.5.1 简约梯度法 设带线性约束的非线性规划问题为 0xbAxxs.t.)(min f (8.65) 其中 A 为 mn 矩阵，b 为 m 维列向量(m0，计算 )()()1(kNkkNkNSxx 若满足边界约束条件NkNNx ) 1

200、(，则转(5)；否则令kk5 . 0，继续计算)1( kNx，直到满足该条件为止。 (5) 解非线性约束方程组(8.96)，求)1( kBx。用牛顿法求解时，给定最大循环次数 I ，按下面步骤进行： (a) 取)1( kBx的初始值)()0)(1(kBkBxx，令 i = 0。 (b) 用式(8.98)或(8.99)计算)1)(1(ikBx，检查它是否满足下降可行条件 BikBBkNkBkNikBffxxxxx)1)(1()()() 1()1)(1(),(),( 和牛顿迭代终止条件 2)1()1)(1(),(kNikBxxh 若两者都满足，则令)1)(1()1(ikBkBxx，置 k = k

201、+1，转(2)；否则转(c)。 (c) 若 i 3 时，就不存在精确解，只能求它的近似解。这时可以用优化设计的方法，把机构的设计运动规律与给定运动规律之间的偏差最小作为设计目标。以 si 表示对应于 i 的摇杆实际输出角，问题可归结为在保证机构正常工作的条件下使目标函数 miisillllf124321)(),( 极小化。 l1 l2 l3 l4 0 0 图 10.1 曲柄摇杆机构简图 在给定机构运动规律的时候，需要用到曲柄的初始输入角 0 与摇杆的初始输出角 0。由图 10.1 知，它们是杆长的函数： 4324232210421242322102)(arccos,)(2)(arccoslll

202、llllllllll (10.2) 另由图 10.2 可知，摇杆的实际输出角 si 与杆长的关系为 )2(,)0(,iiiiiisi (10.3) 式中 iiiiiiiillllllllllcos22arccos,2arccos412421421242322232 图 10.2 曲柄摇杆机构运动学关系 a) 0i b) i 2 l1 l2 l3 l4 i i i i l1 l2 l3 l4 i i i i a) b) 要使机构正常工作，各杆的长度必须满足下面的约束条件： (1) 曲柄能够整周转动的条件 为了实现曲柄整周转动，曲柄 l1必须顺利通过与连杆 l2共线的两个位置，即必须形成图 10.

203、3 中的 ADC1和 ADC2两个三角形。因为三角形的任意两边之和必大于第三边，所以必有 412331244321,llllllllllll (2) 机构能够灵活转动的条件 由图 10.4 知，在曲柄与机架共线的两个位置上，曲柄摇杆机构的传动角 为最小： (,2)(arccos)2/(,2)(arccos322142322322142322minllllllllllll要使机构转动灵活，最小传动角min不得小于其允许值min，因此 l1 l2 l3 l4 A 图 10.3 曲柄整周转动的条件 l1 C1 C2 D l2 min l1 l2 l3 l4 图 10.4 曲柄摇杆机构的最小传动角 C

204、1 C2 D l2 min co2)(,cos2)(322142322min322142322llllllllllll综上所述，使曲柄摇杆机构按给定规律运动问题的数学模型为 )4 , 3 , 2 , 1(, 00cos2)(0cos2)(000s.t.),(minmin322142322min322142322432143214321124321illlllllllllllllllllllllllllllfimiisi (10.4) 目标函数中的 i 与 si 由式(10.1)(10.3)计算。 由于曲柄摇杆机构的 4 根杆长以同样比例放大或缩小，不会改变其运动规律，因此实际使用式(10.4)

205、时，可以把某个杆长取为 1，只计算另外 3 根杆的相对长度。所以式(10.4)的设计变量实际上只有 3个。通常把曲柄 l1 取为 1，把 l2, l3, l4 作为设计变量。 例 10.1 试设计一个曲柄摇杆机构，其机架 l4 是曲柄 l1 的 5 倍，要求曲柄在从初始角 0 转到 0 + 90的区间内，摇杆输出角i与曲柄输入角i之间的对应关系为 )30, 2 , 1 , 0(,60)(32200iiiii 并且要求工作时的最小传动角min不得小于 45。 解 令 l1 = 1，得 l4 = 5。令 TTxxll),(),(2132，并把相关数据代入前述目标函数和约束条件，得本例题的数学模型：

206、 0,045cos236045cos216060404s.t.)(32)(min212122212122212121213012200xxxxxxxxxxxxxxxxfiisix (10.5) 式中 2222101222101025)1 (arccos,)1 (1025)1 (arccosxxxxxx iiiisixxxcos1026cos5arccoscos10262cos1026arccos22122 )30, 2 , 1(,60iii 求解式(10.5)，得 x1 = 4.1286，x2 = 2.3325。因此将曲柄摇杆机构的杆长比例设计为 l1 : l2 : l3 : l4 = 1 :

207、 4.1286 : 2.3325 : 5 就可最大限度地满足设计要求。 10.2 单级直齿圆柱齿轮减速器的优化设计 设计齿轮减速器的已知参数为输入功率 P，输入轴转速 n1以及减速比 i。因为是单级减速器，所以减速比 i 等于齿轮的齿数比 u。现在用优化设计的方法，确定它的其它主要参数，使减速器不仅能够正常工作，而且具有最小的重量。 图 10.5 是单级直齿圆柱齿轮减速器的结构简图。减速器箱体的尺寸由传动零件的尺寸来决定， 只要减速器的一对齿轮和两根轴的支撑轴承间的那部分材料体积减小，整个减速器的重量就会减小。 根据圆柱齿轮的结构图，小齿轮和大齿轮的材料体积分别为 图 10.5 单级直齿圆柱齿

208、轮减速器的结构简图 l2 l1 B1 d1 dz1 B2 d2 dz2 C D2 D3 2321222222221211125. 04)(25. 0)(25. 0)(25. 0DDCBddBVddBVzz 设齿轮的模数 为 m，小齿轮齿数为 z1，齿宽为 b，则其结构尺寸之间具有如下关系： 13122312122112114 . 05 . 225. 0)(25. 0,6 . 11)128(,25. 0) 3 . 02 . 0(, 8)105(,zmumzDDDdDumzmdDbBCbBbbBumzdmzd (10.6) 设两根轴的支撑长度分别为 l1 , l2，它们这部分的体积为 )(25.

209、02222113zzdldlV 因此，问题可归结为在减速器正常工作的条件下把目标函数 2222112212221212132122111)4 . 05 . 225. 0(1592. 0)(25. 0)()8(25. 01),(zzzzzzzdldldmumzbumdumzbdmzbVVVdldlmzbf 极小化。 为了使减速器能够正常工作，各设计参数之间必须满足下面的约束条件。 一、强度条件 (1) 轮齿的接触应力H不得大于齿面较软者的许用值H ： ) 1(25 . 2211HEHubdTuKZ 其中 ZE 是弹性区域系数，对于一对钢齿轮，ZE = 189.8 (MPa)0.5；K 是载荷系数

210、，通常情况下 K = 1.31.6，可取 K = 2；T1 是输入轴的转矩，T1 = 9550000P / n1 。 (2) 轮齿的弯曲应力F不得大于其许用值F： 2,2211222111111FSaFaFFSaFaFmbdTYKYmbdTYKY 其中 YFa1, YFa2 是两个齿轮的齿形系数，YSa1, YSa2 是两个齿轮的应力校正系数。 (3) 两根轴的弯曲应力z不得大于其许用值z。 由于齿轮在轴上对称布置，两根轴上的最大弯矩分别为 2cos24,20cos242212122111111dlTlFMdlTlFlFMnnn并考虑到轴上开一个键槽后对轴的直径应增大 3%，因此 97. 03

211、23297. 0)()(297. 0323297. 0)()(321322122222222311312121121211zzzzzzdTduTMWTMddTdTMWTM 其中Fn是作用在齿轮齿廓上的法向力； 是应力循环特性折合系数， = 0.6。 二、刚度条件 两根轴的最大挠度 yz不得大于其许用值yz： 20cos3464484820cos3432644848242132142322322411311413141311311zzznnzzznznnzyddlETdElFEJlFyddlETEdlFdElFEJlFy 其中E 是材料的弹性模量，对于钢，E = 2.06105 MPa；J 是轴

212、截面的惯性矩；对于安装齿轮的轴，其挠度的许用值yz= (0.010.03)m，可取yz1= yz2= 0.02m。 三、其它设计要求 (1) 齿宽系数 = b/d1 应在 0.9 与 1.4 之间，因此有 0.9d1b 1.4d1 (2) 轴的支撑距离 l1 , l2为 l1= B1+2 +0.5d z1 ， l2= B1+2 +0.5d z2 其中 为齿轮与箱体内壁之间的距离，该距离不得太小。根据设计规范，有 0.025a+1= 0.0250.5(u+1) d1+1 式中 a 为两齿轮的中心距，因此轴的支撑距离 l1 , l2 必须满足如下条件： l1b+0.025(u+1) d1+10 +

213、0.5d z1 l2b+0.025(u+1) d1+10 +0.5d z2 (3) 为了传动平稳，小齿轮齿数一般大于 19，但也不宜太大，取 20 z1 45 使用上述目标函数和约束条件，就可针对具体的设计任务调用优化计算程序求出减速器的主要设计参数。 例 10.2 试设计一个单级直齿圆柱齿轮减速器， 其输入功率 P= 40 kW，输入轴转速 n1= 960 r/min，齿数比 u = 4。 解 减速器的两个齿轮都用45号钢调质， 其许用接触应力 H = 530MPa，许用弯曲应力 F1 = F2 = 180MPa。两根轴都用 35 号钢正火，其许用弯曲应力 z1 = z2 = -1 = 48

214、MPa。把这些数据代入前述目标函数和约束条件，并令设计变量为 TTzzxxxxxxxdldlmzb),(),(765432122111 得本例题的数学模型： )7 , 2 , 1(, 004502005 . 010125. 005 . 010125. 004 . 109 . 0002. 0744849. 1002. 0744849. 1097. 04876. 528312. 04053146097. 04836. 028312. 04053146018015916670180159166705305 .669294s.t.)4 . 05 . 2(692. 0)(4)()8()(min22673

215、214532132113247232364523234372326352324232122232111132276254263321262321252321ixxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxYYxxxYYxxxxxxxxxxxxxxxxxxxxfiSaFaSaFax (10.7) 其中齿轮的齿形系数 YFa1, YFa2和应力校正系数 YSa1, YSa2是齿数 x3的函数，须查表才能获得其数值。 最好是将这两个系数的数表用曲线拟合的方法总结成解析式， 把它们作为等式约束加入到式(10.7)中之后，再进行优化计算。下面使用一个变通方法来解决这两个系数

216、的查取问题。 分析齿形系数和应力校正系数的数表，知道其乘积大概在 4.0 左右，因此先取 YFa1 YSa1 = YFa2 YSa2 = 4.0 对问题(10.7)求解，得 TTTzzxxxxxxxdldlmzb)44.61,75.164,85.43,96.155,69. 2,00.45,91.108(),(),(765432122111 这时的齿轮模数为 2.69，与标准模数 3 接近，因此在式(10.7)中令 x3 =3.0 后再求解，得 TTTzzxxxxxxxdldlmzb)44.61,75.164,85.43,96.155,00. 3,34.40,91.108(),(),(76543

217、2122111 这时的小齿轮齿数为 40.34， 应取为 40。 根据 z1 = 40 和 z2 = 160， 查得齿形系数 YFa1 = 2.40, YFa2 =2.136 和应力校正系数 YSa1= 1.67, YSa2= 1.837。在式(10.7)中令 x3 =3.0、x2 = 40，并使用查得的齿形系数和应力校正系数值进一步求解，得 TTTzzxxxxxxxdldlmzb)47.61,48.166,01.44,75.157, 0 . 3,40,74.110(),(),(765432122111 据此便得齿轮减速器的主要设计参数： TTzzdldlmzb)62,167,44,158,

218、0 . 3,40,111(),(22111 把这些数值代入式(10.6)可获得其它设计参数。 10.3 普通圆柱螺旋弹簧的优化设计 普通圆柱螺旋弹簧是应用最广泛的一种弹簧。设计普通圆柱螺旋弹簧， 通常先根据载荷、 变形量和结构上的要求， 确定弹簧丝直径 d、弹簧的中径 D 和工作圈数 n 等基本参数，然后根据基本参数计算弹簧的节距 p、螺旋角 和高度 H0 等其它参数。 为了保证弹簧能正常工作并有足够寿命，需要兼顾弹簧在刚度、强度、稳定性和共振性等各方面的设计要求。传统设计方法只有通过反复试算才能满足这些要求， 优化设计方法则能较快地获得最优设计方案。对于一般用途的普通圆柱螺旋弹簧，通常以弹簧

219、丝的体积最小作为设计目标。弹簧丝的体积约为 0.252d2Dn，因此可用 f(d, D, n) = d2Dn 作为其目标函数。设计约束有下面几项： 一、弹簧的刚度条件 弹簧的实际刚度为 nDGdkF348 其中 G 为弹簧丝的剪切弹性模量。对于式中的 d、D、n，国家标准有推荐的尺寸系列，因此 kF 只能取一些离散的数值。弹簧的期望刚度为 Fk 其中F是载荷改变量，是相应的变形改变量。k 的取值是连续的，所以 kF 与 k 之间可能存在误差，设计时需要把它的相对误差控制在某个范围之内： kkkF 式中 是刚度相对误差的允许值。因此 kknDGdkk348 二、弹簧的强度条件 (1) 静应力强度

220、条件 最大工作载荷 Fmax 在弹簧丝内侧产生的最大切应力 max 不得超过其允许值： 83maxmaxdDFK 式中 K 是弹簧的曲度系数，根据弹簧的旋绕比 C = D/d 计算： (365. 0615. 025. 0615. 044142DDDdDDddDdDCCCK因此 0)(615. 0365. 08322maxdDddDdDF 其中许用切应力 视弹簧材料及受载情况而定， 对于用碳素弹簧钢丝制成的类弹簧， = 0.3B；类弹簧， = 0.4B；类弹簧， = 0.5B。这些钢丝的拉伸强度极限 B 与钢丝直径 d 有关。 (2) 疲劳强度条件 对于循环次数较多、在变应力下工作的重要弹簧，其

221、疲劳强度的安全系数计算值 Sca不得超过其允许值： FcaSSmaxmin075. 0 075. 00minmaxFS 其中0是弹簧材料的脉动循环剪切疲劳极限，与变载荷作用次数有关；min是安装载荷F1在弹簧丝内侧产生的应力： 31min8dDFK 因此 06801max3FSFdKDF 把曲度系数 K 用 D、d 表示，弹簧疲劳强度条件被整理为 0)(615. 0365. 06803221maxdDddDdDFSFF 许用安全系数 SF 在计算与材料力学性能数据的精确性高时，取 1.31.7；否则取 1.82.2。 三、弹簧的稳定性条件 对于压缩弹簧，为了避免失稳现象的发生，弹簧的长细比 b

222、 不得超过其允许值： 5 . 10bDdnpDHb 其中 p 为节距，可取 p = 0.35D，因此 05 . 135. 0bDdDn 长细比的允许值与弹簧的支撑方式有关。两端固定时，b = 5.3；一端固定一端可转动时，b = 3.7；两端都可转动时，b = 2.6。 四、弹簧的共振性条件 对于承受高速交变载荷的弹簧，为了避免共振现象的发生，弹簧的自振频率 fb 不得低于其工作频率 fw 的 1520 倍： wsFbfGgnDdDndgnDGdmkf)2015(2248212122234 即 0221)2015(2nDdGgfw 其中 ms是弹簧的质量，g 为重力加速度，g = 9800N/

223、s2； 为弹簧材料的重度，对于钢丝， = 7.6410-5N/mm3。 五、弹簧的其它约束条件 (1) 旋绕比 C 的限制 弹簧的旋绕比不能过小或过大，以免弹簧太软或卷绕困难： 4C14 即 4dD14d (2) 最大变形量max的限制 压缩弹簧的最大变形量 max 不应超过轴向总间距： )(84max3maxmaxdpnnGdnFDkFF 其中 是弹簧的轴向间距。因此有 035. 084max3DdGdFD (3) 安装空间对弹簧外径 D2 的限制 D2 = D +d D2max (4) 设计任务对 d、D、n 取值范围的限制 dmin d dmax ，Dmin D Dmax ，nmin n

224、 nmax 根据弹簧特性选用适当的约束条件就可进行普通圆柱螺旋弹簧的优化设计。 例 10.3 试设计一个一端固定一端可转动的普通圆柱螺旋弹簧，其安装载荷 F1 = 200 N，最大工作载荷 Fmax = 480 N，工作行程为10mm，循环次数为 104，外径 D2 不超过 30mm。 解 因循环次数为 104故属类载荷，选用 C 级碳素弹簧钢丝，G = 8104 MPa， 其许用切应力 = 0.4B， 剪切疲劳极限 0 = 0.45B。拉伸强度极限B与钢丝直径 d 有关，须查表才能获得其数值。表中数据可用曲线拟合的方法总结成如下公式： B = 2142 -115d 由载荷和工作行程，得弹簧的

225、期望刚度 N/mm28102004801maxFFFk 现取允许刚度误差 = 0.01， 疲劳强度安全系数 SF = 2.0， 允许长细比 b = 3.7，并令 (d, D, n)T = (x1, x2, x3) T 把相关数据代入上述约束条件，得本例的数学模型： 0,30035. 0048. 00401407 . 35 . 135. 01152142(45. 0)()615. 0365. 0(19.1818)1152142(4 . 0)()615. 0365. 0(31.122201000072.27028.2810000s.t.)(min32112214132211221321123121

226、21221123121212233241332413221xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxf x (10.8) 由于加载频率不高，所以模型中没有出现振动稳定性条件。求解式(10.8)，得 (d, D, n)T = (x1, x2, x3) T= (4.75, 35.25, 4.10) T 据此令 d = 4.5mm，重新求解式(10.8)，得 (d, D, n)T = (x1, x2, x3) T= (4.5, 29.90, 5.43) T 令 d = 4.5mm，n = 5.5 圈，再求解式(10.8)，得 (d, D, n)T = (x1

227、, x2, x3) T= (4.5, 29.76, 5.5) T 因此该弹簧的主要设计参数 d = 4.5mm，D = 29.8 mm，n = 5.5 圈，并由此得 p = 0.35D =10.4 mm， 3 . 6arctanDp H0 = pn+1.5d = 64.0 mm。 需要指出，机械优化设计中的一些设计变量不是连续变量而是离散变量，例如上面例题中的齿轮齿数和模数、弹簧的钢丝直径等。严格地说， 处理这类问题时采用针对离散变量的优化方法在理论上更为合理。但是对设计变量增加离散型限制将使数学模型变得很复杂，问题的求解难度随之增加。从工程实用的角度看，对这些离散变量采用本章的处理办法，即先把它们处理成连续变量，求得结果后，再取结果附近的离散值，是完全可行的。 关于离散变量的优化方法，可参阅参考文献3、4。