趙凱歌 薛創(chuàng) 王立鋒3) 葉文華3) 吳俊峰 丁永坤3) 張維巖3) 賀賢土3)

1)(中國工程物理研究院研究生院,北京 100088)

2)(北京應用物理與計算數(shù)學研究所,北京 100094)

3)(北京大學應用物理與技術研究中心,高能量密度物理數(shù)值模擬教育部重點實驗室,北京 100871)

1 引 言

當兩層流體的界面兩側具有連續(xù)的法向速度與壓力,而其他物理量(如:密度、內(nèi)能)不連續(xù)時,稱為接觸間斷[1].在有勢力場中流體密度梯度與壓力梯度不平行時,那么接觸間斷面上的擾動幅度將被激發(fā)并增長,稱為瑞利-泰勒不穩(wěn)定性(Rayleigh-Taylor instability,RTI)[2?4].RTI是慣性約束聚變(inertial con finement fusion,ICF)研究的重要問題之一[5?8],也是許多天體物理現(xiàn)象的重要過程[9?11],例如:超新星爆炸和星系演化、天體射流等現(xiàn)象[12].在ICF的內(nèi)爆加速和減速階段中均引起RTI,RTI的發(fā)展影響點火熱斑的形成、聚變?nèi)紵湍芰吭鲆鎇13?22].超新星爆炸與ICF在動力學上十分相似,存在加速度引起的RTI.因此,研究RTI的物理機理對于天體物理中各類不穩(wěn)定流動的演化具有十分重要的意義[12,23].

通過對界面不穩(wěn)定性的理論研究[13?16,24?28]、數(shù)值模擬[17,18,29,30]及實驗分析[31?34],將RTI的演化過程劃分為線性階段、變形階段、規(guī)則非線性、不規(guī)則非線性及湍流混合階段.理論研究方法大致分兩類.一類是采用速度勢及界面函數(shù)的微擾展開方法分析RTI的線性和弱非線性階段[4,35,36].在線性階段,理論分析與數(shù)值模擬和實驗結果符合得很好[4,24,28,31].Jacobs和Catton[31]提出真空和流體界面(Atwood數(shù)等于1)的弱非線性(weakly non-linear,WNL)模型;Wang等[13,35,36]將WNL模型推廣適用于任意Atwood數(shù)的情況,很好地描述了弱非線性階段氣泡和尖釘?shù)漠a(chǎn)生過程.由于高次諧波的激發(fā),導致小擾動展至高階時將面對繁瑣的計算.另一類是研究非線性階段的氣泡模型[37?43];Layzer[38]提出了對于真空和流體界面情況在氣泡頂點附近采用勢流理論的模型,描述了氣泡在整個擾動過程中的發(fā)展行為.基于Layzer模型,Zhang[39]對任意時刻氣泡和尖釘?shù)倪\動規(guī)律進行了研究,Goncharov[40]和Sohn[41]將Layzer模型推廣應用于任意Atwood數(shù)下的速度勢模型.對于氣泡和尖釘?shù)拇蟪叨认喔蛇\動,Abarzhi等[42]找到了多重調(diào)和解.Mikaelian[43]給定初始幅值,可以描述任意Atwood數(shù)的界面幅值關系式.陶燁晟等[44]將Layzer的氣泡模型推廣,描述了氣泡由初始線性階段到氣泡以漸近速度增長的非線性階段的發(fā)展過程.然而,Layzer模型僅能描述氣泡和尖釘?shù)捻旤c在擾動演變過程中的運動.

1972年,Ott[45]首先提出非線性RTI的薄層理論,通過分析薄層質(zhì)點的受力,解釋氣泡和尖釘?shù)男纬蛇^程[46,47].然而,Ott理論忽略了薄層兩側流體的分布.此外,當演變時間超過時,在波谷處出現(xiàn)界面卷曲現(xiàn)象,這里k=2π/λ表示波數(shù),λ代表擾動波長,η0是界面擾動的初始幅值,g為重力加速度.Wang等[48]采用WNL模型對具有一定厚度的流體薄層中的RTI進行研究,詳細分析薄層上下界面分別具有初始擾動狀態(tài)時的演變規(guī)律.研究發(fā)現(xiàn)薄層厚度在擾動演變過程中起重要的作用.然而,在實際的流體不穩(wěn)定性中,界面兩側的流體分布范圍比較寬,近似經(jīng)典RTI,而發(fā)生于界面附近的擾動按擾動波長成e指數(shù)形式衰減.當流體層厚度大于數(shù)倍的擾動波長時,在邊界處的擾動已經(jīng)衰減得非常小,因此,邊界位置對界面處的擾動造成的影響可以忽略不計.本文將Ott的薄層模型應用于經(jīng)典RTI界面變形演化過程.

本文在第2部分分析薄層模型的理論分析,給出界面的運動方程組;在第3部分首先數(shù)值求解運動方程組,其次分析薄層模型在RTI中的線性和弱非線性段的增長規(guī)律,再將模型的界面形變與數(shù)值模擬結果進行對比,最后將薄層模型擴展應用于分析擾動分布呈三角波和方波的情形;第4部分為本文的總結.

2 薄層RTI的理論模型

理想二維流體處于重力場(?gey)中,在接觸間斷的界面兩側為有限厚度層流體,其上下兩層的流體密度分別為ρ2和ρ1(ρ2>ρ1).在界面附近選擇一定厚度(h0)的流體薄層,平衡狀態(tài)時,界面處的壓強用p0表示,如圖1(a).已知在重力場中壓力近似由靜壓條件表述,其中p0表示平衡狀態(tài)初始位置處的壓強.當薄層受到擾動后,尖釘處壓力增大,氣泡處壓力減小,流體薄層上下端面的壓強分別用p2,p1表示,如圖1(b).

圖1 流體薄層不穩(wěn)定性的示意圖 (a)平衡狀態(tài);(b)擾動狀態(tài)Fig.1.The RTI of the thin layer:(a)Planar interface in equilibrium;(b)perturbed interface.

平衡狀態(tài)時,選取流體薄層中的一點x=ξ,y=0,ξ是拉格朗日坐標;在t時刻,該點位于r(ξ,t)=x(ξ,t)ex+y(ξ,t)ey.考慮 平衡狀態(tài)時另一點x=ξ+dξ,y=0,該點在t時刻位于r+(?r/?ξ)dξ處.研究分析ξ和ξ+dξ之間流體微團的運動規(guī)律,則流體微團的質(zhì)量dm=(ρ1+ρ2)h0dξ/2=(ρ1+ρ2)h|dr|/2,其中h和|dr|表示擾動過程中流體微團的厚度和寬度.假設微團端面的壓強不受擾動的影響,則上端面和下端面的壓強分別為p2=p0?ρ2g(y+hey·n/2)和p1=p0?ρ1g(y?hey·n/2),如圖1(b),則流體微團的受力關系:

其中,法向量則ey·方程(1)中?gdmey表示流體微團在重力場中的力,(p1?p2)|dr|n表示流場的壓力場作用.由于微團的受力項滿足dF=dm?2r/?t2,則運動方程的分量形式:

AT=(ρ1?ρ2)/(ρ1+ρ2)稱為Atwood數(shù),h0表示界面兩側流體層的初始厚度.方程(2)描述了流體界面RTI的形變規(guī)律,即形變的擾動范圍由線性階段到非線性階段的演化過程.

在擾動初始階段,x和y方向的位置坐標x(ξ,t)=ξ+ δx(ξ,t)和y(ξ,t)=δy(ξ,t), 其中δx(ξ,t)和δy(ξ,t)表示擾動小量.將位置坐標代入方程(2)中保留一階小量,則y方向的擾動演化關系表示為

考慮擾動隨時間呈exp(γt)形式,則增長率γ=

3 結果與討論

3.1 數(shù)值計算方法

方程組(2)為非線性方程組,直接求解比較困難.因此,對方程組采用數(shù)值求解,將時間導數(shù)和空間導數(shù)分別采用隱式中心差分格式,?t和?ξ分別表示時間步長和空間步長.差分后的形式為

則方程組(4)差分后簡為寫如下的方程組:

式中,A表示系數(shù)矩陣,表示未知量矩陣,B為常數(shù)矩陣.上角標n和下角標j分別表示時間變量和空間變量.給定初始狀態(tài),采用Picard迭代求解方程組的結果.

3.2 線性及弱非線性階段

首先,選擇初始擾動的形式y(tǒng)=η0cos(kx),其中η0為初始擾動振幅,k=2π/λ為波數(shù),λ為擾動波長.固定初始擾動幅值η0=0.02μm和加速度g=1μm/ns2,選擇不同的薄層厚度h0,分析薄層模型的線性增長率關于Atwood數(shù)AT和波數(shù)k的變化關系,如圖2所示.已知經(jīng)典RTI在線性段的幅值增長稱為線性增長率.

在薄層模型中分別選擇薄層厚度kh0=1,2,4的情況分析界面不穩(wěn)定性的演變規(guī)律,如圖2所示,將薄層模型的線性增長率與經(jīng)典RTI的結果進行比較.首先,在圖2(a)中固定界面擾動的波長λ=1μm,對比薄層模型和經(jīng)典RTI的線性增長率隨AT數(shù)的變化規(guī)律.其次,在圖2(b)中保持AT=0.6不變,對比薄層模型和經(jīng)典RTI的線性增長率隨擾動波數(shù)k的變化規(guī)律.通過圖2(a)和圖2(b)可以看出,薄層模型的線性增長率依賴于薄層厚度,已知界面處的擾動隨距離界面的寬度呈指數(shù)衰減.當薄層厚度比較薄(kh0=1)時,薄層模型的線性增長率明顯大于經(jīng)典RTI的結果,導致薄層的界面演化加快,擾動增長趨于強烈;當薄層厚度比較厚(kh0=4)時,薄層模型的線性增長率明顯小于經(jīng)典RTI的結果,以致界面演化緩慢,擾動增長趨于舒緩;當且僅當薄層厚度滿足條件kh0=2時,薄層模型的線性增長率隨AT和k的變化規(guī)律與經(jīng)典RTI的結果保持一致.因此,采用薄層模型描述界面的擾動發(fā)展規(guī)律時,應選擇適當?shù)谋雍穸扰c經(jīng)典RTI的結果做比較.這與線性化分析的結果一致,在方程(3)中給出的增長率形式當薄層厚度滿足h0=2/k時,增長率與經(jīng)典RTI結果一致.

圖2 選取不同的薄層厚度,比較薄層理論和經(jīng)典RTI的線性增長率 (a)增長率關于AT數(shù)的關系;(b)增長率關于波數(shù)k的關系Fig.2.Comparision between linear growth rates obtained by the thin layer model and classical RTI at different thin layer thicknesses:(a)The function relation between the linear growth rate and the Atwood number;(b)the function relation between the linear growth rate and the wavenumber.

通過對比線性階段的增長特征,可知在合適的薄層厚度(如h0=2/k)下,可以將模型的界面幅值與WNL模型做比較.已知三階WNL模型給出關于RTI的弱非線性階段界面為

這里η1(x,t),η2(x,t)和η3(x,t)分別表示基模、二次諧波和三次諧波界面,其形式如下:

式中ηL=η0cosh(γct)代表基模的線性增長幅值,η0表示擾動的初始幅值,表示經(jīng)典RTI的線性增長率.

選取初始擾動幅值η0=0.02μm和加速度g=1μm/ns2,分別選擇三組不同的Atwood數(shù)和擾動波數(shù),將薄層模型的界面幅值之和ηb+s=ηb+ηs與WNL模型做對比,如圖3所示,這里ηb,ηs分別表示氣泡和尖釘?shù)姆?圖3(a)中,固定波長λ=5μm,選取AT=0.1,0.4,0.9;圖3(b)中,固定AT=0.8,選取k=0.1,0.2,0.4.

圖3 固定η0=0.02μm和g=1μm/ns2,比較薄層模型和WNL模型的幅值關系 (a)波長λ=5μm,選取AT=0.1,0.5和0.9;(b)AT=0.8,選取k=0.1,0.2和0.4Fig.3.Comparision between the amplitudes obtained by the thin layer model and WNL model for RTI with fixed η0=0.02 μm and g=1 μm/ns2:(a)Choosing AT=0.1,0.5,and 0.9 with the wavelength λ=5 μm;(b)choosing k=0.1,0.2,and 0.4 with the Atwood number AT=0.8.

在圖3中,將薄層模型的幅值之和與WNL模型的結果做比較,圖3(a)中選取不同的AT數(shù),圖中曲線和圓點分別代表WNL模型和薄層模型的結果;圖3(b)中選取不同的波數(shù)k,圖中曲線和方塊分別代表WNL模型和薄層模型的結果.可以看出,薄層模型與WNL模型在線性段及弱非線性區(qū)域內(nèi)符合得相當好,說明薄層模型可以準確描述不穩(wěn)定性的弱非線性段.然而,二者也存在一定的差別.在圖3(a)中,當AT數(shù)比較小(AT=0.1)時,薄層模型僅能夠描述初始線性段的擾動,薄層模型描述至t=6.4 ns時刻,隨后的擾動演變規(guī)律不能準確給出,這是由于AT比較小,數(shù)值計算中導致誤差較大;當AT數(shù)比較大(例如AT=0.9)時,WNL模型可以準確描述至幅值ηb+s=1.1μm的情況,而薄層模型的描述范圍遠大于WNL模型.在圖3(b)可明顯看出,薄層模型在描述不穩(wěn)定性的演變規(guī)律時其適用性大于WNL模型.

3.3 變形與非線性演化

為進一步說明薄層模型的理論,將界面的變形和非線性演化過程與數(shù)值模擬的結果做對比.在數(shù)值模擬中,其上下層的流體密度分別為擾動波長λ=20μm,初始壓強p0=5 MPa,加速度g=17μm/ns2,氣體絕熱指數(shù)γh=5/3.在下面的模擬中,這些參數(shù)保持不變.

首先,對比不穩(wěn)定性在線性階段的界面分布規(guī)律.界面初始幅值η0=0.075μm,圖4給出了不穩(wěn)定性在線性階段的分布,圖4(a)—(d)和(a′)—(d′)分別表示薄層模型和數(shù)值模擬所描述的在不同時刻的不穩(wěn)定性界面.通過對比發(fā)現(xiàn),薄層模型在描述擾動演化的線性階段,界面呈現(xiàn)規(guī)則的余弦形狀;隨擾動的演化,界面的波峰向上增長,波谷向下發(fā)展,發(fā)現(xiàn)薄層模型的不穩(wěn)定界面與數(shù)值模擬的結果在線性階段的界面增長基本一致.

圖4 (a)—(d)薄層模型和(a′)—(d′)數(shù)值模擬所描述的在線性階段的不穩(wěn)定界面 (a),(a′)t=0.0 ns;(b),(b′)t=0.4 ns;(c),(c′)t=0.8 ns;(d),(d′)t=1.0 nsFig.4.Perturbed interfaces obtained by the(a)–(d)thin layer model and(a′)–(d′)numerical simulation in the linear stage:(a),(a′)t=0.0 ns;(b),(b′)t=0.4 ns;(c),(c′)t=0.8 ns;(d),(d′)t=1.0 ns.

圖5 (a)—(d)薄層模型和(a′)—(d′)數(shù)值模擬所描述的在非線性階段的不穩(wěn)定界面 (a),(a′)t=0.0 ns;(b),(b′)t=0.4 ns;(c),(c′)t=0.6 ns;(d),(d′)t=0.82 nsFig.5.Perturbed interfaces obtained by the(a)–(d)thin layer model and(a′)–(d′)numerical simulation in the nonlinear stage:(a),(a′)t=0.0 ns;(b),(b′)t=0.4 ns;(c),(c′)t=0.6 ns;(d),(d′)t=0.82 ns.

其次,對比不穩(wěn)定性在非線性階段的界面分布規(guī)律.假設界面的初始擾動幅值η0=2.0μm,圖5(a)—(d)和(a′)—(d′)分別表示薄層模型和數(shù)值模擬所描述的非線性階段的不穩(wěn)定性界面.通過對比,薄層模型與數(shù)值模擬在描述擾動發(fā)展中的界面幅值的總長度基本保持一致;然而,在薄層模型的發(fā)展后期,波峰幅值明顯大于波谷幅值,這是由于薄層模型中微團上下端面的面積在整個演變過程中是相等的,導致波峰端面受壓力比實際情況要小,波谷端面受力比實際要大,使得氣泡增長加快,尖釘增長變緩,說明薄層模型可以用于描述初始為大擾動幅值的不穩(wěn)定性.另外,薄層模型在擾動發(fā)展的后期(t=0.82 ns),波谷頂端出現(xiàn)類“蘑菇”形結構,可以解釋非線性階段“蘑菇”成型的原因.

3.4 擴展應用

圖6 薄層模型應用于任意初始波形的不穩(wěn)定界面 (a)初始三角波分布;(b)初始方波分布Fig.6.Thin Layer model is applied to obtained the arbitrary perturbed distributions at the initial moment:(a)The perturbed interface with triangular wave;(b)the perturbed interface with square wave.

由于薄層模型是基于對流體微團的受力分析,因此,薄層模型不僅可以應用于初始小擾動的不穩(wěn)定性,而且可以描述初始大擾動的情況.同時,可以描述任意初始波形的不穩(wěn)定界面,例如三角波、方波等.

圖6(a)和圖6(b)分別展示了初始擾動界面為三角波和方波的界面演變分布,其中擾動波長λ=20μm,重力加速度g=1μm/ns2及Atwood數(shù)AT=0.8相同.如圖6所示:在初始階段三角波和方波的擾動界面均表現(xiàn)出規(guī)則的形變規(guī)律,波峰和波谷的增長幅度基本一致;隨著擾動的發(fā)展,波峰的增長幅度稍大于波谷;針對三角波界面,在t=2.5 ns時刻,波谷的頂端處于水平狀態(tài),自此之后,薄層模型描述的波谷位置將會出現(xiàn)凸起現(xiàn)象,薄層模型已經(jīng)不能準確描述擾動的發(fā)展;針對方波界面,在t=1.6 ns時刻,波谷處的界面出現(xiàn)不規(guī)則形變,這是由于薄層模型在描述不穩(wěn)定性時,界面微元向波谷處聚集,導致微元之間相互碰撞.由于薄層模型是基于界面微元的受力分析,因此,初始時界面處于直線狀態(tài)(例如:水平、斜線)時,擾動在演變過程中對此區(qū)域的界面不會帶來形變影響,界面始終保持直線狀態(tài).

4 結 論

本文在Ott薄層理論的基礎上,發(fā)展了能夠描述任意Atwood數(shù)的薄層模型.該模型的數(shù)值解分別與WNL模型、數(shù)值模擬結果比較,薄層模型可以準確描述線性和弱非線性的界面演化和形狀;在非線性階段中“蘑菇”形狀的發(fā)展過程中,界面形狀比較接近數(shù)值模擬的結果.通過研究發(fā)現(xiàn),薄層模型可以描述初始大擾動幅值的不穩(wěn)定性.此外,對于初始時任意分布的擾動界面,該模型均可以很好地描述其演化過程.該模型可以推廣至三維球形幾何,有助于人們理解內(nèi)爆過程中薄層界面的形變規(guī)律和演化機理.

參考文獻

[1]Wang J H 1994Nonstationary Flow and Shock for Two-Dimensional(Beijing:Science Press)p10(in Chinese)[王繼海 1994二維非定常流和激波 (北京:科學出版社)第10頁]

[2]Rayleigh L 1893Proc.R.Math.Soc.14 170

[3]Taylor G I 1950Proc.R.Soc.London:Ser.A201 192

[4]Chandrasekhar S 1961Hydrodynamic and HydromagneticStability(London: Oxford University Press)pp429–514

[5]Nuckolls J H,Wood J,Thiessen A,Zimmerman G 1972Nature239 139

[6]Lindl J D,Amendt P,Berger R L,Glendinning S G,Glenzer S H,Haan S W,Kau ff man R L,Landen O L,Suter L J 2004Phys.Plasmas11 339

[7]Atzeni S,Meyer-ter-Vehn J 2004The physics of Inertial Fusion:Beam Plasma Interaction,Hydrodynamics,Hot Dense Matter(Oxford:Oxford University Press)

[8]He X T,Zhang W Y 2007Eur.Phys.J.D44 227

[9]Remington B A,Drake R P,Ryutov D D 2006Rev.Mod.Phys.78 755

[10]Remington B A,Arnett D,Drake R P,Takabe H 1999Science284 1488

[11]Committee on High Energy Density Plasma Physics Plasma Science Committee Board on Physics and Astronomy Division on Engineering and Physical Science 2001Frontiers in High Density Physics(Washington,DC:Academic Press)

[12]Vlemmings W H T,Diamond P J,Imai H 2006Nature440 58

[13]Wang L F,Ye W H,Li Y J 2010Phys.Plasmas17 052305

[14]Liu W H,Wang L F,Ye W H,He X T 2012Phys.Plasmas19 042705

[15]Wang L F,Wu J F,Fan Z F,Ye W H,He X T,Zhang W Y,Dai Z S,Gu J F,Xue C 2012Phys.Plasmas19 112706

[16]Wang L F,Ye W H,Sheng Z M,Don W S,Li Y J,He X T 2010Phys.Plasmas17 122706

[17]Ye W H,Wang L F,He X T 2010Phys.Plasmas17 122704

[18]Wang L F,Ye W H,He X T,Zhang W Y,Sheng Z M,Yu M Y 2012Phys.Plasmas19 100701

[19]Wang L F,Ye W H,Wu J F,Liu J,Zhang W Y,He X T 2016Phys.Plasmas23 052713

[20]Wang L F,Ye W H,Wu J F,Liu J,Zhang W Y,He X T 2016Phys.Plasmas23 122702

[21]Wang L F,Ye W H,He X T,Wu J F,Fan Z F,Xue C,Guo H Y,Miao W Y,Yuan Y T,Dong J Q,Jia G,Zhang J,Li Y J,Liu J,Wang M,Ding Y K,Zhang W Y 2017Sci.China:Phys.Mech.Astron.60 055201

[22]Zhang W Y,Ye W H,Wu J F,Miao W Y,Fan Z F,Wang L F,Gu J F,Dai Z S,Cao Z Y,Xu X W,Yuan Y T,Kang D G,Li Y S,Yu X J,Liu C L,Xue C,Zheng W D,Wang M,Pei W B,Zhu S P,Jiang S E,Liu S Y,Ding Y K,He X T 2014Sci.Sin.:Phys.Mech.Astron.44 1(in Chinese)[張維巖,葉文華,吳俊峰,繆文勇,范征鋒,王立鋒,谷建法,戴振生,曹柱榮,徐小文,袁永騰,康洞國,李永升,郁曉瑾,劉長禮,薛創(chuàng),鄭無敵,王敏,裴文兵,朱少平,江少恩,劉慎業(yè),丁永坤,賀賢土 2014中國科學:物理學力學 天文學44 1]

[23]Reipurth B,Bally J 2001Annu.Rev.Astron.Astrophys.39 403

[24]Jacobs J W,Catton I 1988J.Fluid Mech.187 353

[25]Kull H J 1991Phys.Rep.206 197

[26]Sanz J 1994Phys.Rev.Lett.73 2700

[27]Garnier J,Raviart P A,Cher fils-Clérouin C,Masse L 2003Phys.Rev.Lett.90 185003

[28]Haan S W 1991Phys.Fluids B3 2349

[29]Youngs D L 1984Physica D12 32

[30]Zhang Y,Drake R P,Glimm J 2007Phys.Plasmas14 062703

[31]Jacobs J W,Catton I 1988J.Fluid Mech.187 329

[32]Waddell J T,Niederhaus C E,Jacobs J W 2001Phys.Fluids13 1263

[33]Wilkinson J P,Jacobs J W 2007Phys.Fluids19 124102

[34]Olson D H,Jacobs J W 2009Phys.Fluids21 034103

[35]Wang L F,Ye W H,Li Y J 2010Chin.Phys.Lett.27 025203

[36]Wang L F,Wu J F,Ye W H,Zhang W Y,He X T 2013Phys.Plasmas20 042708

[37]Davies R M,Taylor G I 1950Proc.Roy.Soc.A200 375

[38]Layzer D 1955Astrophys.J.122 1

[39]Zhang Q 1998Phys.Rev.Lett.81 3391

[40]Goncharov V N 2002Phys.Rev.Lett.88 134502

[41]Sohn S 2003Phys.Rev.E67 026301

[42]Abarzhi S I,Nishihara K,Glimm J 2003Phys.Lett.A317 470

[43]Mikaelian K O 2003Phys.Rev.E67 026319

[44]Tao Y S,Wang L F,Ye W H,Zhang G C,Zhang J C,Li Y J 2012Acta Phys.Sin.61 075207(in Chinese)[陶燁晟,王立鋒,葉文華,張廣財,張建成,李英駿 2012物理學報61 075207]

[45]Ott E 1972Phys.Rev.Lett.29 1429

[46]Manheimer W,Colombant D,Ott E 1984Phys.Fluids27 2164

[47]Colombant D,Manheimer W,Ott E 1984Phys.Rev.Lett.53 446

[48]Wang L F,Guo H Y,Wu J F,Ye W H,Liu J,Zhang W Y,He X T 2014Phys.Plasmas21 122710