磁場對激波沖擊R22重氣柱作用過程影響的數(shù)值模擬＊

林震亞，張煥好，陳志華，劉 迎

（南京理工大學(xué)瞬態(tài)物理國家重點(diǎn)實(shí)驗(yàn)室，江蘇南京210094）

磁場對激波沖擊R22重氣柱作用過程影響的數(shù)值模擬＊

林震亞，張煥好，陳志華，劉 迎

（南京理工大學(xué)瞬態(tài)物理國家重點(diǎn)實(shí)驗(yàn)室，江蘇南京210094）

為研究平面入射激波與磁化R22重質(zhì)圓形氣柱的作用過程，首先通過數(shù)值方法得到了不同初始條件下激波誘導(dǎo)R22氣柱的Kelvin－Helmholtz（KH）及Richtmyer－Meshkov（RM）不穩(wěn)定性導(dǎo)致的重氣柱變形過程，并詳細(xì)討論了不同情況下透射激波在氣柱內(nèi)聚焦誘導(dǎo)射流的過程；然后在加入磁場的情況下，采用CTU＋CT算法進(jìn)行數(shù)值模擬，以保證數(shù)值結(jié)果滿足任意時(shí)刻磁場的散度為零。計(jì)算結(jié)果表明：磁場對激波誘導(dǎo)R22氣柱不穩(wěn)定性具有抑制作用；法向磁場和流向磁場都可以很好地抑制RM不穩(wěn)定性；對于KH不穩(wěn)定性，法向磁場的控制效果更好，不僅可以抑制界面上渦串的卷起，還可以阻止主渦的發(fā)展，而流向磁場做不到后者；磁場對射流影響不大，射流處的磁能量可以一定程度上抑制射流的衰減，同時(shí)法向磁場可以減小聚焦時(shí)壓力及速度峰值。

磁化氣體；磁流體動力學(xué)；磁場；不穩(wěn)定性；激波

激波與兩種不同密度氣體界面相互作用時(shí)，會在界面附近發(fā)生反射、折射及繞射等現(xiàn)象，同時(shí)在界面處還會產(chǎn)生Richtmyer－Meshkov（RM）和Kelvin－Helmholtz（KH）不穩(wěn)定性，導(dǎo)致界面扭曲變形。該現(xiàn)象廣泛存在于超燃沖壓發(fā)動機(jī)中的燃料混合、慣性約束核聚變、水中炸藥爆炸及天體物理中的超新星爆炸等領(lǐng)域，因而具有重要的研究價(jià)值。

目前，對激波與氣體界面相互作用的研究主要集中以下兩方面：一是氣體界面在RM不穩(wěn)定性與KH不穩(wěn)定性作用下的演變過程及其相關(guān)復(fù)雜流動機(jī)理；二是界面演化后期的失穩(wěn)、破碎及湍流混合問題。K．A．Meyer等［1］對激波與氣體界面作用時(shí)的RM不穩(wěn)定性進(jìn)行了數(shù)值模擬，得到了與Meshkov實(shí)驗(yàn)定性一致的結(jié)果；Q．Zhang等［2］研究了二維和三維初始擾動的發(fā)展過程，得到了氣泡、尖釘速度及界面發(fā)展的總體增長率的計(jì)算公式；N．N．Anuchina等［3］對初始擾動幅度下相同入射強(qiáng)度激波與正弦擾動氣體界面的作用過程進(jìn)行了數(shù)值模擬，結(jié)果表明三維數(shù)值計(jì)算得到的擾動增長率高于二維結(jié)果；G．A．Ruev等［4］對馬赫數(shù)為3．5的平面入射激波與正弦擾動氣體界面的作用過程進(jìn)行了數(shù)值模擬，給出了氣體界面RM不穩(wěn)定性的發(fā)展過程及兩種氣體混合區(qū)的寬度變化；J．H．J．Niederhaus等［5］對不同強(qiáng)度平面入射激波與多種不同介質(zhì)的氣柱、氣泡相互作用過程進(jìn)行了大量的數(shù)值研究，重點(diǎn)分析了作用過程中流場環(huán)量的變化和兩種氣體的混合機(jī)理；同時(shí)，B．Thornber等［6］對低馬赫數(shù)的平面入射激波與R22方柱作用過程的數(shù)值模擬，結(jié)果與Holder等的實(shí)驗(yàn)結(jié)果完全相吻；B．Hejazialhosseini等［7］對不同馬赫數(shù)下平面入射激波與氦氣氣泡作用過程進(jìn)行了數(shù)值研究，得到了氣泡的變形過程及氣泡表面在KH不穩(wěn)定作用下次級渦的產(chǎn)生過程；O．Schilling等［8］采用數(shù)值方法研究了平面入射激波與Air／R22氣體界面的相互作用過程，并分析了反射激波與分界面作用前后流場結(jié)構(gòu)和混合層演變；B．Tian等［9］利用ALE（arbitrary Lagrangian－Eulerian）和GALE（global arbitrary Lagrangian－Eulerian）方法對多種激波與氣體界面作用問題進(jìn)行了數(shù)值模擬；S．K．Shankar等［10］采用數(shù)值方法模擬了平面入射激波與重氣柱作用過程，結(jié)果表明流場初始壓力和密度梯度是影響兩個(gè)主渦產(chǎn)生的關(guān)鍵因素；C．Bailie等［11］對不同強(qiáng)度平面入射激波與Air／R22正弦小擾動界面的作用過程進(jìn)行了數(shù)值模擬，并分析入射激波強(qiáng)度對作用過程的影響。

目前，對界面RM不穩(wěn)定性控制方法的研究表明，合適的磁場對氣體界面RM不穩(wěn)定性的產(chǎn)生與發(fā)展具有一定的抑制作用，但仍處在理論研究階段。S．Chandrasekhar等［12］和V．Wheatley等［13－14］采用線性理論方法推導(dǎo)了垂直或平行于界面的磁場對不可壓縮流體界面不穩(wěn)定性的影響，結(jié)果表明：平行于界面的磁場能夠抑制RT和RM不穩(wěn)定性；但是在線性理論下，平面流場中的RM不穩(wěn)定性與垂直于平面的磁場無關(guān)，與實(shí)際情況不符。為進(jìn)一步研究磁場效應(yīng)，M．Khan等［15］在勢流理論的基礎(chǔ)上建立了新模型，發(fā)現(xiàn)考慮非線性效應(yīng)后，磁場能夠影響界面不穩(wěn)定性的發(fā)展。J．Cao等［16］同時(shí)考慮了初始平行于界面的磁場和剪切流對界面不穩(wěn)定性的影響，發(fā)現(xiàn)磁場能夠抑制界面不穩(wěn)定性，而剪切流將會促進(jìn)界面不穩(wěn)定性的發(fā)展，兩者之間存在競爭機(jī)制。M．S．Shin等［17］對不同方向磁場作用下，平面激波與球形氣泡相互作用過程進(jìn)行了數(shù)值模擬，發(fā)現(xiàn)在發(fā)展初期，磁場效應(yīng)很小，而在發(fā)展后期，即使是弱磁場，其對界面形狀的影響也很明顯。此外，磁場越強(qiáng)，界面內(nèi)外物質(zhì)混合率越低，并越能抑制后期湍流的發(fā)展。李源等［18］基于理論分析的方法，綜合考慮黏性和表面張力的影響，考察了磁場中非理想流體中RT不穩(wěn)定性氣泡的演化過程，并推導(dǎo)出描述二維非理想磁流體RT不穩(wěn)定性氣泡運(yùn)動的控制方程組，得到了流體黏性、表面張力和磁場對氣泡界面不穩(wěn)定性發(fā)展的影響。

綜上所述，目前對磁場控制不穩(wěn)定性的研究主要集中在基于不可壓理想磁流場的解析推導(dǎo)及磁場影響下所產(chǎn)生現(xiàn)象的數(shù)值模擬兩個(gè)方面。本文中，基于磁流體動力學(xué)（magneto－h(huán)ydro－dynamics，MHD）方程，采用CTU＋CT算法對激波沖擊磁化R22氣柱進(jìn)行數(shù)值模擬，并以磁壓力為切入點(diǎn)研究磁場對不穩(wěn)定性及射流的影響。

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

1．1 數(shù)值方法

在滿足流體力學(xué)方程的條件下，帶電荷的流體團(tuán)更有利于形成團(tuán)塊，主要原因如下：（1）帶電粒子間的碰撞頻率遠(yuǎn)大于中性粒子之間的碰撞頻率，因此其碰撞平均自由程較短；（2）有磁場存在時(shí)，磁場對帶電粒子有側(cè)向約束。洛倫茲力使帶電粒子在垂直于磁場方向只能保持在回轉(zhuǎn)半徑范圍內(nèi)運(yùn)動。等離子體通常比中性粒子有更低的密度和更高的溫度，由于整個(gè)過程在短時(shí)間內(nèi)完成，本文中考慮高磁雷諾數(shù)的情況，使用MHD方程來研究激波沖擊R22氣柱中出現(xiàn)的不穩(wěn)定性：

式中：ρ為質(zhì)量密度，v為速度矢量，p為壓力，E為能量，μ為磁導(dǎo)率，B為磁感應(yīng)強(qiáng)度，＾b為B的單位矢量，Π為包含各向同性與各向異性部分的黏性應(yīng)力張量，分別由系數(shù)ν0和ν‖控制，γ為比熱比，T為溫度，I為單位張量，符號“∶”表示雙點(diǎn)積運(yùn)算。在各向異性的情況下，黏性通量限制在平行于磁場線方向。Q為熱通量，同樣包含各向同性與各向異性的情況，分別由系數(shù)κ0和κ‖控制。在各向異性的情況下，熱通量同樣限制在平行于磁場線方向。同時(shí)，感應(yīng)方程包含多種非理想MHD效應(yīng)，包括歐姆損耗（由電阻率η控制），霍爾效應(yīng)（由系數(shù)ηH控制）和雙極性擴(kuò)散（由系數(shù)ηAD控制）。

本文中采用CTU＋CT算法對MHD方程組進(jìn)行求解，其中CTU算法是一種用于求解守恒雙曲系統(tǒng)的不分離的二維有限體積算法，其三維變形則由J．Saltzman［19］提出。另外，由于傳統(tǒng)的12－solve CTU算法忽略了MHD方程中始終存在的約束·B＝0，不能作為求解MHD方程的有效算法，因此需要構(gòu)建一個(gè)Godunov通量中網(wǎng)格單元邊緣的平均電場算法，通常稱為CT算法，該算法可以描述為以Godunov電場為預(yù)測值，以CT電場為校正值的預(yù)測過程。鑒于在文獻(xiàn)［20］中εcCT已顯示的最佳性能，本文中采用εcCT算法。該算法包含逆風(fēng)偏差（根據(jù)接觸模式），可以將網(wǎng)格對齊的平行平面流動降低到正確的Godunov電動勢。

1．2 計(jì)算模型

圖1為激波與R22氣柱相互作用的計(jì)算模型，其中R22氣泡直徑為0．05m，氣柱中心與左邊界的距離為0．05m，與右邊界的距離為0．25m，計(jì)算域高度為0．089m。初始時(shí)，氣柱內(nèi)采用R22作為重質(zhì)氣體，外部介質(zhì)為空氣，氣體參數(shù)如表1所示。氣柱內(nèi)外壓力均為101 325Pa，入射激波從左向右傳播，來流馬赫數(shù)Ma＝1．22，上、下邊界設(shè)為固體反射邊界，而右邊界設(shè)為出口條件。計(jì)算域網(wǎng)格采用均勻分布的笛卡爾網(wǎng)格，經(jīng)網(wǎng)格收斂性測試后，區(qū)域網(wǎng)格總數(shù)設(shè)置為3 000×890。本文中，設(shè)初始磁場強(qiáng)度為0．01T（β≈2 500），令平行于來流的方向?yàn)閤方向，垂直于來流的方向?yàn)閥方向。

圖1 計(jì)算模型Fig．1 Computational model

表1 氣體參數(shù)Table 1 Gas parameters used in this paper

為了使氣體受磁場影響，首先通過高溫氣體電離，溫度為12 000K。由于含等離子體，電導(dǎo)率、熱導(dǎo)率、黏性系數(shù)、霍爾系數(shù)以及雙極擴(kuò)散系數(shù)由以下方法確定。電導(dǎo)率采用基于玻爾茲曼方程得到的Spitzer公式進(jìn)行計(jì)算，即由玻爾茲曼方程的碰撞積分出發(fā)，忽略短程碰撞的影響，考慮剩余的碰撞項(xiàng)，得到弱電場情況下由離子與電子所組成的等離子體電導(dǎo)率（σ）為：

式中：j＝m／（2kBT），其中m為質(zhì)量，kB為Boltzmann常數(shù)；γ0＝0．1Z＋0．49；e為電子電荷量；Z為離子電荷數(shù)。

如果碰撞積分到德拜長度為止，則庫侖對數(shù)項(xiàng)為：

計(jì)算可得，初始電導(dǎo)率約為107S／m。

當(dāng)氣體中相當(dāng)部分的原子被電離后，靜電力將對粒子間碰撞起主要作用，若考慮遠(yuǎn)程作用，則離子碰撞將不能視為兩體方式。通過求解玻爾茲曼方程得到一階近似電離氣體熱導(dǎo)率為：

根據(jù)S．Chapman等［22］的理論可以推導(dǎo)出kc、ν的二階近似值。盡管二階近似值與一階近似值有顯著的差異，但其后各階近似所增加的修正量都是很微小的。假定黏度均由重離子引起，電子黏度很小，由于電子的熱運(yùn)動速度遠(yuǎn)大于離子的熱運(yùn)動速度，因此熱擴(kuò)散效應(yīng)主要是由電子引起。經(jīng)過估算：本文中黏性系數(shù)約為3．72×10－5Pa·s，熱導(dǎo)率約為1．4W／（m·K）。

由于準(zhǔn)電中性條件，電子和離子在等離子體內(nèi)不可能獨(dú)立運(yùn)動，當(dāng)電子快速離開某等離子體單元時(shí)，單元內(nèi)將形成由正電荷產(chǎn)生的電場，該電場阻礙電子繼續(xù)離開所在體積單元，并促使離子更快地離去，這種擴(kuò)散狀態(tài)稱為雙極型狀態(tài)。由于電子的擴(kuò)散系數(shù)遠(yuǎn)大于離子，因此：

式中：ηi為正離子的熱擴(kuò)散系數(shù)，Te、Ti分別為電子、離子溫度。由此可見，雙極擴(kuò)散系數(shù)恒小于電子自由擴(kuò)散系數(shù)，而大于離子擴(kuò)散系數(shù)，因此雙極性電場極大影響了電子的定向速度。本文中取雙極擴(kuò)散系數(shù)為（10－5／ρ）m2·Pa·s－1。

Hall效應(yīng)可在磁流體通道中產(chǎn)生橫向洛倫茲力，而且體力的不均勻分布可導(dǎo)致二次流的產(chǎn)生，對流體的運(yùn)動產(chǎn)生重要影響，因此需要考慮霍爾效應(yīng)。此處，霍爾系數(shù)ηH取為10－7m3／C。

2 結(jié)果與討論

為了驗(yàn)證數(shù)值方法與計(jì)算模型的可靠性，圖2給出了不加磁場時(shí)入射激波與R22氣柱相互作用過程計(jì)算結(jié)果與文獻(xiàn)［21］中實(shí)驗(yàn)結(jié)果的對比，兩者具有相同的初始條件。由圖2可知：計(jì)算結(jié)果中激波的反射、繞射及發(fā)展過程均與實(shí)驗(yàn)結(jié)果完全相符，兩者氣柱的變形過程也相似；但是，由于初始時(shí)兩種氣體接觸界面采用間斷面形式，因此計(jì)算結(jié)果中接觸界面上的KH不穩(wěn)定性更明顯，從而導(dǎo)致界面上、下兩側(cè)的湍流轉(zhuǎn)捩過程（見圖2（e））與實(shí)驗(yàn)結(jié)果稍有不同。

圖2無磁場情況下，激波與R22氣柱相互作用過程的計(jì)算紋影（上）與文獻(xiàn)［21］實(shí)驗(yàn)結(jié)果（下）的對比Fig．2 Numerical results of the interaction of shock wave with the R22air column in this paper（upper row），and the experimental result in reference［21］（lower row）．

圖3 為后期激波與R22氣柱作用過程的計(jì)算紋影。從圖3可以看出：由于激波對不同密度氣體界面的作用，界面產(chǎn)生了RM不穩(wěn)定性；在作用初期，隨機(jī)擾動在氣體界面發(fā)生，擾動振幅在界面線性增長；隨著變形的加劇，振幅經(jīng)歷非線性增長，“氣泡”與“尖釘”結(jié)構(gòu)出現(xiàn)；隨著運(yùn)動的繼續(xù)進(jìn)行，在KH不穩(wěn)定性的影響下，“尖釘”破碎，“氣泡”被拉長，最終形成湍流混合狀態(tài)。

圖3后期激波與R22氣柱作用過程計(jì)算紋影Fig．3 Numerical results of the interaction of shock wave with the R22air column in later period．

圖4 為加入均勻磁場后，激波與R22氣柱相互作用過程的密度紋影，其中上方子圖的磁場方向垂直于來流方向（即磁場方向沿y軸正向），而下方子圖中磁場方向平行于來流方向。由圖4可知，與無磁場情況（見圖2）相比，磁作用力能明顯抑制界面處的KH不穩(wěn)定性，且垂直于來流方向的磁作用力對不穩(wěn)定性的控制較平行于來流方向的磁作用力更明顯。為了更清晰地描述此現(xiàn)象，將動量方程中與磁場相關(guān)的項(xiàng)拆分為磁壓力項(xiàng)（B·B）／（2μ）與磁張力項(xiàng)·（BB）／μ兩部分，這兩項(xiàng)均為二階張量。然而，由于計(jì)算結(jié)果中磁壓力值遠(yuǎn)大于磁張力值，下文只對磁壓力的變化情況進(jìn)行討論。

由式（4）可知，垂直于速度方向的磁場分量梯度較大，磁場變化較快，該效應(yīng)會引起磁場在某些位置（如垂直于磁場方向的界面處）的集中。因此，隨著時(shí)間的推移，氣泡界面處的磁場強(qiáng)度不斷增大，妨礙輕質(zhì)氣體和重質(zhì)氣體的相互混合。由圖4可知：當(dāng)初始磁場垂直于來流方向時(shí)，磁壓力作用在x軸方向，可以有效阻止氣柱兩側(cè)的卷起；在t＝450μs之后，由于斜壓作用加劇，氣柱尾部界面開處始形成渦結(jié)構(gòu)；然而，沿x軸方向的磁壓力抑制了不穩(wěn)定性的進(jìn)一步發(fā)展，最終在氣柱尾部界面附近形成了上、下兩個(gè)突起，如圖4（f）所示。當(dāng)磁場平行于來流方向時(shí)，磁壓力作用在y軸方向，因上側(cè)界面處x軸向速度的差異，在界面處卷起渦串，如圖4（c）和圖4（d）所示；同時(shí)，由于斜壓效應(yīng)，氣柱上、下兩側(cè)卷起兩個(gè)渦環(huán)，如圖4（e）所示；隨后，沿y軸方向的磁壓力阻止了主渦的進(jìn)一步發(fā)展，最終在主渦位置形成“尖釘”狀結(jié)構(gòu)，見圖4（f）。由于磁化氣柱的熱導(dǎo)率遠(yuǎn)大于空氣的，界面處氣體混合速度較無磁場時(shí)更快。

圖4 B＝0．01T時(shí)不同磁場方向?qū)げㄅcR22氣柱作用過程的影響Fig．4 Numerical results of the interaction process of shock wave and R22air column when B＝0．01T

圖5和圖6分別給出了磁場方向垂直或平行于來流方向時(shí)氣柱對稱軸上壓力（p）與x方向速度（u）的分布。由圖5和圖6可知，與無磁場的情況［23］相比，初始壓力與速度變化不大。當(dāng)t＝200μs時(shí)，透射激波在x軸上發(fā)生聚焦，并且在匯聚處產(chǎn)生局部高壓區(qū)：初始磁場沿y軸正方向時(shí)，壓力曲線峰值pm＝700kPa，初始磁場沿x軸正方向時(shí)，pm＝760kPa，均小于無磁場情況。此時(shí)，x方向速度同樣處于峰值，但施加沿x軸方向磁場時(shí)的速度（u＝153m／s）遠(yuǎn)大于無磁場（u＝100m／s）及施加沿y軸方向磁場時(shí)的速度（u＝125m／s）。隨著相交點(diǎn)向下游移動，聚焦點(diǎn)的壓力峰值不斷下降，速度逐漸增大，并在t＝250μs時(shí)達(dá)到最大值（u≈170m／s）。隨后，在流體阻力作用下，速度不斷衰減。綜上所述，磁場的加入對壓力影響不大，但卻能提升透射激波聚焦時(shí)速度峰值，而磁場沿x軸方向時(shí)影響更劇烈。此外，在激波衰減過程中，沿y方向磁場較沿x方向磁場更能減緩射流的衰減。

圖5 初始磁場沿y軸正方向時(shí)入射激波誘導(dǎo)射流產(chǎn)生過程Fig．5 Generation process of jet induced by incident shock wave when initial magnetic field is along yaxis

圖6 初始磁場沿x軸正方向時(shí)入射激波誘導(dǎo)射流產(chǎn)生過程Fig．6 Generation process of jet induced by incident shock wave when initial magnetic field is along xaxis

圖7 給出了磁場沿y軸正方向或x軸正方向時(shí)激波與R22氣柱作用過程的計(jì)算紋影。由圖7可知，磁場能夠較好抑制RM不穩(wěn)定性。當(dāng)初始磁場垂直于來流方向時(shí)，界面處的橫向磁壓力對氣體混合的抑制效果較好，同時(shí)氣柱截面被壓扁拉長，兩個(gè)主渦在磁場作用下受到抑制，最終被拉成長條狀；當(dāng)初始磁場平行于來流方向時(shí)，后期氣柱寬度小于磁場垂直于來流方向的情形，雖然形成兩個(gè)主渦，但磁場抑制了氣柱的破碎及次級渦的生成。由此可見，加入磁場后，氣柱界面處沒有在KH不穩(wěn)定性的作用下持續(xù)失穩(wěn)并卷起“珠狀”次級渦，氣柱頂端也未出現(xiàn)RM不穩(wěn)定性，因此，兩種方向的初始磁場對不穩(wěn)定性都有較好的抑制作用。

為更清楚地看到磁場對氣柱界面的影響，圖8給出了初始磁場沿y軸正方向時(shí)不同時(shí)刻磁壓力和磁能量的分布。由圖8可知，初始時(shí)刻，透射激波在氣柱內(nèi)部向下游傳播，圓弧反射激波向上游傳播，磁壓力較大處集中在氣柱界面處，磁能量較大處集中在界面處及透射激波尾端。當(dāng)t＝165μs時(shí)，x方向磁壓力沿x軸對稱分布，而y方向磁壓力沿x軸反對稱分布；同時(shí)x方向磁壓力對氣體界面內(nèi)側(cè)產(chǎn)生向右拉力，對氣體界面外側(cè)則產(chǎn)生向左拉力，而y方向磁壓力對氣體界面內(nèi)側(cè)產(chǎn)生向氣柱中心的拉力，而對氣體界面外側(cè)產(chǎn)生向外的拉力。另外，當(dāng)入射激波與氣柱作用時(shí)，因氣柱內(nèi)外氣體聲阻抗的不同，氣柱內(nèi)外會產(chǎn)生速度差。氣柱外部界面處，空氣x方向速度大于界面速度，而氣柱內(nèi)部界面處，空氣x方向速度小于當(dāng)前界面速度。由于磁場變化率與速度差密切相關(guān)，因此這兩處的速度梯度會帶來磁場的變化，由此產(chǎn)生的磁壓力進(jìn)一步抑制氣柱界面受內(nèi)外剪切力的影響（同理，對空氣y方向速度，也有類似結(jié)論）。由圖7（c）和圖7（d）可知，發(fā)展后期，磁能量與磁壓力分布類似，包裹在界面處。綜上所述，當(dāng)初始磁場垂直于來流方向時(shí)，磁壓力始終包裹著界面，阻止界面受周圍環(huán)境的進(jìn)一步影響，從而防止界面上渦量的生成和擴(kuò)散，抑制界面的不穩(wěn)定性。

圖7 后期激波與R22氣柱作用過程密度紋影Fig．7 Numerical results of interaction process of shock wave with R22air column in later period

圖8 初始磁場沿y軸正方向時(shí)磁壓力及磁能量的分布（每幅子圖中，左上為磁壓力，左下為磁能量，中間和右側(cè)分別為磁壓力在x、y方向的分量）Fig．8 Magnetic pressure and magnetic energy when initial magnetic field is along yaxis（upper left：magnetic pressure；lower left：magnetic energy；middle：magnetic pressure along xaxis；right：magnetic pressure along yaxis）

圖9為初始磁場沿y軸正方向時(shí)磁壓力（pm）和磁能量（Em）沿橫向（x方向）及縱向（y方向）的分布情況。由圖9可知：當(dāng)透射激波在氣柱內(nèi)傳播時(shí)，磁能量在兩個(gè)位置取值較大，一處是氣柱左側(cè)界面處，另一處位于透射激波尾部，與之前分析相符；界面處的磁壓力遠(yuǎn)大于其他位置。由圖9可知，氣柱上、下界面處的磁壓力及磁能量遠(yuǎn)大于左、右界面處。圖9（c）顯示：當(dāng)透射及入射激波先后聚焦后，會在氣柱尾端射流位置形成局部高磁壓力區(qū)域，同時(shí)，該區(qū)域磁能量很低，接近于零。由于圖9（c）和圖9（d）中射流位置的磁壓力及磁能量很低，因此磁場對射流并無明顯影響；但是射流前、后的磁能量可以一定程度上抑制射流的衰減。此外，圖9中不同時(shí)刻磁壓力和磁能量分布曲線的形狀無明顯變化，但上、下界面處磁壓力和磁能量的峰值隨時(shí)間不斷增加，這意味著界面處磁場強(qiáng)度隨時(shí)間不斷增大。

圖10為初始磁場沿x軸正方向時(shí)磁壓力及磁能量的分布。由圖10可知，磁壓力與磁能量的分布類似，其峰值均集中在界面附近。與初始磁場沿y軸正方向的情況相似，x方向磁壓力沿x軸對稱分布，y方向磁壓力沿x軸反對稱分布，同時(shí)磁壓力與磁能量的作用區(qū)域更緊貼界面。由于磁壓力較弱且分布不均勻，圖11給出了邊界處磁壓力矢量分布。由圖11（a）可知，當(dāng)磁場平行于來流方向時(shí)，初始磁壓力分布于界面處，且與剪切力及外界對其壓力的合力方向相反。隨后，在剪切力的作用下，氣柱表面開始卷起渦，磁壓力的方向隨著界面變形而改變，如圖11（b）所示。在主渦生成時(shí)，磁壓力較大的位置依在界面附近，此對渦中心部分影響不大，如圖11（c）和圖11（d）所示。由此可見：當(dāng)渦量大于零時(shí)（渦逆時(shí)針方向），磁壓力指向渦動中心；而當(dāng)渦量小于零時(shí)，磁壓力指向外部。

圖9 初始磁場沿y軸正方向時(shí)的磁壓力及磁能量分布Fig．9 Distribution of magnetic pressure and magnetic energy when initial magnetic field is along yaxis

圖10 初始磁場沿x軸正方向時(shí)磁壓力及磁能量分布（每幅子圖中，左上為磁壓力，左下為磁能量，中間和右側(cè)分別為磁壓力在x、y方向的分量）Fig．10 Magnetic pressure and magnetic energy when initial magnetic field is along xaxis（upper left：magnetic pressure；lower left：magnetic energy；middle：magnetic pressure along xaxis；right：magnetic pressure along yaxis）

圖12給出了初始磁場沿x軸正方向時(shí)磁壓力及磁能量沿橫向及縱向的分布。由圖12（a）和圖12（b）可知，初始時(shí)刻，橫軸和縱軸上磁壓力及磁能量均小于初始磁場沿y軸正方向的情況。當(dāng)t＝165μs時(shí)，橫軸上的磁壓力集中在透射激波尾部，磁能量在頂部界面后迅速減小，接近透射激波時(shí)快速增大并達(dá)到峰值；而縱軸上磁能量和磁壓力的峰值均在上、下界面處。當(dāng)t＝200μs時(shí)，橫軸上的磁壓力集中在氣柱尾端，磁能量仍然在透射激波尾部，縱軸上的磁壓力不斷增大；此時(shí)由于渦在界面處開始生成，縱軸上磁壓力的分布不對稱。當(dāng)t＝370μs時(shí)，透射激波聚集，在尾部形成射流，圓弧形二次激波向外傳播，并最終碰到上、下界面發(fā)生反射；此時(shí)，橫軸上磁壓力與磁能量最大值位置相同，均在氣柱尾部界面處達(dá)到峰值，而縱軸上的磁壓力峰值則繼續(xù)增大，如圖12（c）所示。當(dāng)兩個(gè)主渦形成后，縱軸處主渦位置的磁能量較為明顯，但此時(shí)磁能量及磁壓力最大位置仍在界面附近；橫軸上磁壓力峰值減小，磁能量繼續(xù)增加，并集中在氣柱尾部，如圖12（d）所示。

圖11 初始磁場沿x軸正方向時(shí)邊界處磁壓力矢量圖Fig．11 Vector diagram of magnetic pressure on the interface when initial magnetic field is along xaxis

圖12 初始磁場為x軸正方向時(shí)磁壓力及磁能量分布Fig．12 Distribution of magnetic pressure and magnetic energy when initial magnetic field is along xaxis

3 結(jié) 論

基于CTU＋CT算法，對磁場環(huán)境下激波沖擊重質(zhì)氣體進(jìn)行數(shù)值模擬，并與不加磁場的情況進(jìn)行對比，所得結(jié)論如下。

（1）R22氣柱與入射激波作用后，由于KH不穩(wěn)定性，氣柱上、下兩側(cè)不斷卷起最終形成兩個(gè)主渦。同時(shí)，由于RM不穩(wěn)定性，后期氣柱左端界面處產(chǎn)生了“尖釘”與“氣泡”狀結(jié)構(gòu)。磁場能夠抑制界面不穩(wěn)定性，平行磁場通過減少界面上渦量的生成減緩不穩(wěn)定性的增長；垂直磁場則通過包裹界面的磁力線壓縮界面，抑制剪切力對界面的作用，使界面處不再卷起，從而抑制界面不穩(wěn)定性。

（2）垂直磁場和平行磁場都可以很好地抑制RM不穩(wěn)定性，對于KH不穩(wěn)定性，垂直磁場的控制效果更好，不僅可以抑制界面上渦串的卷起，還可以阻止主渦的發(fā)展，而平行磁場做不到后者。磁場對射流影響不大，射流處的磁能量可以在一定程度上抑制射流的衰減，同時(shí)，垂直磁場可以減小射流聚焦時(shí)壓力及速度峰值。

（3）當(dāng)磁場強(qiáng)度較?。é隆? 500）時(shí)，初始磁場對激波與R22氣柱作用過程沒有顯著影響，隨后磁壓力及磁能量不斷增加，從而抑制不穩(wěn)定性的發(fā)展。磁壓力及磁能量對界面影響很大，在上、下界面處達(dá)到峰值。同時(shí)，磁壓力對界面沒有壓縮作用，而是產(chǎn)生一個(gè)反向作用力減小外界對界面的影響。當(dāng)渦形成后，磁場對渦量影響不大，但可以通過抑制外界的剪切力阻止其進(jìn)一步的發(fā)展。

［1］ Meyer K A，Blewett P J．Numerical investigation of the stability of a shock－accelerated interface between two fluids［J］．Physics of Fluids，1972，15（5）：753－759．

［2］ Zhang Q，Sohn S I．An analytical nonlinear theory of Richtmyer－Meshkov instability［J］．Physics Letters A，1996，212（3）：149－155．

［3］ Anuchina N N，Volkov V I，Gordeychuk V A，et al．Numerical simulations of Rayleigh－Taylor and Richtmyer－Meshkov instability using MAH－3code［J］．Journal of Computational and Applied Mathematics，2004，168（1／2）：11－20．

［4］ Ruev G A，F(xiàn)edorov A V，F(xiàn)omin V M．Development of the Richtmyer－Meshkov instability upon interaction of a diffusion mixing layer of two gases with shock waves［J］．Journal of Applied Mechanics and Technical Physics，2005，46（3）：307－314．

［5］ Niederhaus J H J，Greenough J A，Oakley J G，et al．A computational parameter study for the three－dimensional shock－bubble interaction［J］．Journal of Fluid Mechanics，2008，594：85－124．

［6］ Thornber B，Drikakis D，Youngs D．Large－eddy simulation of multi－component compressible turbulent flows using high resolution methods［J］．Computers and Fluids，2008，37（7）：867－876．

［7］ Hejazialhosseini B，Rossinelli D，Bergdorf M，et al．High order finite volume methods on wavelet－adapted grids with local time－stepping on multicore architectures for the simulation of shock－bubble interactions．［J］．Journal of Computational Physics，2010，229（22）：8364－8383．

［8］ Schilling O，Latini M．High－order weno simulations of three－dimensional reshocked Richtmyer－Meshkov instability to late times：Dynamics，dependence on initial conditions，and comparisons to experimental data［J］．Acta Mathematica Scientia，2010，30（2）：595－620．

［9］ Tian B，Shen W，Jiang S，et al．A global arbitrary Lagrangian－Eulerian method for stratified Richtmyer－Meshkov instability［J］．Computers and Fluids，2011，46（1）：113－121．

［10］ Shankar S K，Kawai S，Lele S K．Two－dimensional viscous flow simulation of a shock accelerated heavy gas cylinder［J］．Physics of Fluids，2011，23（2）：131．

［11］ Bailie C，Mcfarland J A，Greenough J A，et al．Effect of incident shock wave strength on the decay of Richtmyer－Meshkov instability－introduced perturbations in the refracted shock wave［J］．Shock Waves，2012，22（6）：511－519．

［12］ Chandrasekhar S．Hydrodynamic and hydromagnetic stability［M］．Dover Publications，1961．

［13］ Wheatley V，Pullin D I，Samtaney R．Stability of an impulsively accelerated density interface in magnetohydrodynamics［J］．Physical Review Letters，2005，95（12）：125002．

［14］ Wheatley V，Samtaney R，Pullin D I．The magnetohydrodynamic Richtmyer－Meshkov instability：The transverse field case［C］∥18th Australasian Fluid Mechanics Conference．Australasian Fluid Mechanics Society，2012．

［15］ Khan M，Mandal L，Banerjee R，et al．Development of Richtmyer－Meshkov and Rayleigh－Taylor instability in presence of magnetic field［J］．Nuclear Instruments ＆Methods in Physics Research A，2011，653（1）：2－6．

［16］ Cao J，Wu Z，Ren H，et al．Effects of shear flow and transverse magnetic field on Richtmyer－Meshkov instability［J］．Physics of Plasmas，2008，15（4）：445－514．

［17］ Shin M S，Stone J M，Snyder G F．The magnetohydrodynamics of shock－cloud interaction in three dimensions［J］．Astrophysical Journal，2008，680（1）：336－348．

［18］ 李源，羅喜勝．黏性、表面張力和磁場對Rayleigh－Taylor不穩(wěn)定性氣泡演化影響的理論分析［J］．物理學(xué)報(bào)，2014，2（8）：277－285．Li Yuan，Luo Xisheng．Theoretical analysis of effects of viscosity，surface tension，and magnetic field on the bubble evolution of Rayleigh－Taylor instability［J］．Acta Physica Sinica，2014，2（8）：277－285．

［19］ Saltzman J．An unsplit 3Dupwind method for hyperbolic conservation laws［J］．Journal of Computational Physics，1994，115（1）：153－168．

［20］ Gardiner T A，Stone J M．An unsplit Godunov method for ideal MHD via constrained transport［J］．Journal of Computational Physics，2005，205（2）：509－539．

［21］ Haas J，Sturtevant B．Interaction of weak shock waves with cylindrical and spherical gas inhomogeneities［J］．Journal of Fluid Mechanics，1987，181（1）：41－76．

［22］ Chapman S，Cowling T G．The mathematical theory of non－uniform gases［M］．London：Cambridge University Press，1970．

［23］ 沙莎，陳志華，薛大文．激波沖擊R22重氣柱所導(dǎo)致的射流與混合研究［J］．物理學(xué)報(bào)，2013，62（14）：144701．Sha Sha，Chen Zhihua，Xue Dawen．The generation of jet and mixing induced by the interaction of shock wave with R22cylinder［J］．Acta Physica Sinica，2013，62（14）：144701．

Influence of magnetic field on interaction of shock wave
with R22heavy gas column

Lin Zhenya，Zhang Huanhao，Chen Zhihua，Liu Ying

（National Key Laboratory of Transient Physics，Nanjing University of Science and Technology，Nanjing210094，Jiangsu，China）

To study the interaction process of the plane incident shock wave with the magnetized R22 heavy gas column，we numerically simulated the deformation process of the shock－wave－induced R22 heavy gas column resulting from Kelvin－Helmholtz（KH）and Richtmyer－Meshkov（RM）instabilities under different initial conditions，and analyzed the jet focusing and inducing process by the transmitted shock wave．When the magnetic field was taken into consideration，the CTU＋CT algorithm satisfying the divergence equation of the magnetic field at any time was adopted in the numerical simulation．The results show that the magnetic field is capable of restraining the instability of the shockwave－induced R22gas column．Both the normal magnetic field（vertical to the flow direction）and the tangential magnetic field（parallel to the flow direction）can inhibit the RM instability．However，the restraining of the normal magnetic field is more effective than that of the tangential one with regard to the KH instability，as it can not only inhibit the vortex train rolling up on the interface but also prevent the bound vortex from developing．Besides，it is found that the magnetic field has little influence on the jet，and the magnetic energy at the jet point can suppress the jet attenuation to some extent while the normal magnetic field can reduce the peak pressure and velocity when the transmitted shock wave is focused．

magnetization gas；magneto－h(huán)ydro－dynamics；magnetic field；instability；shock wave

O354．5國標(biāo)學(xué)科代碼：1303520

A

10．11883／1001－1455（2017）04－0748－11

（責(zé)任編輯 王玉鋒）

2016－03－24；

2016－08－23

國家自然科學(xué)基金項(xiàng)目（11272156）；中國博士后科學(xué)基金項(xiàng)目（2015M571757）

林震亞（1990－ ），男，博士；通信作者：張煥好，122488989＠qq．com。