基于既往优化知识的稳健性优化方法

吕国臣，程远胜，易家祥，刘均

华中科技大学 船舶与海洋工程学院，湖北 武汉 430074

0 引 言

产品设计会存在很多不确定性因素[1]，例如船舶结构的材料参数、尺寸、载荷。这些不确定性因素的波动可能使产品性能参数发生变化，导致产品性能的波动量超过允许范围而使产品不满足设计要求。稳健性设计优化问题最早由Taguchi[2]提出，目的是在存在不确定性的情况下，使目标/约束函数的值保持在可接受的范围内，同时获得尽可能好的设计问题的解。学者们也对不同问题做 了 工 程 应 用 研 究[3-4]。

不确定性问题主要包括模型不确定性、参数不确定性以及预测不确定性等[5]。本文针对设计变量或参数在区间上存在的不确定性进行研究。传统解决非线性区间不确定性稳健优化问题的方法通常为嵌套优化方法，这会导致大量的资源消耗。Siddiqui 等[6]提出了一种修改的benders 分解方法来降低计算量，即通过离散化不确定性区间来解耦嵌套优化。但由于该方法对不确定性区间进行了离散，并不能保证优化解的稳健性。Zhou等[7]提出了一种区间不确定性鲁棒优化的序列二次规划方法，该方法兼具目标鲁棒性和可行性鲁棒性。Hu 等[8]提出了一种基于系统性能对不确定性的灵敏度雅可比矩阵的有界约束问题稳健性指标计算方法，雅可比矩阵可以定量地测量系统允许的最大变化幅度。辛鹏飞等[9]运用基于切比雪夫多项式的区间扩张函数，将含区间参数的微分代数方程转化为切比雪夫多项式插值点处确定参数的动力学方程，减少了评估稳健性的仿真工作量。王琼等[10]利用区间可能度处理不确定性约束，将稳健性优化模型转换为确定性多目标优化模型。这些方法虽然能够评估稳健性，但都需要利用函数的梯度信息，对于需要采用数值仿真方法计算的实际工程问题，很难利用其梯度信息，所以基于梯度的求解方法不适用。Gunawan 等[11-12]根据逆向思维和最坏情况分析，提出了反映设计方案稳健性的指标，对其施加约束，并作为约束条件加入到优化数学模型中，构成稳健性设计优化问题；但该方法采用了逆向思维来反向映射自变量和函数区间，难以保证稳健优化解的精度。所以对于变化趋势不确定的黑箱工程优化问题，基于嵌套进化算法能够有效解决此问题，但会消耗大量的计算资源。Cheng 等[13]将局部搜索算法应用到嵌套优化中加速其收敛，但所需的计算资源仍然巨大。

根据前人文献的算例经验，若直接利用基于嵌套进化算法来求解稳健性优化问题，通常会调用百万次甚至千万次的原函数，在进行产品设计如船舶结构优化时会消耗大量的计算资源和时间。故本文提出一种基于既往优化知识的稳健性优化方法，以大幅度减少计算资源，提高优化效率。

1 问题描述及求解

1.1 稳健性优化数学模型

不失一般性地，区间稳健性设计优化问题的数学模型可以由式(1)表示：

其中：

式中：x为D维设计变量，xd为 确定性变量；xu为不确定性变量；x和x分别为确定性变量的下界和上界；x和x分别为不确定性变量的下界和上界； δxu为 不确定性变量的波动值； Δxu为不确定性变量的最大波动范围； ηf为目标稳健性指标； ηg为约束稳健性指标； Δf0为目标函数最大允许波动量；f为目标函数；gm为小于等于形式的约束函数。

需要指出的是，对于确定性优化，其约束条件可以统一写成gm≤0 的 形式，当 δxu=0时，式(1)中约束稳健性指标 ηg≤0就表示确定性优化的约束条件。目标稳健性指标为不确定性变量变化下的目标函数的最大波动量与最大允许波动量的比值，其约束稳健性指标为不确定性变量变化下的所有约束函数的最大值。

1.2 基于差分进化的区间不确定性稳健性优化问题求解

1.2.1 差分进化算法

差分进化算法（DE）是工程优化领域应用较多的一种智能进化算法，它由Storn 等[14]提出。该算法能够较好地对函数进行全局搜索，被广泛应用到工程优化问题中。差分进化算法存在多种进化格式，本文采用DE-rand/1/bin 格式。其主要步骤如下：

通常，对于一个种群数为N的种群，第G代的第i个体可以表示为=(p,···,p)，i=1，2，···，N，其下一代对应的变异个体由式(2)产生：

式中：v为（G+1）代变异个体；F为放大因子，是定 义 在 [0,2]范 围 内 的 常 数 因 子；k1，k2和k3为 从{1,2,···,N}中随机产生的3 个整数，且要满足k1≠k2≠k3≠i。

然后，对种群个体进行交叉操作，遵循式（3）的机制将变异个体与当前代的个体交叉产生下一代的新个体。

式中：CR为 在 [0,1]范 围内的交叉因子；rand(j)为[ 0,1]区 间内产生的随机数；jRand为从{1,2,···,N}中随机产生的一个整数。

最后，选择x与x中目标函数较好（或其他准则）的一个进化到下一代。

1.2.2 求解方法

针对式（1）表征的区间不确定性稳健性优化数学模型，一种常规的求解方法就是直接应用差分进化算法进行嵌套优化求解（DE-RO）。外层优化过程是基于差分进化求解最优目标函数值，内层优化是差分进化求解2 个稳健性指标时需要进行2 个无约束的优化，即目标稳健性指标和约束稳健性指标。本文内外层优化均采用相同的进化格式。

2 基 于 既 往 优 化 知 识 的 稳 健 优 化方法

若应用差分进化算法求解方法，在嵌套优化过程中会进行内层优化，从而大量增加目标函数和约束函数的调用次数。本文提出基于既往优化知识的稳健优化方法（DE-BPOK-RO），该方法涉及稳健性指标的近似评估以及基于稳健性指标重新评估的临界距离调整2 个核心策略。

2.1 稳健性指标的近似评估

假设一个稳健性优化问题的不确定性变量的数目为2，画出不确定性变量变化区间的示意图，如图1 所示。图中，p1和p2为2 个外层个体，蓝色点和绿色点分别为p1和p2进行内层优化生成的个体。

由图1 可见，若p1和p2的空间距离较近，两者变化区间的重叠区域较大，则其内层优化生成的个体之间会出现距离极近甚至重合的情况。在内层优化过程中，若能够利用已经精确计算过的个体响应值的信息来近似预测其他个体响应值，则可在很大程度上降低原函数的调用次数。

图1 不确定性变量变化区间示意图Fig. 1 Schematic diagram of the variation interval of uncertain variables

建立精确计算过样本点库Dp={p1,p2,···}及对应的响应值库Dt={t(p1),t(p2),···}。近似评估稳健性指标策略的核心是近似评估内层样本点的响应值，需要定义2 个参数，分别为近似距离阈值R及上限近似点数Ns，Ns的 引入是为了降低R的取值对近似的影响。图2 为近似评估内层个体响应值的示意图，图中蓝色点表示精确计算过的样本点，A，B为2 个个体，其响应值计算公式为：

图2 近似评估内层个体响应值示意图Fig. 2 Schematic diagram of approximate evaluation of inner individual response values

式中：tˆ(p)为 被评估个体的响应值；q为近似距离阈值内精确计算的样本点数；pl为近似距离阈值内样本点库中经距离排序后的样本点。需要指出的是，用来近似评估内层个体响应值的目标函数库和约束函数库均通过外层和内层在既往优化过程中的信息进行更新，是一个动态增加过程。

2.2 临界距离调整策略

在近似评估稳健性指标时，若进化前期样本点库的信息不足，使内层样本点近似得到的响应值误差较大，这可能导致实际满足稳健性约束条件的部分外层个体被误评估为不稳健，从而误导了整体进化过程。而随着代数增加，计算的样本点增加后，重新近似评估其稳健性指标会在一定程度上修正误判解，基于此提出了一种基于稳健性指标重新评估的临界距离调整策略。

在进行第G代外层优化时，建立第G代修正库D={p|f(p)＜f∧gj(p)≤0 ∧(ηf(p)＞1‖ηg(p)＞0)，,p∈WG}， 其中f为第G-1 代的最优目标值，WG为截止到G代每代外层个体的集合。在判断当前代最优之前，随机选取修正库中N个点重新近似评估稳健性指标，定义参数修正率阈值RC0，重新评价后稳健性变化（即由不稳健变为稳健）的个数为Nc，将其与当前种群中min(Nc,N·RC0)（取整）数目的随机个体进行替换。若RC=＞RC0，则近似距离阈值R=R/2 。需要指出的是，用于临界距离调整的修正库在优化的过程中不断变化，其建立和更新是基于外层既往优化过程中可能搜索到的真正最优解附近的信息。

2.3 改进的可行性准则

在应用基于既往优化知识的稳健性优化方法求解稳健性优化问题时，新种群个体的选择策略采用改进的可行性准则，用新生成的个体代替上一代个体需要符合以下几种情况：

1） 新生成的个体是可行解，而原个体是非可行解；

2） 新生成的个体与原个体均是非可行解，而新生成个体约束违反的总体程度更小；

3） 新生成的个体与原个体均是可行解，且新生成的个体满足稳健性约束条件而原个体不满足；

4） 新生成的个体与原个体均是可行解，且新生成的个体和原个体均满足稳健性约束条件而新生成的个体的目标函数值更小；

5） 新生成的个体与原个体均是可行解，且新生成的个体和原个体均不满足稳健性约束条件而新生成的个体的稳健性指标的和更小。

2.4 算法求解框架及流程

算法的基本框架是基于差分进化求解策略的双层循环模式，其中内层优化用于求解外层稳健性指标约束。具体流程图如图3 所示。

图3 基于既往优化知识的稳健性优化方法流程图Fig. 3 Flow chart of robust optimization method based on previous optimization knowledge

步骤1：初始化外层优化的种群。

步骤2：调用原函数精确计算初始外层个体的目标函数f(x)和 约束函数gj(x)。

步骤3：评估初始外层个体的目标稳健性指标 ηf和约束稳健性指标 ηg，即用差分进化算法求解外层个体在不确定区间内目标函数的最大变化量和最大约束函数值（本质是分别以 ηf和 ηg为内层目标函数的无约束优化问题）。

步骤4：用所有精确计算过的点（设计方案）建立初始目标函数库和约束函数库。

步骤5：根据DE-rand/1/bin 的差分进化格式生成新一代外层个体。

步骤6：调用原函数精确计算新一代外层个体的目标函数f(x)和 约束函数gj(x)。

步骤7：根据2.1 节提出的策略近似评估新一代外层个体的目标稳健性指标 ηf和约束稳健性指标 ηg。

步骤8：更新目标函数库和约束函数库。

步骤9：结合2.2 节提出的临界距离调整策略对当代外层个体进行选择并判断当前最优。

步骤10：精确评估当前最优个体的稳健性指标。

步骤11：更新目标函数库和约束函数库。

步骤12：若满足停止条件，则结束外层优化；否则，返回步骤5。

需要指出的是：

1） 只需要对满足gj≤0 的个体进行稳健性指标评价，这在一定程度上减少了计算资源的消耗。

2） 精确评估稳健性指标就是在内层优化中计算所有内层个体的响应值。

3） 对于外层优化过程，其停止准则为进化代数达到了最大进化代数；对于内层优化过程，其停止准则为进化代数达到了最大进化代数或最优值在 μ （一般取20～30）代的差值小于等于10e-3。

3 算例测试结果及分析

采用2 个数学算例和1 个工程算例验证DEBPOK-RO 算法的可行性与有效性。

3.1 算例测试结果

3.1.1 算例1

算例1 是一个目标函数比较复杂的测试函数，其数学表达式为

式中：x1为 不确定性变量，其变化区间为[-0.05,0.05]；x2为确定性变量；目标函数允许的最大波动值 Δf取为0.02。

分别运用DE-BPOK-RO 和DE-RO 算法求解以上算例，并与Cheng 等[13]提出的混合进化方法（DE-SQP-RO）得到的结果比较。对提出的DEBPOK-RO 算法中3 个关键参数进行研究发现，采用表1 所示的参数能得到最佳结果。

表1 DE-BPOK-RO 和DE-RO 算法参数Table 1 Parameters of DE-BPOK-RO and DE-RO algorithms

其中，R0=Δx1/5=0.01。为了降低智能优化算法的不确定性，DE-BPOK-RO 和DE-RO 算法都独立运行了10 次，统计结果如表2 所示，其中FE为所有原函数调用的次数。根据表2，在目标函数最优的结果中，DE-BPOK-RO 算法求解得到的目标函数值与DE-RO 算法得到的目标函数值的误差约为0.072%，与DE-SQP-RO 算法得到的误差约为0.05%，而调用的函数次数仅为另外2 种算法的2.393% 和13.385%；10 次运行目标函数的中位数值要优于DE-RO 算法得到的结果，而原函数的平均调用次数仅为差分进化求解策略的2.620%；最差一次目标函数值虽然劣于DERO 算法得到的结果，但2 种算法得到的最优解均与全局的稳健最优解差距很大，可能由初始种群的随机性导致。

表2 测试函数1 统计结果Table 2 Statistical results of test function-1

单次优化过程的收敛情况如图4 所示。由图可见，2 种算法都收敛到最优值附近。图4 还展示了算法收敛到最优的局部情况。次优化中DEBPOK-RO 算法在67 代时收敛，最优值为-5.952 9，此时原函数调用次数为191 484；而DE-RO 算法要在78 代才能达到相同程度，原函数调用次数为5 259 690，直观展示了DE-BPOK-RO 算法的优势。

图4 算例1 的DE-BPOK-RO 和DE-RO 算法收敛曲线Fig. 4 Convergence curves of DE-BPOK-RO and DE-RO algorithms of numerical example-1

通过以上结果分析，基于既往知识的稳健性优化方法在最优目标误差小于2.4%或更优的情况下，大大减少了真实函数调用的次数，节省了97%以上的计算资源。

3.1.2 算例2

算例2 是一个约束函数比较复杂的测试函数，其数学表达式为

式中：不确定性变量x1和x2的变化区间均为[-0.4,0.4]； 目标函数允许的最大波动值 Δf=2.5。

采用与算例1 相同的参数，其中，R0=min(Δx1,Δx2)/5=0.08。采用DE-RO 算法和DE-BPOK-RO算法对此问题进行10 次独立的优化求解，得到的统计结果如表3 所示，单次优化过程的收敛情况如图5 所示。

表3 测试函数2 的统计结果Table 3 Statistical results of test function-2

图5 算例2 的DE-BPOK-RO 和DE-RO 算法收敛曲线Fig. 5 Convergence curves of DE-RO-DE and DE-RO algorithms of numerical example-2

根据表3，在目标函数最优的结果中, DEBPOK-RO 算法求解得到的目标函数值与DE-RO算法得到的误差约为0.243%，而调用函数次数仅为DE-RO 算法的0.956%；10 次运行的目标函数的中位数值的误差约为2.273%，平均调用原函数的次数仅为DE-RO 算法的0.870%；最差一次目标函数值虽然劣于用DE-RO 算法得到的结果，但2 种算法得到的最优解均与全局的稳健最优解差距很大，表明可能是初始种群的分布对进化的影响较大。

从图5 中可以看出，2 种算法都收敛到最优值附近。根据收敛曲线的局部图，DE-BPOK-RO算法在28 代时收敛，最优值为-1.768 5，原函数调用次数为128 609；而DE-RO 算法要在38 代才能收敛到相同程度，而原函数调用次数为4 759 410，同样证明了DE-BPOK-RO 算法的优势。

由以上分析可知，基于既往知识的稳健性优化方法在最优目标误差小于2.3%的情况下，节省了99%以上的计算资源。

3.2 压力容器设计

选用的工程算例来源于经典的压力容器设计（Kannan 等[15]），如图6 所示。该压力容器也可作为水下结构物内置液舱，获得对设计变量波动不敏感的设计，即稳健性设计对确保这种内置液舱的安全性是非常重要的。压力容器设计经过修改后成为一个稳健性优化问题，本文将比较差分进化求解策略和基于既往优化知识算法对此问题的优化结果，并结合混合优化方法进行对比分析。

图6 压力容器的中心和端部示意图Fig. 6 Schematic diagram of the center and end of pressure vessel

原问题的目标是最小化包括圆柱容器材料、成型和焊接成本的总成本，1 个设计变量分别是x1(Ts，壳体厚度)、x2(Th， 头部厚度)、x3(R，内径)和x4(L，容器圆柱形部分的长度，不包括头部)。修改后的优化问题如式(8)所示：

式中：不确定性变量x4为容器圆柱形部分的长度，其变化区间为 [-0.05,0.05]；其他变量为确定性变量，目标函数（总体设计成本）允许的最大波动值Δf=50。

采用与数学算例相同的算法参数，其中R0=Δx4/5=0.01。用差分进化求解策略和基于既往优化知识算法对此问题进行10 次独立的优化计算，得到的统计结果如表4 所示。

表4 工程算例统计结果Table 4 Statistical results of engineering example

根据表4，最好的优化结果中，DE-BPOK-RO算法求解得到的目标函数值优于DE-RO 算法，与DE-SQP-RO 算法得到的最优值误差仅为0.008%，而调用的函数次数仅为另外2 种算法的5.655%和8.229%；10 次运行的目标函数的中位数值也要优于用DE-RO 算法得到的结果，原函数的平均调用次数也仅仅为DE-RO 算法的5.936%；最差一次目标函数值虽然劣于用DE-RO 算法得到的结果，但出现的概率较低，且误差约为4.193%，在工程优化可接受范围内。通过以上分析可知，基于既往知识的稳健性优化方法在优于或者极少程度牺牲最优解的条件下，大幅减少了真实函数调用次数，对于工程优化问题具有较大意义。

在实际的船舶优化设计中，不确定性因素的存在可能会导致船舶性能达不到预定的设计目标。直接应用嵌套优化算法求解基于区间不确定性的船舶结构稳健性优化问题显然不可行，本文提出的基于既往知识的稳健性优化方法虽大幅度节省了计算资源，但几十万次的有限元计算对船舶结构优化来说仍存在一定压力。

刘婧等[16]将动态代理模型应用到船舶结构可靠性优化中，在未来工作中将把代理模型及本文提出的方法结合进一步降低计算成本，研究提出船舶结构强度或稳定性稳健优化设计方法。

4 结 语

本文针对区间不确定性稳健性优化问题提出了一种基于既往优化知识的稳健性优化方法。该方法基于嵌套差分进化算法的格式，利用临界距离内已经精确计算过的个体响应值近似评估内层个体响应值，以此近似评估个体稳健性指标来大幅度地减少原函数的调用次数，降低计算成本。同时利用进化过程中逐步扩充的精确个体响应值的信息，选择性地重新评估已往个体的稳健性，并根据评估误判率自适应地调整临界距离来降低稳健性误判率。

本文通过2 个数学算例和1 个工程算例验证了算法的有效性。通过与DE-RO 及DE-SQP-RO算法的综合统计结果对比得出，DE-BPOK-RO 算法在最优目标误差为2.5% 以下的情况下，节省了94%以上的计算资源，对工程优化问题也取得了比较好的效果。因此，提出的基于既往优化知识的稳健性优化方法，为解决区间不确定性的产品稳健性优化提供了新思路，大幅节省了基于双层嵌套产品稳健性优化所需要的计算资源。