李 鑫 ,張 韜 ,李偉偉 ,周 玲 ,季路成
(清華大學航空發(fā)動機研究院1,車輛與運載學院2:北京 100084;3.北京動力機械研究所,北京 100024;4.北京理工大學宇航學院,北京 100081)
現(xiàn)代先進航空發(fā)動機的內(nèi)部流動非線性程度高[1]、設(shè)計參數(shù)眾多,完全依靠設(shè)計師憑經(jīng)驗進行手動設(shè)計已難以滿足日益提高的發(fā)動機性能發(fā)展需求[2]。因此,依托計算機和計算技術(shù)對以眾多參數(shù)描述的葉片通道進行精細設(shè)計與優(yōu)化[3]是發(fā)展現(xiàn)實與趨勢。在計算機優(yōu)化算法中,全局隨機類優(yōu)化算法與常規(guī)梯度類優(yōu)化算法均無法應(yīng)對設(shè)計參數(shù)眾多的精細化設(shè)計需求。
伴隨方法具有計算成本與設(shè)計變量數(shù)目無關(guān)的特性,是最為契合精細化設(shè)計概念的合理選擇。該方法由Lions 于1971 年首次提出[4],2000 年由Liu[5]、Dreyer 等[6]用于葉輪機氣動外形優(yōu)化。歷經(jīng)近20 年的發(fā)展,伴隨優(yōu)化已成為葉輪機精細化設(shè)計與優(yōu)化的重要工具,在該方法基礎(chǔ)上,多點優(yōu)化[7]、多排葉輪機優(yōu)化[8-9]、多學科一體化優(yōu)化[10]、魯棒性設(shè)計[11]及氣動不確定性分析[11]等問題均得到了充分研究。在伴隨優(yōu)化系統(tǒng)流程中,參數(shù)化方法描述了離散的設(shè)計變量取值與連續(xù)的外形幾何間一一對應(yīng)的映射關(guān)系,因而成為連接流動數(shù)值仿真與目標函數(shù)梯度分析、梯度尋優(yōu)的重要橋梁。2017 年Xu 等[12]指出,參數(shù)化方法應(yīng)具有全自動、空間大、無震蕩、便于施加幾何約束等系統(tǒng)特征。歷經(jīng)近20 年的發(fā)展,以擾動控制為主的擾動參數(shù)化方法具有更高的設(shè)計空間與靈活性,更符合葉輪機精細化設(shè)計的需求,也因此被廣泛應(yīng)用,如Hicks-henne 函數(shù)、Cst 函數(shù)等[13]。作為3 維空間的隨機擾動施加方法,自由變形技術(shù)(Free-form Deforma?tion,F(xiàn)fd)[14]自2018 年由John 等[15]首次用于葉輪機伴隨優(yōu)化后,便因其設(shè)計空間廣、3 維連續(xù)性強等優(yōu)勢獲得研究人員的關(guān)注。然而,F(xiàn)fd 應(yīng)用過程卻無法與葉輪機幾何分布特征完全契合。從原理上講,F(xiàn)fd 變形過程受限于長方體控制體,難以與圓柱形流道完全匹配,且盡管其能夠保證控制體內(nèi)部實體的曲線連續(xù),于邊界位置也難以施加約束,難以保證回轉(zhuǎn)幾何的周期性。
2013 年研究團隊研發(fā)了基于簡單H 網(wǎng)格、凍結(jié)粘性法、薄層削尖模型的連續(xù)型伴隨優(yōu)化系統(tǒng)Opti?turbo[16],經(jīng)后續(xù)不斷改進,已具備工程可行的伴隨優(yōu)化應(yīng)用能力[17];為解決由網(wǎng)格及湍流模型簡化引發(fā)的數(shù)值精度問題,2017 年進一步研發(fā)了葉輪機氣動外形離散型伴隨優(yōu)化系統(tǒng)Turbosim_Un[18],實現(xiàn)了流動數(shù)值仿真及伴隨敏感性分析精度的顯著改進。然而,其仍使用基于Hicks-henne 函數(shù)[19-20]的擾動參數(shù)化方法,不僅操作復雜,也嚴重限制了設(shè)計空間。為此,從實用的需求出發(fā),順應(yīng)行業(yè)領(lǐng)域發(fā)展趨勢,本文以拓展自由變形技術(shù)(Extended Free-form Deformation,EFFD)為基礎(chǔ)開發(fā)了具有葉片-端壁協(xié)調(diào)控制能力的全新參數(shù)化方法,并對跨聲速壓氣機轉(zhuǎn)子Rotor 67 進行優(yōu)化,在驗證TurboSim_Un 優(yōu)化能力的同時保證參數(shù)化方法的有效性。
拓展自由變形技術(shù)是Coquillart[21]于1900 年對FFD 方法的改進拓展,其主要改進在于對長方體控制體進行變形以構(gòu)造任意形狀控制體,并結(jié)合局部變形概念改善控制體邊界附近曲線的連續(xù)性。文獻[20]已對EFFD 理論進行了詳細描述,本文僅以Rotor 67轉(zhuǎn)子通道為例展示其實施過程:
(1)以通道形狀為控制體,繪制均勻分布的控制點。與流域空間不對等的長方體控制體是FFD 技術(shù)無法實現(xiàn)葉片-端壁一體化設(shè)計的主要限制,而EFFD 方法中通過控制體變形可構(gòu)造任意形狀控制體。為此,直接以通道流域為控制體進行變形處理,并均勻分布控制體內(nèi)部控制點,如圖1(A)所示。
圖1 拓展自由變形過程
(2)對于任意實體位置,應(yīng)用局部變形理念確定其所在的局部控制體位置。為解決控制體邊界附近的曲線變形連續(xù)性等問題,EFFD 方法需使用局部變形,即僅實體位置附近的控制點干擾方能改變實體位置。局部控制體位置計算需利用Bezier 超曲面性質(zhì),以判斷包含該實體位置的具體子控制體微元。
(3)對控制點施加擾動,實現(xiàn)實體變形。施加某隨機擾動后的邊界控制點位置如圖1(B)所示。
(4)自動計算局部控制體內(nèi)虛擬控制點位置。為保證不同子控制體間曲線變形的斜率連續(xù),應(yīng)在子控制體內(nèi)施加虛擬控制點,其不參與參數(shù)化變形,而根據(jù)真實控制點的變化自動調(diào)整位置,如圖1(C)所示。
(5)使用牛頓迭代,計算實體位置在局部控制體內(nèi)的相對坐標。
(6)對子控制體進行FFD 變形,計算公式見文獻[14],最終變形結(jié)果如圖1(d)所示。
Rotor 67 轉(zhuǎn)子葉片是20 世紀70 年代由美國Nasa設(shè)計的低展弦比跨聲速風扇/壓氣機轉(zhuǎn)子葉片,1980年Anthony 等使用激光測速等方法獲得了其詳細的內(nèi)部流場數(shù)據(jù)[22],包括特性線、截面馬赫數(shù)云圖及流向速度分布等,該轉(zhuǎn)子葉片已成為葉輪機領(lǐng)域廣為應(yīng)用的Cfd 數(shù)值校驗、流動機理分析和優(yōu)化設(shè)計的典型案例。Rotor 67轉(zhuǎn)子部分設(shè)計及試驗參數(shù)見表1。
表1 Rotor 67轉(zhuǎn)子部分設(shè)計及試驗參數(shù)
優(yōu)化過程Rotor 67 轉(zhuǎn)子計算流域進口與葉根前緣的距離為0.10 M,流域出口與葉根尾緣的距離為0.12 M,子午面如圖2 所示。采用Autogrid 5 繪制網(wǎng)格,基于Sa(Spalart And Allmaras)模型的網(wǎng)格無關(guān)性校驗結(jié)果如圖3所示。
圖2 Rotor 67子午面
圖3 網(wǎng)格無關(guān)性校驗
經(jīng)過選擇,后續(xù)Rotor 67轉(zhuǎn)子優(yōu)化計算使用58萬網(wǎng)格方案,其壁面y+值如圖4 所示,基于此網(wǎng)格的Ro?tor 67 轉(zhuǎn)子特性線如圖5 所示,其中下標choke 表示堵塞流量,本文后續(xù)堵塞流量均以圖5為基準。
圖5 Rotor 67轉(zhuǎn)子特性線
優(yōu)化過程以效率為目標函數(shù),并采用罰函數(shù)法對流量及壓比進行約束
式中:η為絕熱效率;ω1為流量約束權(quán)值,取值約為1000;ω2為壓比約束權(quán)值,取值為200;M為流量;Pr為總壓比;下標0表示原型的參數(shù)。
在采用EFFD 參數(shù)化方法的優(yōu)化過程中,軸向控制點數(shù)為10,周向控制點數(shù)為7,徑向控制點數(shù)為8,除保證曲線連續(xù)型的若干位置控制點、排除葉尖復雜流動干擾的尖部控制點不變外,總變形控制點數(shù)為462,總設(shè)計參數(shù)的數(shù)量為1386。
在上述設(shè)置下,本文Rotor 67 轉(zhuǎn)子伴隨優(yōu)化過程中相關(guān)參數(shù)變化如圖6所示。
圖6 優(yōu)化過程中性能參數(shù)變化
圖中,橫軸ITER 表示優(yōu)化后所獲得的幾何模型次序編號,左縱軸為流量相對偏差δm及總壓比相對偏差δpr,右縱軸為效率絕對偏差ζη
由此可知,經(jīng)16 次伴隨優(yōu)化,Rotor 67 轉(zhuǎn)子效率提升了0.74%,而流量及總壓比偏差較小,分別為0.42%及0.11%。
為檢驗優(yōu)化前后Rotor 67 轉(zhuǎn)子的整體性能,分別展示了原型及優(yōu)化后結(jié)果在70%、90%及100%轉(zhuǎn)速下的特性線分布,優(yōu)化前后特性線變化如圖7 所示。從圖中可見,本次伴隨優(yōu)化已實現(xiàn)轉(zhuǎn)子氣動性能的顯著改進,并未發(fā)生工況漂移問題,且100%轉(zhuǎn)速狀態(tài)由于近峰值效率點的工況約束程度較高,堵塞流量僅增大約0.24%,而總壓比特性線與原型的仍保持較高的吻合度。具體而言,在100%轉(zhuǎn)速下轉(zhuǎn)子等熵效率提升達到近0.8%;此外,在90%和70%轉(zhuǎn)速下的變工況特性中,峰值效率及近堵塞工況的效率指標同樣有所提升,且隨著轉(zhuǎn)速的降低,流量及壓比偏差減小,但效率提升也減小。
圖7 優(yōu)化前后特性線變化
為描述伴隨優(yōu)化前后Rotor 67 轉(zhuǎn)子葉片發(fā)生的幾何變動,展示了若干展向位置的葉型截面,如圖8 所示。從圖中可見,葉片變形集中在前緣及前緣吸力面附近,而葉片后方變形極小,確保在多排一體化應(yīng)用時優(yōu)化結(jié)果不會對下游工況產(chǎn)生過大干擾。具體而言,前緣均向壓力面?zhèn)绕D(zhuǎn),且往往伴隨著前緣厚度的顯著減小,越靠近根部則絕對葉片厚度減小越明顯。
圖8 葉型變化
厚度是保證優(yōu)化結(jié)果結(jié)構(gòu)強度的重要指標,則進一步提取了相應(yīng)展向位置的葉片厚度及中弧線彎角分布,如圖9 所示。從圖中可見,最大厚度位置前的葉片厚度減小是各截面葉型的共同變化特征,但最大厚度及其位置仍與原型的保持較高的吻合度,基本保證結(jié)構(gòu)強度不會大幅降低;另一方面,中弧線在保持彎角分布規(guī)律一致的前提下,在前緣后方局部區(qū)域內(nèi)彎角減小而后整體增大,且隨展高的增大變化更為明顯。
圖9 葉型厚度與彎角分布
優(yōu)化結(jié)果葉表變形量如圖10 所示,更形象地展示了葉片變化位置。從圖中可見,盡管葉根附近變形量絕對值更大,但考慮輪轂徑向變化的參與,其葉型變化并不十分明顯(圖8);除葉尖幾何未參與氣動優(yōu)化外,80%展高附近葉片變化最為劇烈,且變化以前緣后方位置為主;尾緣附近變形均整體偏小,與圖8變化對應(yīng)。
圖10 優(yōu)化結(jié)果葉表變形量分布
從圖10 中還可見,葉根附近較大的變形量表明輪轂已發(fā)生顯著改變,作為葉片-端壁一體化優(yōu)化的另一主要調(diào)節(jié)對象,輪轂面的徑向變化如圖11 所示。從圖中可見,經(jīng)伴隨優(yōu)化,輪轂面具有典型“波浪形”分布特征:在前緣前、尾緣前葉表兩側(cè)存在3 組較為明顯的通道下壓區(qū),而在吸力面前緣附近、尾緣前通道中部則明顯上抬;整體而言,通道變形形式較為復雜,在任意軸向位置處輪轂壁面型線不再保持等徑或單調(diào)分布,S3截面部分葉片及輪轂形狀如圖12所示。
圖11 模型6輪轂面變形量
圖12 S3截面部分葉片及輪轂形狀
圖中,SS表示吸力面,PS表示壓力面。從圖中可見,經(jīng)伴隨優(yōu)化,除前緣附近外,葉片中后部區(qū)域內(nèi)原本近乎90°的葉片及端壁夾角轉(zhuǎn)化為鈍角,且尤以吸力面角區(qū)更為明顯;此外,端壁曲線呈現(xiàn)典型的“下壓-上抬-下壓”變化:吸力面及端壁交匯區(qū)域以二面角變化為主,通道中部輪轂上抬形成尖峰,而在壓力面與輪轂一體化區(qū)域并未顯著改變二面角,整體表現(xiàn)為平移下沉。
為判斷流動性能提升的內(nèi)在機理,提取了100%轉(zhuǎn)速峰值效率點位置的出口截面內(nèi)若干參數(shù)展向分布,如圖13所示。其中,各相對高度下絕熱效率為
圖13 出口截面若干氣動性能參數(shù)展向分布
式中:Pt為總壓;Tt為總溫;γ為比熱比;下標in表示進口,out表示出口。
經(jīng)伴隨優(yōu)化,近乎在全展范圍內(nèi)效率指標均顯著提升,而60%展高以下區(qū)域總壓均顯著提高,在壓比約束下尖部總壓略有降低。
以20%及80%展高為例,從S1截面的視角展示了優(yōu)化前后流動的變化,如圖14 所示。從圖中可見,盡管20%展高尚未形成激波結(jié)構(gòu),但吸力面加速、尾跡均是流道內(nèi)主要流動現(xiàn)象。在葉根附近,經(jīng)伴隨優(yōu)化,前緣吸力面?zhèn)雀咚賲^(qū)范圍明顯縮減,并結(jié)合葉片前半段的厚度調(diào)整,吸力面加速區(qū)范圍同樣顯著減?。辉诹魉僮兓绊懴?,尾跡范圍明顯減小,表征尾跡損失顯著改善。而在80%展高截面內(nèi)盡管尾跡、激波結(jié)構(gòu)變化并不明顯,但借助厚度調(diào)整,激波前高馬赫數(shù)的范圍和數(shù)值均略有減小,同樣表征激波損失的改善。
圖14 若干展向位置相對馬赫數(shù)分布
優(yōu)化前后葉表吸力面靜壓分布如圖15 所示。從圖中可見,優(yōu)化結(jié)果的低壓區(qū)范圍明顯縮減,對應(yīng)葉表高等熵馬赫數(shù)區(qū)域減小,與前述結(jié)論一致。
圖15 優(yōu)化前后葉表靜壓分布
角區(qū)和輪轂形狀的變化及端區(qū)附近流動性能的顯著提升,表明輪轂附近流動已大幅改善,為此,展示了圖12中截取位置2的周向S3截面馬赫數(shù)分布,如圖16 所示。從圖中可見,相較于原型,優(yōu)化結(jié)果的主流區(qū)流速變化較小,而端區(qū)流速略有提高;以相同的速度等值線與端壁型線間相對面積作為衡量附面層損失的近似判據(jù)判斷其符合二面角原理,增大的二面角分布對改善角區(qū)附面層交匯損失具有正向增益作用。
圖16 截取位置2的S3截面相對馬赫數(shù)分布
優(yōu)化結(jié)果的輪轂面型線具有典型的“下壓-上抬-下壓”分布特征,其中部上抬的輪轂可有效抑制角區(qū)橫向二次流傳播,進而改善吸力面分離渦結(jié)構(gòu),吸力面?zhèn)任簿壐浇谋诿鏄O限流線如圖17 所示。從圖中可見,經(jīng)伴隨優(yōu)化,尾緣后方的渦結(jié)構(gòu)甚至基本消失,吸力面后方的渦尺寸同樣略有減小,表明端區(qū)性能有極大改善。優(yōu)化前后通道出口截面的總壓分布如圖18 所示。從圖中可見,根部區(qū)域的低總壓區(qū)范圍明顯減小,再次表明根部區(qū)域流動性能顯著改善。
圖17 吸力面?zhèn)任簿壐浇谋诿鏄O限流線
圖18 優(yōu)化前后通道出口截面總壓分布
然而,考慮壁面極限流線(圖17)為基于SA 湍流模型的計算結(jié)果,其在復雜湍流狀態(tài)下往往誤差較大。為此,以CFX 商業(yè)軟件為工具,使用剪切應(yīng)力輸運(Shear Stress Transport, SST)湍流模型在相同邊界條件下對原型及優(yōu)化結(jié)果進行氣動仿真,并提取其相應(yīng)位置的流線分布,分離渦結(jié)構(gòu)計算結(jié)果如圖19所示。從圖中可見,經(jīng)伴隨優(yōu)化,尾緣后方渦結(jié)構(gòu)獲得了極大改善,整體表現(xiàn)為展向尺寸及軸向范圍大幅縮減。
圖19 分離渦結(jié)構(gòu)SST模型CFX計算結(jié)果
(1)經(jīng)伴隨優(yōu)化,在100%轉(zhuǎn)速下Rotor 67 轉(zhuǎn)子近峰值效率點效率提升0.74%,且工況約束較高,特性線同樣得到全面改善。
(2)葉片與端壁間二面角的調(diào)整,對調(diào)節(jié)端區(qū)附面層交匯進而調(diào)節(jié)角區(qū)分離具有重要影響,是改善端區(qū)氣動性能時不可忽視的重要設(shè)計空間。
(3)拓展自由變形技術(shù)具備葉片通道的全3維調(diào)節(jié)能力,結(jié)合伴隨方法,可有效探索葉片通道的全3維設(shè)計空間。