面向气动性能的高速列车头型多目标优化设计

初彦彬， 邱广宇，2，陈秉智

(1．大连交通大学 机车车辆工程学院，辽宁 大连 116028；2．中车齐齐哈尔车辆有限公司，黑龙江 齐齐哈尔 161002)

随着列车技术的高速发展，其运行速度不断提高，列车的气动特性对其运行的稳定性、安全性、经济性影响也日趋显著.而列车的气动特性主要评判标准为气动阻力与气动升力[1-2].对于列车气动阻力而言，影响其大小的主要因素为列车运行速度.除此之外，列车头尾外形也是影响气动阻力大小的重要因素.而对于气动升力而言，无论正升力过大，还是负升力过大，均会改变其稳定的动力学性能，使列车运行安全性大大降低[3].而影响列车气动升力大小的主要因素也是列车的流线型外形.故列车头部几何外形优化设计俨然成为改善高速列车的气动性能的主要手段之一.多学科、多目标优化设计可有效缩短新头型设计周期并提升设计效率，故其在高速列车的设计研发过程中应用较为广泛.

因此，国内外学者在高速列车CFD分析和优化设计方面进行了大量的研究.张在中[4]通过对不同列车几何模型进行风洞实验，模拟分析得出列车头部几何外形与阻力位置及大小的关系.其采用的方法实质是优选方案，由于该方法的样本有限，使其不易达到最优解.随着计算机软硬件及优化的方法不断发展，李明等人[5]将外部形状的特征曲线用多项式函数表示，找到函数更优解对应的列车头部最优外形.汪怡平等人[6]则采用网格自由变形技术，直接对车身几何曲面进行形状优化.周家林等人[7]使用响应面代理模型，对列车几何外形进行气动阻力优化.翟建平[8]则将列车的横风特性设定为优化目标对列车头型进行优化.刘加利等人[9]运用几何参数化构建技术自动改变车头外形，以降低列车气动阻力.秦淼[10]和邱利伟等人[11]将气动升力绝对值最小作为目标，对动车组列车性能进行模拟分析并优化.

本文利用建模软件实现高速列车几何外形曲线的参数化定义并采用网格变形软件Sculptor生成计算样本，用拉丁超立方抽样法取得样本中的设计变量，通过Isight软件平台集成流场数值仿真软件Fluent完成样本气动特性的计算.通过构建响应面代理模型，得到优化目标的响应函数.最后采用NSGA-Ⅱ多目标优化算法得到自变量在各自取值范围内的最优解，进而得到气动特性更优的头车外形.

1 多目标优化设计理论

由于本文有两个优化目标，故采用多目标的优化设计方法，其数学模型可以表示为[12-13]：

MinimizeY=(f1(x),f2(x),…,fm(x))

Subject togi(x)≤0,i=1,2,…,k

li≤xi≤ui,i=1,2…,n

(1)

对于多目标优化问题，通常不能找到唯一解，使得目标函数均取得极小值.如果能得到一系列解，其中任何一组解都能满足没有其他解在所有目标函数上都能得到更小值，则称这组解为Pareto最优解集.本文的目标就是要寻找气动升力趋近于0，气动阻力降到最低的Pareto前沿上的一组解.

2 高速列车气动模型

2.1 高速列车几何模型建立

由于高速列车头车表面流线型较为复杂，故采用控制型线来描述头车曲面.控制型线又分基本控制型线和辅助控制型线.基本控制型线是包括纵向的对称面控制型线、横向及水平方向的最大剖面控制型线的三个基本坐标方向上剖切模型的几何外形曲面所得的基本轮廓曲线.辅助控制型线是指用来保证模型过渡连接位置的光顺性连接的除基本控制型线之外的其他三维描述曲线.由于本文的优化目标是高速列车气动特性，通过大量试验性计算及文献的考察，发现列车在高速运行的环境下，其外形中转向架和受电弓的气动特性结果数值占比很小，故可以对列车头部外形做简化处理[7].本文计算模型采用三节编组，即头车、一节中间车厢和尾车，如图1所示.

图1 三节编组列车模型示意图

2.2 优化设计变量

2.2.1 网格变形法

本文采用网格自由变形法对列车的网格模型进行更改，以得到计算样本.该方法的变形原理参考并仿照受外力的弹性物体在发生相应变形时的物理规律.图2所示为列车表面晶格分布示意图.将高速列车头部模型外部几何轮廓由若干个晶格分割固定，在对晶格节点的位置进行移动后，其复杂外形的曲面也随之发生变化.网格自由变形技术相较于传统数模软件可保证变形后物体曲面的几何连续性，同时在对物体进行变形时具有更好的交互性和可控性.

图2 列车整体表面晶格示意图

2.2.2 设计变量的选择

在Sculptor软件中运用网格变形法，并在列车头车选择头尾垂向方向6组晶格节点、分别为dy1、dy2、dy3、dy4、dy5、dy6，纵向方向1组晶格节点dx7和横向方向1组晶格节点dz8，共8组晶格节点作为优化设计自变量，且自变量的控制区间范围依照相应几何位置而不尽相同.而且，明线上运行的列车头车尾车几何外形相同，故本文采用对称处理的方法移动和改变头车与尾车的晶格控制点.头部外形设计变量选择如图3所示.

图3 头部设计变量示意图

由于设计变量数目较多，且变量间无过多交互影响，故设计变量选点应均匀随机.所以选取拉丁超立方设计方法，其相较于正交试验、中心复合等其他试验设计方法在空间上更具均匀覆盖能力，可以充分保证实验样本的均布性.本文的样本数设置为200.

在得到的众多样本中，通过对比可得出，各个方向设计变量对列车头部外形的影响.其中，纵向与横向由于设计变量单一且无交互响应，故可较为明显观察得出：垂向设计变量dy1～dy6的调整直接影响列车的流线型曲线外形；纵向设计变量dx7的调整直接影响列车的鼻间长度；横向设计变量dz8的调整直接影响列车的排障器宽度.图4为几种设计变量组合的列车外形图，可以更加直观看出样本的同异性.

(a) 调整垂向以改变流线型

(b) 调整纵向以改变鼻间长度

(c) 调整横向以改变排障器宽度图4 设计变量组合的列车外形图

2.3 高速列车计算模型建立

2.3.1 空气动力学基本方程

本文选择的在明线上运行的高速列车计算模型其马赫数不大于0.3，认为气体不可压缩.故将其周围流场按照不可压缩粘性流动模型考虑[14].而在此模型问题中，若不计其温度变化，能量方程并非是必要基本方程.故本文进行仿真模拟的流动可用以下质量守恒方程和动量守恒方程结合描述[15-17]：

质量守恒方程：

(3)

质量守恒方程也称连续方程，表述为：单位时间内流体微元体质量的增加等于同一时间间隔内流入该微元体的质量.对于不可压缩流体而言，密度为常数，故连续方程可用以下形式表达：

(4)

各个方向一致的动量守恒方程在x方向如以下表达：

(5)

式中：p为微元体压力；τxx为粘性应力；τxy、τxz为切应力；Fx为体积力.其表述为：微元体中流体动量对时间的变化率等于外界作用在该微元体上的合力.再对其的应力采用张量的表达形式并代入粘性应力和切应力的表达式推导可得：

(6)

该表达式称为Navier-Stokes方程.式中：ui或uj为流场速度，表示三个坐标方向的速度分量；xi或xj同样表示三个坐标方向分量；μ为空气动力粘度；p为热力学压强；ρ为空气密度；δij为克罗内克符号，当i=j时，δij取1；当i≠j时，δij取0.

高速列车周围流场为湍流流动类型流场，本文选用在列车空气动力学工程湍流问题应用最广泛的雷诺时均方程模拟法中涡粘模型下的k-ε两方程湍流模型[18-19].

2.3.2 高速列车外流场网格划分

为保证模拟精度且节约计算成本，故进行仿真模拟时流域选择足够大即可，其长×宽×高分别为299.5 m×92.4 m×90.2 m.车头距流域速度入口处31.2 m,车底距地面高度预留0.525 m.由于气流运动相对规律的列车上方流域较大，而其对列车气动分析结果的影响较小，故对其采用有对边界良好的适应性优点的非结构性网格从列车表面至流域边界大小过渡填充.以确保在精度足够的情况下网格数量较少，计算结果易收敛，提高计算效率.流域网格划分效果如图5所示.其中流场区域共有3 628 270个单元.列车模型共418 748个壳单元.

(a) 整体纵向对称面剖视图

(b) 局部纵向对称面剖视图图5 流场的纵向对称面剖视图

2.3.3 流场边界条件设置

结合列车明线运行实际情况，对边界条件进行以下设置：

(1)入口：速度入口(velocity-inlet)，且大小为320 km/h.

(2)出口：压力出口(pressure-outlet)，压力与大气压相同，即压力出口的表压强为0.

(3)地面：滑移壁面(moving-wall)，按光滑壁面处理，为模拟地面效应，选择滑移边界条件.且速度值取列车运行速度320 km/h.

(4)流域上表面及两侧边界：对称面(symmetry).以模拟流域上方无边界的区域，保证边界不会产生摩擦等干扰条件，以影响计算结果的精度.

(5)列车表面：固定壁面(wall).

图6 车头表面压强分布

3 多目标气动优化

3.1 响应面代理模型

在Isight平台上，集成网格前处理Sculptor软件、流场计算Fluent软件以及自编写批处理文件完成自动计算.运用多项式响应面法(Response Surface Method)来构造气动特性的代理模型，从众多气动阻力响应关系中，选取较易观察出变量间影响关系的部分响应面，如图7所示.并可得出：dy1和dy3对气动阻力的影响大于dy2和dy4这一结论.选取其他10组样本对响应面代理模型进行拟合精度分析，经计算代理模型的R2=0.988.R2接近1，说明所选取的响应面模型与样本符合，构造的响应面模型准确有效.

(a) F对dy1、dy2的响应

(b)F对dy3、dy4的响应图7 优化设计变量与目标气动阻力F的部分响应面

构造的气动升力响应面代理模型如图8所示，拟合精度R2=0.995.也接近1，说明构造的气动升力的响应面模型也准确有效.

图8 优化设计变量与气动升力的部分响应面

3.2 优化目标与设计变量的相关性分析

在得到了优化目标与各设计变量的响应关系后，可以逐步分析得到设计变量对优化目标的影响大小及规律.但分析响应面不如通过Pareto图更加直接准确.图9是各设计变量相对优化目标即气动阻力的贡献率图，即Pareto图.dx7,dy6表示负向效应，其他表示正向效应，且柱形越长表示影响越大，柱形越短则影响越小.由图可知，在影响率的大小方面，对气动阻力影响最大为横向设计变量dz8，其贡献率高达21%；其次是dx7、dy6；而垂向方向的几个设计变量的影响率由大到小依次为dy1、dy3、dy2、dy4、dy5；而在影响率的正负性方面，只有设计变量dx7、dy6对气动阻力的影响为负向效应，即优化目标随着设计变量在可变化范围内的变大而减小.而造成这样影响大小及正负效应的原因是由于在流线型的原模型中，当列车头尾垂向设计变量在各自设计范围内取值变化时，并不能改变主要影响列车气动特性的头部尖锐程度.因而反应横向宽窄的设计变量dz8对于气动阻力影响最大.同时得到dy1、dy3、dy2、dy4、dy5对气动阻力的影响依次降低，这是因为dy2、dy3和dy4处于列车流线型表面几何突变位置，因而此处外形对气流流动影响较大.而dy5晶格控制点位置距离车体顶部较近，对车头尖锐程度只有较小影响.

图9 各设计变量贡献率图

3.3 多目标优化计算方法

在分别得到了优化目标与各设计变量间的响应关系与相关性后，采用NSGA-II多目标遗传算法在其基础上得到近似全局最优解.采用带有精英策略非劣排序的该方法具有种群多样、探索性良好、计算效率高等优点.

3.4 多目标优化结果

本次多目标优化过程的优化目标为列车的气动阻力和气动升力.即让高速列车的气动阻力尽可能地小，同时使气动升力尽可能趋近于零.由于该设计为多目标且各目标限制条件不同，而使优化周期加长且寻得较多对试验造成影响Pareto前沿解.在优化过程中，亦可以依据目标的优化优先级来配比权重以得到实际问题更以应用的优化结果.寻优过程如图10所示.

(a)气动阻力的优化历程

(b)气动升力的优化历程图10 优化历程图

从图10可以看出，星点为求得的Pareto解集中最理想的解，其气动升力趋近于0，且阻力降到了最低.得到的优化头车外形中，横向与纵向设计变量dz8和dx7取值分别为-0.247 4、0.427 0 m，并各自接近可取值区间的最小值和最大值；垂向设计变量dy1值为-0.492 0 m、dy2值为-0.189 2 m、dy3值为-0.234 4 m、dy4值为-0.148 6 m、dy5值为-0.084 6 m、dy6值为0.391 7 m，以上三个方向的优化结果使列车流线型头尾纵向更长、横向更窄、垂向更扁，这与空气动力学的基本理论及工程经验相吻合.优化前，列车所受总气动阻力大小为21 531.29 N，气动升力为23.86 N.优化后列车所受总气动阻力大小为18 589.63 N，降阻比率约为13.66%.气动升力为1.46 N.

4 结论

本文针对某高速列车头型进行多目标优化设计，通过拉丁超立方法结合网格变形法采集200组设计变量组合样本，并通过构造设计变量与优化目标响应面进行分析，最终采用NSGA-II多目标遗传算法设置气动阻力与升力目标，寻得气动特性更优的头车外形全局最优解.在优化过程中构造列车气动阻力、升力与各设计变量间的响应面关系.且各控制变量对不同优化目标的贡献度也不尽相同.对于气动阻力而言，在影响度大小方面：dz8、dx7、dy6、dy1、dy3、dy2、 dy4、dy5对改变气动阻力的贡献依次降低；而在影响的正负性方面，只有设计变量dx7、dy6对气动阻力的影响为负向效应，其余设计变量的影响均为正向效应.同时在优化历程图中也可明显看出，在最优解附近有众多前沿解，为实际应用问题提供大量素材，优化结果得到的列车流线型头尾纵向更长、横向更窄、垂向更扁.优化最优解和优化前相比，列车阻力降低13.66%，且气动升力有效减小.