基于代理模型與參數(shù)敏感性分析的扇形氣膜孔優(yōu)化

管俊俊， 張 祎， 張也平， 鮑阿美， 戴 韌

(1.上海理工大學(xué) 能源與動(dòng)力工程學(xué)院，上海 200093；2.上海電氣電站集團(tuán)本部，上海 201199)

氣膜孔的幾何形狀是影響氣膜冷卻效果的重要因素，為了改善冷卻效果，研究人員提出了各種孔型。雖然冷卻效率很高的孔型很多，但大部分孔型結(jié)構(gòu)復(fù)雜，無法進(jìn)行大規(guī)模工業(yè)加工，而扇形孔不僅冷卻效率高，還具有結(jié)構(gòu)簡(jiǎn)單、加工成本低等優(yōu)點(diǎn)，是氣膜孔優(yōu)化研究中的重點(diǎn)[1-2]。

Gritsch等[3]研究發(fā)現(xiàn)，與圓柱孔相比，扇形孔能夠增大氣膜流向和側(cè)向的覆蓋面積，并在高吹風(fēng)比下仍保持較高的冷卻效率。扇形孔參數(shù)會(huì)影響其自身的冷卻效果，包括傾斜角[4]、前傾角[5]、側(cè)擴(kuò)角[6]、長(zhǎng)徑比[7]和孔間距[8]等。其中，側(cè)擴(kuò)角和前傾角對(duì)扇形孔冷卻效率的影響較大。

基于CFD數(shù)值方法對(duì)扇形孔幾何形狀進(jìn)行優(yōu)化的研究很多，包括多保真度模型[9]和深度學(xué)習(xí)預(yù)測(cè)模型[10]等，但目前成熟的方法仍是基于代理模型的優(yōu)化方法[11-12]。此外，實(shí)驗(yàn)設(shè)計(jì)的優(yōu)化方法也很多，包括拉丁超立方設(shè)計(jì)(LHD)[13]、均勻設(shè)計(jì)[14]和正交設(shè)計(jì)[15]等，但這些設(shè)計(jì)方法也存在不足，其生成的樣本只能均勻分布在樣本空間內(nèi)。對(duì)于氣膜冷卻優(yōu)化問題，不同變量對(duì)于目標(biāo)的影響程度通常不同[16]?；诖?，筆者提出兩步法的扇形孔優(yōu)化方法，對(duì)扇形孔的不同參數(shù)進(jìn)行敏感性分析，確定目標(biāo)的優(yōu)化變量，以減少優(yōu)化變量個(gè)數(shù)，保證代理模型的“凸”性，從而提高優(yōu)化效率和收斂性。

1 優(yōu)化方法

1.1 優(yōu)化變量和目標(biāo)函數(shù)

本文采用新的成型方法對(duì)扇形孔進(jìn)行建模，成型過程見圖1，其中S表示厚度。圓柱段形狀由直徑D、傾斜角α和圓柱段長(zhǎng)度Lm確定。從A向視圖可以看出，擴(kuò)張段形狀由出口前緣長(zhǎng)度L1、尾緣長(zhǎng)度L2和出口寬度H確定，具體過程是用樣條曲線連接圓柱段出口和氣膜孔出口，通過掃掠形成擴(kuò)張段。參考文獻(xiàn)[4]，設(shè)置傾斜角為35°，圓柱孔徑為6 mm，橫向孔間距為13D。優(yōu)化變量為前緣長(zhǎng)度L1、尾緣長(zhǎng)度L2、出口寬度H和圓柱段長(zhǎng)度Lm，其余參數(shù)保持不變。優(yōu)化變量的取值范圍見表1。

表1 優(yōu)化變量范圍Tab.1 Range of optimization variables

圖1 扇形孔成型過程示意圖Fig.1 Schematic diagram of fan-shaped hole forming process

(1)

式中：ηLat為橫向平均冷卻效率。

(2)

(3)

式中：η為絕熱冷卻效率；T∞、Tc和Taw分別為主流溫度、冷氣溫度和絕熱壁溫。

1.2 敏感性分析

在敏感性分析過程中，選擇主效應(yīng)和貢獻(xiàn)率作為比較優(yōu)化變量對(duì)目標(biāo)函數(shù)敏感性大小的標(biāo)準(zhǔn)。其中，主效應(yīng)表示因子在某個(gè)水平時(shí)所有響應(yīng)的平均值，貢獻(xiàn)率表示因子對(duì)響應(yīng)影響程度的百分比，按因子對(duì)響應(yīng)的貢獻(xiàn)排列可得到Pareto圖。

1.3 Kriging代理模型及遺傳算法

Kriging代理模型(KSM)[13]是一種利用隨機(jī)變量估計(jì)非設(shè)計(jì)點(diǎn)上數(shù)值的統(tǒng)計(jì)方法。未知函數(shù)y(x)的近似表達(dá)式為：

y(x)=f(x)+z(x)

(4)

式中：f(x)為用常數(shù)項(xiàng)表示的全局函數(shù)；z(x)為均值為0、帶有協(xié)方差的高斯過程函數(shù)。

識(shí)別2個(gè)或多個(gè)隨機(jī)變量相關(guān)性的z(x)的協(xié)方差可以用總體方差σ2、相關(guān)矩陣R以及相關(guān)函數(shù)R(xi,xj)來表示。

Cov(xi,xj)=σ2R[R(xi,xj)]

(5)

(6)

式中：m為訓(xùn)練樣本個(gè)數(shù);θk為相關(guān)參數(shù);xi、xj均為相關(guān)變量。

如果確定相關(guān)函數(shù)，可得到：

(7)

(8)

(9)

(10)

θk≥0(k=1,2,…,m)

由θ構(gòu)成的Kriging代理模型為擬合精度最優(yōu)的代理模型。

多島遺傳算法(Multi-island GA)通過模仿生物進(jìn)化過程中的遺傳繁殖機(jī)制，對(duì)優(yōu)化問題解空間的個(gè)體進(jìn)行編碼,然后對(duì)編碼后的個(gè)體種群進(jìn)行遺傳操作，通過迭代從新種群里尋找含有最優(yōu)解的組合。流程如下：(1)初始化種群；(2)計(jì)算種群中每個(gè)個(gè)體的適應(yīng)度；(3)按由個(gè)體種群適應(yīng)度所決定的規(guī)則選擇進(jìn)入下一代的個(gè)體；(4)按交叉概率Pc進(jìn)行交叉操作；(5)按變異概率Pm進(jìn)行突變操作；(6)如果沒有滿足停止條件，轉(zhuǎn)到步驟(2)，否則轉(zhuǎn)到步驟(7)；(7)輸出種群里適應(yīng)度最優(yōu)的染色體作為問題的滿意解或最優(yōu)解。

1.4 優(yōu)化流程

圖2為優(yōu)化流程示意圖。根據(jù)表1中的優(yōu)化變量范圍，利用LHD方法生成樣本，基于計(jì)算流體動(dòng)力學(xué)(CFD)數(shù)值方法確定目標(biāo)函數(shù)集，通過Kriging代理模型建立氣膜冷卻的預(yù)測(cè)模型，采用多島遺傳算法尋優(yōu)，尋優(yōu)過程中子群數(shù)量為10個(gè)，島的數(shù)量為10個(gè)，迭代過程為10代。

圖2 優(yōu)化流程示意圖Fig.2 Optimization process diagram

為比較基于敏感性分析的優(yōu)化效果，進(jìn)行兩組優(yōu)化。第1組優(yōu)化一次性生成22個(gè)樣本，不進(jìn)行敏感性分析，使用Kriging代理模型進(jìn)行優(yōu)化,即一步法，優(yōu)化孔命名為優(yōu)化孔A。第2組先生成11個(gè)樣本，進(jìn)行敏感性分析，縮小樣本空間后再生成11個(gè)樣本，總共22個(gè)樣本，使用Kriging代理模型進(jìn)行優(yōu)化,即兩步法,優(yōu)化孔命名為優(yōu)化孔B。

2 CFD數(shù)值方法

2.1 計(jì)算域和邊界條件

圖3給出了模型的幾何尺寸，該模型由主流通道、氣膜孔通道和供氣室組成。主流通道長(zhǎng)、寬、高分別為170D、13D和17D，供氣室長(zhǎng)寬高分別10D、13D和12D。

圖3 計(jì)算模型及邊界條件Fig.3 Calculation model and boundary conditions

主流采用速度入口和壓力出口，主流通道兩側(cè)為周期性邊界條件，冷氣入口為質(zhì)量流量入口，其余邊界為絕熱無滑移壁面，通過改變冷氣質(zhì)量流量qm,c來調(diào)整吹風(fēng)比M。

(11)

式中：ρ∞、ρc分別為主流和冷氣密度；U∞為主流入口速度；Ac為氣膜孔圓柱段的截面積;Uc為氣膜孔圓柱段截面的面平均速度。

主流入口速度U∞為20 m/s，主流溫度T∞為317.35 K，冷氣溫度Tc為299.25 K，主流和冷氣入口湍流度均為1%。選擇Ansys Fluent為數(shù)值模擬的求解器，湍流模型選用Realisablek-ε模型。

吹風(fēng)比對(duì)氣膜優(yōu)化結(jié)果有一定影響[17]。一般吹風(fēng)比范圍為0.5～2.5，本文取1.5，僅進(jìn)行單工況優(yōu)化，將其他吹風(fēng)比的冷卻效果作為后續(xù)檢驗(yàn)。

2.2 網(wǎng)格無關(guān)性驗(yàn)證和數(shù)值驗(yàn)證

選擇Pointwise17.3作為網(wǎng)格劃分軟件。圖4給出了模型網(wǎng)格拓?fù)浣Y(jié)構(gòu)，孔出口和氣膜覆蓋區(qū)域由大密度六面體網(wǎng)格填充，氣膜孔通道由四面體網(wǎng)格填充，其余部分由三棱柱網(wǎng)格填充。采用3套網(wǎng)格對(duì)氣膜橫向平均冷卻效率進(jìn)行無關(guān)性驗(yàn)證，網(wǎng)格數(shù)分別為140萬、310萬和480萬，見圖5。網(wǎng)格數(shù)分別為310萬和480萬時(shí)，橫向平均冷卻效率變化不大，因此選擇網(wǎng)格數(shù)為310萬，第1層網(wǎng)格高度為0.01 mm，拉伸比為1.1，y+<1。

圖4 模型網(wǎng)格拓?fù)浣Y(jié)構(gòu)Fig.4 Topology structure of model grid

圖5 網(wǎng)格無關(guān)性驗(yàn)證Fig.5 Grid independence verification

為驗(yàn)證CFD數(shù)值方法的準(zhǔn)確性，選擇扇形孔作為驗(yàn)證模型，并與文獻(xiàn)[4]中的實(shí)驗(yàn)結(jié)果進(jìn)行對(duì)比。圖6給出了M=1.5時(shí)橫向平均冷卻效率計(jì)算值與實(shí)驗(yàn)值的對(duì)比結(jié)果。從圖6可以看出，基于Realisablek-ε模型的橫向平均冷卻效率計(jì)算值與實(shí)驗(yàn)值更接近，因此選擇Realisablek-ε作為計(jì)算模型。

圖6 橫向平均冷卻效率計(jì)算值與實(shí)驗(yàn)值的對(duì)比

3 結(jié)果與分析

3.1 敏感性分析結(jié)果

表2給出了兩步法第1輪樣本計(jì)算結(jié)果。圖7為通過對(duì)優(yōu)化變量的敏感性進(jìn)行分析得到的Pareto圖。從圖7可以看出，Lm對(duì)目標(biāo)的貢獻(xiàn)率最大，說明Lm發(fā)生變化會(huì)對(duì)面平均冷卻效率造成較大影響；L1和L2的貢獻(xiàn)率次之，H的貢獻(xiàn)率最小，敏感度也最小，其對(duì)面平均冷卻效率的影響很有限。

表2 兩步法第1輪樣本計(jì)算結(jié)果

圖7 Pareto圖Fig.7 Pareto diagram

從圖7可以看出，Lm對(duì)面平均冷卻效率產(chǎn)生的是負(fù)效應(yīng)，面平均冷卻效率隨Lm的增大而減小；而L1、L2和H對(duì)面平均冷卻效率產(chǎn)生的是正效應(yīng)，增大其數(shù)值有利于提高面平均冷卻效率。

為便于比較，采用式(12)對(duì)各優(yōu)化變量進(jìn)行無量綱化，其中X表示各變量，Xmax、Xmin分別為變量上下限，Xi為第i個(gè)變量取值。

(12)

圖8 優(yōu)化變量對(duì)目標(biāo)函數(shù)的影響Fig.8 Influence of optimization variables on the objective function

綜上，Lm、L1和L2的敏感度較大，因此重新確定這3個(gè)優(yōu)化變量的取值范圍來生成第2輪樣本，找出不同優(yōu)化變量對(duì)應(yīng)面平均冷卻效率較高的區(qū)域，以該區(qū)域作為第2輪優(yōu)化變量的取值范圍。具體過程可表示為：

(13)

其中CX是常數(shù)，表示優(yōu)化變量X對(duì)應(yīng)的面平均冷卻效率下限，該下限需要人為確定。第2輪樣本各優(yōu)化變量的取值范圍見表3，基于LHD生成的兩步法第2輪樣本計(jì)算結(jié)果見表4。

表3 敏感性分析后優(yōu)化變量的取值范圍

表4 兩步法第2輪樣本計(jì)算結(jié)果

3.2 優(yōu)化結(jié)果對(duì)比

表5給出了M=1.5時(shí)原型孔和優(yōu)化孔的優(yōu)化變量和面平均冷卻效率的預(yù)測(cè)結(jié)果。2種優(yōu)化孔的面平均冷卻效率均大于原型孔。其中，優(yōu)化孔B的面平均冷卻效率計(jì)算值較優(yōu)化孔A更大，且預(yù)測(cè)誤差更小。這說明基于敏感性分析的優(yōu)化結(jié)果精度更高，預(yù)測(cè)效果更好。

表5 原型孔和優(yōu)化孔的幾何參數(shù)及預(yù)測(cè)結(jié)果對(duì)比

如圖9所示,從氣膜的覆蓋范圍看，尋優(yōu)方向主要是拓展氣膜側(cè)向?qū)挾?，即L1和L2的長(zhǎng)度，氣膜側(cè)向?qū)挾仍龃?，氣膜覆蓋面積增大。如圖10所示，優(yōu)化孔B的氣膜側(cè)向覆蓋范圍大于優(yōu)化孔A的氣膜側(cè)向覆蓋范圍，其更能體現(xiàn)L1和L2對(duì)目標(biāo)的敏感度，因此冷卻效果更好。

(a) 原型孔

(b) 優(yōu)化孔A

(c) 優(yōu)化孔B圖9 M=1.5時(shí)氣膜有效度的對(duì)比Fig.9 Comparison of film effectiveness at M=1.5

圖10 M=1.5時(shí)距離孔出口尾緣2D處絕熱冷卻效率的對(duì)比

圖11給出了吹風(fēng)比為1.5時(shí)3種孔型橫向平均冷卻效率的對(duì)比。與原型孔相比，2種優(yōu)化孔均能夠提高橫向平均冷卻效率，其中優(yōu)化孔B的橫向平均冷卻效率更高，說明基于敏感性分析的優(yōu)化方法能夠提高優(yōu)化精度，優(yōu)化效果更好。

圖11 M=1.5時(shí)不同孔型橫向平均冷卻效率的對(duì)比

3.3 吹風(fēng)比的影響

圖12給出了吹風(fēng)比對(duì)不同孔型氣膜有效度的影響。從圖12可以看出，2種優(yōu)化孔的氣膜有效度分布具有2個(gè)相似的特點(diǎn)。第一，在相同吹風(fēng)比下，相比原型孔，2種優(yōu)化孔在氣膜孔出口區(qū)域的絕熱冷卻效率有較大提升，沿著流向氣膜能延伸到下游更遠(yuǎn)的區(qū)域；優(yōu)化孔A和優(yōu)化孔B均沿橫向增大了氣膜覆蓋面積，使得冷氣能在更大范圍內(nèi)鋪展開，從而提升了冷卻效果。第二，當(dāng)吹風(fēng)比增大時(shí)，原型孔的冷卻效果提升不明顯，其絕熱冷卻效率和氣膜覆蓋面積變化不大，在M=2.5時(shí)原型孔氣膜有效度分布甚至?xí)霈F(xiàn)分叉現(xiàn)象，冷卻效果降低；優(yōu)化孔A和優(yōu)化孔B的絕熱冷卻效率和氣膜覆蓋面積均隨吹風(fēng)比的增大而增大，在大吹風(fēng)比下絕熱冷卻效果最好。與優(yōu)化孔A相比，優(yōu)化孔B的氣膜沿橫向能夠拓寬更大的距離，在相同吹風(fēng)比下其氣膜覆蓋面積更大，隨著吹風(fēng)比的增大,這種差異會(huì)逐漸體現(xiàn)出來。

(a) 原型孔、M=0.5

(d) 原型孔、M=1.0

(g) 原型孔、M=1.5

(j) 原型孔、M=2.0

(m) 原型孔、M=2.5

圖13給出了不同孔型面平均冷卻效率的對(duì)比。在M=2.0時(shí)原型孔存在最大面平均冷卻效率，而優(yōu)化孔B的面平均冷卻效率隨吹風(fēng)比的增大而提高，在M=2.5時(shí)面平均冷卻效率最大。雖然優(yōu)化孔A和優(yōu)化孔B的面平均冷卻效率的變化趨勢(shì)相似，但在相同吹風(fēng)比下優(yōu)化孔B的面平均冷卻效率更高，尤其在M=2.5時(shí)，優(yōu)化孔B的面平均冷卻效率比優(yōu)化孔A高出0.03。因此，基于敏感性分析的優(yōu)化方法不僅精度高，其預(yù)測(cè)結(jié)果也更準(zhǔn)確。

圖13 不同孔型面平均冷卻效率的對(duì)比

4 結(jié) 論

(1) 氣膜孔結(jié)構(gòu)參數(shù)對(duì)冷卻效果的影響存在差異，其中扇形孔出口前緣長(zhǎng)度L1、尾緣長(zhǎng)度L2和圓柱段長(zhǎng)度Lm對(duì)面平均冷卻效率的影響較大，敏感度較強(qiáng)，而出口寬度H的影響最小，敏感度最弱；Lm對(duì)面平均冷卻效率產(chǎn)生的是負(fù)效應(yīng)，而L1、L2和H對(duì)面平均冷卻效率產(chǎn)生的是正效應(yīng)。

(2) 經(jīng)過敏感性分析后的樣本分布更加合理，代理模型的預(yù)測(cè)精度和尋優(yōu)計(jì)算的收斂性均得到提高，結(jié)果均優(yōu)于無敏感性分析。