張 強,楊 永
(西北工業大學 翼型、葉柵空氣動力學國防科技重點實驗室,陜西 西安 710072)
湍流是自然界中普遍存在的一種流動形態。在宏觀上,湍流表現為流體微團的不規則運動,湍流運動產生的質量和能量的輸運遠大于分子熱運動所產生的輸運,湍流運動引起的質量和能量的平均擴散是不能忽略的。因此湍流研究一直是流體力學研究中的一個重要方向。
國內外針對平板、翼型、機翼等外流標準算例,開展的數值模擬研究已經很多,但針對噴管內流的湍流模型研究主要以Spalart-Allmaras(SA)一方程模型或者改進的k-ε兩方程模型為主,其它湍流模型的驗證研究還比較少見。現今兩方程湍流模型的研究按研究角度主要分為三類:構造k-ε和k-ω的混合模型;構造更為精確的剪應力構成模型;構造新的更準確的耗散方程。其代表模型分別是:SSTk-ω兩方程模型,EASMk-ω模型和基于新的耗散方程的k-ζ模型。
國外對射流推力矢量噴管已進行了系統的試驗研究[1],形成了較為完整的數據庫,為CFD數值模擬研究提供了良好的條件。因此本文的研究目的就是針對激波控制射流推力矢量噴管的內外流場、噴管性能參數等,分別采用SSTk-ω,EASMk-ω和k-ζ模型,進行系統的驗證分析,并與SA模型的計算結果和試驗數據進行對比,考察不同湍流模型對噴管內外流場的流場結構、噴管性能參數以及主噴管內壁面的壓力分布等的預測能力。之所以選擇激波控制射流推力矢量噴管作為驗證算例,主要是因為其流場中由于二次射流的干擾,會在主噴流中誘導產生斜激波,存在較強的激波/附面層干擾,這樣的流動對目前湍流模型的研究也是一種極大的挑戰,因此其數值模擬結果也更具說服力。
研究激波控制推力矢量噴管時,所采用的模型是NASA LaRC試驗用的噴管,設計主流壓比為NPRD=8.78,其幾何形狀及尺寸見參考文獻[2],單縫引射對應外形1,雙縫引射對應外形5。
單縫引射的推力矢量噴管,計算網格如圖1所示,圖中生成的計算網格共包括5塊,包圍主噴管的箱體上下各有一塊網格,主噴管和次流噴管各生成一塊網格,主噴管的出口及箱體上下兩塊網格的下游,點對點生成一塊網格,以計算噴管的尾跡流場。雙縫引射噴管外部網格拓撲結構與單縫引射噴管相同,只是在主噴管擴張段上多了另外一個小噴口的網格,如圖2所示,并且兩個噴口總的開口面積和單縫時的開口面積相同。計算中參考長度取為從主噴管出口到其喉道的長度c,由于噴管總長度為2c,因此后面的噴管無量綱化的長度為2,這樣處理以方便與試驗值的比較。

圖1 單縫引射噴管的計算網格及噴口局部圖Fig.1 Computational mesh and close-up of one injection port

圖2 雙縫引射噴管的噴口局部圖Fig.2 Close-up of two injection ports
針對射流推力矢量噴管的數值模擬,采用的是課題組自有的研究型軟件。假定射流推力矢量噴管內外流動是定常的,通過求解不同湍流模型封閉的RANS方程,模擬射流推力矢量噴管的內外流場。RANS流動控制方程的無粘項采用Roe-FDS格式,粘性項采用中心差分,湍流模型的空間項離散也采用二階迎風格式,流動主控方程和湍流模型的控制方程均采用近似因子隱式時間推進,以提高數值計算效率和穩定性,關于算法的詳細介紹可以參考文獻[3-5]。
計算條件取為:單縫引射時外流Ma=0.05,雙縫引射時Ma=0.1,雷諾數Re=1.65×105,分別采用SA、SSTk-ω、EASMk-ω和k-ζ四種湍流模型求解了主流壓比NPR=4.6,7.0,8.78,10.0和SPR=0.7、1.0時的射流推力矢量噴管的內外流場。靜壓固定為標準海平面大氣壓,總溫固定為300K。
利用文獻[2]的試驗結果首先考察了外流馬赫數對數值模擬結果的影響,在兩個外流馬赫數Ma=0.05及Ma=0.1下,采用SA湍流模型,模擬了NPR=4.6、10.0,SPR=0.7時單縫引射的噴管內外流場。計算的主噴管內壁面壓力分布與在靜止大氣中進行的試驗值的比較如圖3(a、b)所示,總的來看,在兩個馬赫數下的數值模擬結果與試驗值比較吻合。在兩個較小的外流馬赫數下,只在NPR=4.6時,噴管出口截面處的壓強略有微小差別,并且我們的數值驗證也表明,在這兩個馬赫數下,噴管性能參數的差別也可以忽略不計,考慮到自由來流馬赫數比較小時(例如Ma=0.05或0.01),流動控制方程不易收斂,為了提高數值計算的穩定性和效率,自由來流馬赫數也可以取0.1,對數值模擬結果的影響不大,可以認為它們是相同的。
Ma=0.05,NPR=4.6,SPR=0.7,SA和SSTk-ω模型給出的噴管流場及二次射流局部流場如圖4所示,由于二次射流的干擾,在次流噴口之前誘導出現斜激波,在次流噴口上游和下游區域均會出現較大的分離區,兩種湍流模型都給出了與試驗很接近的激波結構(圖5),從EASMk-ω和k-ζ模型的計算結果也可以得出相同的結論。四種模型計算給出的主噴管上壁面,次流引起的斜激波位置都比試驗值稍微靠后,約在x/c=1.62處,而試驗測量則是在x/c=1.543處。

圖3 SPR=0.7,Ma=0.05和0.1時,SA模型計算的主噴管內壁面壓力分布與試驗值的比較Fig.3 Experimental and computational pressure distribution at Ma=0.05,0.1with SPR=0.7,SA model

圖4 SA和SST k-ω模型計算馬赫云圖和次噴口局部流場(NPR=4.6,SPR=0.7,Ma=0.05)Fig.4 Mach counters and flow pattern near the injection port for NPR=4.6,SPR=0.7,Ma=0.05with SA and SST model

圖5 NPR=4.6,SPR=0.7時試驗陰影照片Fig.5 Experimental shadowgraph for one injection port,NPR=4.6,SPR=0.7
Ma=0.05,SPR=0.7,NPR分別為4.6、7.0、8.78、10.0時,主噴管內壁面的壓力分布與試驗值的對比如圖6(a~d)所示。在主噴管下壁面,當NPR=4.6時,有一段斜激波在噴管內部,SA模型計算的激波位置最靠近出口,其它三種模型給出的激波位置差別不大,但都優于SA模型的計算結果,更接近實驗值,從喉道到主噴管出口方向依次是k-ζ,SSTk-ω,EASMk-ω模型,在其他NPR下,下壁面不存在斜激波,四種湍流模型計算給出的壓力分布相同,與試驗值吻合良好。
在主噴管上壁面,由于二次射流的影響,在次流噴口的上游和下游位置都會出現較大的分離區,二次噴流誘導產生的斜激波位于上游分離區之前,更接近喉道處,在此激波之前主噴管的上壁面,不同NPR下各湍流模型給出的壓力分布相同,與試驗值完全吻合。在斜激波和次流噴口之間區域,不同NPR下,四種湍流模型給出的激波位置都比試驗值略靠后,但壓力分布差別不大,僅略低于試驗值。在次流噴口以后的區域,受外部回流分離區的影響,不同湍流模型給出的壓力分布略有不同,NPR越大,回流分離區越小,各湍流模型的預測越準確。在圖6(c)中,當NPR=8.78時,在主噴管上壁面次流噴口以后的區域,EASMk-ω的結果與試驗值差別較大,甚至還不如SA模型,這可能與我們沒有考慮流動的非定常效應有關。可見,在有較大分離的區域,四種湍流模型的預測能力還是有差別的。從總的壓力分布來看,相比較而言,SSTk-ω模型最 準確,其次是k-ζ和EASMk-ω,SA模型的誤差最大。
在SPR=0.7時,分別采用SA、SSTk-ω、EASMk-ω、k-ζ四種湍流模型計算的射流推力矢量噴管性能參數隨NPR的變化如圖7所示,不同湍流模型計算給出的噴管性能參數差別不大,隨NPR的增加,流量系數幾乎保持一個常數,推力系數是先增加后減小,矢量偏角和矢量效率都是一直減小的,與試驗得出的規律一致[2]。矢量偏角和矢量效率隨NPR的增加減小的原因是:由于二次射流與主噴管流動相互干擾,誘導生成斜激波,隨著NPR的增加,受主噴流的影響,斜激波末端逐漸外移,遠離主噴管出口,并且強度變弱,導致產生的矢量偏角變小和推力矢量效率降低。

圖6 不同NPR下各湍流模型計算的主噴管壓力分布(SPR=0.7)Fig.6 Pressure distribution of the primary nozzle wall for different NPR(SPR=0.7)

圖7 4種湍流模型計算的噴管性能參數(SPR=0.7)Fig.7 Computational nozzle performance with four turbulence models(SPR=0.7)
圖8是SSTk-ω模型在NPR=4.6,SPR=0.7和1.0時計算給出的主噴管上壁面的壓力分布,可以看出,在相同NPR下,SPR越大,次流噴口前的激波強度越強,位置越靠近喉道,即隨SPR的增加,激波前移,強度增加。不同湍流模型的計算結果均表明,從次噴口進入主噴口的流量與NPR無關,只隨SPR的增加而增加。在不同NPR下,SPR=0.7時,四種湍流模型計算的次噴流流量約是主噴流的4.2%,試驗是4%,SPR=1.0時,次噴流流量約是主噴流的5.7%,試驗是6%。

圖8 不同SPRSST k-ω模型計算的上壁面壓力分布(NPR=4.6)Fig.8 Pressure distribution of the primary nozzle upper wall with SST k-ω model,NPR=4.6,SPR=0.7,1.0
此外,四種湍流模型在SPR=0.7時,數值計算給出的推力矢量偏角隨NPR的變化與試驗值的比較如圖9所示,各湍流模型的計算結果與試驗值吻合良好。在NPR大于4之后,隨著NPR的增加,推力矢量偏角變小,與試驗得出的結論一致。

圖9 推力矢量偏角隨NPR變化與試驗值的比較(SPR=0.7)Fig.9 SPR=0.7,Computational thrust vector angle and experiments varying with NPR
雙縫引射的噴管,在計算過程中取外流馬赫數為0.1,正如我們在前文中所指出的,當馬赫數不大于0.1時,外流馬赫數對噴管性能的影響很小,可以忽略不計。由于試驗中兩個次流注入縫是從同一氣源引氣,因此在數值計算過程中兩個次流噴口的入口給定總溫、總壓邊界條件設置相同(即SPR相同),其它邊界條件的設定與單縫引射時完全一樣。
NPR=4.6,SPR=0.7,不同湍流模型給出的噴管及尾跡流場如圖10所示,雙縫引射會在主噴管中誘導產生兩道斜激波,第一道斜激波位于噴口2(圖2)之前,第二道斜激波在兩個噴口之間,兩個噴口左右兩側都會出現分離區,但噴口1兩側的分離區明顯大于噴口2,各湍流模型都捕捉到了試驗所觀測到的激波結構(圖11)。四種模型計算給出的主噴管內第一道斜激波的位置都比試驗值靠后,約x/c=1.4處,而試驗測量則是在x/c=1.24處,第二道激波位置與試驗相比,k-ζ模型計算給出的激波位置偏靠前,SSTk-ω和EASMk-ω模型與試驗值頗為吻合,SA模型給出的激波位置偏靠后,下壁面出口附近的斜激波試驗測量在x/c=1.89處,三種兩方程湍流模型的計算結果比較接近,都約在x/c=1.9處,SA最靠后,約在x/c=1.94處。
雙縫引射的推力矢量噴管,SPR=0.7時,在一系列NPR下,不同湍流模型計算的主噴管內壁面壓力分布與試驗值的比較如圖12(a~d)所示,總的來看均與試驗值吻合良好。雙縫引射時第一道斜激波的壓力曲線比試驗值更陡,激波位置比試驗值偏靠后。第二道斜激波處的壓力分布比試驗值要低,但對激波位置的捕捉要好于第一道斜激波。總體來說,kζ模型計算給出的激波位置最為靠前,SSTk-ω模型次之,之后是EASMk-ω模型,SA模型最為靠后。

圖10 不同模型計算的馬赫云圖(NPR=4.6,SPR=0.7)Fig.10 Mach counters for two injection ports nozzle at NPR=4.6,SPR=0.7,with different turbulence models

圖11 NPR=4.6,SPR=0.7時試驗陰影照片Fig.11 Experimental shadowgraph for two injection ports,NPR=4.6,SPR=0.7

圖12 雙縫引射噴管上下壁面壓力分布與試驗值對比Fig.12 Pressure distribution of the primary nozzle wall with two injection ports at different NPR
雙縫引射時,不同湍流模型計算的流量系數,推力系數隨NPR、SPR的變化規律與單縫引射相同,在此不再贅述。
圖13是SPR=0.7,不同湍流模型計算的雙縫引射的推力矢量偏角隨NPR的增加與試驗值的比較,從圖中可以看出,各湍流模型都是在NPR=4.6時給出了最大的推力矢量偏角,然后隨NPR增加,矢量偏角減小。在NPR=7.0時,SSTk-ω和k-ζ模型計算的推力矢量偏角稍偏大,主要是由于這兩種湍流模型過高估計了噴口1后回流分離區的壓力,如圖12(b)所示,可能是由于此狀態下,回流分離區的非定常效應較強,采用定常方法計算,這兩種模型無法給出準確的結果。圖14是NPR=4.6,不同SPR下,SSTk-ω模型給出的噴管上壁面壓力分布,從圖中可以看出,與單縫引射類似,隨SPR的增加,兩道斜激波的強度增加,位置前移。

圖13 推力矢量偏角隨NPR的變化(SPR=0.7)Fig.13 SPR=0.7,Computational thrust vector angel and experiments varying with NPR

圖14 不同SPR下SST k-ω模型計算的上壁面壓力分布Fig.14 Pressure distribution of the primary nozzle upper wall with two injection ports,SST k-ω model
當SPR=0.7時,SA湍流模型計算的單/雙縫引射噴管性能參數對比如圖15所示。雙縫引射時,盡管兩個次流注入縫總的開口面積與單縫引射噴管相同,SPR也相同,仍然會導致噴管的流量系數降低約1%,在相同的NPR下,推力系數比單縫引射增加了大約1%,利用其它湍流模型的計算結果也可以得出同樣的結論。
盡管兩個次流注入縫的入口及主噴管入口邊界條件設置完全相同(NPR、SPR相同),但由于次流注入縫在主噴管的位置不同,進入主噴管的流量是不同的,靠近主噴管喉道的噴口(噴口2)要比靠近主噴管出口的噴口(噴口1)的流量稍小一些,進入主噴管的總流量也比單縫引射噴管少,減小了約0.7%。
NPR=4.6,SPR=0.7時,不同湍流模型給出了在6.4°~6.8°范圍內的矢量偏角(試驗是6.5°),小于單次流噴口時不同湍流模型給出的介于7.3°~7.6°之間的推力矢量偏角(試驗是7.0°),但由于雙縫引射時進入主噴管的流量較小,每1%次流注入產生的推力矢量偏角即矢量效率反而是雙縫引射大于單縫引射的噴管。在其它NPR下,雙縫引射產生的推力矢量效率卻比單縫引射小,如圖15(d)所示。

圖15 SA模型計算的單/雙縫引射噴管性能參數比較(SPR=0.7)Fig.15 Computational nozzle performance with SA turbulence model(SPR=0.7)
分別采用四種不同湍流模型,通過對激波控制單/雙縫引射推力矢量噴管的數值模擬研究分析,可以得出以下結論:
(1)四種湍流模型都能比較準確的預測出單/雙縫射流推力矢量噴管的性能參數,并且計算結果差別不大,與試驗得出的規律一致。
(2)各湍流模型計算給出的單/雙縫主噴管內壁面壓力分布均與試驗值吻合良好。在激波位置捕捉方面,SSTk-ω、k-ζ和 EASMk-ω模型較好,SA 模型稍差。針對雙縫引射噴管,在回流分離區內,SSTkω、k-ζ模型的計算結果并不理想,會過高地估計壓力分布。
(3)無論是單縫還是雙縫引射,各湍流模型的計算結果均表明,相同NPR,如果SPR增加,進入主噴管的次流流量增加,次流誘導的斜激波前移,主噴管流量系數損失增加,推力系數損失也增加。
(4)雙縫引射的推力矢量效率并不一定優于單縫引射。如文中噴管構型,NPR=4.6時,雙縫引射的矢量效率較高,NPR大于4.6時,單縫引射的推力矢量效率較高。
[1]DEERE K A.Summary of fluidic thrust vectoring research conducted at NASA langley research center[R].AIAA-2003-3800,2003.
[2]KENRIC A W,DEERE K A.Experimental and computational investigation of multiple injection ports in a convergent-divergent nozzle for fluidic thrust vectoring[R].AIAA-2003-3802,2003.
[3]ZHANG Q,YANG Y。On Exploring predicting precision of an EASMk-ωtwo-equation turbulence model[J].JournalofNorthwesternPolytechnical University,2005,23(1):89-92
[4]ZHANG Q,YANG Y.Applying robinson's new twoequationturbulencek-ζmodel to numerically simulating turbulent flow[J].JournalofNorthwesternPolytechnicalUniversity,2009,27(4):466-470.
[5]ZHANG Q,YANG Y.Numerical analysis of low reynolds number flow over airfoil[J].ActaAerodynamic Sinica,2006,24(4):482-486.