第5章 非线性规划.ppt

第五章非线性规划 NonlinearProgramming 5 1非线性规划的基本概念和定理 一 什么是非线性规划 前面我们已经讨论过线性规划问题 其目标函数和约束条件都是自变量的线性函数 还有一类比线性规划问题范围更广的数学规划 这就是非线性规划问题 什么是非线性规划呢 我们先来看一个例子 例1 已知立式圆柱型油罐的容积为V 设油罐各层圈板厚度与底板和顶板厚度相同 确定油罐的高度与直径使耗钢量最少 非线性规划的基本概念和定理 解 因各圈壁板与顶板 底板厚度相等 故油罐的耗钢量与其总表面积成正比 设油罐的底半径为x 高为y 则V x2y 如果近似把顶做为平板计算 则油罐的表面积为 S 2 x2 2 xy该问题的数学模型为 非线性规划的基本概念和定理 分析该数学模型 可以发现 目标函数不是线性函数 而是决策变量的二次函数 约束条件也不是线性的 而是三次函数 这是一个非线性规划问题 非线性规划的定义 只要目标函数或约束条件中含有一个非线性函数 则称这种规划为非线性规划问题 非线性规划问题的一般形式为 非线性规划的基本概念和定理 其中X xl x2 xn T为n维欧式空间中的向量 f X gi X 为向量函数 gi X 0为由m个非线性等式或不等式组成的约束条件 有三层含义 或 或 对于gi X 0 只要在其两端同乘 1即可化为 gi X 0 任何优化问题都可以分为有约束和无约束优化问题两大类 非线性规划不例外 如在上面的例子中 由约束条件 x2y V可得 非线性规划的基本概念和定理 将y代入目标函数 则原问题变为 显然这是一个无约束非线性规划问题 可用函数求极值的方法求解 结果为y 2x 即罐高等于直径时耗钢量最少 非线性规划的基本概念和定理 由于非线性函数的复杂性和广泛性 一般来说 非线性规划问题的求解比线性规划问题困难得多 而且也不象线性规划那样有单纯形法这一通用方法 在现实世界中 非线性规划的例子很多 下面我们再来讨论几个非线性规划的例子 目前虽己提出了许多非线性规划问题的算法 但每种方法只适用于某一类的问题 因此 当求解一个具体的非线性规划问题时 应根据其特点选择合适的算法 限于时间 这里只介绍几种常见的算法 非线性规划的基本概念和定理 例2 热油管道优化运行问题 某一热油管道 只有一个热泵站 输量一定 确定出站油温和出站压力使总能耗费最小 热油管道的能耗费用包括燃料费和动力费两部分 它们分别是出站油温TR和出站压力Pd的函数 设S为总能耗费 则该问题的数学模型为 非线性规划的基本概念和定理 在目标函数中 f1 TR f2 Pd 一般为非线性函数 约束条件中亦存在不少非线性函数 显然是一个NLP问题 非线性规划的基本概念和定理 例3 最小二乘问题 该问题大量存在于工业生产和科学实验的数据处理中 例如原油的粘度可以表示为 通过实验可测得不同温度Ti下的粘度值 i 如下表 现要确定 0 u T0使实验点构成的曲线与理论曲线误差的平方和最小 这是一个最小二乘问题 决策变量为 0 u T0 其数学模型为 非线性规划的基本概念和定理 最小二乘问题数学模型的一般形式为 二 最优解的概念 对于一般的非线性规划问题 最优解分为局部最优解和全局最优解 简单地说 局部最优解就是目标函数在某一点附近的某个范围内的最优解 点 非线性规划的基本概念和定理 如图所示为某一元函数的曲线 A C E点为局部极大点 B D为局部极小点 对于所研究的整个区域来说 它们不是最优点 非线性规划的基本概念和定理 全局最优解是目标函数在整个可行域内的最优解 它反映目标函数的整体特性 如图中的F点在整个可行域内它是最小的 G点在整个可行域中是最大的 所以它们是全局最优解 点 下面从数学上给出局部极小和全局极小的定义 对于极大点 可以类似定义 非线性规划的基本概念和定理 一般来说 f X 在区域D上可能有若干个局部最优解 而我们总是希望求得全局最优解 然而 目前已有的方法只能求得局部最优解 为了求得全局最优解 通常是找出多个局部最优解 比较求得全局最优解 非线性规划的基本概念和定理 三 无约束问题的最优性条件 函数极值存在的条件 1 一元函数 必要条件 f x 在x 处取得极值的必要条件是f x 0 充分条件 若f x 0 则x 为极大点 若f x 0 则x 为极小点 2 多元函数 必要条件 f X 在D域内存在极值点X 的必要条件为 即f X 在X 处的所有一阶偏导数等于0 非线性规划的基本概念和定理 充分条件 X 为f X 的极小点的充分条件是X 点的海森矩阵A为半正定 若A为正定 则X 为严格局部极小点 f X 在X 点的海森矩阵是指f X 在X 点的二阶偏导数组成的矩阵 若矩阵A左上角的各阶主子式均大于或等于0 则称矩阵A为半正定的 若各阶主子式均大于0 则称矩阵A为正定的 非线性规划的基本概念和定理 四 凸函数与凸规划 对于一般函数 其局部最优解并不一定是全局最优解 然而存在这样一类具有特殊性质的函数 其局部最优解就是全局最优解 这类函数就是凸函数 1 凸函数的定义 则称f X 为一凸函数 若将上式中的不等号反向 则称f X 为凹函数 若将式中的等号去掉则称f X 为D上的严格凸函数或凹函数 非线性规划的基本概念和定理 根据定义 线性函数既是凸函数 又是凹函数 凸函数的几何意义 对于一元函数f x 若函数曲线上任意两点之间的连线永远不在曲线的下方 则f x 为凸函数 参见右图 非线性规划的基本概念和定理 对于二元函数f x1 x2 若函数曲面上任意两点之间的连线永远不在曲面的下方 则f x1 x2 为凸函数 参见右图 非线性规划的基本概念和定理 对于一元函数f x 若函数曲线上任意两点之间的连线永远不在曲线的上方 则f x 为凹函数 参见右图 凹函数的几何意义 非线性规划的基本概念和定理 对于二元函数f x1 x2 若函数曲面上任意两点之间的连线永远不在曲面的上方 则f x1 x2 为凹函数 参见右图 非线性规划的基本概念和定理 2 凸函数的性质 性质1 设f X 是定义在凸集D上的凸函数 则对于任意实数 0 函数 f X 为定义在D上的凸函数 性质2 有限个凸函数的非负组合仍为凸函数 即 性质3 设f X 是定义在凸集D上的凸函数 则对于任意实数 集合 非线性规划的基本概念和定理 3 凸函数的判定 定理1 一阶条件 即用一阶导数判断 设f X 在开凸集D上可微 则f X 为D上凸函数的充要条件是 对于任意X 1 X 2 恒有 即每一点的切线在图形的下方 即梯度向量 开凸集是指不包括边界的凸集 非线性规划的基本概念和定理 定理2 二阶条件 用二阶导数判断 设f X 在开凸集D上二阶可微 则f X 为D上凸函数的充要条件是 对于所有 其海森矩阵为半正定 若为正定 则为严格凸函数 4 凸函数的极值对于凸函数的极值 有以下两个重要的结论 凸函数的任一局部极小值等于全局极小值 若是严格凸函数 则局部极小点是唯一的 凸函数的极小点形成一个凸集 非线性规划的基本概念和定理 根据凸函数极值的性质 如果我们能预先证明目标函数f X 在可行域D内是凸函数 则可以肯定所求得的局部极小点就是全局极小点 而不用求出其所有的极小点并进行比较 5 凸规划 在非线性规划模型 中 如果f X 和 gi X 都是凸函数 或者说f X 为凸函数 gi X 为凹函数 则称这种规划为凸规划 非线性规划的基本概念和定理 凸规划具有如下性质 可行解集为凸集 最优解集为凸集 若最优解存在的话 任何局部最优解均为其全局最优解 若目标函数为严格凸函数 且最优解存在 则最优解必唯一 根据凸规划的定义 线性规划是一种凸规划 非线性规划的基本概念和定理 例 判断下列非线性规划是否为凸规划 非线性规划的基本概念和定理 解 先验证约束条件gi X 为凹函数 g1 X g3 X g4 X 均为自变量的线性函数 为凹函数 g2 X 的海森阵为 为半负定 故g2 X 为凹函数 f X 的海森阵为 为正定 故f X 为凸函数 上述规划为凸规划 有唯一极小点X 0 58 1 34 T f X 3 8 5 2无约束非线性规划寻优方法概述 无约束非线性规划问题的寻优方法大致可分为两类 一 间接寻优方法 也称为解析法 利用求导数寻求函数极值的方法即古典的微分法就属于这一类 这类方法要求把一个非线性规划问题用数学方程式描述出来 然后按照函数极值的必要条件用数学分析的方法求出其解析解 再按照充分条件或者问题的实际物理意义间接地确定是否最优解 因此称为间接法 这类间接寻优方法适用于目标函数具有简单而明确的数学形式的非线性规划问题 而对于目标函数比较复杂 或甚至无明确的数学表达式的情况 间接法就无能为力了 这时需要采用直接寻优方法 非线性规划寻优方法概述 二 直接寻优方法 也称搜索法 这是一种数值方法 利用函数在某一局部区域的性质或一些已知的数值 来确定下一步计算的点 这样一步步搜索逼近 最后达到最优点 这种方法一般称为下降迭代算法 对于最小化问题 下降迭代算法的基本思想是 从某一初始点分X 0 出发 按照一定规则 即算法 转换到某一个比X 0 更好的点X 1 对极小化问题 要求f X 1 f X 0 故称为下降迭代算法 再转换到X 2 如此继续下去就产生了一个点列 X k X 0 X 1 X k 相应的目标函数值满足 非线性规划寻优方法概述 f X 0 f X 1 f X 2 f X k 若该点列收敛于f X 的最优点X 即 则称算法是收敛的 对于某一算法来说 我们要求其产生的点列 X k 中的某一点本身就是问题的精确最优点 或者该点列的极限点是问题的精确最优点 但该极限点是不能达到的 在这种情况下 经过有限次迭代不可能得到问题的精确最优解 而只能得到一个近似最优解 对于一般的非线性函数 大多数迭代算法都只能求得近似最优点 非线性规划寻优方法概述 下降迭代算法的迭代过程可分为4步 选择初始点X 0 令k 0 确定搜索方向 从迭代点X k 出发 确定一个搜索方向P k 要求沿这个方向能找到使目标函数下降的点 确定步长 沿P k 方向前进一步 得到新的迭代点X k 1 检验X k 1 是否满足计算终止准则 如满足则迭代结束 此时X k 1 就是要求的极小点或近似极小点 否则 令k k 1 转 在以上步骤中 选定搜索方向是最关键的一步 它是区分各种算法的主要标志 非线性规划寻优方法概述 在许多算法中 搜索方向P k 确定之后 步长是按使目标函数值下降最大的原则选定的 即求 k 使得 这实际上是求 的一元函数f X k P k 的极小点 k 故常把按上述原则确定 k的过程称为一维搜索 由此确定的步长称为最佳步长 非线性规划寻优方法概述 对于线性规划问题 我们可以利用检验数方便地判断一个基本可行解是否为最优解 但对于非线性规划问题则要复杂得多 由于精确极值点X 在求出之前并不知道 因此 在实用上只能根据所得到的迭代点的信息来判断是否应该终止迭代过程 常用的计算终止准则有以下几种 根据相邻两个迭代点的绝对误差判断 非线性规划寻优方法概述 根据相邻两个迭代点的相对误差判断 根据函数梯度的模判断 非线性规划寻优方法概述 式中 1 5是事先给定的足够小的正数 视问题的精度要求确定 需要指出的是 上述判敛准则没有一种是对任意函数都普遍有效的 例如当函数值随X变化不剧烈时 只用函数值的差值判敛可能会出现错误 虽然满足了 但离最优解X 可能还很远 同样对于函数值随X变化剧烈的情况 只用判断也是不可靠的 最好的办法是将上述判敛准则联合使用 非线性规划寻优方法概述 几乎所有的非线性规划寻优方法求得的结果都是局部最优解 但若非线性规划的目标函数是凸函数 则局部最优解就是全局最优解 在实际工作中 当我们处理一个具体的非线性规划问题时 目标函数是否为凸函数 有时可以验证 有时并不容易验证 而一般求出的极小值往往是局部极小值 这种情况下 对于很多具有迭代性的方法 为了求出全局极小值 可从多个初始点出发进行迭代 求出多个局部极小值 比较确定全局极小值 5 3几种常用的一维搜索算法 一维搜索算法用于求单变量函数极值 单变量极值问题是最简单的极值问题 它不但存在于实际问题中 而且是某些多变量函数最优化算法的基础 许多多元函数无约束最优化搜索算法实质上是由一系列一维搜索构成的 而一维搜索的成功与否对整个算法的效率影响很大 如上节提到的求X k 1 X k P k 中的步长因子 k的单变量极值问题 因此一维搜索算法是最优化理论中很重要的一种算法 几种常用的一维搜索算法 在所讨论的算法中 大多数属于序列搜索法 所谓序列搜索法就是利用过去的信息来生成以后的搜索点的方法 它是一类直接利用目标函数值信息的一维搜索方法 它通过比较一系列试点的目标函数值而逐渐缩小搜索区间 其基本原理是 在搜索区间 a b 内任取两点xl x2且xl x2 计算函数值f xl f x2 将f xl 与f x2 比较 有三种可能情况 f xl f x2 几种常用的一维搜索算法 去掉 x2 b 最小点x 必在 a x2 内 x1为保留点 几种常用的一维搜索算法 f xl f x2 去掉 a x1 最小点x 必在 x1 b 内 x2为保留点 几种常用的一维搜索算法 f xl f x2 同时去掉 a x1 和 x2 b 最小点x 必在 x1 x2 内 这种情况很少遇到 而且在计算机中 浮点数是不存在相等关系的 几种常用的一维搜索算法 然后在缩小的区间内只取一个新点 并将该点的函数值与保留点的函数值比较 按上述原则缩小区间 这样继续下去 当区间缩短到足够小时 即得到了满足给定精度要求的近似解 一 穷举搜索法 我们先讨论一种最简单的搜索方法 穷举搜索法 它是通过将搜索区间划分成若干份 然后计算各个点的函数值并比较而求得最优解 其计算方法和搜索步骤如下 1 确定区间缩短的精度或允许误差 即令要考察的试点xi之间的间隔为 几种常用的一维搜索算法 2 将搜索区间 a b 划分成宽为 的k个区间 k b a 得k 1个试点 3 计算k 1个试点的函数值f xi 4 取k 1个函数值中最小者所对应的xi作为最优解x 则最大误差为 几种常用的一维搜索算法 二 对分搜索法 对分法的搜索步骤为 1 在区间 a b 的中点x a b 2两侧各 2处取两个点 对分搜索法属于序列搜索法 它要求f x 是单峰的 否则只能求得局部最优解 对于min问题 要求f x 是下单峰的 即在x 的左侧f x 严格下降 在x 的右侧f x 严格上升 几种常用的一维搜索算法 并计算其函数值f x1 f x2 其中 为x 的最大允许误差 若f x1 f x2 则删除区间 a x1 并令a x1 2 比较f x1 和f x2 几种常用的一维搜索算法 3 若b a 则x a b 2 f x f x1 f x2 2 算法终止 否则转1 对分法计算时 每次删除的区间约为当前区间的一半 但每次需计算两个点 它虽比穷举法效率高 但仍不太满意 那么有没有一次只计算一个点的更有效的方法呢 回答是肯定的 几种常用的一维搜索算法 三 0 618法 也叫黄金分割法 0 618法是序列搜索法的一种 由于其计算次数较少 且计算过程又比较简单 是一种常用的方法 0 618法的取点原则 对称取点 即所取的两个试点关于区间中点对称 如图所示 在区间内取两个试点xl和x2 由于消去的区间事先不能预测 因此消去区间 a x1 和 x2 b 都是可能的 所以最有利的方案应是要求它们一样长 即 几种常用的一维搜索算法 在迭代过程中 始终保持两个试点的相对位置 即距离的比例关系 不变 也就是说在舍去某一区间后 保留的点仍处于留下区间的相应位置上 即如果消去的是区间 x2 b 留下了 a x2 则保留点x1在 a x2 中的位置应与x2在 a b 中的位置相仿 即比值相同 几种常用的一维搜索算法 由式 2 得 将x2代入式 1 得 由式 3 和 4 得 即 解上述方程得 几种常用的一维搜索算法 这就是0 618法的来历 其中 0 618称为区间缩短率 表示每次缩短后的区间与上一次搜索区间的长度之比 0 618是一个很奇特的数字 在古代有人认为这种比例在人体造型中具有美学价值 故叫黄金分割 如工艺美术和日用品的长宽设计用此比例容易引起美感 大家考察一下是不是这样 事实上人体中的许多比例也符合0 618 如人体肚脐以下部分与身高之比为0 618 上肢与下肢之比 膝关节以下部分与整条腿长之比 肘关节到指尖的长度与整条胳膊长度之比等均为0 618 如果人体的这些部位的比例偏离0 618太多 则人的体型是不美观 不匀称的 另外 股市的涨跌幅预测中也常用黄金分割比例 因此 黄金分割比例广泛存在与自然界中 几种常用的一维搜索算法 由上面的推导可知 只要一个点取在原区间的0 618处 第二个点在它的对称位置 就能保证无论经过多少次舍去 保留点始终在新区间的0 618处 要进一步缩短区间 在其保留点的对称位置再作一次计算就可以了 0 618法的搜索步骤 设给定初始搜索区间 a b 及误差限 0 1 在区间 a b 中选两个点x1和x2 使 几种常用的一维搜索算法 3 若y1 y2 则删去 a x1 并令 a x1 x1 x2 x2 a b x1 y1 y2 y2 f x2 转2 若y1 y2 则删去 x2 b 并令 b x2 x2 x1 x1 a b x2 y2 y1 y1 f x1 转2 几种常用的一维搜索算法 0 618法的FORTRAN语言程序 使用该程序时不需要输入最优解存在的区间 a b 只需要输入迭代的初始点 确定最优解存在区间时使用的步长和迭代的误差限 程序会自动确定最优解存在的区间 程序中目标函数是通过语句函数定义的 必须在第一条可执行语句前定义 该函数定义好后 可以在程序中用F X 的形式直接调用 当目标函数不能用一个简单的函数表示时 可以通过子程序 包括函数子程序和子例行子程序 定义 下面的程序是通过主程序 主要完成输入输出功能 调用子例行子程序 子例行子程序中又调用函数子程序实现的 几种常用的一维搜索算法 C黄金分割法计算机程序 对于不同的问题 请修改函数子程序C该程序先用外推法求最优解存在的区间 然后用黄金分割法求最优解 最小值 C若要求极大点 请将目标函数转换为求极小值 目标函数乘以 1 CMAINPROGRAMDOUBLEPRECISIONEPS X DXWRITE 11 11FORMAT 1X 请输入初始点X 最优解区间搜索步长DX和误差限EPS READ X DX EPSCALLA12 X DX EPS WRITE 最优解为 10WRITE 14 X DX14FORMAT 5X X G15 7 5X S G15 7 STOPEND 几种常用的一维搜索算法 SUBROUTINEA12 X DX EPS DOUBLEPRECISIONFO F1 F2 FA FB X X1 X2 A B T EPS DXDX0 DX2FO F X IF F X DX FO LE 0 GOTO1X2 XDX DX1X1 X DXF1 F X1 IF F1 FO GT 0 GOTO4CDX 2 DXX2 XX X1FO F1IF ABS F1 GT 1E20 OR ABS X1 GT 1E20 THENWRITE 从所输入的初始点出发找不到最优解 1 请重新输入初始点X 几种常用的一维搜索算法 READ XDX DX0GOTO2ENDIFGOTO14A X2B X2IF X2 GT X1 A X1IF X2 LT X1 B X1T SQRT 5 0 1 0 5FA F A FB F B X1 B T A B F1 F X1 X2 A T B A F2 F X2 5IF F1 F2 GE 0 GOTO7B X2FB F2 几种常用的一维搜索算法 X2 X1F2 F1X1 B T A B F1 F X1 GOTO97A X1FA F1X1 X2F1 F2X2 A T B A F2 F X2 9IF ABS B LT EPS GOTO15IF B A GT EPS GOTO5GOTO1615IF ABS B A GT EPS GOTO516IF ABS F1 LT EPS OR ABS F2 LT EPS GOTO17IF ABS F1 F2 F2 GT EPS GOTO5GOTO18 几种常用的一维搜索算法 17IF ABS F1 F2 GT EPS GOTO518IF F1 GT F2 GOTO47X X1DX F1RETURN47X X2DX F2RETURNENDCC定义目标函数的函数子程序FUNCTIONF X DOUBLEPRECISIONF XF X 2 10 X 36RETURNEND 几种常用的一维搜索算法 四 分数法 分数法是与0 618法很相似的一种方法 它们的取点原则相同 只是二者的区间缩短率不同而己 0 618法的区间缩短率为常数0 618 而分数法的区间缩短率为变数Fn 1 Fn 其中 Fn 为Fibonacci数列 它满足 根据这两个关系式可以算出每步的区间缩短率 见下表 几种常用的一维搜索算法 几种常用的一维搜索算法 分数法的计算步骤与0 618法相似 只要将0 618法修改一下即可 当n 时 Fn 1 Fn 0 618 由上表可见 当n 7时 分数法的区间缩短率已接近0 618 因此0 618法是分数法的一个很好的近似 但其计算步骤比分数法简单得多 便于上机计算 因此一般多用0 618法 而不采用分数法 但分数法也有一个优点 就是对于给定的计算精度 可以预先知道迭代的次数 避免出现死循环 几种常用的一维搜索算法 1 确定试点的个数n 设原区间为 a b 将其缩短为原来的 倍 即 查表由Fn确定n 2 令k 1 在 a b 内取两个试点x1和x2 满足 3 若k n 转4 否则 令k k 1 继续下面的计算 如果y1 y2 则删去 x2 b 令 几种常用的一维搜索算法 如果y1 y2 则删去 a x1 令 计算结束 几种常用的一维搜索算法 五 抛物线插值法 1 外推法求极值存在的区间在前面讨论的方法中 都假设最优解存在的区间是已知的 但对于许多实际问题 并不知道最优解存在的区间 这时 无论用前面讨论的方法还是插值法都必须首先确定最优解存在的范围 通常采用的方法是外推法 为了加快寻找区间的速度 通常采用成倍加大步长的措施 几种常用的一维搜索算法 计算步骤 设从某点x1出发 原始步长为h0 令x2 x1 h0 计算f x1 f x2 令k 2 如果f x1 f x2 转 否则转 k k 1 xk xk 1 2k 2h0 计算f xk 若f xk f xk 1 转 否则转 本步开始处 将步长加倍 求x3 x2 2h0 x4 x3 4h0 等点的函数值 直到函数值增加为止 几种常用的一维搜索算法 k k 1 xk xk 1 2k 3h0 若k 3 xk xk 2 2k 3h0 计算f xk 若f xk f xk 1 转 否则转 本步开始处 将步长加倍 求x3 x1 h0 x4 x3 2h0 x5 x4 4h0 等点的函数值 直到函数值增加为止 几种常用的一维搜索算法 令xk 1 xk 1 xk 2 计算f xk 1 若f x1 f x2 则d 2k 3h0 否则d 2k 4h0 令f xb min f xk 2 f xk 1 f xk f xk 1 xa xb d xc xb d 则 xa xc 为最优解存在的区间 该步的含义是 在xk 1和xk中点取点xk 1 得4个等间距的点xk 2 xk 1 xk 1 xk 其间距为d 比较这4个点的函数值 将其中最小者命名为xb 该点两侧各为d的范围即为最优点存在的区间 几种常用的一维搜索算法 2 抛物线插值法 用外推法求出最优解存在的区间后 可以采用前面讨论的任何一种方法求得最优解 为了能充分利用外推法的计算结果 加速寻优过程 可采用抛物线插值法求最优解 外推法确定的最优解存在区间为 xa xc 且xb为该区间内的一点 其函数值满足 f xa f xb f xc f xb 可通过三点 xa f xa xb f xb xc f xc 作一条二次抛物线P x 来逼近函数f x 设 P x a0 a1x a2x2 几种常用的一维搜索算法 将三点的坐标代入可得三个线性方程 联立求解可得 P x 的极小点即为抛物线的顶点 几种常用的一维搜索算法 若xa xc xb三点等间距 上式简化为 迭代步骤 设最优解存在区间为 xa xc 在 xa xc 中取一点xb 计算f xa f xb f xc 如果最优解存在区间是由外推法求得的 可直接利用其xa xb xc f xa f xb f xc 值 按式 1 计算近似最优解xm 若三点等间距 可按式 2 计算 几种常用的一维搜索算法 若 xb xm 为给定的允许误差 0 则x xm 停止计算 否则转 若f xm f xb 转 否则转 若xm xb 令xa xm f xa f xm 即保留xb xc不变 删掉 xa xm 见图1 转 几种常用的一维搜索算法 若xm xb 令xc xm f xc f xm 即保留xa xb不变 删掉 xm xc 见图2 转 几种常用的一维搜索算法 若xm xb 令xa xb f xa f xb xb xm f xb f xm 即保留xc不变 删掉 xa xb 见图3 转 几种常用的一维搜索算法 若xm xb 令xc xb f xc f xb xb xm f xb f xm 即保留xa不变 删掉 xb xc 见图4 转 几种常用的一维搜索算法 抛物线插值法的收敛速度在很大程度上取决于目标函数的性质 当目标函数比较接近二次抛物线时 这种方法可能比分数法或0 618法收敛快 然而 当相差较大时 该方法的收敛速度可能很慢 甚至不收敛或收敛点与f x 的极值点相差很远 因此 在对函数性质缺乏了解时 采用0 618法更可靠 几种常用的一维搜索算法 3 外推内插法的FORTRAN语言程序 程序中目标函数采用语句函数定义 欲求解其它目标函数时 可修改程序第一句中的目标函数定义 亦可采用子程序进行定义 使用本程序时 应输入下列数据 x0 初始点 E 精度 另外 本程序中 给定的初始步长h0取值为1 0 由H0 1 0赋值 迭代时步长缩短倍数P取值为0 1 由P 0 1赋值 如果这两个量选择其它数值 修改上述两个赋值语句即可 几种常用的一维搜索算法 CCCC外推内插法求解一元函数的最小值F X X 3 10 X X 36H0 1 0P 0 1WRITE 请输入初始点X0 精度E READ X0 EX1 X0F1 F X1 40N 0H H0X2 X1 HF2 F X2 IF F1 GT F2 GOTO105H HX3 X1F3 F180X1 X2 几种常用的一维搜索算法 F1 F2X2 X3F2 F3GOTO115105H H HN N 1115X3 X2 HF3 F X3 IF F2 GT F3 GOTO150IF N GT 0 GOTO165135C1 F3 F1 X3 X1 C2 F2 F1 X2 X1 C1 X2 X3 GOTO220150H H HN N 1GOTO80165X4 0 5 X2 X3 F4 F X4 IF F4 GT F2 GOTO205 几种常用的一维搜索算法 X1 X2F1 F2X2 X4F2 F4GOTO135205X3 X4F3 F4GOTO135220IF ABS C2 LT E GOTO285X4 0 5 X1 X3 C1 C2 F4 F X4 IF F2 LT 1 0 GOTO250F5 F2GOTO255250F5 1 0 几种常用的一维搜索算法 255IF ABS F4 F2 F5 LT E GOTO300IF F4 GT F2 GOTO285X1 X4F1 F4275H0 P H0GOTO40285X1 X2F1 F2GOTO275300WRITE 305 F4 X4305FORMAT 1X 最优目标函数值 G15 6 1X 最优解X G15 6 STOPEND 多元函数无约束最优化问题的直接搜索法 多元函数无约束最优化问题一般采用直接搜索法 工程上许多问题的数学模型往往难以用解析法求出目标函数的极值 只能采用逐步逼近的直接法 在最优化方法中有许多直接搜索法 可以分成两大类 即 1 不求导的算法 不必计算目标函数的导数 只靠计算目标函数值来搜索 这种方法发展较早 其中有的算法收敛较慢 主要有坐标轮换法 步长加速法 模矢搜索法 方向加速法 共扼方向法 Powell法 单纯形法和随机搜索法等 2 以求目标函数的导数为基础的算法 这种方法提出较晚但发展较快 也较为有效 主要有最速下降法 一阶梯度法 共扼梯度法 拟牛顿法 变尺度法等 本章主要给大家介绍坐标轮换法 步长加速法 方向加速法单纯形加速法和变尺度法等 5 4坐标轮换法 坐标轮换法也叫变量轮换法 这是一种最古老的直接优化方法 它的思路非常简单 即依次在各个坐标方向上作一维搜索 且每次搜索都以前一次搜索的最好点作为起点 下面先以二元函数为例来进行分析 然后给出n元函数的一般算法 设某二元函数S f X f x1 x2 其等高线见下图 设极值存在的区间为a1 x1 b1 a2 x2 b2 用分别表示x1 x2坐标方向 一 坐标轮换法的算法 坐标轮换法 搜索过程如下 1 设从X 0 x1 0 x2 0 出发 先固定x1 x1 0 不变 求以x2为单变量的目标函数f x2 f x2 0 f 的极值 可利用上一节介绍的一维搜索方法求以为单变量的目标函数的极值 从而可定出最佳步长因子 得 坐标轮换法 2 然后固定x2 x2 1 x2 0 不变 变动x1 求以x1为单变量的目标函数f x1 f x1 0 f 的极值 亦即求以为单变量的目标函数的极值 从而可定出最佳步长因子 得 3 因为S 2 优于S 1 因此在进行了一轮坐标轮换后进一步搜索时 搜索区间已缩短为 x1 1 x1 b1 x2 1 x2 b2 在新的区间中重复以上步骤 每次搜索不仅可使目标函数值得到改善 而且可以把搜索区间进一步缩小 直到达到给定得精度为止 坐标轮换法 对n元函数 坐标轮换法的一般步骤为 记n个坐标方向为 则在第j个坐标方向上搜索即为单变量极值问题 坐标轮换法 其中X j 1 为在方向上搜索得到的最好点 所谓在方向上搜索是指固定其它n 1个变量不变 求以xj为单变量的目标函数的极值 算法如下 作一维搜索 确定 若j n 已完成一轮坐标轮换 转 若j n 令j j 1 转 坐标轮换法 否则 令X 0 X n j 1 转 上述算法中 X 1 X 2 X n 为进行一轮坐标轮换后得到的最好点集 坐标轮换法的收敛速度很慢 只适用于变量数目少 且目标函数中变量间的交互作用弱的情况 如果目标函数中没有自变量的乘积项 即一个变量的变化与其它变量无关 则称无交互作用 这时坐标轮换法最有效 如果目标函数中有自变量的乘积项 则称为弱交互作用 这时坐标轮换法的效率很低 坐标轮换法 如果目标函数中出现了脊线 见下图 称为强交互作用 此时如果某一步搜索正好落在脊线上 则无论向哪个方向搜索 目标函数值均增大 搜索过程终止于脊线上 坐标轮换法无效 坐标轮换法 二 坐标轮换法的FORTRAN语言程序 本程序中的目标函数由语句函数定义 当目标函较复杂时 目标函数可改为由子程序计算 最佳步长采用0 618法确定 欲求解其它目标函数时 应修改本程序第二句中的目标函数定义 若所求解的目标函数不是二元函数 还需修改程序中调用该函数的语句 坐标轮换法 使用本程序时 要输入下列数据 N 目标函数的维数 E 精度 A I 第I维变量的初始区间的起点 B I 第I维变量的初始区间的终点 X I 初始点的第I维坐标 I 1 2 N 坐标轮换法 CCCC坐标轮换法求解多元函数无约束极值的计算程序DIMENSIONA 20 B 20 X 20 F X1 X2 X1 X1 X2 X2 X1 X2 10 X1 4 X2 60OPEN 21 FILE LUN IN READ 21 N EDO100I 1 NREAD 21 A I B I X I 100CONTINUER 0 1M 1155DO385I 1 NA1 A I B1 B I X0 X I Y2 A1 0 381966 B1 A1 X I Y2F2 F X 1 X 2 195Y1 A1 0 618034 B1 A1 X I Y1F1 F X 1 X 2 坐标轮换法 215IF F1 LT F2 GOTO235G F2X I Y2GOTO245235G F1X I Y1245IF ABS Y1 Y2 GT R GOTO265IF ABS G LT 1 0 GOTO260IF ABS F1 F2 G LT E GOTO330GOTO265260IF ABS F1 F2 LT E GOTO330265IF F1 LT F2 GOTO310B1 Y1Y1 Y2F1 F2Y2 A1 0 381966 B1 A1 X I Y2F2 F X 1 X 2 GOTO215310A1 Y2Y2 Y1F2 F1 坐标轮换法 GOTO195330IF M I EQ 1 GOTO380IF ABS G LT 1 0 GOTO355IF ABS G G0 G LT E GOTO400GOTO360355IF ABS G G0 LT E GOTO400360IF X0 LT X I GOTO375B I X0GOTO380375A I X0380G0 G385CONTINUEM M 1GOTO155400OPEN 22 FILE LUN OUT WRITE 22 401 G401FORMAT 1X 最优目标函数 E10 4 WRITE 22 最优解为 WRITE 22 416 I X I I 1 N 416FORMAT 5X X I2 E10 4 CLOSE 21 CLOSE 22 STOPEND 5 5步长加速法 步长加速法也叫模矢搜索法PatternSearch 该法是由胡克 Hooke 和基夫斯 Jeeves 于1961年提出的 是一种特别有效的直接搜索法 该方法的特点是将多维搜索分成探测性移动和模矢移动进行 前者的目的是寻求函数下降 对min问题 的有利方向 后者的目的是在前者的基础上 沿着有利的方向以加速的步伐寻求较好的点 两种移动交替进行 一 步长加速法的基本思想 步长加速法 将探测性移动 用S表示 的起点称为参考点 用R0 R1 表示 其终点称为基点 用B0 B1 表示 模矢移动 用P表示 以基点为起点 以参考点为终点 整个搜索过程从某一点出发 用B0 R0表示 搜索过程总是以探测性移动开始 搜索过程可表示为 二 步长加速法算法 下面以n元实函数f X 的极小值为例说明两种移动的过程 步长加速法 1 任选一初始点X0 B0 R0 求出该点的目标函数值f X0 2 对于每个独立变量xi选择一个探索移动步长 i 令 即为第i个分量为 i 其余分量为0的向量 3 计算参考点R0处的目标函数值f R0 并在R0处作局部探测如下 与坐标轮换法类同 先对x1分量增加一个 1 即先沿xl方向探测 记作 计算处的目标函数值 步长加速法 即 步长加速法 然后在R01处对第二个独立变量x2作同样的探测性试验 直到进行到第n个独立变量xn为止 对于第j个独立变量xj 步长加速法 4 在对n个变量进行探测后 将最后一个临时参考点R0 n记作基点B1 R0 n 由基点B0到基点B1的方向就构成了n维空间的第一个模矢 这就是对基点B0来说使目标函数值得以改善的最有利的方向 沿这一方向前进 目标函数值下降最快 第一个模矢移动得到一个新的参考点R1 B0 2 B1 B0 2B1 B0 即把B0到B1的距离沿模矢方向延伸一倍 5 若f R1 f B1 即无改善 转8 否则 在R1点进行局部探测 方法同步骤3 对n个变量进行局部探测后 得到第三个基点B2 R1 n 步长加速法 6 像前述那样 建立新的模矢 得新的参考点R2 2B2 B1 一般情况为 Ri Bi 1 2 Bi Bi 1 2Bi Bi 1 7 如果对第i个参考点Ri进行所有的局部探测均不能改善结果 但f Ri f Bi 则将 i缩小一半再进行局部探测 8 若f Ri f Bi 则退至Bi 即令Bi 1 Bi Ri 1 依此作为参考点进行新的局部探测 如果得到一个好点 则进行模矢移动 否则缩短步长 即令 i i 2 9 计算终止准则 i缩小到 i 时计算终止 为要求的计算精度 步长加速法 该方法求得的是局部极小点 当目标函数的交互作用很强时 搜索容易失败 模矢图如下 步长加速法可用于油库 计量站 转油站和联合站等的最佳位置的确定 步长加速法 例 拟在下图所示的场地上建一小型油库 铁路岔道 公路 上下水管道 煤气管道及电力线均需由图示的枢纽点接入 要求铁路支线由垂直方向引进 单位长度道路 铁路及管道等的修建费见下表 工程地质勘察表明 距铁路以下5米的硬粘土层是由铁路向公路方向倾斜的 坡度为1 100 其上被一层软沉积土所覆盖 故建造油罐和泵房时需打支撑桩150根 试确定最佳的库址 解 如上图所示 取横向为xl轴 纵向为x2轴 油库中心位置可表示为 xl x2 上述问题实际上是一个二元目标函数求极小值的问题 各种管道 道路 铁路岔道和电力线的长度分别为 铁路岔道长度为 x2千米 cl 150元 米 15万元 千米 上水管和煤气管长度为 x12 x2 2 2 1 2千米 c2 5万元 千米 下水管长度为 x1 0 2 2 5 6 x2 2 1 2千米 c3 4万元 千米 电力线长度为 5 x1 2 x22 1 2千米 c4 3万元 千米 道路长度为 3 x1 2 4 8 x2 2 1 2千米 c5 12万元 千米 随库址不同 每根支撑桩需增加的长度为 0 01x2千米 c6 1 5万元 千米 步长加速法 设总的建造费用 只考虑随库址而变的费用 相同部分不考虑 S f x1 x2 则问题的数学模型为 步长加速法 求解过程如下 1 将初始参考点取在接近场地中心的位置 即 f R0 f 2 5 2 5 110 164万元 2 先沿xl轴作探测性移动 步长加速法 3 再沿x2轴作探测性移动 4 得第一次探测性移动的基点B1 R02 2 4 2 4 f B1 109 400进行模矢移动 R1 2B1 B0 4 8 2 5 4 8 2 5 2 3 2 3 f R1 f 2 3 2 3 108 713 f B1 5 进行第二次探测性移动 先沿xl轴作探测性移动 步长加速法 再沿x2轴作探测性移动 得第二次探测性移动的基点B2 R12 2 2 2 2 f B2 108 100 6 进行第二次模矢移动 R2 2B2 B1 4 4 2 4 4 4 2 4 2 0 2 0 f R2 f 2 0 2 0 107 095 f B2 步长加速法 7 按上述同样的步骤进行下去 到第7步时 B7 1 9 0 8 f B7 104 791 R7 2 1 0 7 f R7 f 2 1 0 7 104 767 f B7 以R7为参考点进行探测性移动 步长加速法 8 作模矢移动 R8 2B8 B7 4 0 1 9 1 6 0 8 2 1 0 8 f R8 f 2 1 0 8 104 767 f B8 故应退回到B8 即令 R8 B8 2 0 0 8 f R8 104 759 以R8为参考点进行探测性移动 步长加速法 此时无论向哪个方向移动 目标函数值均增大 所以当采用步长0 1km时 其最优点为 2 0 0 8 满足约束条件 这时xl x2的最大误差均为0 1km 如果给定的允许最大误差为0 05km 则应以 2 0 0 8 为初始点 将步长缩短一半继续重复上面的步骤 这里就不再做了 大家可在课下练习一下 5 6方向加速法 方向加速法也叫共扼方向法 该方法是由M J D Powell在1964年首先提出的 所以又叫Powell法 它可在有限步内找出n元二次函数的极小点 被认为是直接法中比较有效的一种方法 一 方向加速法的基本思想 从选定的初始点X1 0 出发 先依次沿n个坐标方向搜索 得到X1 n 坐标轮换法 然后沿X1 0 到X1 n 的连线方向 X1 n X1 0 再搜索一次 得新点X1 n 1 把从X1 0 到X1 n 1 所经历的n 1次一维搜索称为一轮迭代 方向加速法 方向加速法 上述性质表明 对于二次函数f X 从某点X1 0 出发 沿共扼方向搜索 可以很快得到其极值点 且对于n元正定二次函数 沿n个关于A共额的方向相继进行一轮一维搜索 最后得到的点Xn n 1 即为f X 的极值点 也就是说至多需要n轮迭代即可求得函数的极值点 方向加速法 正定二次函数的上述性质表明 方向加速法具有二次收敛性 即对于正定二次函数 算法可在有限步内收敛 因为一般函数的极小点附近可用二次函数近似 因此该方法对于一般函数也具有较好的收敛性 但是 对于一般函数按上述迭代步骤得到的n个方向可能是线性相关的 从而只能得到函数在某一子空间上的极小点 为此 Powell对上述方法进行了改进 每轮迭代不一定换掉一个坐标方向 而且也不一定按顺序换掉坐标方向 得到了改进的Powell法 方向加速法 二 方向加速法的算法 我们首先以二元函数为例说明方向加速法的具体做法 然后再归纳出n元函数的一般步骤 方向加速法 返回 方向加速法 Powell法的关键在于如何确定应替换掉的坐标方向 其原则是 应替换掉使目标函数值有最大下降的搜索方向 方向加速法 方向加速法 根据上面对二元函数的分析 可以得出n元函数的Powell算法 方向加速法 方向加速法 方向加速法 Powell法可用于确定热油管道的最优加热温度 对于正定二次函数 上述算法中第 步的判敛准则是多余的 因为前面已经讲过进行n轮迭代后 得到的优点就是最优点X 例 求目标函数f X 60 10 x1 4x2 x12 x22 xlx2的最小值及其最优点 解 第一阶段迭代 取X0 0 0 f0 f X0 60 方向加速法 方向加速法 方向加速法 方向加速法 第二阶段迭代 方向加速法 方向加速法 方向加速法 由于正定二次函数的二次收敛性 第二阶段的最后一个点即为最优点 即 X 8 6 f X 8 5 7变尺度法 1959年戴维敦Davidon 提出变尺度法 1963年付莱丘 Fletcher 和鲍威尔 Powell 改进这一方法 简称DFP法 是一种基于目标函数导数的求解无约束极值问题的有效方法 与其它同类算法相比 它不需要计算目标函数的二阶导数矩阵及矩阵求逆 收敛速度快 特别是对于多维问题具有显著的优越性 被认为是求解多维无约束极值问题的有效算法之一 变尺度法 一 DFP法的基本原理 假定无约束极值问题的目标函数f X 二阶连续可微 Xk为极小点的某一近似值 在这个点附近f X 的二阶泰勒多项式逼近为 其梯度为 该近似函数的极小点满足 从而得到 变尺度法 其中H Xk 为f X 在Xk点的海森矩阵 X X Xk 如果f X 为二次函数 则H X 为常数阵 这时逼近式 1 是准确的 在这种情况下 从任一点Xk出发 用式 3 只要一步即可求出f X 的极小点 当f