李逸翔 汪 球 羅 凱 李進(jìn)平 趙 偉
* (中國科學(xué)院力學(xué)研究所高溫氣體動(dòng)力學(xué)國家重點(diǎn)實(shí)驗(yàn)室,北京 100190)
? (中國科學(xué)院大學(xué)工程科學(xué)學(xué)院,北京 100049)
新一輪的高超聲速技術(shù)對(duì)飛行器在“極端”環(huán)境和“極端”動(dòng)力條件下的飛行提出了新的挑戰(zhàn),它也面臨降熱、減阻、控制等一系列的難題.當(dāng)飛行器以高超聲速飛行時(shí),頭部弓形會(huì)由于激波壓縮和黏性阻滯減速而產(chǎn)生高溫,其溫度高達(dá)幾千度甚至上萬度,高溫下氣體化學(xué)反應(yīng)產(chǎn)生的離子和電子在飛行器周圍形成一個(gè)等離子體鞘層.高超聲速M(fèi)HD流動(dòng)控制是利用飛行器內(nèi)部布置的磁場與電離氣體的相互作用來達(dá)到流動(dòng)控制的目的.它作為一種非接觸的主動(dòng)流動(dòng)控制方法,在飛行器氣動(dòng)力操控、熱環(huán)境管理、“電磁窗口”減弱“黑障”、“電磁舵面”飛行姿態(tài)控制、磁流體發(fā)電等領(lǐng)域均有廣闊的應(yīng)用前景,受到國內(nèi)外的廣泛關(guān)注[1-3].
實(shí)驗(yàn)研究、數(shù)值模擬和理論分析是研究科學(xué)問題的常用手段,三者相輔相成.高超聲速M(fèi)HD 流動(dòng)控制的研究對(duì)實(shí)驗(yàn)條件要求較高,需要高超聲速設(shè)備能夠模擬高速飛行狀態(tài)下的來流以及產(chǎn)生高磁場強(qiáng)度的磁場發(fā)生裝置;同時(shí),電磁干擾以及脈沖設(shè)備較短的實(shí)驗(yàn)時(shí)間也會(huì)帶來測(cè)量方面的困難[4].由于當(dāng)前實(shí)驗(yàn)可進(jìn)行測(cè)量的物理量較為有限,所以研究主要是針對(duì)激波脫體距離[5-8]、模型氣動(dòng)力(阻力/磁阻力等)[8-12]以及模型表面熱流[13]等方向來開展,但總體來講,隨著高超聲速地面設(shè)備模擬能力及電磁技術(shù)的發(fā)展,越來越多的學(xué)者考慮通過地面實(shí)驗(yàn)來研究磁場流動(dòng)控制的效果及規(guī)律.
數(shù)值模擬相較于實(shí)驗(yàn)研究而言,所需要的研究成本大為降低,通過數(shù)值模擬可以得到全流場信息,有助于直觀地了解流場規(guī)律.經(jīng)過近二三十年的發(fā)展,尤其是對(duì)計(jì)算格式的不斷修正[14-18],MHD 數(shù)值模擬發(fā)揮著越來越重要的作用.然而,MHD 涉及多種復(fù)雜物理現(xiàn)象,包括產(chǎn)生和維持帶電粒子的過程、帶電粒子與磁場的作用、非平衡能量輸運(yùn)機(jī)理以及電磁場中的霍爾效應(yīng)等,每一種效應(yīng)都包含復(fù)雜的規(guī)律,即使是全MHD 的數(shù)值模擬也是基于一些簡化假設(shè)(如電磁介質(zhì)各向同性、無磁極化、忽略粒子滑移等).因此,MHD 數(shù)值模擬也亟待可靠的MHD 流動(dòng)控制實(shí)驗(yàn)數(shù)據(jù)來對(duì)計(jì)算方法及結(jié)果進(jìn)行校驗(yàn).同時(shí),在相同的計(jì)算條件下,它往往要比CFD模擬耗費(fèi)更大的計(jì)算量;除了CFD 模擬存在的問題外,MHD 模擬還存在一些如磁場散度、全MHD 奇異性、導(dǎo)電壁面邊界等特有的問題[19].
基于高超聲速M(fèi)HD 的理論研究雖然起步非常早,但相比于實(shí)驗(yàn)研究和數(shù)值模擬,其研究工作仍然非常有限,且成果主要集中于二十世紀(jì)五六十年代.較為經(jīng)典的是Bush[20]的研究工作,他在激波層內(nèi)常密度、常電導(dǎo)率、無黏、球形激波等假設(shè)下,采用分離變量的方法求解流動(dòng)方程和電磁方程,將復(fù)雜的偏微分方程組化簡為常微分方程組,以轉(zhuǎn)化為一個(gè)初值問題進(jìn)行離散求解.Smith 和Wu[21]、Porter和Cambel[22]沿用Bush 的方法,分別研究了黏性系數(shù)和Hall 效應(yīng)的影響.該領(lǐng)域其他方面的理論研究工作也在同時(shí)期展開,Ludford 和Murray[23]研究了不可壓縮無黏導(dǎo)電流體流經(jīng)磁化球的穩(wěn)態(tài)流動(dòng)問題.Meyer[24]采用Newton-Plus-Centrifugal 近似方法,發(fā)展了MHD回轉(zhuǎn)體穩(wěn)態(tài)流動(dòng)理論.Lykoudis[25]引入Newton-Busemann 壓力定律,通過對(duì)無黏動(dòng)量方程進(jìn)行積分得到激波脫體距離的表達(dá)式.Eriscon 和Maciulaitis[26]研究了高度?速度圖上的電離流與磁場相互作用的分布范圍,并進(jìn)一步評(píng)估了較強(qiáng)相互作用流態(tài)下的磁流體控制作用力.Hooks 和Lewis[27]基于激波坐標(biāo)系和物體坐標(biāo)系之間的關(guān)系,得到了駐點(diǎn)處的無黏速度梯度、物體周圍的靜壓分布和邊界層外沿的速度分布等.激波脫體距離作為高超聲速飛行器飛行問題中的一個(gè)重要參數(shù),是氣動(dòng)領(lǐng)域研究的關(guān)注重點(diǎn)之一.然而,上述理論工作對(duì)于磁場影響下激波脫體距離的研究相對(duì)較少,相關(guān)結(jié)果也較為有限.有學(xué)者采用的方法雖然可以得到駐點(diǎn)線附近的流場信息,但其求解過程和結(jié)果不能直觀快速的反映磁控效果,有一定局限性[20-22];另一方面,Lykoudis[25]、Hooks 和Lewis[27]對(duì)激波脫體距離解析形式的理論推導(dǎo)是基于激波層內(nèi)磁場方向處處垂直于流速方向的條件下進(jìn)行的,該類型磁場可與電離流充分作用,然而其分布形式較為理想,且積分形式的激波脫體距離表達(dá)式需要進(jìn)一步數(shù)值求解.因此,高超聲速M(fèi)HD 激波脫體距離的快速求解方法仍有待進(jìn)一步發(fā)展和完善,以達(dá)到從理論上快速評(píng)估磁控效果以及能夠反映MHD 流動(dòng)控制本質(zhì)的目標(biāo).
本文通過理論分析的方法得到高超聲速M(fèi)HD無量綱激波脫體距離的解析表達(dá)式.求解過程主要分為兩部分,首先對(duì)連續(xù)方程進(jìn)行積分,將無量綱激波脫體距離作為待求量包含在內(nèi);進(jìn)一步,引入切向速度切向偏導(dǎo)數(shù)沿駐點(diǎn)線上的線性近似分布,利用徑向動(dòng)量方程建立磁場與激波脫體距離的聯(lián)系,尋找球頭駐點(diǎn)處切向速度沿切向的偏導(dǎo)數(shù)值.最后,通過求解關(guān)于激波脫體距離的一元二次方程,獲得包含磁場影響的解析表達(dá)式.
磁場添加后對(duì)流動(dòng)影響最直觀的體現(xiàn)是脫體激波距離或激波結(jié)構(gòu)的變化,學(xué)者在多種設(shè)備中均通過光學(xué)方法觀測(cè)到了磁場對(duì)脫體激波距離的影響,一般而言,磁場越強(qiáng),脫體激波距離越大.
本文的理論研究是基于以下條件:
(1)低磁雷諾數(shù)假設(shè),對(duì)于高超聲速飛行條件,激波層內(nèi)氣體電離度一般較低,磁雷諾數(shù)遠(yuǎn)小于1,此時(shí)誘導(dǎo)磁場對(duì)于外加磁場的影響可以忽略不記[1,28-29].Porter 和Cambel[29]詳細(xì)研究了不同磁雷諾數(shù)條件下的磁場分布,其結(jié)果也表明在低磁雷諾數(shù)條件下,電流導(dǎo)致的誘導(dǎo)磁場對(duì)外加磁場的影響很小,在激波層內(nèi)只考慮外加磁場即可;
(2)激波層內(nèi)流體不可壓縮,盡管自由來流速度是高超聲速,但強(qiáng)激波后激波層內(nèi)駐點(diǎn)線附近馬赫數(shù)較小,該區(qū)域內(nèi)的流體可當(dāng)不可壓縮流體近似處理,Bush[20]、Smith 和Wu[21]在求解高超聲速M(fèi)HD問題時(shí)也采用了該假設(shè);
(3)忽略黏性影響,高超聲速來流條件下,流體雷諾數(shù)較大,在關(guān)注脫體激波距離問題時(shí)可以忽略黏性項(xiàng),將問題簡化為無黏流動(dòng)問題[20,25];
(4)激波層內(nèi)電導(dǎo)率為常值,電導(dǎo)率一般與溫度相關(guān),在無黏假設(shè)條件下,激波層內(nèi)不考慮邊界層,沿駐 點(diǎn)線溫度變化較小,可將其以常值近似處理[17,20,25].
對(duì)連續(xù)方程的處理借鑒了Olivier[30]的方法,他通過對(duì)連續(xù)方程積分的方式獲得了包含激波脫體距離的方程,本文在其基礎(chǔ)上進(jìn)一步考慮添加磁場對(duì)脫體激波距離的影響,對(duì)該方法的具體應(yīng)用過程如下.球坐標(biāo)系下的連續(xù)方程為
其中 ρ∞和u∞為來流密度和速度,ρs為激波層內(nèi)密度,rs和rb分別為球頭半徑和激波半徑,Δ 為激波脫體距離,各變量示意圖如圖1 所示.
圖1 鈍頭體和激波結(jié)構(gòu)示意圖Fig.1 Schematic of the blunt body and the shock
將式(2)沿著駐點(diǎn)線由駐點(diǎn)向激波處進(jìn)行積分,積分的目的是將無量綱激波脫體距離加入到方程中,將其作為待求量來進(jìn)一步求解,積分過程具體如下
由于求解過程在駐點(diǎn)線附近,此時(shí)可以通過引入駐點(diǎn)線假設(shè)簡化問題的求解過程,此過程實(shí)質(zhì)上是去除方程中的高階小量,同時(shí)保留低階小量,從而保證了低階量的可求解性.根據(jù)駐點(diǎn)線假設(shè),引入如下近似條件[30]
Olivier[30]近似地給出了/?θ在駐點(diǎn)線上的線性分布表達(dá)式,既可以使/?θ滿足激波處的邊界條件,同時(shí)提供了壁面處 (/?θ)b的待求量,通過 (/?θ)b可以建立激波脫體距離與磁場和來流的聯(lián)系.以下給出/?θ在駐點(diǎn)線附近的線性分布形式
動(dòng)量方程包含了待求問題的特征元素,即電磁作用項(xiàng),同時(shí)需要通過對(duì)動(dòng)量方程的處理,尋找激波脫體距離與磁場影響的聯(lián)系,因此這一部分較為關(guān)鍵.通過處理徑向動(dòng)量方程,得到駐點(diǎn)處/?θ的具體數(shù)值,即 (/?θ)b.此處選取徑向動(dòng)量方程求解,這是考慮到壓力項(xiàng)的量級(jí)在分析的過程中會(huì)被排除掉,有助于求解的簡化.同時(shí),需要額外地說明一下,對(duì)動(dòng)量方程取旋度的方式是不可行的,原因是分離變量后,在駐點(diǎn)處方程兩邊恒為0,對(duì)于問題的求解無作用.對(duì)徑向的動(dòng)量方程進(jìn)行處理具體步驟如下.包含電磁項(xiàng)的徑向動(dòng)量方程為
其中σ為電導(dǎo)率.根據(jù)電磁感應(yīng)定律,變化的磁場會(huì)導(dǎo)致感生電場,但由于本文涉及的是穩(wěn)態(tài)磁場,所以電場需要滿足?×E=0,對(duì)于二維軸對(duì)稱問題,在駐點(diǎn)線上滿足E=0[20],于是式(9)可以進(jìn)一步整理
對(duì)于磁場分布,選取較為常用的模型為均勻磁化球形磁體的磁場分布,該模型的磁場分布簡單,存在解析形式的表達(dá),便于實(shí)驗(yàn)和理論的研究[12,20,22],這種磁場的方向與激波層內(nèi)電離流速度方向趨近于垂直,可以較為充分地與電離流發(fā)生相互作用,其分布形式見圖2.由式(10)知,電流為環(huán)向,垂直于磁場和速度矢量所在的面,于是MHD 動(dòng)量方程中的電磁力項(xiàng)J×B會(huì)受到增益的影響.這里參考Porter給出的偶極子磁場分布形式[22]
圖2 均勻磁化球體磁感線分布示意圖Fig.2 Distribution of magnetic induction lines of uniformly magnetized sphere
根據(jù)量級(jí)分析,式(16) 的壓力項(xiàng)不包含小量O(θ2),所以當(dāng)θ →0 時(shí),等式兩端包含O(θ2)的項(xiàng)都趨近0,而此時(shí)壓力項(xiàng)為0,將該結(jié)果代入式(16),可以得到
由于磁場會(huì)明顯影響流場,切向速度需要包含磁場以及來流所帶來的影響,所以這里選取第二個(gè)解.額外地,由式(16)還會(huì)得到于是可以發(fā)現(xiàn),當(dāng)出現(xiàn)磁場時(shí),駐點(diǎn)附近地切向速度有向駐點(diǎn)方向流動(dòng)的趨勢(shì),而該處的徑向速度可由微分形式的連續(xù)方程得到
此結(jié)果說明,駐點(diǎn)附近徑向速度方向有與來流相反的趨勢(shì).
根據(jù)激波前后質(zhì)量守恒關(guān)系,激波處駐點(diǎn)線上的速度有[19]
根據(jù)圖3,對(duì)于相同的相互作用系數(shù)Q,隨著ε的增加(來流速度降低)無量綱激波脫體距離變大.對(duì)于相同的ε,隨著Q增加激波無量綱脫體距離增加.從圖3 中也可以觀察到隨著ε的減小,無量綱激波脫體距離關(guān)于ε的靈敏度降低,在高超聲速M(fèi)HD實(shí)驗(yàn)或數(shù)值模擬開始之前,可以通過合理估計(jì)的ε值以及由初始條件構(gòu)成的Q來預(yù)測(cè)結(jié)果,這有助于快速評(píng)估磁控效果.
圖3 無量綱激波脫體距離與磁相互作用系數(shù)Q 的關(guān)系曲線Fig.3 The curve of dimensionless shock stand-off distance vs.MHD interaction parameter Q
圖4 對(duì)比了本文結(jié)果與Lefevre 等[31]在膨脹管中的實(shí)驗(yàn)結(jié)果.C2~ C6 代表Lefevre 選取的不同實(shí)驗(yàn)條件,每一個(gè)實(shí)驗(yàn)條件對(duì)應(yīng)著不同的相互作用系數(shù)Q以及激波前后氣流密度比ε,但由于實(shí)驗(yàn)時(shí)間很短,不足以形成穩(wěn)定流場,文章以20 μs 為間隔進(jìn)行數(shù)據(jù)取點(diǎn),Lefevre 同時(shí)對(duì)比了MHD 的CFD 模擬結(jié)果.從圖4 中可以看到,理論結(jié)果、實(shí)驗(yàn)和數(shù)值結(jié)果符合得較好.
圖4 無量綱激波脫體距離比/NM 與磁相互作用系數(shù)Q 的關(guān)系曲線Fig.4 The curve of ratio of shock stand-off distance /NM vs.MHD interaction parameter Q
Lykoudis[25]同樣在低磁雷諾數(shù)以及激波層內(nèi)常值密度、常值電導(dǎo)率、無黏的假設(shè)條件下推導(dǎo)了球頭激波脫體距離表達(dá)式,他采用了較為理想的磁場分布模型:J=σvBφ,即磁場處垂直于流線,這樣的磁場分布有助于問題簡化.在對(duì)動(dòng)量方程的處理過程中,他認(rèn)為壓力分布隨外加磁場的變化很小,于是通過引入Newton-Busemann 壓力分布近似公式進(jìn)一步補(bǔ)充了流場信息,結(jié)合流管上下游的質(zhì)量守恒關(guān)系,最終推導(dǎo)得到了激波脫體距離積分形式的解析表達(dá)式,分析結(jié)果反映了磁場對(duì)氣流流動(dòng)的影響效果.
圖5 給出了本文與Lykoudis 理論結(jié)果的對(duì)比,兩者結(jié)果趨勢(shì)一致,都體現(xiàn)了磁場與導(dǎo)電流體相互作用產(chǎn)生的阻滯效應(yīng),但兩者也存在一定偏差,主要原因如下:Lykoudis 采用的磁場分布模型與本文存在一定的差別,他的推導(dǎo)是基于激波層內(nèi)磁場處處垂直于流速的條件,相較于偶極子磁場而言,Lykoudis使用的磁場可與流場進(jìn)行更加充分的相互作用,這也會(huì)產(chǎn)生更強(qiáng)的阻滯作用,使激波層厚度增益效果更加明顯,從而導(dǎo)致其結(jié)果在一定Q值范圍內(nèi)高于本文所得的結(jié)果.當(dāng)然,由于Lykoudis 采用的磁場過于理想,在現(xiàn)實(shí)中難以用磁場裝置產(chǎn)生該類磁場,在實(shí)驗(yàn)研究中難以采用.另外,需要注意的是,包括Lykoudis 文章在內(nèi)的一些文獻(xiàn)中的相互作用系數(shù)表達(dá)式為Qs=σB2rs/(ρ∞u∞),與本文使用的特征長度不一致,需要進(jìn)行轉(zhuǎn)換.
圖5 本文結(jié)果與Lykoudis 結(jié)果對(duì)比圖[25]Fig.5 Comparison between results of this paper and the Lykoudis’s
通常情況下,理論分析引入的假設(shè)條件會(huì)導(dǎo)致模型存在一定的局限性.根據(jù)線性近似分布假設(shè),/?θ 是僅關(guān)于r?的函數(shù),當(dāng)來流速度越高時(shí),對(duì)于相同磁相互作用系數(shù)Q條件下,慣性力的影響越強(qiáng),/?θ在駐點(diǎn)線上的分布越符合線性近似假設(shè),然而隨著激波脫體距離的增加,該假設(shè)會(huì)出現(xiàn)一定程度的失效,即/?θ 將與r?和兩個(gè)量相關(guān),/?θ的分布與激波層厚度有關(guān)聯(lián),于是式(5)的積分結(jié)果將會(huì)受到影響,導(dǎo)致無法將其化為關(guān)于激波脫體距離的一元二次方程形式.總之,在Q較大時(shí),/?θ分布偏離線性假設(shè)是限制本文模型應(yīng)用的關(guān)鍵點(diǎn).然而,在Q>1 時(shí)就會(huì)出現(xiàn)明顯的MHD 流動(dòng)控制效果,且當(dāng)前實(shí)驗(yàn)條件還難以達(dá)到較大的Q值,因此本文的模型仍然有較強(qiáng)的指導(dǎo)意義.當(dāng)然,后續(xù)會(huì)通過更多的理論和實(shí)驗(yàn)來對(duì)該模型在Q較大時(shí),/?θ分布偏離線性假設(shè)是限制本文模型應(yīng)用值下進(jìn)行改進(jìn)優(yōu)化.
本文的結(jié)果中有一個(gè)值得關(guān)注的地方,當(dāng)磁場存在時(shí),駐點(diǎn)附近的切向速度有指向駐點(diǎn)的趨勢(shì),通過連續(xù)方程得到與來流相反的徑向速度,如圖6.這里認(rèn)為該結(jié)果是低磁條件假設(shè)造成的,原因是磁場的旋度方程中,如果考慮電流對(duì)磁場的影響,該影響將會(huì)是一個(gè)小量,這個(gè)小量會(huì)影響徑向動(dòng)量方程的量級(jí)分析過程,從而改變求解結(jié)果.然而這種流動(dòng)趨勢(shì)是實(shí)際存在的,由于電離流與磁場相互作用的影響,流體會(huì)受到與流動(dòng)相反方向的洛倫茲力,從而使激波層內(nèi)的流速降低,從而導(dǎo)致了激波脫體距離的增加.
圖6 駐點(diǎn)處流動(dòng)示意圖Fig.6 Schematic of flow near stagnation point
本文利用連續(xù)方程建立了關(guān)于激波脫體距離的待求方程,通過徑向動(dòng)量方程尋找 (?v/?θ)b的具體數(shù)值,從而使磁場影響與待求方程得以關(guān)聯(lián),將問題化為求解一個(gè)關(guān)于激波脫體距離的一元二次方程,得到了較為直觀簡潔的解析表達(dá)式.本文主要結(jié)論如下:
(1)理論分析結(jié)果表明,當(dāng)磁相互作用系數(shù)增加時(shí),激波脫體距離增加,即磁控效果增強(qiáng).將本文結(jié)果與Lefevre 的實(shí)驗(yàn)數(shù)據(jù)、Lykoudis 的結(jié)果進(jìn)行了對(duì)比,三者吻合較好.本文理論模型對(duì)高超聲速M(fèi)HD實(shí)驗(yàn)設(shè)計(jì)和結(jié)果分析具有較好的指導(dǎo)意義;
(2)本文模型的建立是基于特定的假設(shè)條件,存在一定的適用范圍,在Q較大時(shí),/?θ分布偏離線性假設(shè)會(huì)限制本文模型的應(yīng)用.通常情況下,當(dāng)Q>1 時(shí)即會(huì)出現(xiàn)明顯的MHD 流動(dòng)控制效果,本文結(jié)果的適用性可以達(dá)到并遠(yuǎn)大于該范圍.