黨進(jìn)鋒，郜 冶，劉加寧

(哈爾濱工程大學(xué) 航天與建筑工程學(xué)院，哈爾濱 150001)

?

雙向渦數(shù)值模擬雷諾應(yīng)力各向異性研究與分析

黨進(jìn)鋒，郜 冶，劉加寧

(哈爾濱工程大學(xué) 航天與建筑工程學(xué)院，哈爾濱 150001)

以實(shí)驗(yàn)數(shù)據(jù)和理論結(jié)果作為參考，探討基于布茨涅克各向同性假設(shè)的RNGk-ε模型、SST模型和各向異性RSM湍流模型在雙向渦流數(shù)值計(jì)算中的適用性，并對(duì)不同入射速度下的雙向渦流結(jié)構(gòu)進(jìn)行了比較分析。計(jì)算結(jié)果表明，基于布茨涅克各向同性假設(shè)的RNG模型與SST模型在不同入射速度下均無(wú)法得到內(nèi)部受迫渦、外部自由渦的雙向渦流結(jié)構(gòu)，而RSM模型則可以準(zhǔn)確得到雙向渦流結(jié)構(gòu)。內(nèi)、外渦分界面處雷諾正應(yīng)力各向異性特征顯著，且隨入射速度的提高而增強(qiáng)，最大差值達(dá)到了198%，軸線附近流場(chǎng)雷諾正應(yīng)力各向異性特征相對(duì)較弱，但是最大差值依然達(dá)到46%，因此基于布茨涅克各向同性假設(shè)的湍流模型不適用于渦流冷壁發(fā)動(dòng)機(jī)數(shù)值模擬研究。

雙向渦流；渦流冷壁發(fā)動(dòng)機(jī)；雷諾應(yīng)力；各向異性

0 引言

渦流是自然界中普遍存在的一種流動(dòng)形式，大到天體星云的運(yùn)動(dòng)，小到湍流中的粘性耗散渦，均是不同形式的渦流。由于渦流具有豐富多樣的結(jié)構(gòu)形式，近年來(lái)被越來(lái)越多的應(yīng)用在流體機(jī)械中，如渦流燃燒器、渦流固液組合發(fā)動(dòng)機(jī)、渦流分離器等。2002年，美國(guó)Orbital科技公司Chiaverini等提出了一種新概念液體火箭發(fā)動(dòng)機(jī)形式[1]，氧化劑從發(fā)動(dòng)機(jī)噴管上方沿切向噴入發(fā)動(dòng)機(jī)內(nèi)形成外部渦流，外部渦流向發(fā)動(dòng)機(jī)頂端旋轉(zhuǎn)流動(dòng)，在發(fā)動(dòng)機(jī)頂端改變流動(dòng)方向，并與頂端噴入的燃料混合形成內(nèi)部燃燒渦流，燃燒的內(nèi)部渦流從噴管流出產(chǎn)生推力，該發(fā)動(dòng)機(jī)被稱為渦流冷壁發(fā)動(dòng)機(jī)(Vortex Combustion Cold-Wall Chamber VCCWC)。VCCWC的優(yōu)點(diǎn)在于氧化劑所形成的外部渦流對(duì)于發(fā)動(dòng)機(jī)壁面起到了冷卻作用，大大降低了發(fā)動(dòng)機(jī)內(nèi)壁熱防護(hù)層厚度，核心流內(nèi)氧化劑與燃料旋轉(zhuǎn)燃燒，增加了反應(yīng)時(shí)間從而縮短了反應(yīng)距離，最終可減小發(fā)動(dòng)機(jī)長(zhǎng)度，該兩方面優(yōu)點(diǎn)均有助于降低發(fā)動(dòng)機(jī)成本，因此美國(guó)NASA將VCCWC確定為未來(lái)廉價(jià)發(fā)動(dòng)機(jī)的候選方案[2-3]。由于VCCWC內(nèi)部渦與外部渦流動(dòng)方向相反，因此該組合渦稱為雙向渦流(Bidirectional Vortex)。

從VCCWC中雙向渦流問(wèn)題的提出至今十余年時(shí)間里，有限研究工作主要分為兩部分，一部分為Chiaverini等所做VCCWC冷流與燃燒實(shí)驗(yàn)[4]，實(shí)驗(yàn)驗(yàn)證了雙向渦流的存在并得到內(nèi)渦旋轉(zhuǎn)燃燒圖像。另一部分為Majdalani等開展的雙向渦流理論研究，文獻(xiàn)[5]中推導(dǎo)得到不可壓無(wú)粘流雙向渦流結(jié)構(gòu)解析式，文獻(xiàn)[6]中推導(dǎo)得到不可壓粘性核心雙向渦流結(jié)構(gòu)解析式，文獻(xiàn)[7]得到定常剪切應(yīng)力核心流模型，文獻(xiàn)[8-9]中分別通過(guò)Bragg-Hawthorne流函數(shù)法和Beltramian流場(chǎng)近似法得到可壓縮雙向渦流模型。

與之相比，雙向渦流CFD研究進(jìn)展較慢。由于雙向渦流三維效應(yīng)顯著，且涉及到VCCWC化學(xué)反應(yīng)數(shù)值計(jì)算，計(jì)算量巨大，因此LES、DNS等方法的數(shù)值研究還未開展。Fang在文獻(xiàn)[10]中結(jié)合Standardk-ε模型和RSM模型數(shù)值計(jì)算模擬了VCCWC冷流與燃燒過(guò)程，從其冷流計(jì)算結(jié)果來(lái)看速度曲線與理論曲線趨勢(shì)相似，能夠得到內(nèi)、外雙向渦流結(jié)構(gòu)。文獻(xiàn)[11]中通過(guò)數(shù)值計(jì)算指出外部渦流厚度約占發(fā)動(dòng)機(jī)半徑的25%～29%，并提出冷卻劑噴嘴傾斜一定的角度可增強(qiáng)內(nèi)部渦流強(qiáng)度、改善推進(jìn)劑摻混程度提高燃燒效率。文獻(xiàn)[12]中數(shù)值計(jì)算研究了以GH2/GO2為燃料的VCCWC燃燒過(guò)程。

文獻(xiàn)[13]中通過(guò)二維物理模型比較了RNG 模型和RSM模型在VCCWC數(shù)值模擬中的準(zhǔn)確性，指出RSM模型更適合于強(qiáng)旋流場(chǎng)。本文在文獻(xiàn)[13]工作基礎(chǔ)上，對(duì)VCCWC內(nèi)流場(chǎng)進(jìn)行了三維完整物理模型、3種不同湍流模型、3種不同入射速度下的數(shù)值計(jì)算研究，定量探究RSM模型更適合VCCWC數(shù)值計(jì)算的內(nèi)在原因。計(jì)算發(fā)現(xiàn)雙向渦流結(jié)構(gòu)對(duì)于湍流模型的選擇性非常大，不合適的湍流模型計(jì)算所得結(jié)果與實(shí)驗(yàn)結(jié)果差別非常明顯，甚至可以說(shuō)是完全錯(cuò)誤的，分析原因認(rèn)為：雙向渦流與普通流動(dòng)的區(qū)別在于，雙向渦流在遠(yuǎn)離壁面的核心流區(qū)域內(nèi)湍流脈動(dòng)雷諾正應(yīng)力各向異性特性依然異常顯著，而基于布茨涅克各向同性假設(shè)(Boussinesq Isotropic Hypothesis)的湍流模型則認(rèn)為核心流區(qū)域各向雷諾正應(yīng)力均相等[14]，這就導(dǎo)致了傳統(tǒng)基于布茨涅克各向同性假設(shè)的湍流模型無(wú)法得到正確結(jié)果。

本文以文獻(xiàn)[15]中實(shí)驗(yàn)數(shù)據(jù)及文獻(xiàn)[5，16]中理論結(jié)果作為參考，探討基于布茨涅克各向同性假設(shè)的RNGk-ε模型、SST模型和各向異性RSM湍流模型在雙向渦流數(shù)值計(jì)算中的適用性，并對(duì)不同入射速度下的雙向渦流結(jié)構(gòu)進(jìn)行了比較分析，結(jié)論可用于指導(dǎo)雙向渦流數(shù)值模擬湍流模型的選擇。

1 物理模型及數(shù)值計(jì)算方法

圖1(a)、(b)分別為雙向渦流物理模型示意圖及入射噴口分配示意圖。物理模型直徑d=100 mm，長(zhǎng)度L=100 mm，即長(zhǎng)徑比α=1，模型頂部為平面。入射噴口橫截面為5 mm×4 mm矩形，模型底部12個(gè)入射噴口周向均勻排列，采用切向速度入口形式。本文中物理模型與文獻(xiàn)[10，13]中類似，因此采用與文獻(xiàn)[10，13]相同的Outflow出口邊界條件，出口直徑為30 mm，其余均為無(wú)滑移壁面。

(a)VCCWC物理模型 (b)入射噴口分配示意圖

物理模型網(wǎng)格示意圖如圖2所示。計(jì)算中，3種湍流模型(RNGk-ε模型、SST模型及RSM模型)均采用SIMPLE算法求解，梯度差分采用G-S節(jié)點(diǎn)基格式，壓力差分采用二階格式，其余通量差分均為二階迎風(fēng)格式。Majdalani在文獻(xiàn)[5]中指出雙向渦流在不可壓縮與可壓縮流動(dòng)中均存在，因此計(jì)算中采用低速不可壓縮流形式，流體介質(zhì)為理想不可壓空氣，3種入口速度分別為10、15、20 m/s。計(jì)算工作在64位服務(wù)器上完成，耗時(shí)320 h，所使用的CFD軟件為ANSYS FLUENT11.0。

圖2 物理模型網(wǎng)格示意圖Fig.2 Geometric mesh of physical model

2 計(jì)算可信性驗(yàn)證

圖3為入射速度vin=20 m/s時(shí)L=50 mm處無(wú)量綱壓力差徑向分布曲線，3條曲線分別為RNG、SST以及RSM模型計(jì)算所得，并與文獻(xiàn)[15]中實(shí)驗(yàn)數(shù)據(jù)進(jìn)行了比較，文獻(xiàn)[15]中只給出了0.00

圖4所示為本文中RSM模型計(jì)算所得速度矢量與Majdalani在文獻(xiàn)[5]中所得理論速度矢量比較。

圖3 L=50 mm無(wú)量綱壓力徑向分布與文獻(xiàn)[15]實(shí)驗(yàn)數(shù)據(jù)Fig.3 Non-dimensional pressure difference distribution inradial direction at L=50 mm and experimentalresults in Reference[15]

可看到在不考慮出口對(duì)流場(chǎng)的影響時(shí)，CFD結(jié)果與理論結(jié)果非常一致，Majdalani在文獻(xiàn)[5]推導(dǎo)過(guò)程中忽略了出口直徑對(duì)于流場(chǎng)的影響，假設(shè)出口直徑與內(nèi)渦半徑相等，因此內(nèi)、外渦流動(dòng)方向相反，不存在回流區(qū)。當(dāng)出口直徑小于內(nèi)渦半徑時(shí)就會(huì)存在如本文中數(shù)值計(jì)算中所示的回流區(qū)。

綜合以上壓力與速度矢量圖比較可知，本文中RSM計(jì)算結(jié)果與實(shí)驗(yàn)數(shù)據(jù)及理論結(jié)果非常符合，說(shuō)明本文中所采用物理模型、計(jì)算網(wǎng)格、邊界條件及數(shù)值求解方法具有足夠的可信度。

3 計(jì)算結(jié)果及討論

3.1 不同湍流模型速度及雙向渦流線分析

圖5為不同入射速度下L=50 mm處切向速度vt徑向分布曲線(圖5(a)～(c)所示)。

vin分別為10、15、20 m/s，以及文獻(xiàn)[16]雙向渦vt理論曲線圖如5(d)所示，理論曲線使用模型半徑和入口速度對(duì)橫、縱軸進(jìn)行了無(wú)量綱化。RSM模型所得速度曲線：vin=10 m/s時(shí)vt，max位于r≈±0.045 m處，vin=15 m/s時(shí)vt，max向軸線位置處移動(dòng)位于r≈±0.017 m處，vin=20 m/s時(shí)vt，max位于r≈±0.015 m處。同時(shí)vt，max隨vin的提高而增大，說(shuō)明內(nèi)部受迫渦強(qiáng)度隨vin提高而增大。RNG模型和SST模型在不同入射速度情況下vt，max始終位于r≈±0.045 m處。

圖6所示為vin=20 m/s時(shí)RNGk-ε模型、SST模型和RSM模型所得雙向渦流結(jié)構(gòu)流線圖，在物理模型、入射速度及其他邊界條件均相同的情況下，3種湍流模型所得流線圖差別非常顯著。RNG模型和SST模型也可得到雙向旋轉(zhuǎn)流動(dòng)結(jié)構(gòu)，結(jié)合圖3、圖5曲線可知此時(shí)并沒(méi)有形成內(nèi)部受迫渦，整個(gè)流場(chǎng)均為自由渦，所以雙向流動(dòng)并沒(méi)有延伸至整個(gè)內(nèi)部流場(chǎng)，旋流長(zhǎng)度僅為流場(chǎng)長(zhǎng)度的2/3，外部渦流半徑也僅為流場(chǎng)半徑的2/3。RSM模型所得流線則可以形成完整的雙向渦流結(jié)構(gòu)。

從引言中可知，基于雙向渦流的VCCWC優(yōu)點(diǎn)在于外部渦流能夠?qū)τ诎l(fā)動(dòng)機(jī)起到冷壁作用，所以外部渦流半徑應(yīng)接近于發(fā)動(dòng)機(jī)半徑，內(nèi)部渦流高速旋轉(zhuǎn)能夠縮短化學(xué)反應(yīng)距離。因此，vt，max應(yīng)位于內(nèi)部受迫渦范圍內(nèi)，這均與RSM計(jì)算結(jié)果相一致。

綜合圖3中壓力曲線、圖5中速度曲線及圖6流線圖可知，SST模型和RNG模型不能得到正確的雙向渦流場(chǎng)分布，RSM模型能夠準(zhǔn)確計(jì)算得到雙向渦流結(jié)構(gòu)。

(b)vin=15 m/s

(c)vin=20 m/s

(d)文獻(xiàn)[16]理論結(jié)果

3.2 RSM模型雙向渦結(jié)構(gòu)分析

圖7為不同入射速度RSM模型計(jì)算所得L=50 mm處軸向速度va徑向分布曲線。雙向渦流內(nèi)、外渦軸向速度方向相反，因此va曲線與0.00 m/s速度曲線相交點(diǎn)為內(nèi)、外渦分界面位置，軸線到該分界面距離即為渦幔半徑rm(Mantle Radius)。vin=10m/s時(shí)，rm≈0.042 m，vin=15m/s時(shí)，rm≈0.035 m，vin=20 m/s時(shí)rm與vin=15 m/s時(shí)基本相同。

(a)RNG模型

(b)SST模型

(c)RSM模型

圖8所示為vin=20 m/s時(shí)RSM模型所得不同軸向位置處va分布曲線，軸向位置分別為L(zhǎng)=25 mm、L=50 mm和L=75 mm。從圖8看到，渦幔半徑從頂端到出口處略有變化，但是基本保持在0.035 m左右，這與現(xiàn)有文獻(xiàn)結(jié)論非常吻合[5，11，17]。

圖7 RSM 模型計(jì)算所得L=50mm處軸向速度徑向分布Fig.7 Axial velocity distribution in radial direction atL=50mm using RSM model

圖8 vin=20 m/s時(shí)RSM模型計(jì)算所得不同位置軸向速度徑向分布Fig.8 Axial velocity distribution in radial direction in inlet velocity of vin=20 m/s using RSM model

3.3 雷諾應(yīng)力各向異性分析

從圖3～圖6看到，基于布茨涅克各向同性假設(shè)的RNG模型與SST模型均無(wú)法得到正確的雙向渦流結(jié)構(gòu)，而RSM模型則能夠準(zhǔn)確得到完整的雙向渦流結(jié)構(gòu)，壓力差分布曲線(圖3)與實(shí)驗(yàn)數(shù)值吻合的非常好，對(duì)稱面速度矢量分布(圖4)與理論結(jié)果符合較好。從國(guó)內(nèi)外學(xué)者關(guān)于雙向渦流數(shù)值計(jì)算文章中可看到，部分文章依然采用基于布茨涅克各向同性假設(shè)的湍流模型進(jìn)行數(shù)值計(jì)算研究。為了證明基于布茨涅克各向同性假設(shè)湍流模型的不足，圖9～圖11將RSM模型計(jì)算所得不同入射速度情況下各軸向位置處雷諾正應(yīng)力分布曲線進(jìn)行比較分析。

(a)vin=10 m/s (b)vin=15 m/s (c)vin=20 m/s

(a)vin=10 m/s (b)vin=15 m/s (c)vin=20 m/s

(a)vin=10 m/s (b)vin=15 m/s (c)vin=20 m/s

比較3.1節(jié)與本節(jié)雷諾正應(yīng)力分析可知，RNG模型與SST模型無(wú)法得到正確雙向渦流的原因在于雙向渦流結(jié)構(gòu)中雷諾正應(yīng)力各向異性特性顯著，內(nèi)、外渦分界面附近雷諾正應(yīng)力最大差值達(dá)到198%，且隨內(nèi)部受迫渦強(qiáng)度增大還會(huì)進(jìn)一步提高，此時(shí)布茨涅克各向同性假設(shè)已經(jīng)無(wú)法滿足雙向渦流物理特性，因此RNG模型與SST模型無(wú)法得到正確的雙向渦流，推論可知基于布茨涅克各向同性假設(shè)的其他湍流模型，如Standardk-ε模型，Standardk-ω模型等也不適用于雙向渦流數(shù)值計(jì)算。雷諾應(yīng)力各向異性決定了湍流脈動(dòng)的各向異性，化學(xué)反應(yīng)數(shù)值計(jì)算準(zhǔn)確性很大程度上決定于湍流場(chǎng)計(jì)算的準(zhǔn)確性，因此未來(lái)無(wú)論進(jìn)行VCCWC冷流或燃燒數(shù)值模擬時(shí)均應(yīng)摒棄基于布茨涅克各向同性假設(shè)的湍流模型。

4 結(jié)論

(1)基于布茨涅克各向同性假設(shè)的RNG模型與SST模型在入射速度從vin=10 m/s增大到vin=20 m/s時(shí)均無(wú)法準(zhǔn)確得到內(nèi)部受迫渦、外部自由渦的雙向渦流結(jié)構(gòu)，RSM模型能夠準(zhǔn)確得到與實(shí)驗(yàn)與理論結(jié)果相符的雙向渦流結(jié)構(gòu)。

(2)對(duì)于本文中所使用的物理模型尺寸，當(dāng)入射速度vin=10 m/s時(shí)無(wú)法形成穩(wěn)定的內(nèi)部受迫渦，當(dāng)vin≥15 m/s時(shí)形成穩(wěn)定的雙向渦結(jié)構(gòu)，在一定的速度變化范圍內(nèi)渦幔半徑保持在rm=0.035 m左右。隨入射速度提高內(nèi)部受迫渦強(qiáng)度增大，因此最大切向速度vt，max增大，其位置隨入射速度提高而靠近軸線位置。

(3)當(dāng)雙向渦流形成時(shí)，內(nèi)、外渦分界面附近雷諾正應(yīng)力各向異性特征顯著，且隨入射速度的提高而增強(qiáng)，最大差值達(dá)到了198%，軸線附近流場(chǎng)雷諾正應(yīng)力各向異性特征相對(duì)較弱，但是最大差值依然達(dá)到46%，因此基于布茨涅克各向同性假設(shè)的湍流模型不適用于VCCWC雙向渦流數(shù)值模擬研究。

[1] Chiaverini,M J,Malecki,M J,Sauer,J A,et al.Vortex combustion chamber development for future liquid rocket engine applications[R].AIAA2002-4149.

[2] Turner A E.Aquarius low cost launch main engine study[C]//3rd Responsive Space Conference April 25-28,2005 Los Angeles,CA.

[3] Munson S M,Sauer J A,Rocholl J D,et al.Development of a low-cost vortex-cooled thrust chamber using hybrid fabrication techniques[C]//47th AIAA/ASME/SA.E/ASEE Joint Propulsion Conference & Exhibit,31 July-03 August 2011,San Diego,California.

[4] Chiaverini M J,Malecki M J,Sauer A,et al.Vortex thrust chamber testing and analysis for O2-H2propulsion applications[R].AIAA 2003-4473.

[5] Vyas A B,Majdalani J.The bidirectional vortex.part 1:an exact inviscid solution[R].AIAA2003-5052.

[6] Vyas A B,Majdalani J.The bidirectional vortex.Part 2:viscous core corrections[R].AIAA2003-5053.

[7] Maicke B A,Majdalani J.A constant shear stress core flow model of the bidirectional vortex[J].Proceedings of the Royal Society,2009,465:915-935.

[8] Maicke B A,Majdalani J.On the compressible bidirectional vortex.Part 1:A bragg-hawthorne stream function formulation[C]//50th AIAA Aerospace Sciences Meeting including the New Horizons Forum and Aerospace Exposition,09-12 January 2012,Nashville,Tennessee.

[9] Maicke B A,Majdalani J.On the compressible bidirectional Vortex.Part 2:A beltramian flowfield approximation[C]//50th AIAA Aerospace Sciences Meeting including the New Horizons Forum and Aerospace Exposition,09-12 January 2012,Nashville,Tennessee.

[10] Fang D,Majdalani J.Simulation of the cold-wall swirl driven combustion chamber[R].AIAA 2003-5055.

[11] 唐飛,李家文,?？擞?渦流冷卻推力室中渦流結(jié)構(gòu)的分析與優(yōu)化[J].推進(jìn)技術(shù),2010,31(2):165-169.

[12] 吳東波,李家文,?？擞?GH2/GO2渦流冷卻推力室設(shè)計(jì)與數(shù)值計(jì)算[J].火箭推進(jìn),2010,36(5):17-22.

[13] 孫得川,白榮博,劉上.渦流燃燒發(fā)動(dòng)機(jī)燃燒室數(shù)值模擬[J].彈箭與制導(dǎo)學(xué)報(bào),2011,31(2):111-114.

[14] Schmitt F G.About boussinesq's turbulent viscosity hypothesis:historical remarks and a direct evaluation of its validity[J].Comptes,Rendus,Mecanique,2007,9-10,335:617-627.

[15] Maicke B A,Majdalani J.Characterization of the Bidirectional vortex using particle image velocimetry[D].PhD.Giovanna Cavazzini (Ed.),ISBN:978-953-51-0625-8,InTech.

[16] Majdalani J,Halpenny E K.The bidirectional vortex with sidewall injection[C]//44th AIAA/ASME/SAE/ASEE Joint Propulsion Conference & Exhibit,21-23 July 2008,Hartford,CT.

[17] Smith,J L.An analysis of the vortex flow in the cyclone separator[J].Journal of Basic Engineering Transactions of the ASME,1962:609-618.

(編輯：呂耀輝)

Research and analysis of reynolds stress anisotropic in bidirectional vortex simulation

DANG Jin-feng，GAO Ye，LIU Jia-ning

(College of Aerospace and Civil Engineering,Harbin Engineering University,Harbin 150001,China)

In this paper,simulation results were compared with experimental data and theoretical results.The simulation validity of RNGk-εmodel and SST model based on Boussinesq isotropic hypothesis and RSM model in bidirectional vortex were discussed and bidirectional vortex properties in different inlet velocities were analyzed as well.It is pointed out that the RNG model and SST model could not achieve the bidirectional vortex which frees in outer area and is forced in inner space under any inlet velocity based on Boussinesq isotropic hypothesis,while the RSM model simulation results fit the experiment data and theoretical solutions well.In mantle surface of inner vortex and outer vortex,Reynolds stress anisotropy characteristics is significant and is strengthened with the increase of inlet velocity.The biggest difference between normal stresses inx,y,zdirection can reach up to 198%.In core flow area close to axial,Reynolds stress anisotropy is relatively weak,but the greatest difference can still reach as large as 46%.Therefore,it is believed that turbulent models based on Boussinesq isotropic hypothesis are invalid in VCCWC simulation.

bidirectional vortex；vortex combustion cold-wall chamber；reynolds stress；anisotropic

2014-10-30；

：2014-12-11。

自然科學(xué)基金面上項(xiàng)目(11372079)。

黨進(jìn)鋒(1987—)，男，博士生，研究方向?yàn)榘l(fā)動(dòng)機(jī)內(nèi)流場(chǎng)數(shù)值仿真。E-mail：xiandjf@163.com

V430

A

1006-2793(2015)05-0640-06

10.7673/j.issn.1006-2793.2015.05.007