初彥彬， 邱廣宇，2，陳秉智

(1．大連交通大學 機車車輛工程學院，遼寧 大連 116028；2．中車齊齊哈爾車輛有限公司，黑龍江 齊齊哈爾 161002)

隨著列車技術(shù)的高速發(fā)展，其運行速度不斷提高，列車的氣動特性對其運行的穩(wěn)定性、安全性、經(jīng)濟性影響也日趨顯著.而列車的氣動特性主要評判標準為氣動阻力與氣動升力[1-2].對于列車氣動阻力而言，影響其大小的主要因素為列車運行速度.除此之外，列車頭尾外形也是影響氣動阻力大小的重要因素.而對于氣動升力而言，無論正升力過大，還是負升力過大，均會改變其穩(wěn)定的動力學性能，使列車運行安全性大大降低[3].而影響列車氣動升力大小的主要因素也是列車的流線型外形.故列車頭部幾何外形優(yōu)化設(shè)計儼然成為改善高速列車的氣動性能的主要手段之一.多學科、多目標優(yōu)化設(shè)計可有效縮短新頭型設(shè)計周期并提升設(shè)計效率，故其在高速列車的設(shè)計研發(fā)過程中應(yīng)用較為廣泛.

因此，國內(nèi)外學者在高速列車CFD分析和優(yōu)化設(shè)計方面進行了大量的研究.張在中[4]通過對不同列車幾何模型進行風洞實驗，模擬分析得出列車頭部幾何外形與阻力位置及大小的關(guān)系.其采用的方法實質(zhì)是優(yōu)選方案，由于該方法的樣本有限，使其不易達到最優(yōu)解.隨著計算機軟硬件及優(yōu)化的方法不斷發(fā)展，李明等人[5]將外部形狀的特征曲線用多項式函數(shù)表示，找到函數(shù)更優(yōu)解對應(yīng)的列車頭部最優(yōu)外形.汪怡平等人[6]則采用網(wǎng)格自由變形技術(shù)，直接對車身幾何曲面進行形狀優(yōu)化.周家林等人[7]使用響應(yīng)面代理模型，對列車幾何外形進行氣動阻力優(yōu)化.翟建平[8]則將列車的橫風特性設(shè)定為優(yōu)化目標對列車頭型進行優(yōu)化.劉加利等人[9]運用幾何參數(shù)化構(gòu)建技術(shù)自動改變車頭外形，以降低列車氣動阻力.秦淼[10]和邱利偉等人[11]將氣動升力絕對值最小作為目標，對動車組列車性能進行模擬分析并優(yōu)化.

本文利用建模軟件實現(xiàn)高速列車幾何外形曲線的參數(shù)化定義并采用網(wǎng)格變形軟件Sculptor生成計算樣本，用拉丁超立方抽樣法取得樣本中的設(shè)計變量，通過Isight軟件平臺集成流場數(shù)值仿真軟件Fluent完成樣本氣動特性的計算.通過構(gòu)建響應(yīng)面代理模型，得到優(yōu)化目標的響應(yīng)函數(shù).最后采用NSGA-Ⅱ多目標優(yōu)化算法得到自變量在各自取值范圍內(nèi)的最優(yōu)解，進而得到氣動特性更優(yōu)的頭車外形.

1 多目標優(yōu)化設(shè)計理論

由于本文有兩個優(yōu)化目標，故采用多目標的優(yōu)化設(shè)計方法，其數(shù)學模型可以表示為[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)

對于多目標優(yōu)化問題，通常不能找到唯一解，使得目標函數(shù)均取得極小值.如果能得到一系列解，其中任何一組解都能滿足沒有其他解在所有目標函數(shù)上都能得到更小值，則稱這組解為Pareto最優(yōu)解集.本文的目標就是要尋找氣動升力趨近于0，氣動阻力降到最低的Pareto前沿上的一組解.

2 高速列車氣動模型

2.1 高速列車幾何模型建立

由于高速列車頭車表面流線型較為復雜，故采用控制型線來描述頭車曲面.控制型線又分基本控制型線和輔助控制型線.基本控制型線是包括縱向的對稱面控制型線、橫向及水平方向的最大剖面控制型線的三個基本坐標方向上剖切模型的幾何外形曲面所得的基本輪廓曲線.輔助控制型線是指用來保證模型過渡連接位置的光順性連接的除基本控制型線之外的其他三維描述曲線.由于本文的優(yōu)化目標是高速列車氣動特性，通過大量試驗性計算及文獻的考察，發(fā)現(xiàn)列車在高速運行的環(huán)境下，其外形中轉(zhuǎn)向架和受電弓的氣動特性結(jié)果數(shù)值占比很小，故可以對列車頭部外形做簡化處理[7].本文計算模型采用三節(jié)編組，即頭車、一節(jié)中間車廂和尾車，如圖1所示.

圖1 三節(jié)編組列車模型示意圖

2.2 優(yōu)化設(shè)計變量

2.2.1 網(wǎng)格變形法

本文采用網(wǎng)格自由變形法對列車的網(wǎng)格模型進行更改，以得到計算樣本.該方法的變形原理參考并仿照受外力的彈性物體在發(fā)生相應(yīng)變形時的物理規(guī)律.圖2所示為列車表面晶格分布示意圖.將高速列車頭部模型外部幾何輪廓由若干個晶格分割固定，在對晶格節(jié)點的位置進行移動后，其復雜外形的曲面也隨之發(fā)生變化.網(wǎng)格自由變形技術(shù)相較于傳統(tǒng)數(shù)模軟件可保證變形后物體曲面的幾何連續(xù)性，同時在對物體進行變形時具有更好的交互性和可控性.

圖2 列車整體表面晶格示意圖

2.2.2 設(shè)計變量的選擇

在Sculptor軟件中運用網(wǎng)格變形法，并在列車頭車選擇頭尾垂向方向6組晶格節(jié)點、分別為dy1、dy2、dy3、dy4、dy5、dy6，縱向方向1組晶格節(jié)點dx7和橫向方向1組晶格節(jié)點dz8，共8組晶格節(jié)點作為優(yōu)化設(shè)計自變量，且自變量的控制區(qū)間范圍依照相應(yīng)幾何位置而不盡相同.而且，明線上運行的列車頭車尾車幾何外形相同，故本文采用對稱處理的方法移動和改變頭車與尾車的晶格控制點.頭部外形設(shè)計變量選擇如圖3所示.

圖3 頭部設(shè)計變量示意圖

由于設(shè)計變量數(shù)目較多，且變量間無過多交互影響，故設(shè)計變量選點應(yīng)均勻隨機.所以選取拉丁超立方設(shè)計方法，其相較于正交試驗、中心復合等其他試驗設(shè)計方法在空間上更具均勻覆蓋能力，可以充分保證實驗樣本的均布性.本文的樣本數(shù)設(shè)置為200.

在得到的眾多樣本中，通過對比可得出，各個方向設(shè)計變量對列車頭部外形的影響.其中，縱向與橫向由于設(shè)計變量單一且無交互響應(yīng)，故可較為明顯觀察得出：垂向設(shè)計變量dy1～dy6的調(diào)整直接影響列車的流線型曲線外形；縱向設(shè)計變量dx7的調(diào)整直接影響列車的鼻間長度；橫向設(shè)計變量dz8的調(diào)整直接影響列車的排障器寬度.圖4為幾種設(shè)計變量組合的列車外形圖，可以更加直觀看出樣本的同異性.

(a) 調(diào)整垂向以改變流線型

(b) 調(diào)整縱向以改變鼻間長度

(c) 調(diào)整橫向以改變排障器寬度圖4 設(shè)計變量組合的列車外形圖

2.3 高速列車計算模型建立

2.3.1 空氣動力學基本方程

本文選擇的在明線上運行的高速列車計算模型其馬赫數(shù)不大于0.3，認為氣體不可壓縮.故將其周圍流場按照不可壓縮粘性流動模型考慮[14].而在此模型問題中，若不計其溫度變化，能量方程并非是必要基本方程.故本文進行仿真模擬的流動可用以下質(zhì)量守恒方程和動量守恒方程結(jié)合描述[15-17]：

質(zhì)量守恒方程：

(3)

質(zhì)量守恒方程也稱連續(xù)方程，表述為：單位時間內(nèi)流體微元體質(zhì)量的增加等于同一時間間隔內(nèi)流入該微元體的質(zhì)量.對于不可壓縮流體而言，密度為常數(shù)，故連續(xù)方程可用以下形式表達：

(4)

各個方向一致的動量守恒方程在x方向如以下表達：

(5)

式中：p為微元體壓力；τxx為粘性應(yīng)力；τxy、τxz為切應(yīng)力；Fx為體積力.其表述為：微元體中流體動量對時間的變化率等于外界作用在該微元體上的合力.再對其的應(yīng)力采用張量的表達形式并代入粘性應(yīng)力和切應(yīng)力的表達式推導可得：

(6)

該表達式稱為Navier-Stokes方程.式中：ui或uj為流場速度，表示三個坐標方向的速度分量；xi或xj同樣表示三個坐標方向分量；μ為空氣動力粘度；p為熱力學壓強；ρ為空氣密度；δij為克羅內(nèi)克符號，當i=j時，δij取1；當i≠j時，δij取0.

高速列車周圍流場為湍流流動類型流場，本文選用在列車空氣動力學工程湍流問題應(yīng)用最廣泛的雷諾時均方程模擬法中渦粘模型下的k-ε兩方程湍流模型[18-19].

2.3.2 高速列車外流場網(wǎng)格劃分

為保證模擬精度且節(jié)約計算成本，故進行仿真模擬時流域選擇足夠大即可，其長×寬×高分別為299.5 m×92.4 m×90.2 m.車頭距流域速度入口處31.2 m,車底距地面高度預留0.525 m.由于氣流運動相對規(guī)律的列車上方流域較大，而其對列車氣動分析結(jié)果的影響較小，故對其采用有對邊界良好的適應(yīng)性優(yōu)點的非結(jié)構(gòu)性網(wǎng)格從列車表面至流域邊界大小過渡填充.以確保在精度足夠的情況下網(wǎng)格數(shù)量較少，計算結(jié)果易收斂，提高計算效率.流域網(wǎng)格劃分效果如圖5所示.其中流場區(qū)域共有3 628 270個單元.列車模型共418 748個殼單元.

(a) 整體縱向?qū)ΨQ面剖視圖

(b) 局部縱向?qū)ΨQ面剖視圖圖5 流場的縱向?qū)ΨQ面剖視圖

2.3.3 流場邊界條件設(shè)置

結(jié)合列車明線運行實際情況，對邊界條件進行以下設(shè)置：

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

(2)出口：壓力出口(pressure-outlet)，壓力與大氣壓相同，即壓力出口的表壓強為0.

(3)地面：滑移壁面(moving-wall)，按光滑壁面處理，為模擬地面效應(yīng)，選擇滑移邊界條件.且速度值取列車運行速度320 km/h.

(4)流域上表面及兩側(cè)邊界：對稱面(symmetry).以模擬流域上方無邊界的區(qū)域，保證邊界不會產(chǎn)生摩擦等干擾條件，以影響計算結(jié)果的精度.

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

圖6 車頭表面壓強分布

3 多目標氣動優(yōu)化

3.1 響應(yīng)面代理模型

在Isight平臺上，集成網(wǎng)格前處理Sculptor軟件、流場計算Fluent軟件以及自編寫批處理文件完成自動計算.運用多項式響應(yīng)面法(Response Surface Method)來構(gòu)造氣動特性的代理模型，從眾多氣動阻力響應(yīng)關(guān)系中，選取較易觀察出變量間影響關(guān)系的部分響應(yīng)面，如圖7所示.并可得出：dy1和dy3對氣動阻力的影響大于dy2和dy4這一結(jié)論.選取其他10組樣本對響應(yīng)面代理模型進行擬合精度分析，經(jīng)計算代理模型的R2=0.988.R2接近1，說明所選取的響應(yīng)面模型與樣本符合，構(gòu)造的響應(yīng)面模型準確有效.

(a) F對dy1、dy2的響應(yīng)

(b)F對dy3、dy4的響應(yīng)圖7 優(yōu)化設(shè)計變量與目標氣動阻力F的部分響應(yīng)面

構(gòu)造的氣動升力響應(yīng)面代理模型如圖8所示，擬合精度R2=0.995.也接近1，說明構(gòu)造的氣動升力的響應(yīng)面模型也準確有效.

圖8 優(yōu)化設(shè)計變量與氣動升力的部分響應(yīng)面

3.2 優(yōu)化目標與設(shè)計變量的相關(guān)性分析

在得到了優(yōu)化目標與各設(shè)計變量的響應(yīng)關(guān)系后，可以逐步分析得到設(shè)計變量對優(yōu)化目標的影響大小及規(guī)律.但分析響應(yīng)面不如通過Pareto圖更加直接準確.圖9是各設(shè)計變量相對優(yōu)化目標即氣動阻力的貢獻率圖，即Pareto圖.dx7,dy6表示負向效應(yīng)，其他表示正向效應(yīng)，且柱形越長表示影響越大，柱形越短則影響越小.由圖可知，在影響率的大小方面，對氣動阻力影響最大為橫向設(shè)計變量dz8，其貢獻率高達21%；其次是dx7、dy6；而垂向方向的幾個設(shè)計變量的影響率由大到小依次為dy1、dy3、dy2、dy4、dy5；而在影響率的正負性方面，只有設(shè)計變量dx7、dy6對氣動阻力的影響為負向效應(yīng)，即優(yōu)化目標隨著設(shè)計變量在可變化范圍內(nèi)的變大而減小.而造成這樣影響大小及正負效應(yīng)的原因是由于在流線型的原模型中，當列車頭尾垂向設(shè)計變量在各自設(shè)計范圍內(nèi)取值變化時，并不能改變主要影響列車氣動特性的頭部尖銳程度.因而反應(yīng)橫向?qū)捳脑O(shè)計變量dz8對于氣動阻力影響最大.同時得到dy1、dy3、dy2、dy4、dy5對氣動阻力的影響依次降低，這是因為dy2、dy3和dy4處于列車流線型表面幾何突變位置，因而此處外形對氣流流動影響較大.而dy5晶格控制點位置距離車體頂部較近，對車頭尖銳程度只有較小影響.

圖9 各設(shè)計變量貢獻率圖

3.3 多目標優(yōu)化計算方法

在分別得到了優(yōu)化目標與各設(shè)計變量間的響應(yīng)關(guān)系與相關(guān)性后，采用NSGA-II多目標遺傳算法在其基礎(chǔ)上得到近似全局最優(yōu)解.采用帶有精英策略非劣排序的該方法具有種群多樣、探索性良好、計算效率高等優(yōu)點.

3.4 多目標優(yōu)化結(jié)果

本次多目標優(yōu)化過程的優(yōu)化目標為列車的氣動阻力和氣動升力.即讓高速列車的氣動阻力盡可能地小，同時使氣動升力盡可能趨近于零.由于該設(shè)計為多目標且各目標限制條件不同，而使優(yōu)化周期加長且尋得較多對試驗造成影響Pareto前沿解.在優(yōu)化過程中，亦可以依據(jù)目標的優(yōu)化優(yōu)先級來配比權(quán)重以得到實際問題更以應(yīng)用的優(yōu)化結(jié)果.尋優(yōu)過程如圖10所示.

(a)氣動阻力的優(yōu)化歷程

(b)氣動升力的優(yōu)化歷程圖10 優(yōu)化歷程圖

從圖10可以看出，星點為求得的Pareto解集中最理想的解，其氣動升力趨近于0，且阻力降到了最低.得到的優(yōu)化頭車外形中，橫向與縱向設(shè)計變量dz8和dx7取值分別為-0.247 4、0.427 0 m，并各自接近可取值區(qū)間的最小值和最大值；垂向設(shè)計變量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，以上三個方向的優(yōu)化結(jié)果使列車流線型頭尾縱向更長、橫向更窄、垂向更扁，這與空氣動力學的基本理論及工程經(jīng)驗相吻合.優(yōu)化前，列車所受總氣動阻力大小為21 531.29 N，氣動升力為23.86 N.優(yōu)化后列車所受總氣動阻力大小為18 589.63 N，降阻比率約為13.66%.氣動升力為1.46 N.

4 結(jié)論

本文針對某高速列車頭型進行多目標優(yōu)化設(shè)計，通過拉丁超立方法結(jié)合網(wǎng)格變形法采集200組設(shè)計變量組合樣本，并通過構(gòu)造設(shè)計變量與優(yōu)化目標響應(yīng)面進行分析，最終采用NSGA-II多目標遺傳算法設(shè)置氣動阻力與升力目標，尋得氣動特性更優(yōu)的頭車外形全局最優(yōu)解.在優(yōu)化過程中構(gòu)造列車氣動阻力、升力與各設(shè)計變量間的響應(yīng)面關(guān)系.且各控制變量對不同優(yōu)化目標的貢獻度也不盡相同.對于氣動阻力而言，在影響度大小方面：dz8、dx7、dy6、dy1、dy3、dy2、 dy4、dy5對改變氣動阻力的貢獻依次降低；而在影響的正負性方面，只有設(shè)計變量dx7、dy6對氣動阻力的影響為負向效應(yīng)，其余設(shè)計變量的影響均為正向效應(yīng).同時在優(yōu)化歷程圖中也可明顯看出，在最優(yōu)解附近有眾多前沿解，為實際應(yīng)用問題提供大量素材，優(yōu)化結(jié)果得到的列車流線型頭尾縱向更長、橫向更窄、垂向更扁.優(yōu)化最優(yōu)解和優(yōu)化前相比，列車阻力降低13.66%，且氣動升力有效減小.