理想磁流體中激波與矩形密度界面相互作用的數(shù)值研究

李 源， 羅喜勝

（中國科學(xué)技術(shù)大學(xué)近代力學(xué)系先進(jìn)推進(jìn)實(shí)驗(yàn)室，合肥 230026）

理想磁流體中激波與矩形密度界面相互作用的數(shù)值研究

李 源， 羅喜勝

（中國科學(xué)技術(shù)大學(xué)近代力學(xué)系先進(jìn)推進(jìn)實(shí)驗(yàn)室，合肥 230026）

發(fā)展一套采用三階WENO格式和混合GLM方法的理想磁流體數(shù)值方法，并對(duì)激波與矩形密度界面相互作用進(jìn)行數(shù)值研究.通過圓極化阿爾芬波和旋轉(zhuǎn)激波管問題對(duì)數(shù)值方法的穩(wěn)定性和可靠性進(jìn)行驗(yàn)證.在入射激波馬赫數(shù)為10，界面內(nèi)外氣體密度比為10的情況，對(duì)比不同磁場(chǎng)中矩形密度界面的演變過程.結(jié)果表明，磁場(chǎng)能夠減少界面上渦量的生成從而抑制界面不穩(wěn)定性，并且磁場(chǎng)對(duì)界面的加速過程以及界面內(nèi)外氣體混合率有影響；而界面的存在將會(huì)使波后部分區(qū)域磁場(chǎng)增強(qiáng)；由于尖角的存在，矩形界面的發(fā)展與圓形界面不同.

磁流體；矩形界面；激波；不穩(wěn)定性

0 引言

激波沖擊被擾動(dòng)的物質(zhì)界面時(shí)，會(huì)產(chǎn)生Richtmyer-Meshkov（RM）不穩(wěn)定性［1－2］.在宇宙中，由于熱不穩(wěn)定性，星際介質(zhì)的密度分布并不均勻［3］，會(huì)形成大量的團(tuán)狀介質(zhì)和星云.大量的觀測(cè)表明，在超新星爆發(fā)、恒星風(fēng)和星云碰撞等天體現(xiàn)象中會(huì)產(chǎn)生在星際介質(zhì)中傳播的激波.激波與不均勻的星際介質(zhì)相互作用，會(huì)引發(fā)RM不穩(wěn)定性，使得星際介質(zhì)中能夠形成多相結(jié)構(gòu)，并有可能觸發(fā)恒星的形成［4－5］；另外，激波在不同相的星際介質(zhì)之間的質(zhì)量和能量的交換過程中也起到了重要作用.因此，了解激波與星際介質(zhì)相互作用的機(jī)理對(duì)研究星際介質(zhì)的結(jié)構(gòu)和演變有著重要意義.

激波與不均勻的星際介質(zhì)之間的相互作用一般用激波撞擊物質(zhì)界面的模型來刻畫.由于宇宙中存在磁場(chǎng)，并且部分星際介質(zhì)被宇宙射線電離，因此該模型可用磁流體方程描述.Wheatley和Samtaney等人［6－7］對(duì)不同方向磁場(chǎng)中激波與單模物質(zhì)界面相互作用進(jìn)行了理論和數(shù)值研究，得到磁場(chǎng)能夠抑制RM不穩(wěn)定性的結(jié)論，并且發(fā)現(xiàn)在入射激波和磁場(chǎng)很弱時(shí)，數(shù)值結(jié)果與不可壓線性理論解吻合很好；Samtaney等人［8］還考察了磁流體中激波與斜接觸面的相互作用，發(fā)現(xiàn)在磁場(chǎng)中，密度界面上產(chǎn)生的渦量會(huì)被輸運(yùn)到慢磁聲波或中間波上，阻止了渦量在密度界面上的積累，從而抑制了界面不穩(wěn)定性的發(fā)展；Sano等人［9］對(duì)二維磁流體方程組進(jìn)行了數(shù)值模擬，考察了不同磁場(chǎng)方向下，強(qiáng)激波與單模界面的相互作用，發(fā)現(xiàn)RM不穩(wěn)定性是一種穩(wěn)定的增強(qiáng)星際間磁場(chǎng)的機(jī)制，并且是超新星遺跡中激波位置處強(qiáng)磁場(chǎng)的形成原因；Mac Low等人［10］對(duì)平行磁場(chǎng)（與激波運(yùn)動(dòng)方向相同）中激波與二維氣體云的相互作用進(jìn)行了數(shù)值研究，得出了平行磁場(chǎng)能夠抑制RM不穩(wěn)定性和Kelvin-Helmholtz（KH）不穩(wěn)定性，并且磁場(chǎng)的存在能有效減少了云內(nèi)氣體與環(huán)境氣體之間的混合等重要結(jié)論；Fragile等人［11］則模擬了在輻射狀態(tài)下，激波與二維氣體云的相互作用，發(fā)現(xiàn)氣體云外部磁場(chǎng)會(huì)加強(qiáng)激波對(duì)氣體云的壓縮，并使其輻射變快，磁場(chǎng)變強(qiáng)，而氣體云內(nèi)部磁場(chǎng)則起到相反作用，在磁場(chǎng)很弱時(shí)，內(nèi)部磁場(chǎng)可以完全抑制住低溫氣體（T＜100 K）的輻射冷卻.

前人工作主要集中在磁場(chǎng)中單模界面、斜平面、圓形界面或球面受到激波沖擊后的演變和發(fā)展.這些界面形狀簡(jiǎn)單，適用于理論分析，但沒有考慮到星際介質(zhì)的不均勻性，密度界面形狀往往并不規(guī)則.而矩形界面或多邊形界面既包含有規(guī)則的分界同時(shí)又包含有尖角，是構(gòu)造復(fù)雜形狀界面的基本組成單元之一.范美茹等人［12］考慮了不同形狀界面與激波的相互作用，并得到界面形狀不同會(huì)引起流場(chǎng)波系結(jié)構(gòu)和渦量分布的不同的結(jié)論，但并沒有考慮磁場(chǎng)的影響.本文綜合考慮界面形狀與磁場(chǎng)的影響，對(duì)磁流體中矩形氣體云在強(qiáng)激波沖擊下的發(fā)展進(jìn)行數(shù)值模擬，并對(duì)過程中流場(chǎng)和磁場(chǎng)的變化進(jìn)行探討.

1 方法

忽略粘性，輻射冷卻、熱傳導(dǎo)以及自重等物理現(xiàn)象，則激波與氣體云相互作用過程可用二維理想磁流體方程組（IMHD）描述［13］：

其中ρ為密度，v為速度，B為磁場(chǎng)，pt＝p＋B2／2，ψ為標(biāo)量勢(shì)函數(shù)，E為總能

其中γ＝5／3.C是表示流場(chǎng)中各點(diǎn)含有界面內(nèi)部物質(zhì)多少的變量.初始時(shí)，在界面內(nèi)C＝1，界面外C＝0，隨著界面的演變，兩種物質(zhì)會(huì)發(fā)生混合，則有0≤C≤1.氣體云總質(zhì)量為

將0.1≤C≤0.9的部分認(rèn)為是混合氣體，則混合氣體中氣體云總質(zhì)量為

記氣體混合率fmix＝Mmix／Mcl，用來描述界面內(nèi)氣體與環(huán)境氣體的混合程度.類似［13］，對(duì)任一函數(shù)f，定義其質(zhì)量加權(quán)值為

并定義氣體云在x，y方向上的有效尺寸分別為δx，δy

同時(shí)，為了更好地理解氣體云被激波加速的過程，可引入激波運(yùn)動(dòng)方向上的氣體云平均速度〈u〉，其中u是x方向上的速度分量.

由于自然界中不存在磁單極子，即?·B＝0恒成立.而事實(shí)上，在MHD的數(shù)值模擬過程中，由于離散誤差、數(shù)值求解誤差等原因，該條件不一定能夠很好滿足，因此要對(duì)磁場(chǎng)散度進(jìn)行特殊處理.采取了混合GLM修正法［14］，在方程（4）中引入標(biāo)量勢(shì)函數(shù)ψ來保持?·B＝0成立，ψ滿足

綜上，控制方程寫為

其中，

類似［15］，將每個(gè)時(shí)間步的求解分為兩部分：①令S＝0，對(duì)方程（10）進(jìn)行數(shù)值離散求解；②通過積分，求出S，如（12）所示

S＝0時(shí)，數(shù)值離散方程（10），

其中，Lsi＋1／2和Rsi＋1／2分別為矩陣?F／?U的第s個(gè)左右特征向量

時(shí)間采用三階TVD Runge-Kutta格式

2 程序驗(yàn)證

其中，urefi，j是參考解.

2.1 圓極化阿爾芬波傳播問題

圓極化阿爾芬波傳播問題有理論解，因此可以用來考察程序處理連續(xù)流的能力，并且評(píng)估程序的計(jì)算精度.

設(shè)圓極化阿爾芬波的波矢k＝（kx，ky）＝（2π，4π），則該波的傳播方向與x軸的夾角α＝arctan（ky／kx）≈63.4°，文中計(jì)算區(qū)域?yàn)椋?，1］×［0，0.5］，網(wǎng)格數(shù)為Nx×Nx／2，采用周期邊界條件，γ＝5／3．初始條件為：ρ＝1，p＝0.1，v‖＝0，B‖＝1，v⊥＝B⊥＝0.1sinφ，vz＝Bz＝0.1cosφ，φ＝2π（xcosα

通過圓極化阿爾芬波（Circular Polarized Alfvén波）傳播問題和二維MHD激波管問題［17］對(duì)程序進(jìn)行驗(yàn)證，分別考察數(shù)值方法處理連續(xù)流和捕捉間斷流動(dòng)的能力，同時(shí)考察網(wǎng)格大小對(duì)計(jì)算結(jié)果的影響.對(duì)于Nx× Ny的網(wǎng)格，定義變量u的誤差L1的模為＋ysinα）；下標(biāo)‖表示沿波矢的分量，⊥表示垂直波矢的分量.可知，Alfvén波速vA＝B‖／ρ＝1，波的周期T＝2π／（｜k｜（vA＋v‖））＝5／5，在不同網(wǎng)格下對(duì)該問題進(jìn)行數(shù)值模擬，計(jì)算終止時(shí)間tend＝5T＝5．表1列出了Nx＝16，32，64，128，256時(shí)，磁場(chǎng)誤差值和精度分析；由表1可知，隨著網(wǎng)格增多，磁場(chǎng)誤差逐漸減小，并且程序的計(jì)算精度達(dá)到了預(yù)期的三階.圖1給出了不同網(wǎng)格下，沿波矢方向上的B⊥分布，可以看出，網(wǎng)格分辨率越高，數(shù)值解與理論解的偏差就越小，當(dāng)Nx≥64之后，數(shù)值解與理論解幾乎無差別.

表1 t＝5T時(shí)，圓極化阿爾芬波傳播問題在不同網(wǎng)格下的誤差和精度Table 1 Error and accuracy of circular polarized Alfvén wave propagation problem at t＝5T

圖1 t＝5T時(shí)，不同網(wǎng)格下，沿波矢方向上B⊥分布（數(shù)字表示Nx大小，ref為理論解.）Fig.1 B⊥along direction of wave propagation with different grids at t＝5T（ref means analytical solution.）

2.2 二維MHD激波管問題

二維MHD激波管問題是將一維MHD激波管問題旋轉(zhuǎn)后得到.該問題除了可以考察程序捕捉間斷流能力外，可以用來考察程序能否很好地保持斜方向流動(dòng)與正方向流動(dòng)的數(shù)值解的對(duì)稱性.計(jì)算區(qū)域?yàn)椋郏?，1］×［－0.01，0.01］，網(wǎng)格數(shù)為Nx×Nx／100，令Nx＝600，采用出口邊界條件，γ＝2.初始條件為：x‖＜0時(shí)，（ρ，v‖，v⊥，B‖，B⊥，p）＝（1，0，0，0.75，1，1），x‖≥0時(shí)，（ρ，v‖，v⊥，B‖，B⊥，p）＝（0.125，0，0，0.75，－1，0.1），其中x‖＝xcosα＋ysinα，α＝π／4．圖2給出了t＝2／10時(shí)，在激波運(yùn)動(dòng)方向上的ρ，B‖和B⊥分布，并與參考解對(duì)比.參考解是一維問題在網(wǎng)格數(shù)1 024時(shí)的數(shù)值解.由密度圖可知，MHD波從左到右分別為膨脹波、復(fù)合波、接觸間斷、激波、膨脹波，可以看到，程序準(zhǔn)確地捕捉到了所有的間斷并且與參考解一致；而B‖則有微小震蕩，但最大相對(duì)誤差不超過0.25%.

圖2 t＝2／10時(shí)，沿激波運(yùn)動(dòng)方向上的物理量分布（從左到右ρ，B‖，B⊥，圓符號(hào)為二維解，實(shí)線為一維參考解.）Fig.2 Profiles of ρ，B‖，B⊥along direction of shock propagation at t＝2／10（Circle symbol represents solution of 2D computations whereas solid line is reference solution.）

3 數(shù)值結(jié)果和討論

3.1 問題描述

考慮一道平面激波沖擊一個(gè)矩形氣體云界面.如圖3所示，計(jì)算區(qū)域?yàn)閤∈［－0.5，0.5］，y∈［－0.5，2.5］，激波從左向右傳播，左邊界為超音速來流，上、下及右邊界為出口；入射激波強(qiáng)度Ms＝10；界面中心在原點(diǎn)處，尺寸為a＝0.1，b＝0.15；界面內(nèi)為重氣體，壓強(qiáng)與環(huán)境壓強(qiáng)相等，密度比χ＝10；波前環(huán)境壓強(qiáng)P＝1，聲速Cs＝1；氣體比熱比γ＝5／3．根據(jù)程序驗(yàn)證部分網(wǎng)格收斂性研究的結(jié)論，計(jì)算中網(wǎng)格大小取為800×2 400，可以很好保證計(jì)算結(jié)果的精度.

對(duì)該問題在不同強(qiáng)度的磁場(chǎng)中進(jìn)行了模擬，磁場(chǎng)強(qiáng)度用波前環(huán)境氣體的無量綱值β＝2P／B2描述，選用三種β，分別為1、10和∞；另外，文中還考慮兩種不同的初始磁場(chǎng)方向：①初始磁場(chǎng)與激波運(yùn)動(dòng)方向平行，即x方向磁場(chǎng)，下稱平行磁場(chǎng)（BX）；②初始磁場(chǎng)與激波運(yùn)動(dòng)方向垂直，即y方向磁場(chǎng)，下稱垂直磁場(chǎng)（BY）.所有算例參數(shù)見表2.

圖3 初始設(shè)定（上）和初始?xì)怏w云形狀（下）Fig.3 Initial conditions（top）and initial shape of the cloud（bottom）

表2 算例參數(shù)Table 2 Initial parameters of simulation cases

入射激波與界面接觸后，會(huì)產(chǎn)生透射激波，將透射激波通過氣體云所需的特征時(shí)間記為tcc［18］，則

根據(jù)初始條件，tcc≈0.047 4，設(shè)激波與界面開始接觸的時(shí)刻為0.

3.2 界面變形

圖4列舉了t／tcc＝1.37時(shí)，所有算例的log［（?ρ／?x）2＋（?ρ／?y）2＋1］灰度圖.從圖4（a）中可以看出，激波與界面相互作用后，不穩(wěn)定性會(huì)產(chǎn)生.而由圖5可知，磁場(chǎng)能夠抑制界面不穩(wěn)定性.初始磁場(chǎng)為平行磁場(chǎng)時(shí)，在界面發(fā)展初期，磁場(chǎng)對(duì)界面不穩(wěn)定性的抑制并不明顯，如圖4（b），（d）所示，磁場(chǎng)中的氣體云形狀與無磁場(chǎng)時(shí)差別不大.相對(duì)于平行磁場(chǎng)，垂直磁場(chǎng)對(duì)不穩(wěn)定性的抑制更為強(qiáng)烈，見圖4（c）.而由圖4（e）可知，t＝1．37tcc時(shí)，弱垂直磁場(chǎng)對(duì)界面發(fā)展已經(jīng)有一定的影響.雖然在界面發(fā)展初期，平行磁場(chǎng)對(duì)界面外形影響不大，但到了發(fā)展后期，平行磁場(chǎng)也能夠很大程度地影響界面外形，如圖5所示.

Chandrasekhar的線性理論［19］指出，當(dāng)?shù)谹lfvén波速超過界面兩邊速度差，即β＜2／M2s時(shí)，KH不穩(wěn)定性被抑制；當(dāng)Alfvén波穿越時(shí)間小于加速時(shí)間尺度，即β＜1．2（χ／Ms）2時(shí)，Rayleigh-Taylor（RT）不穩(wěn)定性被抑制.本文中χ＝Ms＝10，可知KH不穩(wěn)定性只有在很強(qiáng)磁場(chǎng)（β＜0．02）中才被抑制；而只要β＜1．2時(shí)，RT不穩(wěn)定性就被抑制.該結(jié)論與圖5的計(jì)算結(jié)果符合.

3.3 有效長(zhǎng)度

為了更好地理解界面的演變過程，對(duì)氣體云的有效尺寸進(jìn)行討論.圖6給出了δx／δx0和δy／δy0隨時(shí)間的演變過程，其中δx0、δy0分別為初始時(shí)刻的縱向和橫向有效長(zhǎng)度.當(dāng)t／tcc＜1時(shí)，氣體云被入射激波壓縮，δx急劇下降到不足初始時(shí)刻的40%；與δx不同，在這個(gè)時(shí)間內(nèi)，橫向有效尺寸δy幾乎沒有改變（不包括強(qiáng)垂直磁場(chǎng)）；之后，氣體云會(huì)再膨脹，在這個(gè)過程中，界面不穩(wěn)定性會(huì)迅速增長(zhǎng)，最終導(dǎo)致氣體云破碎摧毀［18］.再膨脹過程能夠被磁場(chǎng)限制，如圖6所示，雖然弱平行磁場(chǎng)很難限制住氣體云的再膨脹過程，但強(qiáng)平行磁場(chǎng)在3tcc后，開始明顯地限制住氣體云δx、δy的增長(zhǎng)；與平行磁場(chǎng)不同，即使是弱垂直磁場(chǎng)，也能夠明顯地限制氣體云的膨脹.

不同方向的磁場(chǎng)都能夠限制氣體云的再膨脹過程，但兩者機(jī)制不同.激波撞擊氣體云后，界面上會(huì)產(chǎn)生渦量，在沒有磁場(chǎng)存在時(shí)，渦會(huì)帶著界面內(nèi)氣體四處擴(kuò)散，使得氣體云體積增大，平行磁場(chǎng)可以減少界面上渦量的生成，從而限制住氣體云膨脹；相比而言，垂直磁場(chǎng)不僅可以減少界面上渦量的生成，事實(shí)上，由于是二維模擬，磁力線無法繞過氣體云而只能夠包裹住氣體云前沿，隨著這些磁力線的累積，界面將被強(qiáng)烈地壓縮，使得再膨脹過程受限.

圖4 t＝1．37tcc時(shí)，算例（a）N，（b）BX1，（c）BY1，（d）BX10，（e）BY10的log［（?ρ／?x）2＋（?ρ／?y）2＋1］灰度Fig.4 Gray-scale plots of log［（?ρ／?x）2＋（?ρ／?y）2＋1］for（a）N，（b）BX1，（c）BY1，（d）BX10，（e）BY10 at t＝1.37 tcc

圖5 t＝5．06 tcc時(shí)，算例N，BX10，BX1，BY10，BY1（從左往右）的log（ρC＋1）灰度Fig.5 Gray-scale plots of log（ρC＋1）for N，BX10，BX1，BY10，BY1（from left to right）at t＝5.06 tcc

3.4 平均速度和氣體混合率

圖7（左）給出了氣體云〈u〉／Cs隨時(shí)間的變化過程，其中Cs是初始環(huán)境氣體聲速.氣體云與激波相互作用后會(huì)被加速，Klein［18］將其分為兩個(gè)過程：①激波在氣體云內(nèi)傳播時(shí)，界面內(nèi)氣體會(huì)被激波加速；②當(dāng)激波通過氣體云后，氣體云會(huì)被波后氣流“拖著”向前走.在第一個(gè)過程中，由于氣體云受磁場(chǎng)影響很小，加速過程基本一致，如圖7（左）所示；而第二個(gè)過程與氣體云橫向有效尺寸δy相關(guān)，δy越小，說明氣體云能被波后氣體牽引的區(qū)域越小，導(dǎo)致氣體云加速度越小，在發(fā)展后期，氣體云的速度就越小.由圖7（左）可知，平行磁場(chǎng)越強(qiáng)，后期氣體云平均速度越小，與上述討論一致；但在垂直磁場(chǎng)中，界面前沿累積了大量磁力線，由于磁力線上的張力會(huì)給界面內(nèi)氣體一個(gè)加速度，并且磁場(chǎng)越強(qiáng)，這個(gè)加速度就越大，所以盡管垂直磁場(chǎng)中的δy很小，氣體云還是能夠被很快地加速.

圖6 δx／δx0（左）和δy／δy0（右）隨時(shí)間變化Fig.6 Time histories of δx／δx0（left）and δy／δy0（right）

圖7 〈u〉／Cs（左）和fmix（右）隨時(shí)間變化Fig.7 Time histories of〈u〉／Cs（left）and fmix（right）

激波與界面相互作用后，由于斜壓，界面上會(huì)產(chǎn)生渦量并向外擴(kuò)散，使得界面內(nèi)外氣體相混合.故fmix與界面不穩(wěn)定性相關(guān)，界面不穩(wěn)定性越強(qiáng)，氣體混合率越大.圖7（右）列出了所有算例的fmix隨時(shí)間變化.可知，磁場(chǎng)的存在有效地降低了fmix，即磁場(chǎng)一定程度上抑制了界面不穩(wěn)定性，并且磁場(chǎng)強(qiáng)度越強(qiáng)，這種抑制就越強(qiáng)烈；對(duì)同種強(qiáng)度的磁場(chǎng)而言，平行磁場(chǎng)對(duì)界面內(nèi)外氣體混合過程的抑制不如垂直磁場(chǎng)強(qiáng)烈.而由計(jì)算結(jié)果可知，盡管磁場(chǎng)可以降低fmix，但并不能徹底阻止氣體之間的混合.

3.5 磁場(chǎng)變化

流場(chǎng)中沒有界面存在時(shí)，入射激波對(duì)熱力學(xué)壓強(qiáng)的增幅為M2s.初始磁場(chǎng)為平行磁場(chǎng)時(shí)，波后磁場(chǎng)不變.初始磁場(chǎng)為垂直磁場(chǎng)時(shí)，若Ms→∞，波后磁場(chǎng)的增幅為（γ＋1）／（γ－1）；故波后磁壓強(qiáng)遠(yuǎn)小于熱力學(xué)壓強(qiáng)，即β→∞．當(dāng)流場(chǎng)中存在界面時(shí)，在波后某些區(qū)域，磁壓強(qiáng)可以增強(qiáng)到與熱力學(xué)壓強(qiáng)同一個(gè)量級(jí)，即β～1．

圖8給出了t＝5．06 tcc時(shí)，算例BX10和BY10的log β灰度.初始磁場(chǎng)為平行磁場(chǎng)時(shí)，最小的β出現(xiàn)在氣體云下風(fēng)處的磁場(chǎng)通量管［11］處，如圖8（a）所示.激波撞擊到氣體云后，環(huán)境氣體中的激波迅速掃過氣體云上下界面，并在氣體云背面對(duì)稱軸處相交，形成一個(gè)馬赫桿，見圖4.經(jīng)馬赫桿加速后的氣體速度遠(yuǎn)大于周圍氣體速度，隨著該部分氣體向前運(yùn)動(dòng)，氣體原位置處會(huì)形成一個(gè)稀薄空間，周圍氣體只能通過壓縮磁力線來填滿這個(gè)空間，從而形成磁場(chǎng)通量管，并使得通量管中的磁場(chǎng)增強(qiáng).除此之外，在氣體云界面有渦量處，磁場(chǎng)也會(huì)被放大.初始磁場(chǎng)為垂直磁場(chǎng)時(shí)，最小的β發(fā)生在界面前沿處，如圖8（b）所示.正如前面所討論的，磁力線會(huì)在界面前沿處累積，導(dǎo)致磁場(chǎng)增強(qiáng).與平行磁場(chǎng)不同的是，在氣體云背面并沒有形成單一的磁場(chǎng)通量管，而是兩條平行的絲，并且氣體云下風(fēng)處對(duì)稱面上形成了強(qiáng)電流層.

圖8 t＝5．06tcc時(shí)，算例（a）BX10，（b）BY10的log β灰度圖（從β＝0．1（黑）到β＝200（白））Fig.8 Gray-scale plots of log β for（a）BX10 and（b）BY10 at t＝5.06 tcc（from β＝0．1（black）to β＝200（white））

4 總結(jié)

對(duì)二維理想磁流體中強(qiáng)激波與矩形界面相互作用進(jìn)行數(shù)值模擬，在入射激波Ms＝10，界面內(nèi)外氣體密度比χ＝10的前提下，考慮不同的初始磁場(chǎng)強(qiáng)度和方向，并與無磁場(chǎng)情況對(duì)比，得到以下結(jié)果.

1）磁場(chǎng)能夠抑制界面不穩(wěn)定性，磁場(chǎng)強(qiáng)度越強(qiáng)，抑制越明顯.初始磁場(chǎng)方向不同，抑制的機(jī)制也不同.平行磁場(chǎng)通過減少界面上渦量的生成減緩不穩(wěn)定性的增長(zhǎng)；垂直磁場(chǎng)則通過包裹界面的磁力線來壓縮界面，并防止界面上渦量的生成和擴(kuò)散，從而抑制界面不穩(wěn)定性.

2）由于界面上渦量的生成和擴(kuò)散被磁場(chǎng)抑制，使得界面內(nèi)外氣體混合率下降，從而延長(zhǎng)了氣體云破碎時(shí)間.初始磁場(chǎng)越強(qiáng)，氣體混合率越低.磁場(chǎng)會(huì)影響氣體云加速過程，平行磁場(chǎng)中，磁場(chǎng)越強(qiáng)，氣體云加速越慢；垂直磁場(chǎng)中，由于界面前沿磁力線張力為界面內(nèi)氣體提供了額外的加速度，使得磁場(chǎng)越強(qiáng)，氣體云加速越快.

3）界面的存在使得波后某些區(qū)域的磁場(chǎng)增強(qiáng)，平行磁場(chǎng)中，這些區(qū)域主要集中在磁場(chǎng)通量管處；垂直磁場(chǎng)中，這些區(qū)域則集中在界面前沿.

［1］ Richtmyer R D.Taylor instability in a shock acceleration of compressible fluids［J］.Communication on Pure and Applied Mathematics，1960，13（2）：297－319.

［2］ Meshkov E E.Instability of the interface of two gases accelerated by a shock wave［J］.Fluid Dynamics，1969，4（5）：101－104.

［3］ Koyama H，Inutsuka S.An origin of supersonic motions in interstellar clouds［J］.The Astrophysical Journal，2002，564：97－100.

［4］ McKee C F，Ostriker J P.A theory of the interstellar medium：Three components regulated by supernova explosions in an inhomogeneous substrate［J］.The Astrophysical Journal，1977，218：148－169.

［5］ Elmegreen B G，Lada C J.Sequential formation of subgroups in OB associations［J］.The Astrophysical Journal，1977，214：725－741.

［6］ Wheatley V，Pullin D I，Samtaney R.Stability of an impulsively accelerated density interface in magnetohydrodynamics［J］. Physical Review Letter，2005，95：125002.

［7］ Wheatley V，Samtaney R，Pullin D I.The magnetohydrodynamics Richtmyer-Meshkov instability：The transverse field case［C］.Proceedings of the 18th Australian Fluid Mechanics Conference，Launceston，2012.

［8］ Samtaney R.Suppression of the Richtmyer-Meshkov instability in the presence of a magnetic field［J］.Physics of Fluids，2003，15（8）：53－56.

［9］ Sano T，Nishihara K，Matsuoka C，et al.Magnetic field amplification associated with the Richtmyer-Meshkov instability［J］. The Astrophysical Journal，2012，758：126－138.

［10］ 范美茹，翟志剛，司廷，等.激波與不同形狀界面作用的數(shù)值模擬［J］.中國科學(xué)：物理學(xué)，力學(xué)，天文學(xué)，2011，41（7）：862－869.

［11］ Mac Low M M，McKee C F，Klein R I，et al.Shock interactions with magnetized interstellar clouds.I.Steady shocks hitting nonradiative clouds［J］.The Astrophysical Journal，1994，433：757－777.

［12］ Fragile P C，Anninos P，Gustafson K，et al.Magnetohydrodynamic simulation of shock interactions with radiative clouds［J］. The Astrophysical Journal，2005，619：327－339.

［13］ Shin M S，James M S，Gregory F S.The magnetohydrodynamics of shock-cloud interaction in three dimensions［J］.The Astrophysical Journal，2008，680：336－348.

［14］ Dender A，F(xiàn)emm F，Kroner D，et al.Hyperbolic divergence cleaning for the MHD equations［J］.Journal of Computational Physics，2002，175：645－673.

［15］ Mignone A，Tzeferacos P.A second-order unsplit Godunov scheme for cell-centered MHD：The CTU-GLM scheme［J］. Journal of Computational Physics，2010，229（6）：2117－2138.

［16］ Mignone A，Tzeferacos P，Bodo G.High-order conservative finite difference GLM-MHD schemes for cell-centered MHD［J］. Journal of Computational Physics，2010，229（17）：5896－5920.

［17］ Yee H C，Sj?green B.Efficient low dissipative high order schemes for multiscale MHD flows，II：Minimization of?·B numerical error［J］.Journal of Scientific Computing，2006，29（1）：559－575.

［18］ Klein R I，Mckee C F，Colella P.On the hydrodynamic interaction of shock waves with interstellar clouds.I.Nonradiative shocks in small clouds［J］.The Astrophysical Journal，1994，420：213－236.

［19］ Chandrasekhar S.Hydrodynamic and hydromagnetic stability［M］.Oxford：Clarendon Press，1961：428－511.

Numerical Study of Shock Interactions with Rectangular Density Interface in Magnetohydrodynamics

LI Yuan，LUO Xisheng
（Advanced Propulsion Laboratory，University of Science and Technology of China，Hefei 230026，China）

A magnetohydrodynamic simulation method is developed for shock interacting with rectangular density interface in a magnetic field.The method employs 3rd WENO scheme and mixed GLM method.With circular polarized Alfvén wave propagation test and rotated shock tube problem，the method is validated.Under conditions that Mach number of shock is 10 and ratio of density of cloud to ambient gas is 10，evolutions of shocked interface in different initial orientations and strengths of magnetic field are compared. It shows that magnetic field decreases vorticity formed on surface and reduces growth of hydrodynamic instabilities；Field influences acceleration and mixing rate of cloud；And field is greatly amplified in some regions behind shock when cloud is presented.It is also found that，due to sharp corner，evolution of rectangular interface is different from circular one.

magnetohydrodynamics；rectangular interface；shock wave；instability

date： 2013－10－28；Revised date： 2014－02－01

O361.3

A

2013－10－28；

2014－02－01

李源（1990－），男，碩士生，從事磁流體界面不穩(wěn)定性研究，E-mail：yuli＠m(xù)ail.ustc.edu.cn通訊作者：羅喜勝，E-mail：xluo＠ustc.edu.cn

1001-246X（2014）06-0659-09