黄锡龙,杨宝锋,严俊峰,李春乐,许开富
(西安航天动力研究所,陕西 西安 710100)
预压泵能够显著提高主泵的空化性能,减轻推进剂贮箱质量,提高火箭的有效载荷,在新一代大型液体运载火箭中得到了普遍应用[1]。但国内针对预压泵性能优化的研究还较少[2-3]。国内学者分别基于正交优化设计、权矩阵分析、NSGA-III算法和二次回归正交组合设计等对离心轮和诱导轮性能开展了大量优化工作[4-11]。Derakhshan等基于神经网络与人工蜂群算法对离心泵叶轮进行优化设计[12]。
随着对泵精细化设计的要求越来越高,多参数、多目标的优化逐渐引起学者重视。Bellary等基于参数化建模并采用多目标优化算法对离心泵叶轮形状进行了优化[13]。Park等基于人工神经网络建模和多目标粒子群优化算法对轴流泵进行了优化[14]。袁寿其等基于最优拉丁超立方法和多岛遗传算法求得最优的叶轮参数组合[15]。薛城等利用正交试验方法,采用五因素四水平设计了16组试验方法,对离心泵的性能进行了优化[16]。
在多参数优化领域中,响应面法因其独特的优势也得到了广泛的应用。张人会等采用响应面法对叶片载荷分布规律进行了优化研究[17]。张德胜等利用响应面法构建了6个叶轮外形参数与优化目标的数学模型,使泵性能得到明显提升[18]。
本文基于响应面法针对高速预压泵中的诱导轮开展多参数水力优化研究,利用罚函数法对响应面模型进行寻优,获得最优组合,提高诱导轮的水力性能。
响应面法是数据回归方法中的一种,它综合了试验设计和数学建模[19]。该方法一般包括试验、建模、模型检验、最优化和结果验证等几个步骤。通过对设计域内部分具有代表性的点进行试验,回归拟合影响因素与优化目标之间的近似函数关系,基于一定的算法求得各因素最优水平[20]。该方法在机械优化领域内得到了广泛的应用,是近年来高速发展的优化理论方法。
在响应面法中,影响因素的选择直接决定了优化效果,选择合适的影响因素不仅是优化成功的关键,也是减少资源消耗、加快优化进程的重要手段。目前在叶片泵的水力优化中,一般选取1~3个参数进行分析,但由于影响叶片泵性能参数众多,得到的优化模型往往还有进一步优化的空间。但是从优化效率和速度考虑,影响因素又不宜太多。
综上,选取6个结构参数对诱导轮进行优化,包括叶片倾角γ、入口安放角β1、出口安放角β2、分流叶片轴向位置l、入口轮毂直径D1、出口轮毂直径D2等,如图1所示。
图1 诱导轮参数示意图
为解决影响因素不同量纲所带来的不便,方便后续统一处理,将所有参数作如下线性变换:
(1)
式中影响因素Xi的设计范围是[X1i,X2i],i=1,2,3,…,6。
采用旋转性及外推稳健性均较好的中心复合有界设计开展试验方案设计,该试验设计的试验点分布情况见图2。
图2 中心复合有界设计试验点分布
由于共有6个影响因素,故试验包含12个轴向点、64个立方点,另取2个中心点,共计78个试验点,即78个诱导轮方案。
对78个诱导轮方案进行数值仿真,获得诱导轮扬程系数和效率。同时,对计算结果进行极差分析,获得各影响因素对诱导轮扬程系数和效率的影响。
利用ANSYS ICEM CFD软件对流域划分网格。采用非结构网格,并对壁面进行加密。计算域如图3所示。
图3 计算域
为提高计算方法的准确性并降低对计算资源的需求,对网格进行了无关性验证,结果如图4所示。为保证计算结果的准确性,选定总网格单元数为574万。
图4 网格无关性验证
采用ANSYS CFX15.0进行稳态数值计算。将常温水作为工作介质,以压力入口、质量流量出口作为边界条件,并设入口总压为0.073 MPa、质量流量为6.2 kg/s。壁面均为无滑移壁面,采用scalable壁面函数。诱导轮转速设为12 500 r/min。采用k-ε湍流模型,ZGB空化模型,饱和蒸气压为3.574 kPa,残差收敛精度设为10-6。利用ANSYS CFX软件的session功能简化设置过程。
为使计算结果更准确,先计算无空化条件下的流场,并将其作为初场再计算空化条件下的流场。部分计算结果如表1所示(其中000-1模型和000-2模型是原始模型的两次独立的仿真计算,即两个中心点)。扬程系数定义为
表1 部分计算结果
(2)
式中:pout为出口压力;pin为入口压力;ρ为流体密度;vt为泵叶尖速度。
考察各影响因素对诱导轮扬程和效率的影响大小,对64个立方点的计算结果进行极差分析,结果如表2和表3所示。其中Kij是第i个因素、第j个水平所对应的计算结果之和,kij是Kij的平均值,Ri是第i个因素的极差,极差越大反映影响因素对诱导轮扬程或效率的影响越大。
表2 扬程系数极差分析
表3 效率极差分析
从表2可得,影响因素对扬程的影响从大到小依次是:出口安放角、出口轮毂直径、入口安放角、叶片倾角、入口轮毂直径、分流叶片轴向位置。从表3可得,6个影响因素对效率的影响从大到小依次是入口安放角、出口安放角、分流叶片轴向位置、入口轮毂直径、叶片倾角、出口轮毂直径。
根据试验样本建立扬程系数和效率的响应面模型。采用常用的二阶响应面模型。该模型具有较高的拟合精度。
由于包含6个影响因素,二阶响应面模型公式为
(3)
式(3)共包含28个未知参数。根据得到的样本数据,通过线性回归可求解该模型中的所有未知参数。
为方便计算,将式(3)进行如下换算。
将式(3)改写成多元线性模型,由此得到响应面模型的线性表达形式,即
(4)
(5)
式中:Z为78组方案的影响因素水平矩阵;Y为78组方案的扬程系数矩阵或效率矩阵。
经过矩阵计算,分别得到扬程系数响应面模型和效率响应面模型,即
(6)
对响应面模型进行误差分析、方差分析及残差分布分析,验证模型的准确性、可靠性和适应性。
表4给出了响应面模型误差分析和方差分析的结果。结果表明,两个响应模型均具有较高的准确性和可靠性。
表4 评估结果
下面对响应面模型的残差分布进行分析,验证响应面模型的适应性。图5和图6分别给出了残差正态概率分布图和残差与预测值分布图。图中每个点代表一个试验方案,不同颜色代表不同的值,其值从蓝到红逐渐增大。
图5 残差正态概率分布图
图6 残差与预测值分布图
从图5可以看到,残差正态概率总体呈线性分布,主要位于一根直线上,表明两种响应面模型具有良好的适应性。而在图6中,所有的残差均分布在两根红线内,且无明显的分布规律,亦能表明响应面模型具有良好的适应性。
综上分析,扬程系数响应面模型和效率响应面模型具有良好的准确性、可靠性和适应性,能够对优化目标进行有效、准确的预测。
将扬程系数作为优化目标,效率大于60%作为约束条件,采用内点罚函数法在不同范围内对扬程系数进行寻优。优化问题的数学形式为
(7)
式中:F(X)为目标函数;M(X)为约束函数;a、b为寻优范围的上下界。
构成内点罚函数的寻优数学模型为
(8)
式中rk为惩罚因子,rk+1=βrk,β=0.1。
持续迭代,惩罚因子将趋于0,惩罚项也趋于0,从而使罚函数逐渐逼近目标函数。
利用Design-Expert软件分别求解不同寻优范围下扬程系数的最优解。下面给出了4个不同范围下的寻优结果。
1)-1≤x1,x2,x3,x4,x5,x6≤1
在[-1,1]寻优范围内求得最优解为x=[-1,-1,1,-0.271,-1,-1],扬程系数的预测值为0.326 6,效率的预测值为61.53%。
根据最优解给出的诱导轮参数,建立诱导轮三维模型,并进行仿真计算,得到扬程系数的仿真值为0.316 3,效率为61.21%,优化后的诱导轮扬程系数相比原始模型提升了11.06%,与预测值的相对误差为3.26%。将该优化模型命名为优化模型A。
2)-1.2≤x1,x2,x3,x4,x5,x6≤1.2
在[-1.2,1.2]范围内求得最优解为x=[-1.2,-1.2,1.2,-0.527,-1.2,-1.075],扬程系数的预测值为0.334 2,效率的预测值为61.06%。
经仿真计算,得到扬程系数的仿真值为0.323 2,效率为60.67%,优化后的诱导轮扬程系数较原始模型提升了13.48%,与预测值的相对误差为3.40%。将该优化模型命名为优化模型B。
3)-1.5≤x1,x2,x3,x4,x5,x6≤1.5
在[-1.5,1.5]范围内,求得最优解为x=[-1.5,-1.5,1.5,-0.578,-1.5,-0.979],扬程系数的预测值为0.346 0,效率的预测值为60.64%。
经仿真计算,得到扬程系数的仿真值为0.328 7,效率为61.37%,优化后的诱导轮扬程系数较原始模型提升了15.41%,与预测值的相对误差为5.26%。将该优化模型命名为优化模型C。
4)-2≤x1,x2,x3,x4,x5,x6≤2
在[-2,2]范围内,求得最优解为x=[-2,-2,2,-0.844,-2,-0.541],扬程系数的预测值为0.366 0,效率的预测值为60.00%。
经仿真计算,得到扬程系数的仿真值为0.319 5,效率为61.04%,优化后的诱导轮扬程系数较原始模型提升了12.18%,与预测值的相对误差为14.55%。将该优化模型命名为优化模型D。
优化结果汇总如表5所示。
表5 优化结果
从表5可以看到,在一定范围内,随着寻优范围的扩大,最优解的诱导轮扬程系数越大,预测值与仿真值之间的相对误差越大,且相对误差处于较低水平,体现了响应面模型对扬程系数预测的准确性。但是当寻优范围扩至[-2,2]时,两者的相对误差急剧增加到14.55%,响应面模型预测得到的结果明显失真,且扬程系数未继续增大,表明在该寻优范围下,扬程系数响应面模型已不再适用,无法对诱导轮扬程系数作出准确的预测。在4个优化模型中,优化模型C增压能力最高,其扬程系数较原始模型提升了15.41%。
将扬程系数最高的优化模型C与原始模型做进一步对比分析。
图7给出了优化前后的诱导轮流道轴截面压力分布对比。从图7中可以看到,优化模型的流道入口低压区明显减少,叶轮入口条件明显改善,具有更好的进口性能,且叶轮出口处高压区也明显增多,表明优化模型增压能力得到了提升。
图7 压力分布对比
图8给出了90%叶高处湍动能分布对比,湍动能可以反映流道内能量的耗散程度。从图8中可以看到,在分流叶片与长叶片之间的流道的湍动能明显高于其他区域,且优化后湍动能小幅下降,出口处的能量损失减少。此外,优化模型与原始模型相比,入口处的湍动能均处于较低的水平,这可能是因为优化后减弱了进口回流,改善了诱导轮入口条件,使得入口能量损失减少。
图8 90叶高处湍动能分布对比
综上对比分析,优化后的诱导轮做功能力更强,入口条件更好,流场改善明显。
利用基于中心复合有界设计的响应面法,选取叶片倾角γ、入口安放角β1、出口安放角β2、分流叶片轴向位置l、入口轮毂直径D1和出口轮毂直径D2这6个诱导轮结构参数作为影响因素,对诱导轮进行了多参数水力优化,使诱导轮扬程提升15.4%。
通过对64个立方点的计算结果进行极差分析,得到各影响因素对诱导轮扬程和效率的影响大小。其中,出口安放角和出口轮毂直径对扬程的影响较大,分流叶片轴向位置和入口轮毂直径对扬程的影响较小;入口安放角对效率的影响最大,出口轮毂直径对效率的影响最小。