潘志榕， 朱光宇， 伊德景

(福州大学机械工程及自动化学院， 福建 福州 350116)

0 引言

传统的螺旋桨设计， 主要根据实际需求查找相应的螺旋桨图谱， 其设计的桨型只是一组固定数值， 叶型误差较大且设计过程复杂. 已有研究人员将遗传算法等优化算法运用于螺旋桨的设计中. Mirjalili等[1]采用遗传算法设计船用螺旋桨， 改变螺旋桨设计参数， 改善了螺旋桨效率和空泡等问题. Xie[2]使用遗传算法对螺旋桨参数进行优化， 以扭力系数和敞水效率的最大值作为目标函数， 最后得到两组分别有利于扭力和效率的计算值， 根据设计要求及运行条件对二者进行折衷选优求得最优解. 叶礼裕等[3]为提高螺旋桨优化设计的质量和效率, 通过将试验设计方法、 近似模型和粒子群算法相结合, 建立了基于近似模型的螺旋桨优化设计方法. 为提高螺旋桨的水动力性能， 研究人员尝试采用多种优化方法， 但都是基于等效的做法， 即通过修改螺距和剖面的拱度使三维桨各半径处的等效拱线与二维剖面的拱线一致. 这类做法需对剖面拱度和螺距同时进行修改且需反复迭代， 各家采取的等效方式存在差异， 迭代调整的方法不同， 存在很大的经验性和不确定性. 另一方面， 多目标遗传算法多采用Pareto支配[4]与启发式技术来指导算法搜索， 在优化领域取得较好的效果， 但也存在一些不足： 算法参数需要根据经验来进行调整， 另外， 目前的算法很难在Pareto解集的分布性和收敛性两方面同时取得满意的效果. 由于多目标优化的输出结果为一组解集， 因而不仅要考虑Pareto解集与Pareto解的逼近程度， 还要考虑到整个解集在目标空间分布的广泛程度和均匀程度， 解决这些问题的关键是确定Pareto解集的评价指标. 针对评价指标， 在遗传算法中引入超体积准则概念. Schütze等[5]验证了超体积最大化理论， 即当算法所得解集的Pareto前沿越逼近真实Pareto前沿， 并且拥有相对均匀和广泛的分布时， 对应的超体积指标越大.

本研究提出结合超体积准则[6]的多目标遗传算法用于螺旋桨水动力性能优化. 集成基于面元法的螺旋桨性能计算， 大大减少螺旋桨性能的计算时间. 在遗传算法中加入超体积评价准则， 保障Pareto解集的分布性和收敛性， 解决之前优化算法存在的问题. 以降低阻升比、 桨型面积差和压力方差为优化目标， 将无法直接比较的三个子目标结果转变为具体数值， 减少计算次数和时间的同时将流场计算程序与面元法结合， 实现对螺旋桨叶片形状的优化设计.

1 螺旋桨叶片力学性能模型

采用面元法进行螺旋桨的数值计算. 面元法通过在桨型表面分布面涡或面源求解桨型的运动特性， 求解过程中伴有直均流叠加， 它是当前求解绕流问题的主要分析方法.

图1 螺旋桨叶片轮廓Fig.1 Shape of the propeller blades

螺旋桨叶片轮廓如图1所示. 图中剖出9个叶剖面， 是通过与螺旋桨共轴的圆柱面和桨叶相截所得的截面， 从螺旋桨中心往外依次在半径为0.2R、 0.3R、 0.4R、 0.5R、 0.6R、 0.7R、 0.8R、 0.9R(R为螺旋桨半径). 螺旋桨求解剖面绕流问题时， 首先生成代表叶剖面的精确坐标值， 并使之有光顺的函数表达式. 之后， 对每个叶剖面， 沿着曲线方向， 按照确定的面元之间的距离， 得出面元的总数. 面元的划分并不固定, 可以是等分也可以是非均匀的划分， 本研究是对所有叶剖面从前缘到尾缘使用余弦函数进行分割. 如图2所示， 在叶剖面弦线上按下列余弦方程取结点xj：

(1)

其中：C为叶剖面的弦长；m为面元切割数量. 这样所构成的面元分布， 在前缘和尾缘附近比较集中， 正好是性能分析主要部分. 对每个面元上布置相应的涡和偶极子， 或等价的面源， 建立相应的速度求解方程[7].

图2 0.3R叶剖面上的面元分布和控制点坐标示意图Fig.2 Schematic diagram of the surface element distribution and the control point coordinates on the 0.3R wing section

0.3R处叶剖面如图2所示. 每个小段近似直线， 形成m个面元，mi表示第i个面元. 每个面元布置恒定的涡密度为γj(j=1, 2, …,m)， 攻角为α， 水速为VA的来流叠加， 每个面元上的中点为该面元的控制点. 计算螺旋桨升力和阻力时， 添加边界条件， 使尾缘和前缘部分上下两个面元上压力相等， 此条件得到的方程是非线性方程组， 为此增加近似的线性Kutta条件， 该条件需速度为固定的才符合速度势ψ的计算公式， 这样所建立的方程是线性方程[3]. 通过计算控制点处的速度值， 可建立满足边界条件的m个方程组, 另加Kutta条件得到m+1个方程组, 计算过程如下.

任一面元j上的涡密度对面元i的控制点的坐标(xi,yi)作用的流函数[9]为：

(2)

其中：rij表示面元i到面元j的距离； dS表示面积分.

(3)

将面元与来流叠加在一起等于速度势ψ， 可写为：

(4)

若满足叶剖面上法向速度为零的条件， 则ψ0=0， 可建立m个方程：

(5)

其中：ni为面元上单位法向矢量. 求解上述公式还需要一个方程， 在此处使用Kutta条件， 根据经验到尾缘距离为 2%处取 0 点[10]， 使0点处的法向速度为零， 则可写出方程为：

(6)

其中：n0为尾缘的单位法向矢量， 联立式(5)便可求解涡密度γ值. 叶剖面面上任一点速度Vi=∂ψ/∂S的计算式， 可从式(4)得出, 即有：

(7)

根据所建立的方程组， 便可求解压力系数Cp, i=1-(Vi/VA)， 压力系数和弧长乘积的升力系数为：

(8)

由此得到螺旋桨叶剖面的力学性能模型， 为之后的优化建立理论模型.

2 超体积指标与遗传算法

2.1 超体积指标及计算

超体积准则作为一种解集的综合性能评价方法用于多目标遗传算法中， 其主要作用是评价多目标优化的Pareto解集， 引导算法搜索, 实现多目标的兼顾优化， 具有较强的多目标空间搜索能力和清晰的数学定义[11]. 超体积算法可定量评价算法得到的种群质量好坏， 通过超体积指标便可得到最优结果. 超体积指标不偏向任何一个目标函数， 所以也无需设置目标函数的权重， 减少了人为干预的情况[12].

超体积准则的计算方式为： 将多个目标函数构成的目标空间分割成多个小的切块，n个目标函数构成n维空间， 通过多次迭代计算将高维空间中的面逐次映射到低维空间， 再计算其体积. 这种计算方法符合空间扫描解集的特点， 可以有效地实现空间扫描算法， 将多目标问题简化. 设任意一个由n个目标函数值构成的点p:(f1(p),f2(p), …，fn(p))， 一个参考点ref:(f1(ref),f2(ref), …,fn(ref))， 点p与参考点ref在第i维方向上的差值， 即为最大深度， 即：

max depthp=|fi(ref)-fi(p)|

(9)

以三维空间为例计算超体积指标， 如图3所示. 图3(a)在空间中有 4 个点P1～P4， 坐标表示为(0, 7, 8), (2, 3, 6), (4, 5, 2), (8, 2, 4)， 参考点ref:(10, 10, 10). 点坐标值为对应的多目标优化的各目标函数值， 组成集合S={P1,P2,P3,P4}. 图3(b)中， 平行于f2-f3面沿f1轴将参考点ref与点P1～P4所围成的多面体依次切为 4 个小切块. 对于切块1，P1沿坐标轴f1方向上的坐标值f1(P1)为0， 点P2在f1方向上的坐标值f1(P2)为2， 根据式(9), 则切块1的深度为|f1(P2)-f1(P1)|=|2-0|=2. 将切块1投影到f2-f3平面上， 得到如图3(c)相应的投影面， 切块的体积即为其对应切块的深度与投影面的面积之积， 切块1的体积为((10-7)×(10-8))×2=12 mm. 其他切块的体积亦可按照同样的方法计算， 计算出总体积12+56+192+116=376 mm， 可知参考点的设置直接决定了超体积指标的大小.

图3 三维情况下超体积指标计算示例Fig.3 Calculation example of hypervolume index under 3 dimensional cases

2.2 多目标遗传算法及螺旋桨叶片多目标优化模型

遗传算法由模拟生物界自然选择和进化过程而来， 针对单目标优化， 因其适用于求解高度复杂的非线性问题， 并且其不依赖函数原型， 具有较强的通用性， 被广泛地应用于多个学科. 但在工程领域往往存在多目标和多约束条件， 且大部分情况下这些子目标之间存在冲突， 为寻找到合适的多目标优化结果， 通常需要对各子目标进行综合考虑， 即对多个子目标进行折衷， 并且多目标优化结果为一个解集， 即Pareto最优解集. 多目标优化的数学模型一般表达式为:

min{f1(x),f2(x), …,fn(x)}, s.t.gi(x)≤0

(10)

其中:f1(x),f2(x)和fn(x)为目标函数;gi(x)为约束条件.

目标函数的建立主要基于以下3点考虑. ①保证优化后螺旋桨的水动力性能， 即提高螺旋桨叶片效率. 在相同的进给系数下， 优化桨型的升力系数不能低于原桨值. ②新型桨外形曲线需保持平滑. 原始桨和新桨外形差值较小， 保证螺旋桨的光滑性和制造的可行性， 减少桨型变化面积. ③螺旋桨叶片应避免压力曲线突变， 防止出现受力不均现象， 叶片压力方差尽可能小. 根据上述条件， 提出具体的目标函数以提高升力系数, 减小优化前后桨型面积差, 降低压力方差.

为便于统一， 将升力系数的倒数作为目标函数， 使优化过程转换为求三个目标函数的最小值. 为此， 目标函数包括：

minf1(1/CL), minf2(ΔS)， minf3(D(x))

(11)

其中：CL由式(8)得到， ΔS=S-S′为优化前后桨型面积差，D(x)=E{∑[CL(x)-E(CL(x))]2}表示桨面上的压力方差.

将超体积准则作为适应度值引导遗传算法进化， 寻求优化问题的多目标最优解. 算法的具体步骤如下：

步骤1 初始化遗传算法的参数. 设置种群规模M， 算法迭代代数T， 交叉概率Pc， 变异概率Pm. 生成初始化种群.

步骤2 计算超体积指标. 设置参考点值， 将种群中的个体按照上述方法计算超体积指标， 并以此值作为适应度值.

步骤3 选择操作. 选择种群中超体积指标值较大的个体进入下一步操作.

步骤4 交叉操作. 根据交叉概率Pc， 采用多点交叉方式产生新个体.

步骤5 变异操作. 根据变异概率Pm， 使用均匀变异方法保证个体多样性.

步骤6 维护和更新外部档案. 根据超体积指标值的大小对每一代产生的非劣解进行删除和添加， 结合精英保留策略使得种群的质量和多样性得到提升.

步骤7 终止条件. 当迭代代数达到预先设定值或者迭代结果的精度满足条件要求时， 程序运行终止， 并输出优化解； 否则转步骤2.

3 多目标叶片参数优化设计

图4 螺旋桨多目标优化设计流程图Fig.4 Flow chart of the multi-objective optimization of the propeller blade

螺旋桨多目标优化设计的思路如图4所示. 第一节建立的螺旋桨力学性能模型计算公式为连续性方程， 而螺旋桨图谱为离散的坐标点. 为进行力学性能计算同时解决螺旋桨图谱桨型数据误差较大的问题， 对螺旋桨图谱上的坐标点进行参数化拟合. 采用 B 样条曲线对螺旋桨剖面进行参数化表达， 通过图谱上的坐标点确定螺旋桨的几何形状. 根据取值范围产生新的螺旋桨控制点， 以26个控制点来确定桨型， 设置点的上下限范围， 在此范围内随机生成控制点组成一定数量的设计方案， 每组控制点坐标(x1,y1), (x2,y2), (x3,y3)， …， (x26,y26)代表一个解. 计算每个桨型的三个目标函数值， 由于三个目标函数值量级、 量纲不同， 对函数值进行归一化处理. 然后计算超体积指标， 将其作为多目标优化的适应度值实现选择操作. 在遗传算法中对解进行多点交叉及均匀变异操作， 产生新的解. 保留较优的解并剔除劣解. 重复上述过程， 直到迭代数达到预设的最大进化代数.

为验证所提方法的有效性， 在螺旋桨图谱中选取 B 型桨作为优化的原型桨， 在原始桨坐标点基础上， 建立优化桨型坐标. 螺旋桨数据： 螺旋桨半径R为 290 mm， 螺距比P/D=0.8， 盘面比AE/AO=0.35， 攻角为α=10°. 利用超体积准则与遗传算法对螺旋桨8处叶剖面进行多目标优化， 优化图1中的0.2R、 0.3R、 0.4R、 0.5R、 0.6R、 0.7R、 0.8R、 0.9R处的叶剖面. 为使螺旋桨形状分布更均匀， 同时保证螺旋桨的基本特性， 对叶剖面弧长的26个控制点坐标值进行约束， 取值范围为原厚度的[-3%， 3%]之间， 使桨型变化更加合理， 以提高搜索效率.

图5 0.3 R处螺旋桨性能指标Fig.5 Propeller performance index in 0.3 R

多目标遗传算法参数为： 群体规模M=50、 交叉概率Pc=0.3、 变异概率Pm=0.03、 最大进化代数T=50. 算法中个体编码方式使用实数编码. 算法中参考点选取时， 考虑到三个目标值越小越好， 所以选取每个目标函数值中的最大值， 由三个目标中每个目标最大值作为超体积准则的参考点.

0.3R处叶剖面50次迭代后三个目标函数值如图5所示. 四方形红点为原始螺旋桨三目标函数值， 三角形蓝点为超体积指标值最大的解对应的目标函数值. 发现优化后的解更靠近坐标原点， 三个目标函数值都优于原始桨. 表1为优化前后的螺旋桨函数值(性能参数)， 优化后， 升力系数(CL)有明显提高， 压力方差下降明显， 桨型面积变化不大， 符合目标期望.

表1 B型桨0.3R处优化性能结果

六处原始桨与优化桨的叶剖面形状对比如图6所示. 从0.2R到0.6R处优化桨型相较原始桨型最大厚度明显变大， 下表面两边向外侧偏移， 中间部分内移， 拱线曲率变大， 起到增加升力的作用. 而从0.7R到0.9R处桨型变化不大， 但却呈现出整体上移的趋势， 应是受攻角影响.

图6 原始桨与优化桨形状曲线比较Fig.6 Shape of the initial blade and the optimal blade

仿真表明， 0.2R到0.6R处桨型两侧由于形状变化的原因压力曲线也相应出现明显变化， 而0.7R到0.9R处没有明显变化. 选取的2处原始桨与优化桨的叶剖面压力曲线， 如图7所示.

图7 原始桨与优化桨压力曲线比较Fig.7 Pressure curve of the initial blade and the optimal blade

优化后桨型的压力曲线较原先有较大提高， 但压力曲线不如原先曲线平缓. 由于遗传算法的随机性特点， 桨型优化后在尾缘处下表面可能发生突变， 这将影响桨型的其它性能， 并会在尾缘处产生不规律的漩涡. 该情况在今后螺旋桨叶剖面定型中需要考虑， 并在实际的生产中加以改进， 以免漩涡产生对螺旋桨的空泡、 剥蚀等不利影响.

为验证螺旋桨优化效果， 将优化后的螺旋桨模型导入Fluent中进行仿真分析. 螺旋桨优化前后压力云图如图8所示， 压力呈现出由随边到导边逐渐变小的规律. 由于叶片厚度变小， 优化后螺旋桨叶片上受力较原先更小， 达到减小阻力的目的. 优化后螺旋桨桨型较原始桨型产生较大变化， 所以叶片上受力存在不均匀的现象， 靠近随边的位置出现受力突变. 根据得到的推力及力矩计算螺旋桨的敞水效率. 在相同的进给系数J=0.3下， 螺旋桨的敞水效率η从0.425提高到0.540.

图8 优化桨型和原始桨型压力云图Fig.8 Pressure cloudmap of the original and the optimal propeller type

4 结语

本研究在B型桨基础上， 通过面元法建立螺旋桨数学模型， 计算螺旋桨力学性能， 使用超体积准则作为遗传算法的评价指标对桨型进行多目标优化设计. 结果表明： 在面积差变化很小的情况下， 优化后的桨型实现了降低阻升比、 减少压力方差的目标. 在进给系数J=0.3情况下， 螺旋桨的敞水效率η从0.425提高到0.540， 效率提升12%， 表明基于超体积准则多目标遗传算法能够得到更好的设计方案.

参考文献：

[1] MIRJALILI S, LEWIS A, MIRJALILI S A M. Multi-objective optimisation of marine propellers[J]. Procedia Computer Science, 2015(51): 2247-2256.

[2] XIE G M. Optimal preliminary propeller design based on multi-objective optimization approach[J]. Procedia Engineering， 2011(16): 278-283.

[3] 叶礼裕， 王超， 孙文林， 等. 近似模型方法在螺旋桨优化设计中的应用[J]. 上海交通大学学报， 2016， 50(8)： 1173-1179.

[4] COELHO R F. Bi-objective hypervolume-based Pareto optimization[J]. Optimization Letters, 2016, 9(6): 1091-1103.

[5] SCHÜTZE O, HERNNDEZ V A S, TRAUTMANN H,etal. The hypervolume based directed search method for muti-objective optimization problems[J]. Journal of Heuristics, 2016， 22(3)： 1381-1231.

[6] BADER J, ZITZLER E. An algorithm for fast hypervolume-based many-objective optimization[J]. Evolutionary Computation， 2011, 19(1): 45-76.

[7] 王国强， 董世汤. 船舶螺旋桨理论与应用[M]. 2版. 哈尔滨： 哈尔滨工程大学出版社， 2005.

[8] 王献孚. 船用翼理论[M]. 北京： 国防工业出版社， 1998.

[9] 黄树权. 基于遗传算法的螺旋桨性能优化研究[D]. 大连： 大连理工大学, 2009.

[10] 陆志良. 空气动力学[M]. 北京： 北京航空航天大学出版社， 2009.

[11] 郑金华, 李珂, 李密青, 等. 一种基于Hypervolume指标及其在多目标进化算法[J]. 计算机研究与发展, 2012, 49(2): 312-326.

[12] 袁萃. 遗传算法与CFD技术在翼型优化设计中的应用研究[D]. 西安： 西北工业大学， 2007.