第4章确定最小安全系数的最优化方法 4.1概述 4.1.1应用最优化方法确定最小安全系数 边坡稳定分析应包含下面两个步骤 )对滑坡体内某一滑裂面按第2或第3章介绍的方法,确定其抗滑稳定安全系数: 2)在所有可能的滑裂面中,重复上述步骤,找出相应最小安全系数的临界滑裂面 在讨论了计算一个滑裂面的安全系数的方法后,本章介绍边坡稳定分析极限平衡法的 第二步,就是要寻找最小安全系数的方法。 如果滑裂面曲线为y(x)。那么,这个问题具体化为寻找下列泛函的极值 F=F) 岩土工程中的边坡的几何形状各异,材料通常是非均质性,纯解析的变分原理很难进 行极值计算。用最优化方法通过数值方法求解,是一个比较现实可行的途径。 我国学者较早开展应用数值规划方法求解安全系数极值的问题(闫中华,1983;张天宝 1981;周文通,1984)。20世纪80年代初期,孙君实(1983)和 Nguyen(1985)分别提出了使 用复形法和单形法搜索任意形状和圆弧滑裂面的最小安全系数的方法。80年代中期有更多 的学者发表了有关研究工作 Celestino and Duncan,1981; Li. and White,1987)。Chen& Shao(1988)采用单形法和牛顿法进行任意形状滑裂面搜索的研究成果。这篇论文所附的几 道例题,后来在国外多篇论文中被引用,作为检验新的优化方法的考题。借此书的机会, 作者将这些例题和原数据文件在本章介绍。 国内外有关的研究成果和STAB软件十几年在全国推广的实践证明,应用计算机自动 搜索临界滑裂面是可行的。 同时作者也发现,相对于三维或二维斜分条的极限平衡法,垂直条分法的极小值搜索 问题比较简单。采用任何一种优化计算方法配合本章介绍的随机搜索法,即可快速地找到 临界滑裂面。因此,目前在STAB程序中向用户提供的只是单形法这个功能。以后作者还 将把本章中介绍的各种方法补充进程序,供用户选用 4.1.2最优化方法 最优化方法是近代数学规划中十分活跃的一个领域。目前,已有许多十分成熟的计算 方法。这些计算方法总的来看,可以分为以下三大类。 1.枚举法 枚举法的基本思想是,根据一定的模式,比较不同自变量的目标函数,经过筛选,最
第4章 确定最小安全系数的最优化方法 4. 1 概述 4. 1. 1 应用最优化方法确定最小安全系数 边坡稳定分析应包含下面两个步骤 1) 对滑坡体内某一滑裂面按第 2 或第 3 章介绍的方法 确定其抗滑稳定安全系数 2) 在所有可能的滑裂面中 重复上述步骤 找出相应最小安全系数的临界滑裂面 在讨论了计算一个滑裂面的安全系数的方法后 本章介绍边坡稳定分析极限平衡法的 第二步 就是要寻找最小安全系数的方法 如果滑裂面曲线为 y(x) 那么 这个问题具体化为寻找下列泛函的极值 F = F( y) (4.1) 岩土工程中的边坡的几何形状各异 材料通常是非均质性 纯解析的变分原理很难进 行极值计算 用最优化方法通过数值方法求解 是一个比较现实可行的途径 我国学者较早开展应用数值规划方法求解安全系数极值的问题 闫中华, 1983 张天宝, 1981 周文通, 1984 20 世纪 80 年代初期 孙君实(1983)和 Nguyen (1985)分别提出了使 用复形法和单形法搜索任意形状和圆弧滑裂面的最小安全系数的方法 80 年代中期有更多 的学者发表了有关研究工作(Celestino and Duncan, 1981; Li. and White, 1987) Chen & Shao(1988)采用单形法和牛顿法进行任意形状滑裂面搜索的研究成果 这篇论文所附的几 道例题 后来在国外多篇论文中被引用 作为检验新的优化方法的考题 借此书的机会 作者将这些例题和原数据文件在本章介绍 国内外有关的研究成果和 STAB 软件十几年在全国推广的实践证明 应用计算机自动 搜索临界滑裂面是可行的 同时作者也发现 相对于三维或二维斜分条的极限平衡法 垂直条分法的极小值搜索 问题比较简单 采用任何一种优化计算方法配合本章介绍的随机搜索法 即可快速地找到 临界滑裂面 因此 目前在 STAB 程序中向用户提供的只是单形法这个功能 以后作者还 将把本章中介绍的各种方法补充进程序 供用户选用 4. 1. 2 最优化方法 最优化方法是近代数学规划中十分活跃的一个领域 目前 已有许多十分成熟的计算 方法 这些计算方法总的来看 可以分为以下三大类 1. 枚举法 枚举法的基本思想是 根据一定的模式 比较不同自变量的目标函数 经过筛选 最
88土质边坡稳定分析一原理·方法程序 终找到最小值。这是最原始、简单的方法 如图41中示,任一圆弧可用其圆心坐标x,y)和半径r确定。其相应的安全系数F可 表达为 F=f(x0,y0,D3) 式中:D为滑弧深度,即圆弧最低点的坐标,可知 D Vo 显然这是三个自由度的问题。采用枚举法,不断地改变x。y和D的数值,逐一比较 相应的安全系数,最终找到最小的安全系数。在具体操作中,先固定一个D,然后在圆 可能的位置中布置一个网格,见第12章图121。设网格的中心坐标为x和y2。在左、右方 向,各布置了N格,在上、下方向各布置了N格,则共计有(2N2+1)×(2N+1)个网格点,分 别以该网格点为圆心,以D,为滑弧深度计算相应安全系数,找出最小的安全系数。然后改 变一个D值,重复相同的步骤。在这一过程中,有可能出现圆弧和边坡不相交的情况,则 应令程序自动抛弃该圆弧。同样,D也是以一个中心值起算,在其上、下各布置N层,这 样总计计算(2N2+)×(2N+1)x(2N+1)个圆弧。 在枚举法的基础上,用0618法或其他方法,提高搜索最小安全系数工作的效率。这 方面的工作可参阅有关文献。 图4.1圆弧滑裂面 2.数值分析方法 随着计算机的发展,数值分析方法逐步形成一门完整的学科,统称为最优化方法( Method of optimization)。这一方法又可分为两大类 第一类称模式搜索法( Pattern search method)。其基本思想是,根据一定的模式,比较不 同自变量的目标函数,确定最优的搜索方向,最终找到最小值。 另一类称牛顿法。它要通过解析手段寻找使目标函数F对自变量κ的偏导数为零的极 值点(aFOx2=0,i=1,2…,m)。同时,从理论上讲,还需要满足由二阶导数形成的 Hessian矩 阵正定这个达到极小值的充分条件。此类方法中以导数为研究的主要对象,因此,也称为 以导数为基础的方法( Gradient- based- method)。一般认为,当自由度较多时,直接搜索法效 率较低。此时需要考虑牛顿法体系的分析方法
88 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 终找到最小值 这是最原始 简单的方法 如图 4.1 中示 任一圆弧可用其圆心坐标(x0, y0)和半径 r 确定 其相应的安全系数 F 可 表达为 ( , , ) 0 0 Ds F = f x y (4.2) 式中 Ds为滑弧深度 即圆弧最低点的坐标 可知 0 D r y s = + (4.3) 显然这是三个自由度的问题 采用枚举法 不断地改变 x0, y0 和 Ds 的数值 逐一比较 相应的安全系数 最终找到最小的安全系数 在具体操作中 先固定一个 Ds 然后在圆心 可能的位置中布置一个网格 见第 12 章图 12.1 设网格的中心坐标为 xc和 yc 在左 右方 向 各布置了 Nx 格 在上 下方向各布置了 Ny 格 则共计有(2Nx+1)×(2Ny+1)个网格点 分 别以该网格点为圆心 以 Ds 为滑弧深度计算相应安全系数 找出最小的安全系数 然后改 变一个 Ds 值 重复相同的步骤 在这一过程中 有可能出现圆弧和边坡不相交的情况 则 应令程序自动抛弃该圆弧 同样 Ds 也是以一个中心值起算 在其上 下各布置 Nd 层 这 样总计计算(2Nx+1)×(2Ny+1)×(2Nd+1)个圆弧 在枚举法的基础上 用 0.618 法或其他方法 提高搜索最小安全系数工作的效率 这 方面的工作可参阅有关文献 图 4. 1 圆弧滑裂面 2. 数值分析方法 随着计算机的发展 数值分析方法逐步形成一门完整的学科 统称为最优化方法(Method of optimization) 这一方法又可分为两大类 第一类称模式搜索法(Pattern search method) 其基本思想是 根据一定的模式 比较不 同自变量的目标函数 确定最优的搜索方向 最终找到最小值 另一类称牛顿法 它要通过解析手段寻找使目标函数 F 对自变量 zi 的偏导数为零的极 值点(∂F/∂zi = 0, i = 1,2,...,n) 同时 从理论上讲 还需要满足由二阶导数形成的 Hessian 矩 阵正定这个达到极小值的充分条件 此类方法中以导数为研究的主要对象 因此 也称为 以导数为基础的方法(Gradient−based−method) 一般认为 当自由度较多时 直接搜索法效 率较低 此时需要考虑牛顿法体系的分析方法
第4章确定最小安全系数的最优化方法89 本节介绍直接搜索法中的单纯形法和 Powell法以及牛顿法中负梯度法和DFP法 ( Davidon- Fletcher-Powell method)由于这些方法的原理在众多的文献及教科书中都有所介 绍,这里不拟全面阐述其原理,而着重讨论在边坡稳定分析领域需解决的特殊问题。 3.非数值分析方法 第三类方法是在近期计算机发展基础上形成的,称为非数值方法。这一方法被广泛应 用于管理科学、计算机科学、分子物理学及超大规模集成电路设计中,用于解决组合优化 问题。非数值分析利用计算机具有容量大、计算速度快的优点,通过大量随机采样来找到 目标函数的最优值。近期涌现了诸如模拟退火,遗传算法,神经网络和蚂蚁算法等。均属 此类。在边坡稳定分析中,上述算法均有人尝试。本章47节将对此类方法作一简要介绍 4.2任意形状滑裂面的模拟和目标函数的确立 4.2.1任意形状滑裂面的模拟 最优化问题的提法是:对于一个具有n个自变量的向量z=(1,z2,…,m),确定其目标 函数F的最小安全系数Fm,相应的自变量为zn 为此,需要对式(41)中的曲线y(x)用若干参数来模拟。也就是说,需要将任意形状滑 裂面yx)用z来近似表达。将滑裂面曲线用m个点A1,A2,…,Am离散,见图42。也就是将 此m个点用直线或光滑的曲线连起来,以近似模拟此曲线。此m个点的坐标用x(i=1,2,…m) 表示: 次弱夹层 图4.2任意形状滑裂面 旦这种连接的模式确定,安全系数F即可表达成此m个点的坐标x,y2x2y2…xmyn 的函数 F=F(x,yi, x2, y2,. xm,y
第 4 章 确定最小安全系数的最优化方法 89 本节介绍直接搜索法中的单纯形法和 Powell 法以及牛顿法中负梯度法和 DFP 法 (Davidon−Fletcher−Powell method) 由于这些方法的原理在众多的文献及教科书中都有所介 绍 这里不拟全面阐述其原理 而着重讨论在边坡稳定分析领域需解决的特殊问题 3. 非数值分析方法 第三类方法是在近期计算机发展基础上形成的 称为非数值方法 这一方法被广泛应 用于管理科学 计算机科学 分子物理学及超大规模集成电路设计中 用于解决组合优化 问题 非数值分析利用计算机具有容量大 计算速度快的优点 通过大量随机采样来找到 目标函数的最优值 近期涌现了诸如模拟退火 遗传算法 神经网络和蚂蚁算法等 均属 此类 在边坡稳定分析中 上述算法均有人尝试 本章 4.7 节将对此类方法作一简要介绍 4. 2 任意形状滑裂面的模拟和目标函数的确立 4. 2. 1 任意形状滑裂面的模拟 最优化问题的提法是 对于一个具有 n 个自变量的向量 z = (z1, z2, ..., zn) 确定其目标 函数 F 的最小安全系数 Fm 相应的自变量为 zm 为此 需要对式(4.1)中的曲线 y(x)用若干参数来模拟 也就是说 需要将任意形状滑 裂面 y(x)用 z 来近似表达 将滑裂面曲线用 m 个点 A1, A2, …, Am离散 见图 4.2 也就是将 此 m 个点用直线或光滑的曲线连起来 以近似模拟此曲线 此 m 个点的坐标用 zi (i=1,2, ,m) 表示 = i i i y x z (4.4) 图 4. 2 任意形状滑裂面 一旦这种连接的模式确定 安全系数 F 即可表达成此 m 个点的坐标 x1, y1, x2, y2, .... xn, yn 的函数 F = F( ) x1, y1, x2 , y2 ,Lxm , ym (4.5)
90土质边坡德定分析一原理方法程序 在进行最优化搜索过程中,A1,A2…,Am将移到临界滑裂面的位置B,B2,Bn,见图 42,此处m=6,其中端点A,Am原来在边坡线上,有可能移到边坡线外或内,如图42中 的B1,Bn。为此,需要通过一定的处理方式,分别找到它们和边坡线的交点B1和Bn,仍 以B1,B2,…,Ba作为硏究的对象。对均匀的土质边坡,通常希望滑裂面比较光滑,此时,可 采用三次或更高次的样条函数连接这些点。第11章将详细介绍使用样条函数构筑光滑曲线 滑裂面的方法。当然,也可以采用直线和光滑曲线的组合构筑滑裂面。例如图42,A3,A,A3 A用曲线相联;A2,A3用直线相联。 应用样条函数构筑光滑滑裂面,对于减少自由度,提高数值计算效率具有重要意义 例41]折线和光滑滑裂面比较 如图43所示例,如果用四个点模拟滑裂面,那么,折线滑裂面和曲线滑裂面相应的 最小安全系数分别为1489和1.364,差别颇大。换句话说,如果使用折线模式,要得到与 曲线模式相同的精度,就要更多的离散点,这意味着更多的自由度。而数值分析的机时和 收敛难度是随自由度按指数增加的。值得注意的是,在笔者所知的所有这方面的文献中, 还没有其他研究者使用过曲线模拟任意形状滑裂面 在稳定分析中,对土质比较均匀的边坡,常采用圆弧滑裂面。此时,自变量可以是圆 弧的圆心的x,y坐标和半径r,由于自变量仅三个,最优化方法搜索最小安全系数收敛性能 很好,因此将不作为重点在本章中讨论。 图4.3说明用不同模式模拟滑裂面的差别 滑裂面1-直线模式;滑裂面2-曲线模式 4.1.2目标函数的确定 当滑裂面z的离散模型确定后,安全系数便是m个控制点A1,A2,,A的函数。在优 化计算过程中,这m个点中有n个点各沿某一设定方向β向临界滑裂面移动,或者不规定 方向任其自由移动。其余m-n个点由于问题本身的要求可以固定。滑裂面上任意一点的z 可以用相应于一个初始滑裂面的相对坐标来代表,见图42 =+ 式中:=1,2,n,d为第个点沿β移动的距离
90 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 在进行最优化搜索过程中 A1, A2, …, Am将移到临界滑裂面的位置B1,B2,...,Bm ′ ′ 见图 4.2 此处 m = 6 其中端点 A1, Am原来在边坡线上 有可能移到边坡线外或内 如图 4.2 中 的B1, Bm ′ ′ 为此 需要通过一定的处理方式 分别找到它们和边坡线的交点 B1 和 Bm 仍 以B1,B2,...,Bm 作为研究的对象 对均匀的土质边坡 通常希望滑裂面比较光滑 此时 可 采用三次或更高次的样条函数连接这些点 第 11 章将详细介绍使用样条函数构筑光滑曲线 滑裂面的方法 当然 也可以采用直线和光滑曲线的组合构筑滑裂面 例如图 4.2, A3, A4, A5, A6用曲线相联 A2, A3用直线相联 应用样条函数构筑光滑滑裂面 对于减少自由度 提高数值计算效率具有重要意义 [例 4.1] 折线和光滑滑裂面比较 如图 4.3 所示例 如果用四个点模拟滑裂面 那么 折线滑裂面和曲线滑裂面相应的 最小安全系数分别为 1.489 和 1.364 差别颇大 换句话说 如果使用折线模式 要得到与 曲线模式相同的精度 就要更多的离散点 这意味着更多的自由度 而数值分析的机时和 收敛难度是随自由度按指数增加的 值得注意的是 在笔者所知的所有这方面的文献中 还没有其他研究者使用过曲线模拟任意形状滑裂面 在稳定分析中 对土质比较均匀的边坡 常采用圆弧滑裂面 此时 自变量可以是圆 弧的圆心的 x, y 坐标和半径 r 由于自变量仅三个 最优化方法搜索最小安全系数收敛性能 很好 因此将不作为重点在本章中讨论 图 4. 3 说明用不同模式模拟滑裂面的差别 滑裂面 1−直线模式 滑裂面 2−曲线模式 4. 1. 2 目标函数的确定 当滑裂面 z 的离散模型确定后 安全系数便是 m 个控制点 A1, A2, ..., Am的函数 在优 化计算过程中 这 m 个点中有 n 个点各沿某一设定方向 βi向临界滑裂面移动 或者不规定 方向任其自由移动 其余 m−n 个点由于问题本身的要求可以固定 滑裂面上任意一点的 zi 可以用相应于一个初始滑裂面的相对坐标来代表 见图 4.2 = o i o o i i y x z (4.6) = + ∆ = + i i i o t i o i i d β β sin cos z z z z (4.7) 式中 i=1,2,…,n, di为第 i 个点沿βi移动的距离
第4章确定最小安全系数的最优化方法 对于固定的点,其自由度为零,如图42中的A2,A3,当某点需沿一个软弱夹层移动时, 则其自由度为1。如对该点的移动方向无特殊要求,如图42中的其它点,该点的自由度为 2。问题的自由度为各点自由度的总和,搜索最小安全系数的问题具体化为求下列函数的最 小值问题 F=F(d1,d2,…,dn,阝1,阝2,,阝n) 如前所述,滑裂面的一部分可能是光滑曲线,另一部分受软弱夹层的控制为直线。STAB 程序是这样规定这一定义格式的。设滑裂面被m个点A,A2 Am分为m-1段,从上 交点向下交点编号为1,2,…,m-1。在STAB程序中,用LNO代表此m-1段中为直线段的 线段总数,在数组LOO①,I=12,,LNO中存入这些直线段的编号。其它线段则默认为曲 线。当LNO为零时,滑裂面为一光滑曲线;当LNO为一大于m-1的数时,则程序默认为 全部为直线段。这两种情况都无须填写LOO(D 4.3模式搜索法 4.3.1单形法(朱伯芳等1984) 1.建立初始单形 对某一初始向量z,按下面模式构筑n个向量z(i=1,2,,n,组成m+1个顶点 z=[=+p,=2+q,…2+q z2=[=8+q-2+p…,z0 z=[1+g,22+q,……,m+p 其中 (4.10) (n+1) (4.11) 式中:a为初始步长,可根据实际情况进行调试 2.确定搜索方向 计算单纯形n+1个顶点的目标函数值,比较其大小,从中找出目标函数最大点z和次 大点和最小点z,然后按下式确定下一步的搜索方向 式中上标ν表示迭代次数
第 4 章 确定最小安全系数的最优化方法 91 对于固定的点 其自由度为零 如图 4.2 中的 A2, A3 当某点需沿一个软弱夹层移动时 则其自由度为 1 如对该点的移动方向无特殊要求 如图 4.2 中的其它点 该点的自由度为 2 问题的自由度为各点自由度的总和 搜索最小安全系数的问题具体化为求下列函数的最 小值问题 ( , , , , , ,..., ) F = F d1 d2 L dn β1 β2 βn (4.8) 如前所述 滑裂面的一部分可能是光滑曲线 另一部分受软弱夹层的控制为直线 STAB 程序是这样规定这一定义格式的 设滑裂面被 m 个点 A1, A2, ……, Am分为 m−1 段 从上 交点向下交点编号为 1, 2, …, m−1 在 STAB 程序中 用 LNO 代表此 m−1 段中为直线段的 线段总数 在数组 LOO(I), I = 1,2,…,LNO 中存入这些直线段的编号 其它线段则默认为曲 线 当 LNO 为零时 滑裂面为一光滑曲线 当 LNO 为一大于 m−1 的数时 则程序默认为 全部为直线段 这两种情况都无须填写 LOO(I) 4. 3 模式搜索法 4. 3. 1 单形法 朱伯芳等, 1984 1. 建立初始单形 对某一初始向量 z 0 按下面模式构筑 n 个向量 z i (i=1,2,..,n) 组成 n+1 个顶点 T m n T m T m z q z q z p z q z p z q z p z q z q [ , , , ] [ , , , ] [ , , , ] 0 0 2 0 1 0 0 2 0 1 2 0 0 2 0 1 1 = + + + = + + + = + + + LL LL LL LL z z z (4.9) 其中 a n n n p 2 ( +1) + −1 = (4.10) a n n q 2 ( +1) −1 = (4.11) 式中 a 为初始步长 可根据实际情况进行调试 2. 确定搜索方向 计算单纯形 n+1 个顶点的目标函数值 比较其大小 从中找出目标函数最大点 zH 和次 大点 zG 和最小点 zL 然后按下式确定下一步的搜索方向 = ∑ − − + v H v i n i v n n z z z 1 2 1 (4.12) 式中上标 v 表示迭代次数
92土质边坡稳定分析一原理·方汝猩序 显然,从=到=+2,目标函数是下降的。因此以和=n+2两点的联线作为搜索方向。 3.优化计算 (1)反射。沿搜索方向=和=2点向前再走一步,步长为a(zn+2-zH),到达=+3点 二y+3称为反射点,∞>0称为反射系数。计算反射点的目标函数,并根据它的大小来决 定下一步的走法 (2)扩张。如果反射点的函数值小于z1点的函数值,即 F(zn+3)V(zH),则缩小单 纯形的边长,以最好点zL为顶点,其它各顶点向z移近一半距离,即按下式计算 z}=x+0.5(x}-z),i=0.1,2,,n (418) 得到新的单纯形,转入第(4)步重复计算。 4.收敛判断 按照一定的方式通过反射,扩充和收缩,使单形不断更新逼近极值点。收敛准则为
92 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 显然 从 ν zH 到 ν n+2 z 目标函数是下降的 因此以 ν zH 和 ν n+2 z 两点的联线作为搜索方向 3. 优化计算 (1) 反射 沿搜索方向 ν zH 和 ν n+2 z 点向前再走一步 步长为 ( ) 2 v H v n z − z α + 到达 ν n+3 z 点 ( ) 3 2 2 v H v n v n v nz = z + z − z + + α + (4.13) ν n+3 z 称为反射点 α>0 称为反射系数 计算反射点的目标函数 并根据它的大小来决 定下一步的走法 (2) 扩张 如果反射点的函数值小于 zL点的函数值 即 ( ) ( ) v F F L z z v n+3 1 为扩张系数 如果 ( ) ( ) 4 v L v F zn V z + (4.16) 表明反射点走得太远 需要缩回去 按下式计算 ( ) 5 2 2 v n v H v n v n+ = + + − + z z β z z (4.17) 式中 0V z + 则缩小单 纯形的边长 以最好点 Lz 为顶点 其它各顶点向 Lz 移近一半距离 即按下式计算 i n v L v i v i z = z + 0.5(z − z ), = 0,1,2,..., v L (4.18) 得到新的单纯形 转入第(4)步重复计算 4. 收敛判断 按照一定的方式通过反射 扩充和收缩 使单形不断更新逼近极值点 收敛准则为 ∑ [ ] − ≤ ε + + − 2 2 0 ( ) ( ) 1 1 v n v i n i v v n z z (4.19)
第4章确定最小安全系数的最优化方法 式中:ε为要求的计算精度。如式(419)满足,则结束计算,并以z作为极小点。否则, 置ν=ν+1,转至第3步,重复计算 据已有经验,可取ax=1;,04≤β≤0.6,2.0≤y≤3.0。经 Nelder和Mead论证,为使单纯形 适应函数的性态及便于收敛,a不宜比1大很多,而α<1的计算次数比a=1要多,故折衷 取α=1。另外,认为β的数值变化对搜索效率的影响比γ要大。他们推荐采用 28≤y≤3.0 421) [例42]说明单形法计算过程例 为了形象直观地了解单形法在搜索最小安全系数时的工作状况,我们考察图44所示 的一个有两个自由度的例子滑裂面。它由ABC组成。计算时令C点固定不动,A,B两点 沿水平线移动,则该滑裂面的安全系数F由A点的x坐标x1和B点的x坐标x决定。图45 示F相应x1,x2的等值线图。根据枚举法可以发现在x1=92.0,x2=143.0时安全系数获得最小 值1257,相应临界滑裂面如图44中标5的那个滑裂面 图4.4具有两个自由度的算例 1,2,3,4-初始滑裂面;5-临界滑裂面 如果使用单形法,则按式(49)初始生成的单形即三个滑裂面如图45在左下角三角形 所示,最大、次大和最小分别为A、B、C三点。D为按式(412)所得的新的顶点。AD代 表了两点联线的方向,D代表了第一次反射和扩张后达到的点。第一次迭代后,形成了B、 C、D构成的新的单形。开始新的一轮迭代。依次循环,最终达到安全系数的极值E。搜索 过程如图45中折线1所示,最终收敛到Fm=1257相应的z=(92.00,14300y
第 4 章 确定最小安全系数的最优化方法 93 式中 ε为要求的计算精度 如式(4.19)满足 则结束计算 并以 v Lz 作为极小点 否则 置 v=v+1 转至第 3 步 重复计算 据已有经验 可取α=1; 0.4 ≤ β ≤ 0.6; 2.0 ≤ γ ≤ 3.0 经 Nelder 和 Mead 论证 为使单纯形 适应函数的性态及便于收敛 α不宜比 1 大很多 而α < 1 的计算次数比α=1 要多 故折衷 取α=1 另外 认为β的数值变化对搜索效率的影响比γ要大 他们推荐采用 0.4 ≤ β ≤ 0.6 (4.20) 2.8 ≤ γ ≤ 3.0 (4.21) [例 4.2] 说明单形法计算过程例 为了形象直观地了解单形法在搜索最小安全系数时的工作状况 我们考察图 4.4 所示 的一个有两个自由度的例子滑裂面 它由 ABC 组成 计算时令 C 点固定不动 A, B 两点 沿水平线移动 则该滑裂面的安全系数 F 由 A 点的 x 坐标 x1和 B 点的 x 坐标 x2决定 图 4.5 示 F 相应 x1, x2的等值线图 根据枚举法可以发现在 x1=92.0, x2=143.0 时安全系数获得最小 值 1.257 相应临界滑裂面如图 4.4 中标 5 的那个滑裂面 图 4. 4 具有两个自由度的算例 1, 2, 3, 4−初始滑裂面 5−临界滑裂面 如果使用单形法 则按式(4.9)初始生成的单形即三个滑裂面如图 4.5 在左下角三角形 所示 最大 次大和最小分别为 A B C 三点 D 为按式(4.12)所得的新的顶点 AD 代 表了两点联线的方向 D 代表了第一次反射和扩张后达到的点 第一次迭代后 形成了 B C D 构成的新的单形 开始新的一轮迭代 依次循环 最终达到安全系数的极值 E 搜索 过程如图 4.5 中折线 1 所示 最终收敛到 Fm=1.257 相应的 z m=(92.00, 143.00)T
94土质边坡德定分析一原理方法程序 安全系数等值线 初始单 x1,A点x坐标(m) 图4.5对图44示例用单形法计算过程 D (D2,D3) C1(C2,C3) 软弱夹层c=9.8kPa,∮=1 图4.6具有软弱夹层的算例 [例4.3]使用单形法计算最小安全系数 图46具有分层剖面和一个倾斜的软弱夹层。图中示4个不同的初始滑裂面。点A和 B可以不受限制地自由移动,故自由度为2。同时,此两点有意未设在边坡面上,以说明 程序在搜索过程中将进行判断和调整,自动找到滑面与坡面的交点。C和B点位于夹层上 因此,自由度为1。该两点在搜索过程中将沿一定方向β移动,β为该夹层的倾角。图47中 滑裂面3为采用单纯形法获得的临界滑裂面,相应安全系数为1025
94 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 图 4. 5 对图 4.4 示例用单形法计算过程 图 4. 6 具有软弱夹层的算例 [例 4.3] 使用单形法计算最小安全系数 图 4.6 具有分层剖面和一个倾斜的软弱夹层 图中示 4 个不同的初始滑裂面 点 A 和 B 可以不受限制地自由移动 故自由度为 2 同时 此两点有意未设在边坡面上 以说明 程序在搜索过程中将进行判断和调整 自动找到滑面与坡面的交点 C 和 B 点位于夹层上 因此 自由度为 l 该两点在搜索过程中将沿一定方向β移动 β为该夹层的倾角 图 4.7 中 滑裂面 3 为采用单纯形法获得的临界滑裂面 相应安全系数为 1.025
第4章确定最小安全系数的最优化方法 C1(C2,C3,C1) B2,B3) 图4.7对图46示例用最优化方法计算过程 临界滑裂面:1-改良后的DFP法,F=1.025;2-负梯度法,F=1.025 3-单形法,F=1.025;4-DFP法,F=1.025 4.1.2Powe法 1.基本原理 有关此法的详细原理可参阅文献(朱伯芳等,1984)。其基本算法如下。 (1)给定初始点,并令 S=(0.0,,1.,0)7 (422) 式(422)右端的1处于第i位,i=1,2,,n (2)相应迭代步i=1,进行一次一维搜索,第一次搜索方向P=S,使 (乙1-1+axP)=极小 并令 P 然后令i=计1,重复上述过程,直至i=n为止。 (3)取P=p1,=12,,n1,并令 Pn=2n-20 (4.25) 通过这一步,增加zn-,丟弃原来的坐标方向P1。 (4)沿xn-z作一次一维搜索,求xn,使x2+anPn)=极小,并令 然后转向第(2)步。 (5)重复进行上述运算n次,即沿n个互相共轭的方向都进行了一维搜索,然后作如 下判断 kmP
第 4 章 确定最小安全系数的最优化方法 95 图 4. 7 对图 4.6 示例用最优化方法计算过程 临界滑裂面 1−改良后的 DFP 法 Fm=1.025 2−负梯度法 Fm =1.025 3−单形法 Fm =1.025 4−DFP 法 Fm =1.025 4. 1. 2 Powell法 1. 基本原理 有关此法的详细原理可参阅文献 朱伯芳等 1984 其基本算法如下 (1) 给定初始点 z0 并令 i T S = (0,0,...,1,...,0) (4.22) 式(4.22)右端的 1 处于第 i 位 i=1,2,…,n (2) 相应迭代步 i=1 进行一次一维搜索 第一次搜索方向 pi=Si 使 V(zi−1 +αi pi) = 极小 (4.23) 并令 i i i i z = z +α p −1 (4.24) 然后令 i=i+1 重复上述过程 直至 i=n 为止 (3) 取 pi=pi+1, i=1,2,…,n−1 并令 0 p = z − z n n (4.25) 通过这一步 增加 0 z − z n 丢弃原来的坐标方向 1p (4) 沿 zn−z0作一次一维搜索 求αn 使 V(zn+αn pn)=极小 并令 n n pn z0 = z +α (4.26) 然后转向第(2)步 (5) 重复进行上述运算 n 次 即沿 n 个互相共轭的方向都进行了一维搜索 然后作如 下判断 α < ε n n p (4.27)
96土质边坡德定分析一原理方法程序 如上述判据满足则结束计算,此时x0即为所求极小点,否则转向第(1)步。 2.应用 Powell法时所做的改进 在采用 Powel法进行计算时,第一步需要确定n个独立的搜索方向。在大多数情况下, 都采用与n个自变量坐标轴方向一致的单位向量,即由式(422)定义的S。 这种处理方式并不适于边坡稳定分析。参见图48(a)所示一个有两个自由度例子。按 Powell法,第一步固定B点不同,让A 点沿水平方向移动到A,使安全系数达 到极小值,显然能够移动的范围是很有限 的。第二步固定A不动,B在其移动方向 到达B’,再次得一个极小值,其移动范围 更为有限,这样的搜索效率很低,当自由 度较多时,事实上无法实现有效的搜索。 设想如果在A点移动时,让B点也向前 第二次迭代 移动一定值,第一次一维搜索就可以使初 第一次迭代值z向极值zm迈进一大步,如图48(b 因此,在构筑S时,我们的思路是,当某 控制点在移动方向搜索最小值时,要求 其邻近的点也按一定的比例移动,使滑裂 图.8在边坡稳定分析中使用 Powell法作的改进整体而不是局部在优化计算中不断地调 (a)未作改进;(b)改进后 整,以迅速地达到极值点。 为此按下式构筑n个向量S(=12,3,,n) (428) 按下式确定e(=2,3,,n) cs(B/-B)~布 当≤时 cost(B-B.)I-31 当i=1时,总是用(429)的下式;i=n时,总是用式(429)的上式。c的物理意义是当 第i个变量沿其搜索方向移动一个单位值时,邻近的点j都将在其移动方向获得一个增量 e,这个增量值取决于第j点距端点(与j同侧的那个)的水平距离,在=0和i=n时等 于零,ij时等于1,其余点按线性分布原则确定内插。参见图49
96 土质边坡稳定分析 原理 ⋅ 方法 ⋅ 程序 如上述判据满足则结束计算 此时 z0即为所求极小点 否则转向第(1)步 2. 应用 Powell 法时所做的改进 在采用 Powell 法进行计算时 第一步需要确定 n 个独立的搜索方向 在大多数情况下 都采用与 n 个自变量 zi坐标轴方向一致的单位向量 即由式(4.22)定义的 Si 这种处理方式并不适于边坡稳定分析 参见图 4.8(a)所示一个有两个自由度例子 按 图 4. 8 在边坡稳定分析中使用 Powell 法作的改进 (a) 未作改进 (b) 改进后 Powell 法 第一步固定 B 点不同 让 A 点沿水平方向移动到 A′ 使安全系数达 到极小值 显然能够移动的范围是很有限 的 第二步固定 A′不动 B 在其移动方向 到达 B′ 再次得一个极小值 其移动范围 更为有限 这样的搜索效率很低 当自由 度较多时 事实上无法实现有效的搜索 设想如果在 A 点移动时 让 B 点也向前 移动一定值 第一次一维搜索就可以使初 值 z0 向极值 z m 迈进一大步 如图 4.8 (b) 因此 在构筑 Si 时 我们的思路是 当某 一控制点在移动方向搜索最小值时 要求 其邻近的点也按一定的比例移动 使滑裂 面整体而不是局部在优化计算中不断地调 整 以迅速地达到极值点 为此按下式构筑 n 个向量 (i 1,2,3,..., n) i S = Τ = ( , ,..., ) 1 2 i n i i i S e e e (4.28) 按下式确定e (i 1,2,3,..., n) i j = > − − − ≤ − − − = 当 时 当 时 i j x x x x i j x x x x e n i n j j i i j j i i j cos( ) cos( ) 1 1 β β β β (4.29) 当 i=1 时 总是用(4.29)的下式 i=n 时 总是用式(4.29)的上式 i j e 的物理意义是当 第 i 个变量沿其搜索方向移动一个单位值时 邻近的点 j 都将在其移动方向获得一个增量 i j e 这个增量值取决于第 j 点距端点 与 j 同侧的那个 的水平距离 在 i=0 和 i=n 时等 于零 i=j 时等于 1 其余点按线性分布原则确定内插 参见图 4.9