波瓣混合器的大渦模擬

許錦錦，邊一帆，趙馬杰，葉桃紅

（中國科學技術大學熱科學和能源工程系，安徽合肥230027）

波瓣混合器是一種增強同向混合流動的裝置，在航空發(fā)動機的排氣系統中，中心射流與次流的混合效率對航空發(fā)動機的噴射噪聲有重要影響，因而波瓣混合器最初用于降低排氣噪聲（見圖1）[1-4]，隨后研究發(fā)現波瓣混合器能夠增加發(fā)動機推力[5]、減少燃油消耗[6-7]，降低紅外輻射[8]，此外波瓣混合器也應用于增強燃料與空氣在燃燒室的混合作用，用以提高燃燒效率和減少污染物的形成[9]，因而波瓣混合器在航空發(fā)動機中有著廣泛的應用。所以開展對波瓣混合器增強混合機制的研究具有重要意義。

圖1 波瓣混合器

許多研究人員利用實驗對波瓣混合器的流場結構進行了捕捉，用以分析其增強混合作用的機制。Paterson[10]最早利用實驗的方法對波瓣混合器的混合機制進行了研究，研究發(fā)現波瓣混合器內存在大尺度的流向渦、正交渦以及小尺度的渦流結構，但對其中的混合過程并不十分清楚。Werle[11]等通過對二維波瓣結構的實驗測量分析，揭示了流向渦的形成過程，圖2所示為波瓣混合器內存在的幾種渦流結構。McCormick[12]等通過實驗的方法指出由于K-H不穩(wěn)定性形成的正交渦對加強混合具有同等重要的作用，同時發(fā)現正交渦和流向渦的相互作用產生了更小尺度的渦流結構。Glauser[13]等利用熱線測速的方法對二維波瓣混合器進行了研究，研究發(fā)現流向渦導致了正交渦的變形，最終致使其破碎，從而形成了小尺度湍流結構。Hu[14]751等利用實驗分析了正交渦和流向渦隨著流向距離的衰減情況，發(fā)現它們具有同樣的衰減趨勢，進一步證明了正交渦與流向渦之間的相互作用。Belovich[15]總結了波瓣混合器的混合過程主要受到三個因素控制：①由于波瓣混合器的幾何外形而誘導產生的流向渦；②自由射流的剪切層中K-H不穩(wěn)定性而形成的正交渦；③由于波瓣混合器特殊的幾何外形結構增大了射流與次流的接觸面積，加強了剪切作用。但是，實驗方法對波瓣混合器流場的測量受到諸多因素的限制，對波瓣混合器的研究具有局限性，因而開展對波瓣混合器的數值模擬具有重要作用。

圖2 渦流結構示意圖

目前對波瓣混合器的數值模擬主要采用雷諾平均方法。Koutmos[16]利用k-ε雙方程模型模擬計算，并與實驗測量結果進行對比，發(fā)現可以很好地再現實驗中的結果。Cooper[17]等則采用不同的湍流模型模擬了波瓣混合器內的流場結構，驗證了不同的湍流模型對實驗結果的吻合情況，研究發(fā)現k-εRealizable湍流模型相對于其他雙方程模型，對波瓣混合器內最大流向渦和正交渦的值吻合的最好。而更多的數值模擬研究集中在優(yōu)化波瓣混合器的幾何外形，進而提高混合效率。O’Sullivan[18]和Du[19]研究了不同波瓣張角對混合效率的影響，Wright[20]研究了不同波瓣數目以及波瓣面積對波瓣混合器增強同向混合效率的影響，Lei[21]研究發(fā)現流向渦的破碎過程對于波瓣混合器提高混合效率具有重要作用，因此設計出新的波瓣混合器幾何結構用以增強混合效率。

但是，基于雷諾平均方法對波瓣混合器的數值模擬研究只能得到混合器內大尺度的流場特性，而缺少了小尺度的結果和波瓣混合器內的瞬時流場特性，盡管Ooba[22]等采用大渦模擬方法（LES）對波瓣混合器的流場進行了數值模擬，但也只給出了數值的平均結果，對波瓣混合器內小尺度的湍流結構變化缺乏研究，因此波瓣混合器增強同向混合作用的混合機制仍然有待完善。本文采用大渦數值模擬計算，研究了波瓣混合器流場中渦流結構，分析了不同渦流結構的發(fā)展過程。

1 數學物理模型及數值方法

本研究中波瓣混合器的具體尺寸參照Hu[14]759等人的實驗數據，在建立三維模型時，圓管部分的特征直接通過對圓環(huán)進行拉伸得到，而波瓣處的外形特征，由于其復雜的波瓣外形則需要首先通過二維尺寸添加波瓣出口處的外張角和內張角的輔助線，最后通過各個輔助線的位置創(chuàng)建外形曲面，如圖3所示，其主要參數如下：D為波瓣混合器內徑40mm；H為波瓣高度15mm；d為波瓣寬度6mm；θout-in為波瓣外壁內張角22°；θout-out為波瓣外壁外張角28°；θin-in為波瓣內壁內張角17°；θin-out為波瓣內壁外張角9°。

圖3 波瓣混合器三維結構

在計算域的選取方面，對于波瓣混合器的內部計算域，直接采用幾何抽取的方法，抽取計算域如圖4所示，需要注意的是計算域只包含管道內部流動空間，對于管壁厚度不予考慮。對于波瓣混合器的外部計算域，為了消除邊界對于數值計算結果的影響，采用10D的圓柱狀，其流向距離為12.6D，最終得到的波瓣混合器的計算域如圖4所示。

圖4 波瓣混合器計算域

在網格劃分方面，考慮邊界層的影響，在波瓣混合器壁面處采用邊界層網格，第一層網格高度為0.01，網格層數為15層。由于波瓣混合器出口位置射流與周圍流體之間存在很強的剪切作用，因而對波瓣混合器出口截面位置進行局部加密，加密長度為波瓣混合器出口截面沿下游方向6D，最終劃分的網格數量為38 867 454，如圖5所示。

圖5 非結構化網格劃分

初始和邊界條件方面，本文根據Hu[14]760等人實驗測得波瓣混合器出口處的速度U0=20 m/s，通過波瓣混合器幾何外形上出口與入口面積之比為1.36，流體的入口速度為14.5 m/s。流體介質為空氣，密度為1.225 kg/m3，動力黏度為1.789 4 Pa·s。波瓣混合器的外部流場區(qū)域采用壓力出口邊界條件，靜壓為1個大氣壓；波瓣混合器壁面采用無滑移邊界條件。

數值求解的方程包括不可壓縮的連續(xù)方程和動量方程。在大渦模擬中，采用密度加權過濾（Favre平均）后得到的控制方程為

式中：“-”代表網格尺度過濾后的物理量，“～”表示為Favre過濾；上標sgs表示亞格子量（subgrid scale）和分別為過濾后的速度分量和壓力，v為運動黏度，指小尺度脈動與可解尺度之間的動量輸運，本文采用動態(tài)Smagorinsky模型對亞格子應力進行封閉，表達式為

式中：δij為克羅內克算子為應變率張量是亞格子渦黏系數為亞格子湍動能；CD和C1是通過局部動態(tài)模型得到的模型系數；Δ為過濾尺度。

本文基于流體力學計算軟件OpenFOAM開展波瓣混合器的大渦模擬研究，物理模型的流體進口為常溫常壓的空氣，馬赫數較小，忽略其可壓縮性，求解質量守恒方程和動量守恒方程，采用動態(tài)Smagorinsky渦粘模型對亞格子應力進行模化。采用Piso算法求解壓力速度的耦合，采用有限體積法對控制方程進行離散，對流項和擴散項使用二階中心差分格式，而時間導數項使用二階隱式差分。

2 結果討論

定義沿X、Y、Z方向上各截面上的無量綱渦量ωx、ωy、ωz分別為

式中：U0為波瓣混合器出口中心線速度，U0=20 m/s，u，v，w均為瞬時速度分量。

2.1 流向渦分布

流向渦的產生對波瓣混合器提高混合作用具有關鍵作用，流向渦的形成是由于波瓣混合器的特殊波瓣結構。

圖6為不同截面上波峰處ωz分布，在波瓣混合器的圓管與波瓣噴管的連接截面上，部分射流流體開始進入到波瓣噴管當中，進入到波瓣處的流體分別向兩側流動，形成了兩個旋渦，隨著流動距離的增加，由于波瓣張角的存在，使得進入到波瓣內的流體越來越多，從而導致了流向渦越來越大并沿著波瓣截面分布。而后隨著距離波瓣混合器出口截面距離進一步減少，流向渦進一步增大，并始終分布在波瓣兩側?？梢灶A見波瓣張角、高度以及寬度直接影響著進入到波瓣內的流體質量，從而影響著流向渦的大小以及其分布形狀。隨著流向距離的進一步發(fā)展，當射流氣體通過混合器截面以后，在波瓣混合器的出口截面，盡管外部流體的速度為0，但由于流向渦對周圍流體的卷吸作用使得在波瓣混合器的外壁面存在較小的速度，并與波瓣混合器外壁面相互作用，存在較小的渦量。此時，內側涵道產生沿徑向方向向外的二次流動，外側流體則有向內的二次流，內外兩側涵道誘導產生二次流發(fā)展成為流向渦，最終形成了一對反向旋轉的渦對。同時，由于在波瓣壁面邊界層的存在，在流向渦內靠近波瓣壁面處形成了狹長的沿著波瓣壁面分布與流向渦方向相反的渦量。值得注意的是，在流向渦的形成過程中，在波瓣內部徑向向外的壓力梯度的作用下，邊界層沿波瓣表面向波峰方向卷起，形成了一對反向旋轉的渦對，稱為馬蹄渦，這部分渦量并沒有像流向渦一樣隨著距離的增加逐漸增大，而是始終存在于波瓣頂部，并在波瓣出口截面以后迅速消散，而由于波瓣壁面產生的這部分渦量也在下游位置迅速消散。這兩部分渦量的值以及分布相對流向渦較小，且耗散速率較快，因而對提高混合效率的影響較少。

在不同截面處大渦數值模擬得到的流向渦與實驗結果的對比如圖7所示。從圖7中可知，大渦數值模擬很好地反映了流向渦沿流向方向的變化。在波瓣混合器的出口位置，內側涵道產生沿徑向方向向外的二次流動，外側流體則有向內的二次流，內外兩側流體相互作用，使得流向渦發(fā)展成為一對反向旋轉的渦對，以波瓣波峰線為對稱軸，兩側的渦量幾何相等，但符號相反。相對于實驗的測量值，大渦數值模擬得到的流向渦更為狹長，分布在波瓣混合器的波瓣兩側。隨著流動距離的增加，速度不斷下降，流向渦的值并不是變得更小，而是整個渦結構開始破碎，由最初的一對反向旋轉的渦對，逐漸分解，形成多個尺度更小的渦結構，在Z/D=1和1.5截面上，大渦數值模擬以及實驗測量都再現了流向渦的破碎過程。隨著流動距離的進一步增大，可以看到在Z/D=2和3截面上，已經破碎為小尺度結構的流向渦開始逐漸耗散，直至消亡。因而流向渦在波瓣混合器流場中的過程可以分為產生、發(fā)展、破碎以及消亡四個階段。相比于Hu[14]761等人基于雷諾平均方法的計算結果，大渦數值模擬更好地反映了流向渦變化，尤其在流向渦破碎為小尺度的渦流階段。

圖6 不同截面上波峰處ωz分布

2.2 正交渦分布

圖7 不同截面上流向渦分布

在內外兩側涵道中，射流與外部流體存在較大的速度差，從而形成較強的剪切層，導致在圍繞波瓣混合器的出口截面上形成了正交渦量。圖8顯示了正交渦的瞬時特征，大渦數值模擬較好地反映了正交渦的分布。相比于流向渦，正交渦的形成是由于剪切層的作用，因而不同于流向渦，正交渦形成在波瓣混合器的出口截面位置，其形狀特征與波瓣混合器截面的幾何結構一致。同流向渦一樣，正交渦也存在著產生、發(fā)展、破碎以及消亡四個階段。在Z/D=0.25截面上，正交渦連接成環(huán)狀結構，而在波峰與波谷的連接處，正交渦出現了扭曲變形，而隨著流向距離的增加，在Z/D=1處截面上，正交渦開始破碎，破碎的位置發(fā)生在波峰與波谷相連的扭曲變形處，破碎的正交渦分成兩部分，一部分在沿著波峰分布呈現月牙狀，另一部分沿著波谷分布形成花冠結構，隨著流向距離的進一步增大，這兩部分逐漸破碎，正交渦也不再連接成環(huán)狀，各自分別形成更小尺度的渦量，最終耗散消亡。采用RANS的計算方法盡管捕捉到了正交渦最初的破碎過程，但是并沒有反映出之后的衰減情況。值得注意的是，Hu[14]762等人對不同截面上的最大流向渦以及正交渦進行了研究，發(fā)現它們具有相同的衰減規(guī)律，這表明兩種渦結構存在著相互作用。

圖8 不同截面上正交渦分布

2.3 速度分布

圖9（a）為沿波瓣混合器波峰方向橫斷截面的速度分布。由于波瓣混合器的出口截面相對于圓截面的面積減少，因而在波瓣噴管處，流體速度將增加。在出口位置的流向方向，速度呈現“三叉”結構，波峰出口處速度與中心射流速度形成三個明顯的“尖叉”，這使得剪切層不僅存在于外部流體與射流流體之間，在中心射流與波瓣波峰流體之間的區(qū)域也存在剪切層，增強了混合作用。隨著流向距離的增加，速度逐漸降低，尤其在波峰處的流體速度降低的最為明顯，這是由于這部分流體的速度梯度較大，受到剪切層的影響較大，耗散效應更為明顯。在Z/D=2的位置處，從波峰處射流出的流體速度已經基本耗散，這與流向渦和正交渦的消亡階段保持一致。而在沿波瓣混合器波谷方向橫斷截面上如圖9（b）所示，此時波瓣混合器的外形是向內凹進，因而不會出現“三叉”結構，中心射流僅受到外部流體的剪切作用。

圖9 速度分布（t=0.467 009 56）

3 結論

本文基于大渦模擬方法（LES）開展了波瓣混合器渦結構的數值模擬，在波瓣混合器流場中主要存在大尺度的流向渦、正交渦以及小尺度的馬蹄渦，通過對比PIV實驗結果，大渦數值模擬很好地再現了波瓣混合器內不同渦結構的發(fā)展過程。流向渦產生于波瓣混合器的波瓣處，流體進入到波瓣內并產生旋渦，隨著進入到波瓣內的流體越來越多，旋渦逐漸增大并形成流向渦，在波瓣混合器出口處，射流流體卷吸次流流體，使得流向渦進一步發(fā)展增大，流向渦隨著流向距離的增加，破碎成更小尺度的渦量而并不是逐漸衰減變弱。相比于流向渦，正交渦的形成是由于射流流體與次流流體的剪切作用，因而正交渦形成與波瓣混合器出口截面，并呈現波瓣狀，同時在波瓣混合器流場中存在多個剪切作用，一是中心射流與周圍流速較低的流體之間的剪切作用，二是次流流體與射流之間的剪切作用。