基于MULTI2D-Z程序的Z箍縮動態黑腔形成過程模擬?

陳忠旺 寧成

(北京應用物理與計算數學研究所,北京 100088)

基于MULTI2D-Z程序的Z箍縮動態黑腔形成過程模擬?

陳忠旺 寧成?

(北京應用物理與計算數學研究所,北京 100088)

(2016年12月26日收到;2017年4月12日收到修改稿)

對輻射流體力學程序MULTI-2D進行改造,增加磁場演化方程程序模塊,自洽地在運動方程模塊中增加洛倫茲力,在能量方程模塊中增加歐姆加熱,將它改造成輻射磁流體力學程序MULTI2D-Z.驗證了新增磁場程序模塊的可靠性,并發現溫度和密度的增大會抑制磁場的擴散,負徑向速度梯度的流體對流也會抑制磁場的擴散.利用改造好的MULTI2D-Z程序模擬了峰值為8 MA的脈沖電流驅動的鎢絲陣Z箍縮動態黑腔形成過程.得到了X光功率(約30 TW)和能量(約300 kJ)、泡沫輻射溫度(約120 eV)、箍縮軌跡等模擬結果.在動態黑腔形成過程中,磁場主要分布在鎢主體等離子體中;輻射向內傳播,燒蝕泡沫柱而使它膨脹;輻射熱波在被撞擊的泡沫柱中傳播,其傳播速度比物質溫度傳播得快,當輻射熱波傳播到中心軸時泡沫柱中的輻射場變得比較均勻,并且除了沖擊波處外輻射溫度與物質溫度基本上沒有分離.這些模擬結果可增強人們對磁場擴散和對流規律以及動態黑腔形成機制的理解,同時表明了MULTI2D-Z程序可成為Z箍縮及其應用的新的程序模擬工具.

Z箍縮,動態黑腔,輻射磁流體,MULTI-2D程序

1 引 言

Z箍縮(Z-pinch)是由驅動電流產生的強磁場作用于自身載流等離子體(負載)上,使其受洛倫茲力作用而向軸心內爆的物理過程.當等離子體箍縮到心滯止時,產生高溫高密度等離子體,并輻射出很強的X光,形成等離子體強輻射源.這類強X光輻射源在慣性約束聚變、實驗室天體物理、材料輻射特性和輻射效應等研究領域有著廣泛的應用[1].

近二十年來,國內外的Z箍縮實驗和數值模擬研究均取得了很大的進展.國內在“強光一號”和“聚龍一號”(PTS)等脈沖功率裝置上進行了大量的絲陣Z箍縮實驗研究[2?5].在美國Sandia實驗室的Z裝置上進行的絲陣實驗中,產生了峰值功率200 TW、能量2 MJ、脈寬5 ns的X射線脈沖,能量轉換效率高達15%[6?8].在英國帝國理工大學的Magpie裝置上也進行了大量的絲陣Z箍縮基礎物理實驗研究,發展了多種先進的診斷技術,對絲陣消融和內爆的細致物理過程進行了深入的觀察和研究[9?11].

在Z箍縮的諸多應用中,其驅動慣性約束聚變應用研究是最令人矚目的.在Z箍縮間接輻射驅動慣性約束聚變的概念設計中,需要形成一個黑腔為氘氚的靶丸內爆創造一個盡可能均勻、對稱的輻射環境.目前人們提出了三種通過Z箍縮而形成的黑腔,即雙端驅動黑腔、靜態壁黑腔和動態黑腔[12].前兩種黑腔的輻射源和靶丸相對分離,其輻射場的對稱性均比較好,大體上能滿足慣性約束聚變對輻射場對稱性的要求,但其能量耦合效率較低,需要相對較高的驅動電流才能達到聚變點火條件.在動態黑腔的形成中,其核心思想是利用單層或者雙層絲陣Z箍縮內爆等離子體撞擊內置的中心泡沫塑料柱,內爆等離子體撞擊塑料泡沫柱時損失的動能以及形成的沖擊波產生輻射,這些輻射在泡沫柱中傳播和均勻化,并且被外面高Z絲陣等離子體俘獲,而在靶丸的周圍形成一個受約束的輻射環境.因為它是在絲陣內爆動態過程中形成的,故稱其為動態黑腔.在動態黑腔中Z箍縮內爆等離子體、動能-X光能量轉換體(泡沫柱)和聚變靶丸等緊密地耦合在一起,從而使它們之間的能量耦合效率比較高,進而可能在驅動電流水平較低的Z箍縮裝置上實現慣性約束聚變點火.目前動態黑腔成了Z箍縮間接驅動慣性約束聚變的主要技術途徑之一[13].

Z裝置上的動態黑腔實驗[14,15]表明:泡沫柱的密度由15 mg/cc被壓縮到了35 mg/cc,而電子溫度由400 eV因輻射降到了300 eV;受沖擊壓縮的泡沫等離子體產生了大于180 kJ的輻射能量,其中靶丸吸收了超過40 kJ的輻射能量,而此能量已是靶丸點火所需能量的1/4;沖擊波對稱性受內爆等離子體不均性和不穩定性的影響較小,從而有利于中心聚變靶丸實現球對稱輻射燒蝕內爆壓縮.此外在Z裝置上還首次進行了利用Z箍縮動態黑腔驅動充氘氣靶丸出中子的實驗[16],產生了1×1010—5×1010個熱核中子.

美國著名的激光慣性約束聚變模擬程序LASNEX已增加了磁場模塊,已應用于Z箍縮驅動的慣性約束聚變物理研究[17,18].原來LASNEX程序是典型的三溫(電子、離子和輻射溫度)輻射流體力學程序,與國內大部分的激光慣性約束聚變和Z箍縮模擬程序的物理模型相似,所采用的計算方法也主要是拉格朗日流體算法.它在所增加的磁場模塊中考慮了非理想磁場演化因素,如Hall效應和Nernst效應.在空間二維上它對動態黑腔中絲陣內爆到泡沫柱轉換體和嵌套在泡沫柱中的靶丸內爆,以及電極板對輻射的吸收等過程進行了整體數值模擬[18].絲陣內爆時間在100 ns左右,而緊接著的靶丸內爆僅有幾個納秒.模擬得到的中子產額比實驗測量的值大兩倍,沖擊波速度比測量值高30%,靶丸心部的氣體溫度比由光譜推測的值高.模擬結果再現了Z箍縮內爆到泡沫柱和緊接著產生的輻射沖擊波過程,還預言了泡沫柱吸收等離子體早期因歐姆加熱和內爆而產生的輻射所形成的消融沖擊波,從而幫助人們更好地理解動態黑腔的形成過程.

在國內也開展了Z箍縮驅動動態黑腔的實驗和數值模擬研究[5,19?21].目前已發表的主要是在驅動電流比較低的實驗裝置(如強光I號)上獲得的實驗結果,得到了徑向和軸向的X光輻射功率、外邊界的內爆軌跡、泡沫柱的發光圖像等結果.數值模擬主要是在空間一維上再現了動態黑腔形成過程的主要物理特征,如各物質界面的內爆軌跡,密度、電子溫度和輻射溫度的時空演化及其相互關系,沖擊波形成和傳播過程等.

MULTI系列程序是西班牙的Ramis教授等開發的,用于模擬高能量密度等離子體中發生的輻射和物質相互作用過程,如激光慣性約束聚變中的激光-等離子體相互作用過程.最早的MUTLI一維程序是使用Fortran語言編寫的,1994年Ramis教授開發了一種新的程序語言r94,用于編寫MULTI程序.經過二十多年的發展,目前MUTLI系列程序有MULTI-IFE,MULTI-fs,MULTI-1D,MULTI-2D和MULTI-3D等版本.由于它們都是開源的,因此得到了比較廣泛的應用[22?24].相對而言,一維和二維程序已發展得比較成熟.

動態黑腔負載是由柱對稱的多層物質嵌套組合而形成,至少需要空間二維的物理模型才能比較好地描述動態黑腔的形成過程.因此,本文對國外開源的激光慣性約束聚變模擬程序MULTI-2D進行改造,使它們能應用于Z箍縮物理、Z箍縮驅動慣性約束聚變及其聚變能應用方面的物理研究.在物理模型上,它最大的特點是僅有一個物質溫度,輻射溫度是通過等離子體物質的輻射參數和近似輻射能流計算得到,因此,具有非常高的計算效率.它采用任意歐拉-拉格朗日算法和非結構三角形網格,能處理在內爆過程中產生的流體大變形,較好地避免在流體模擬中因網格畸變而無法繼續計算的問題,以及能處理由不同幾何構形部件組合而成的復雜的內爆物理問題.

2 物理模型和程序框架

采用輻射磁流體力學方程對等離子體Z箍縮動力學過程進行模擬[25,26].質量守恒方程為

這里ρ是物質密度,u為流體速度,t為時間.

動量守恒方程為

其中p為總壓強,包括流體壓強和黏性;J×B為導電流體受到的單位體積洛倫茲力,J為電流密度,B為磁感應強度.在柱對稱的Z箍縮中,可以僅考慮磁場的角向分量Bθ,此時,洛倫茲力在二維柱坐標系(r,z坐標系)下的表達式為

這里fr和fz分別表示f在r方向和z方向的分量.能量守恒方程為

式中,ε是流體單位質量內能;?·(ρεu)表示單位時間流入流體單元的內能;p?·u是壓強做功;?·F是熱傳導傳輸的能量;?·S是輻射能流沉積;W是歐姆加熱,其在二維柱坐標系下的表達式為

其中,η⊥為垂直于磁場方向的等離子體電阻率,c為光速.

在輻射輸運中,流體速度遠小于光速,在任一時刻,輻射場都可以認為是準靜態的.在局部熱力學平衡附近,忽略散射和折射,輻射輸運方程為[23]

其中,n是傳播方向;qν是單位體積輻射的沉積能;λ′ν為輻射頻率為ν的輻射平均自由程,與材料的溫度和密度有關[27];Iν是輻射強度;IPν是普朗克分布函數,

這里h和k分別是普朗克常數和玻爾茲曼常數,T為物質溫度.若考慮單群方法處理輻射輸運方程,則對所有頻率積分就得到

這里平均輻射自由程λ定義為

其中I,q和IP=σT4/π分別表示Iν,qν和IPν對所有頻率的積分值;σ為Stefan-Boltzmann常數.在局部熱力學平衡附近,考慮近似n·?I≈n·?IP,可得輻射強度近似表達式

從而可得輻射能流

利用輻射強度,空間某一點的輻射溫度通過以下表達式進行計算:

磁場演化方程為

只考慮θ方向的磁場,在二維柱坐標系下(7)式變成

其中,uz,ur分別為流體沿著z和r方向的速度分量;磁擴散系數vm=η⊥c2/4π.平行磁場或無磁場下Spitzer電阻率η =1.15× 10?14z lnΛT?3/2(s),其中z為平均電離度.

本文利用實用的解析表達式來計算電離度.在Thomas-Fermi理論框架下,More提出了一個平均電離度計算公式[28].設ρ為質量密度,單位g/cm3,Te為電子溫度,單位eV(約11605 K),Z為原子序數,A為相對原子質量.定義如下過程量R=ρ/(ZA),T0=Te/Z4/3,T1=T0/(1+T0),α =14.3139,β =0.6624,x= αQβ,f(x)=則平均電離度z=f(x)Z.其中Q的計算過程如下:

相對于原來MULTI-2D程序的物理模型,為了使它能模擬研究Z箍縮及其相關的物理問題,需要增加磁場演化方程,在動量守恒方程中添加洛倫茲力,在能量守恒方程的源項中添加歐姆加熱項.此外,為了計算電阻率而增加一個電離度的解析計算表達式.這樣就把它改造成了輻射磁流體力學程序,并仍然具有高的計算效率.

MULTI-2D是結構化程序,由若干相對獨立的功能模塊組成.這些模塊包括流體模型建立模塊、循環模塊、數據輸出和繪圖模塊,如圖1所示.在流體模型建立模塊中,包括網格生成、參數輸入和數據初始化.在主體循環模塊中,有外加驅動源、流體運動、能量輸運以及狀態方程等模塊.驅動源是外加能量輸入項,如激光輻射、重離子束注入,及本文增加的電流或磁場.能量模塊包括物質能量方程,以及輻射輸運、熱傳導和做功等的計算.狀態方程采用QEOS程序生成的狀態方程參數[29,30].

圖1 MULTI-2D程序框架結構Fig.1.Programm frame of MULTI-2D.

MULTI-2D程序采用r94語言編寫.該語言于1994年開發,是基于C語言的高級語言,適用于UNIX操作系統.它通過構建大量函數庫而進行編程,使得程序精煉,可使用X窗口進行操作.

本文在MULTI-2D程序中加入磁場演化程序模塊,在流體運動模塊中加入洛倫茲力,在物質能量方程模塊中加入歐姆加熱,將它發展成二維輻射磁流體力學數值模擬程序MULTI2D-Z,并應用于Z箍縮動態黑腔形成過程的模擬研究.

3 模擬結果

3.1 磁場演化模塊驗證與磁場演化分析

為了驗證所添加的磁場演化方程程序模塊的正確性,考慮如圖2所示的等離子體柱.橫坐標軸為柱等離子體的對稱軸(z軸),縱坐標軸為柱等離子體的半徑(r軸).等離子體柱的初始半徑為1.25 cm,半高度為0.5 cm.考慮的是鎢等離子體,密度為8.6×10?4g/cm3,初始溫度為1 eV.由于柱等離子體的對稱性,因此,實際需要的計算區域僅為實線的右半部分.在z方向上加載8 MA的恒定電流,該電流在θ方向產生角向磁場,該磁場將向對稱軸方向擴散演化.為了驗證的方便,給定等離子體速度,不求解流體運動方程.此時,磁場的演化過程由方程(8)描述.磁場的邊界條件為

其中rm為等離子體柱的徑向邊界,L為z方向邊界.

明顯地,該磁場演化過程可以由徑向(r方向)一維磁場演化方程來描述.一維方程為

其邊界條件為

圖2 磁場演化模塊驗證模型Fig.2.Plasma model for veri fi cation of magnetic fi eld evolution module.

一維磁場演化方程程序模塊在柱對稱Z箍縮等離子體的模擬中是普遍使用的,其磁場演化規律清楚,程序模塊可靠.將由新添加的二維磁場演化方程程序模塊模擬出的結果與一維磁場演化程序模擬的結果進行比較,就可以判斷基于非結構網格上的新增磁場程序模塊的正確性.當流體靜止時,僅有純磁擴散,而當存在流體運動速度梯度時,則磁場既有擴散,也有對流.下面就這兩種情況對磁場演化模塊和相應的磁場演化規律進行考察.

圖3(a)是在假定流體靜止時,當8 MA的恒定電流加載到該等離子體柱時,利用一維和二維磁場演化方程對其中磁場演化過程的模擬結果進行比較.恒定電流所產生的角向磁場僅隨r而變化,磁場逐漸向中心軸擴散,并且很快就達到穩定分布.可見一維和二維的模擬結果符合得很好.

圖3(b)給出的是在流體流動、并且存在速度梯度時,利用一維和二維磁場演化方程對該磁場演化過程的模擬結果.這里假定流體徑向速度隨半徑的變小而下降,為某一徑向位置,um=?5×106cm/s.磁場外邊界值Bm=2I/crm隨著流體向內箍縮運動而增大,磁場既隨流體對流,又向內擴散.可見,此時二維和一維模擬結果也很好地重合,說明基于非結構網格上的新增二維磁場演化方程程序模塊是正確的.

下面利用我們新增的基于非結構網格上的二維磁場演化方程程序模塊來研究磁擴散系數和流體對流對磁場演化的影響.在模擬磁擴散系數對磁場演化的影響時,通過改變等離子體溫度或密度而改變電阻率,使磁擴散系數發生變化.這里采用圖2所示的基準等離子體溫度或密度的100倍來改變電阻率(磁擴散系數),觀察(100倍溫度,基準密度),(基準溫度,100倍基準密度)和(基準溫度,基準密度)等三種情況下t=100 ns時的磁場分布情況.在考察流體對流對磁場演化的影響時,使等離子體處于三種運動狀態,即零速度、恒定速度u1=um=?5×106cm/s和速度隨半徑變化,并假定磁場的邊界值為恒定的Bm.為了描述磁場擴散或向軸心演化的快慢,我們定義一個磁場擴散速度.磁場向軸心擴散而逐漸減小,因此,定義磁場為Bm/e的位置r隨時間變化率為磁場擴散速度(自然常數e=2.71828).

圖4(a)是磁場擴散的模擬結果,可見密度變化對磁場分布的影響比較小,密度增加導致電離度減小,從而電阻率減小,磁擴散系數變小,進而抑制磁場的擴散;溫度對磁場分布的影響比較復雜,雖然溫度升高使電離度有增大的趨勢,使得電阻率增大,但電阻率與T?3/2成正比,溫度升高又使得電阻率下降得更快.綜合起來,本模擬結果表明溫度升高抑制了磁場的擴散,磁場更集中于等離子體柱的表面.這表明等離子體的導電性好,則電流更趨膚,是符合物理理解的.圖5(a)是與圖4(a)相同條件下的磁擴散速度隨時間的變化情況.在電流加載的初始階段,磁擴散速度比較大,但隨著時間的增加,擴散速度快速下降,最后趨于零,即磁場分布趨于穩定.但等離子體溫度較高時,磁擴散速度初始時相對較小,下降也緩慢,或者說磁擴散系數減小,磁場擴散速度減慢,需要更長的時間才能到達穩定.

圖3 (網刊彩色)二維磁場演化程序與一維磁場演化程序模擬得到的磁場徑向分布演化的比較 (a)流體靜止;(b)存在流體速度梯度Fig.3. (color online)The comparison of magnetic fi eld pro fi les in radial direction calculated byrtwodimensional and one-dimensional magnetic inductance equations:(a)ur=uz=0;(b)ur(ri)= ium,u=0.rmz

圖4 (網刊彩色)磁場的分布演化 (a)流體靜止時等離子體參數變化對磁場擴散的影響;(b)流體對流對磁場分布演化的影r響,u1表示等離子體以恒定速度um=?5×106cm/s向軸心運動,u2表示等離子體以空間變化速度ur(ri)=ium向軸心運動;邊界磁場固定Bm=128 T rmFig.4.(color online)The evolutions of magnetic fi eld pro fi les:(a)The e ff ects of plasma parameters on the magnetic fi eld evolution when the the plasma is at rest;(b)the e ff ect of fl uid motion on the magnetic fi eld evolution when the plasma is moving towards the axis at a constant velocity um=?r5×106cm/s,which is denoted by u1,as well as the plasma is moving at a convective velocity ur(ri)=ium,which is denoted by u.The magnetic fi eld at boundary interface is fi xed,B=128 T.rm2m

圖5 (網刊彩色)磁場擴散速度隨時間變化 (a)流體靜止;(b)流體運動,uo表示流體靜止;其他與圖4中的情況相對應Fig.5.(color online)The variation of velocity of magnetic fi eld di ff usion with time:(a)The plasma is at rest;(b)the plamsa is moving at a certain velocity.The conditions used in the simulations are the same as that in Fig.4,besides the uowhich denotes the plasma at rest.

圖4(b)是在等離子體運動的情況下磁場擴散的結果.當等離子體以恒定的速度向軸心運動時,磁場在不同時刻的徑向分布是相似的(這是因為在拉氏坐標下當流體速度恒定均勻時,對流項消失,從而磁場擴散與流體靜止時的情況一樣),但隨著等離子體的壓縮,內部的磁場逐漸增大.當等離子體以隨r逐漸減少的速度箍縮時,磁場的擴散效應減弱,磁場的分布更趨膚.表明此時的等離子體對流抑制了磁場的擴散,或者說在實際的Z箍縮中,在箍縮初始階段已擴散進入等離子體中的磁場,當等離子體變得快速向心運動時,磁場將有向外反擴散而再趨膚的過程.圖5(b)是與圖4(b)相同條件下的磁擴散速度隨時間的變化情況.有意思的是,當等離子體靜止時磁擴散速度逐漸趨于零,而當它以恒定速度流動時磁擴散速度等于等離子體運動速度5×106cm/s,即磁場被凍結在等離子體中,隨流體一起運動;當速度不恒定,存在負的速度梯度時,磁擴散速度略小于流體運動速度,這是因為此時的磁擴散被減弱了.

3.2 絲陣Z箍縮驅動的動態黑腔形成過程

新添加的磁場演化方程程序模塊得到驗證后,改造后MULTI-2D程序的驅動源除了原來的激光和粒子束外,又增加了磁場力,從而它就具有了模擬Z箍縮及其相關物理過程的能力.典型而簡單的動態黑腔負載構形就是在絲陣的中心嵌入一根低密度泡沫柱(其組成成分為碳和氫,通常將它簡記為CH),當在強脈沖電流的作用下絲陣內爆到泡沫柱,撞擊泡沫柱而在泡沫柱中形成沖擊波,絲陣和泡沫等離子體被熱化,產生輻射向泡沫柱中心傳播,在絲陣等離子體的部分約束下在泡沫等離子體中形成比較勻化的輻射場.為了能利用動態黑腔來驅動慣性約束聚變,這一整體動力學過程的發展、演化過程需要通過數值模擬(當然也需要實驗)研究、理解清楚.

下面在二維柱坐標系下,利用改造完成的MULTI2D-Z程序對峰值電流為Im=8 MA的脈沖電流驅動的嵌套了中心泡沫柱的絲陣Z箍縮內爆過程進行模擬.

圖7是計算得到的絲陣和泡沫柱的箍縮軌跡、X光功率和能量、絲陣和泡沫柱的收縮比、電流波形等隨內爆過程的變化.整體箍縮到心在150 ns左右,X光功率峰值接近30 TW,總能量在330 kJ左右,絲陣和泡沫的最大壓縮比分別為15和5左右.泡沫柱在被絲陣等離子體碰撞壓縮以前,略微膨脹.

圖6 動態黑腔負載模型Fig.6. The used model in simulation of dynamic hohlraum.

圖7 (網刊彩色)動態黑腔形成過程中一些宏觀參量隨時間的變化 (a)絲陣外半徑R、電流I、X光功率P和徑向輻射能量E;(b)絲陣和泡沫的柱外半徑(rWA,rCH),絲陣和泡沫柱的收縮比(Ratio 2,Ratio 1)Fig.7.(color online)The variations of some macroscopical quantities with time during the formation process of dynamic hohlraum:(a)Outer radius of wire-array(R),current(I),X-ray power(P)and energy(E);(b)outer radii of wire-array(rWA)and foam(rCH),convergence ratioes of wire-array(Ratio 2)and foam(Ratio 1).

圖8—圖10分別是絲陣等離子體內爆碰撞泡沫柱后三個時刻(146,150,152 ns)密度、等離子體溫度和輻射溫度、內爆速度、磁場、輻射場(溫度)等的空間分布情況.t=146 ns時,絲陣內爆到了半徑為0.3 cm的地方,離泡沫柱的初始位置還有0.1 cm.此時在絲陣等離子體中形成了很強的沖擊波(見圖8(a)—(c)中的密度尖峰位置),內爆速度達到了4×107cm/s,物質溫度達到了370 eV左右,泡沫中的最大輻射溫度接近100 eV,邊界磁場達到約520 T.在靠近中心軸處的輻射溫度比物質溫度高,表明輻射熱波向軸心方向傳播比物質溫度快.在箍縮過程中磁場主要分布在主體鎢等離子體中,形成很強的磁壓,使等離子體加速向內箍縮,在壓縮做功和歐姆加熱作用下,主體等離子體溫度較高,密度較大.溫度升高和密度增加使磁擴散系數減小,另外流體壓縮也抑制磁場的擴散.高溫高密度主體絲陣等離子體向泡沫內部和外表面及兩端電極發射輻射,向內輻射加熱燒蝕泡沫柱,使泡沫柱因溫度升高而向外膨脹,此時,泡沫柱的外半徑膨脹到了大約0.24 cm.向外表面(r方向)和電極(z方向)的輻射導致內爆等離子體的能量漏失,使得靠近電極板區域和等離子體外則的輻射溫度均比較低.

到了150 ns(見圖9),絲陣等離子體內爆到全部進入了原泡沫柱所在的區域,到達半徑為0.17 cm的位置,表明泡沫柱已被絲陣等離子體碰撞,泡沫柱得到了一點壓縮.此時在絲陣等離子體中的沖擊波(見圖9(a)—(c)中的密度尖峰位置)得到增強,溫度和密度增高,但內爆速度下降了一點,最大內爆速度下降到3.7×107cm/s,物質溫度升高到460 eV左右,泡沫中輻射溫度在中心區域基本均勻化,其值達約120 eV,并且其中的物質溫度和輻射溫度幾乎沒有分離,達到了平衡,邊界磁場超過了900 T.鎢等離子碰撞低密度泡沫后形成物質沖擊波,速度減緩,動能轉換為內能,溫度迅速升高,形成輻射沖擊波;輻射沖擊波迅速向內傳播,在物質沖擊波到來之前形成了較為均勻的黑腔輻射環境,黑腔中物質溫度和輻射溫度較為重合,溫度達120 eV左右;與此同時,主體等離子體向外輻射出X光的第一個小波峰,功率在7 TW左右(見圖7(a)).在動態黑腔概念設計中,希望在絲陣碰撞壓縮泡沫柱到達靶丸之前,所形成的黑腔輻射場對靶丸的燒蝕壓縮,并使其達到聚變點火條件.因此,碰撞后主體等離子的速度大幅度減小,損失的動能轉化為內能和輻射,大部分輻射被泡沫吸收,這是黑腔輻射場形成的重要機制.

圖8 (網刊彩色)146 ns時刻各參量的隨空間分布 (a)輻射溫度Tr、物質溫度T和密度分布ρ;(b)速度ur和密度ρ分布;(c)磁場B和密度ρ分布;(d)輻射溫度Tr空間分布Fig.8.(color online)The spatial pro fi les of some physical quantities at 146 ns:(a)Radiation,plamsa temperatures,and density pro fi les;(b)velocity and density pro fi les;(c)magnetic fi eld and density pro fi les;(d)spatial distribution of radiation temperature.The part of a line with symbols of circle belongs to foam,and it with symbols of plus belongs to wire-array.

圖9 (網刊彩色)150 ns時各時刻參量隨空間變化 (a)輻射溫度、物質溫度和密度分布;(b)速度和密度分布;(c)磁場和密度分布;(d)輻射溫度空間分布Fig.9.(color online)The spatial pro fi les of some physical quantities at 150 ns:(a)Radiation,plamsa temperatures,and density pro fi les;(b)velocity and density pro fi les;(c)magnetic fi eld and density pro fi les;(d)spatial distribution of radiation temperature.The part of a line with symbols of circle belongs to foam,and it with symbols of plus belongs to wire-array.

t=152 ns以后,絲陣和泡沫等離子體繼續箍縮,泡沫等離子體整體得到進一步的加速,泡沫中的溫度和密度進一步增大,在泡沫中出現小的密度峰,中心附近的輻射和物質溫度稍高.由于磁壓力做功,使得泡沫的溫度和輻射溫度繼續增加,最高達到約170 eV,它們仍然基本上沒有分離.輻射溫度在r方向分布的均勻性變弱(見圖10(a)),但對比圖9(d)可以看到,高溫和比較均勻的輻射區域在z方向有所增大.由于邊界輻射溫度的增大,使得邊界輻射能量以更高的速率損失(輻射能流約為導致徑向輻射溫度相對更低.因此,對于動態黑腔驅動聚變靶丸的實際應用而言,需要根據驅動電流、靶丸大小仔細地設計動態黑腔的負載構形,以明確能夠利用的動態黑腔中輻射場的時間窗口和空間區域,使其靶丸實現對稱壓縮和點火成功.在整體箍縮到心時,等離子體減速滯止,產生強烈的X光輻射,此時X光輻射脈沖的峰值功率達約28 TW.

在徑向一維上密度、物質溫度(電子溫度)和輻射溫度的空間分布和隨時間的演化關系的整體特征,如絲陣等離子體的前期膨脹、輻射燒蝕泡沫柱而膨脹、輻射熱波傳播較快、輻射熱波與輻射沖擊波的空間分離等,與一維拉氏非平衡輻射磁流體力學程序CRMHA的模擬結果[20]以及國外的模擬結果[18]都是相似的.并且預期的黑腔輻射溫度與文獻[20]的結果也接近,都是120 eV左右.說明MULTI2D-Z程序模擬動態黑腔形成過程所得到的物理結果是基本正確、可靠的.但本文中泡沫柱的輻射溫度的均勻化過程要慢一些,即泡沫等離子體的輻射吸收更強烈一些,并且本文能給出輻射場在徑向和軸向方向上的均勻化情況.

4 結 論

從輻射磁流體力學基本方程出發,通過在MULTI-2D程序中增加磁場演化方程程序模塊,同時自洽地在流體動量守恒方程中增加導電流體所受到的洛倫茲力,在能量守恒方程中增加歐姆加熱能量源項,并在相應的程序中進行改動.通過新增程序模塊的驗證,使得MULTI-2D程序改造、發展成了MULTI2D-Z程序,成為基于非結構網格的二維輻射磁流體力學程序,可以應用于Z箍縮驅動的動態黑腔形成動力學過程模擬.

圖10 (網刊彩色)152 ns時各時刻參量隨空間變化 (a)輻射溫度、物質溫度和密度分布;(b)速度和密度分布;(c)磁場和密度分布;(d)輻射溫度空間分布Fig.10.(color online)The spatial pro fi les of some physical quantities at 152 ns:(a)Radiation,plamsa temperatures,and density pro fi les;(b)velocity and density pro fi les;(c)magnetic fi eld and density pro fi les;(d)spatial distribution of radiation temperature.The part of a line with symbols of circle belongs to foam,and it with symbols of plus belongs to wire-array.

利用新研制的磁場演化方程程序模塊,對磁場的擴散和對流演化規律進行了一些研究.發現溫度和密度的增大會抑制磁場的擴散;負徑向速度梯度的流體(如Z箍縮等離子體的箍縮運動)也會抑制磁場的擴散.當等離子體靜止時,磁擴散速度隨著時間的增加,迅速減少,然后趨于零,磁場分布變得穩定;當等離子體以恒速箍縮時,磁擴散速度也是電流加載后快速下降,然后趨同于流體運動速度,表明磁場被凍結等離子體中,隨流體一起運動;當等離子體以隨半徑變小而減少的速度做箍縮運動時,磁場的擴散速度略小于運動速度,表明磁場擴散過程被削弱了,出現磁場的再趨膚現象.通過對磁場演化的分析,有助于理解Z箍縮過程中各種因素對磁場演化的影響趨勢.

利用MULTI2D-Z程序對8 MA電流水平的絲陣Z箍縮驅動的動態黑腔形成動力學過程進行了模擬,再現了該動態黑腔的形成過程.動態黑腔形成過程的整體特征與國內外的模擬結果基本相似,黑腔輻射溫度與相同電流水平下一維CRMHA程序的模擬結果接近,說明了該程序用于動態黑腔形成過程模擬的可靠性.本模擬中發現,在Z箍縮等離子體的高溫高密度和對流的共同作用下,磁場主要集中在鎢等離子體表面,使絲陣等離子體快速內爆,并撞擊泡沫柱;鎢主體等離子體的輻射向泡沫柱傳播,使泡沫燒蝕而膨脹;絲陣等離子體碰撞泡沫柱,形成很強的物質沖擊波和輻射,輻射向泡沫等離子體中傳播,形成輻射熱波;在泡沫等離子體中輻射熱波的傳播比物質溫度的傳輸速度快,當輻射熱波傳播到中心軸后,泡沫中的輻射達到了較好的均勻化(約120 eV),并且除了沖擊波附近外,輻射溫度和物質溫度基本上重合,沒有分離現象;此后,隨著整體箍縮的進行,泡沫等離子體繼續被加速和壓縮,其中的輻射溫度繼續提升,輻射場在軸向方向的均勻性有所增加,但在徑向方向的均勻性有所下降.當泡沫柱中嵌套有聚變靶丸時,需要通過程序模擬謹慎、精確地設計絲陣和泡沫柱的尺寸和質量,以便能得到滿足聚變點火要求的動態黑腔輻射場時間窗口和空間區域.

動態黑腔形成過程的數值模擬對相關實驗的設計和物理過程的理解,以及Z箍縮慣性約束聚變概念設計是非常重要的.今后還需要繼續完善和發展MULTI2D-Z程序,改善其輻射模型,評估非結構網格劃分對模擬結果的影響.

感謝北京應用物理與計算數學研究所的文武老師、薛創老師、鄔吉明研究員和高志明老師以及中國工程物理研究院研究生部的吳福源師兄在研究工作中給予的指導和幫助.感謝中國工程物理研究院聚龍一號(PTS)裝置的Z箍縮實驗團隊對絲陣Z箍縮和動態黑腔實驗結果的交流和分享.

[1]Matzen M K 1997 Phys.Plasmas 4 1519

[2]Qiu A C,Kuai B,Wang L P,Wu G,Cong P T 2008 High Power Laser and Particle Beams 20 1761(in Chinese)[邱愛慈,蒯斌,王亮平,吳剛,叢培天2008強激光與粒子束20 1761]

[3]Huang X B,Zhou S T,Dan J K,Ren X D,Wang K L,Zhang S Q,Li J,Xu Q,Cai H C,Duan S C,Ouyang K,Chen G H,Ji C,Wei B,Feng S P,Wang M,Xie W P,Deng J J,Zhou X W,Yang Y 2015 Phys.Plasmas 22 072707

[4]Dan J K,Ren X D,Huang X B,Zhang S Q,Zhou S T,Duan S C,Ouyang K,Cai H C,Wei B,Ji C,He A,Xia M H,Feng S P,Wang M,Xie W P 2013 Acta Phys.Sin.62 245201(in Chinese)[但加坤,任曉東,黃顯賓,張思群,周少彤,段書超,歐陽凱,蔡紅春,衛兵,計策,何安,夏明鶴,豐樹平,王勐,謝衛平2013物理學報62 245201]

[5]Jiang S Q,Ning J M,Chen F X,Ye F,Xue F B,Li L B,Yang J L,Chen J C,Zhou L,Qin Y,Li Z H,Xu R K,Xu Z P 2013 Acta Phys.Sin.62 155203(in Chinese)[蔣樹慶,甯家敏,陳法新,葉繁,薛飛彪,李林波,楊建倫,陳進川,周林,秦義,李正宏,徐榮昆,許澤平 2013物理學報62 155203]

[6]Deeney C,Nash T J,Spielman R B,Seaman J F,Chandler G C,Struve K W,Porter J L,Stygar W A,McGum J S,Jobe D O,Gilliland T L,Torres J A,Vargas M F,Ruggles L E,Breeze S,Mock R C,Douglas M R,Fehl D L,McDaniel D H,Matzen M K,Peterson D L,MatuskaW,Roderick N F,MacFarlane J J 1997 Phys.Rev.E 56 5945

[7]Spielman R B,Deeney C,Chandler G A,Douglas M R,Fehl D L,Matzen M K,McDaniel D H,Nash T J,PorterJ L,Sanford T W L,Seaman J F,Stygar W A,Struve K W,Breeze S P,McGurn J S,Torres J A,Zagar D M,Gilliland T L,Jobe D O,McKenney J L,Mock R C,Vargas M,Wangoner T,Peterson D L 1998 Phys.Plasmas 5 2105

[8]Stygar W A,Ives H C,Fehl D L,Fehl D L,Cuneo M E,Mazarakis M G,Bailey J E,Bennett G R,Bliss D E,Chandler G A,Leeper R J,Matzen M K,McDaniel D H,McGurn J S,McKenney J L,Mix L P,Muron D J,Porter J L,Ramirez J S,Ruggles L E,Seamen J F,Simpson W W,Speas C S,Spielman R B,Struve K W 2004 Phys.Rev.E 69 046403

[9]Lebedev S V,Mitchell I H,Aliaga-Rossel R,Bland S N,Chittenden J P,Dangor A E,Haines M G 1998 Phys.Rev.Lett.81 4152

[10]Harvey-Thompson A J,Lebedev SV,Patankar S,Bland S N,Burdiak G,Chittenden J P,Colaitis A,de Grouchy P,Doyle H W,Hall G N,Khoory E,Hohenberger M,Pickworth L,Suzuki-Vidal F,Smith R A,Skidmore J,Suttle L,Swadling G F 2012 Phys.Rev.Lett.108 145002

[11]Swadling G F,Lebedev S V,Harvey-Thompson A J,Rozmus W,Burdiak G,Suttle L,Patankar S,Smith R A,Bennett M,Hall G N,Suzuki-Vidal F,Bland S,Yuan J 2015 Phys.Plasmas 22 072706

[12]Leeper R J,Alberts T E,Asay J R,Baca P M,Baker K L,Breeze S P,Chandler G A,Cook D L,Cooper G W,Deeney C,Derzon M S,Douglas M R,Fehl D L,Gilliland T,Hebron D E,Hurst M J,Jobe D O,Kellogg J W,Lash J S,Lazier S E,Matzen M K,McDaniel D H,McGurn J S,Mehlhorn T A,Moats A R,Mock R C,Muron D J,Nash T J,Olson R E,Porter J L,Quintenz J P,Reyes P V,Ruggles L E,Ruiz C L,Sanford T W L,Schmidlapp F A,Seamen J F,Spielman R B,Stark M A,Struve K W,Stygar W A,Tibbetts-Russell D R,Torres J A,Vargas M,Wagoner T C,Wake fi eld C,Hammer J H,Ryutov D D,Tabak M,Wilks S C,Bowers R L,McLenithan K D,Peterson D L 1999 Nucl.Fusion 39 1283

[13]Ning C,Feng Z X,Xue C 2014 Acta Phys.Sin.63 125208(in Chinese)[寧成,豐志興,薛創 2014物理學報63 125208]

[14]Bailey J E,Chandler G A,Slutz S A,Bennett G R,Cooper G,Lash J S,Lazier S,Lemke R,Nash T J,Nielsen D S,Moore T C,Ruiz C L,Schroen D G,Smelser R,Torres J,Vese y R A 2002 Phys.Rev.Lett.89 095004

[15]Rochau G A,Bailey J E,Maron Y,Chandler G A,Dunham G S,Fisher D V,Fisher V I,Lemke R W,Mac-Farlane J J,Peterson K J,Schroen D G,Slutz S A,Stambulchik E 2008 Phys.Rev.Lett.100 125004

[16]Ruiz C L,Cooper G W,Slutz S A,Bailey J E,Chandler G A,Nash T J,Mehlhorn T A,Leeper R J,Fehl D,Nelson A J,Franklin J,Ziegler L 2004 Phys.Rev.Lett.93 015001

[17]Zimmerman G B,Kruer W B 1975 Comments Plasma Phys.Control.Fusion 2 51

[18]Slutz S A,Peterson K J,Vesey R A,Lemke R W,Bailey J E,Varnum W,Ruiz C L,Cooper G W,Chandler GA,Rochau G A,Mehlhorn T A 2006 Phys.Plasmas 13 102701

[19]Meng S J,Huang Z C,Ning J M,Hu Q Y,Ye F,Qin Y,Xu Z P,Xu R K 2016 Acta Phys.Sin.65 075201(in Chinese)[蒙世堅,黃展常,甯家敏,胡青元,葉繁,秦義,許澤平,徐榮昆2016物理學報65 075201]

[20]Xiao D L,Sun S K,Xue C,Zhang Y,Ding N 2015 Acta Phys.Sin.64 235203(in Chinese)[肖德龍,孫順凱,薛創,張揚,丁寧2015物理學報64 235203]

[21]Xiao D L,Ding N,Ye F,Ning J M,Hu Q Y,Chen F X,Qin Y,Xu R K,Li Z H,Sun S K 2014 Phys.Plasmas 21 042704

[22]Ramis R,Ramirez J,Meyer-ter-Vehn J 1988 Comput.Phys.Commun.49 475

[23]Ramis R,Meyer-ter-Vehn J,Ramirez J 2009 Comput.Phys.Commun.180 977

[24]Ramis R,Meyer-ter-Vehn J,Ramirez J 2016 Comput.Phys.Commun.203 226

[25]Liberman M A,De Groot J S,Toor A 2003 Physics of High-Density Z-Pinch Plasmas(Beijing:National Defence Industry Press)(in Chinese)[Liberman M A,De Groot J S,Toor A著(孫承緯譯)2003高密度Z箍縮等離子體物理學(北京:國防工業出版社)

[26]Ning C,Yang Z H,Ding N 2003 Acta Phys.Sin.52 415(in Chinese)[寧成,楊震華,丁寧 2003物理學報52 415]

[27]Murakami M,Meyer-ter-Vehn J,Ramis R 1990 J.X-Ray Sci.Technol.2 127

[28]Salamann D 1998 Atom Physics in Hot Plasmas(New York:Oxford University Press)p28

[29]Kemp A J,Meyer-ter-Vehn J 1998 Nucl.Instrum.Meth.A 415 674

[30]More R,Warren K H,Yong D A,Zimmermann G 1988 Phys.Fluids 31 3059

PACS:52.58.Lq,52.59.Qy,52.30.Cv,52.38.PhDOI:10.7498/aps.66.125202

Simulation of forming process of Z-pinch dynamic
hohlraum based on the program MULTI2D-Z?

Chen Zhong-Wang Ning Cheng?

(Institute of Applied Physics and Computational Mathematics,Beijing 100088,China)

26 December 2016;revised manuscript

12 April 2017)

The radiation hydrodynamics code MULTI-2D,which was developed by Ramis et al.in 2009(2009 Comput.Phys.Commun.180 977)and adopted the single temperature fl uid and unstructured lagrangian mesh,is modi fi ed into a radiation magnetohydrodynamics code MULTI2D-Z by adding the program module of evolution equation of magnetic fi eld,and self-consistently considering the Lorentz force in the module of motion equation and the Ohmic heating in the module of energy equation.The newly developed module for magnetic fi eld was validated to be reliable.The module is used to study the magnetic fi eld di ff usion process,and it is found that the di ff usion is weakened due to the increasing of plasma temperature and density and the fl uid convection,in which there is minus grads of velocity in radial direction.The new code MULTI2D-Z is used to simulate the formation process of dynamic hohlraum driven by tungsten wire-array Z-pinch at an 8 MA current level.The obtained results are that X-ray power and energy are,respectively,～30 TW and～300 kJ,radiation temperature in foam is～120 eV,and the implosion trajectory of wire-array is also obtained.The calculated results reveal that the magnetic fi eld is mainly distributed in the outside of tungsten plasma during the hohlraum formation.The foam expands due to the radiation heating from the shock wave created by the collision between wire-array plasma and the foam.The thermal radiation wave,which is characterized by radiation temperature,spreads towards the central axis faster than the plasma temperature.When the thermal radiation wave spreads to the central axis,the radiation temperature becomes comparatively uniform in space,and is almost equal to the plasma temperature except at the place of the shock wave.These results help the people to better understand the magnetic fi eld di ff usion and convection in Z-pinch,as well as the formation mechanism of dynamic hohlraum driven by wire-array Z-pinch.It is also indicated that the newly developed code MULTI2D-Z can be considered as a new tool for simulating Z-pinch and its applications,such as inertial con fi nement fusion and magnetically accelerated fl yer plates.

Z-pinch,dynamic hohlraum,radiation magnetohydrodynamics,MULTI-2D code

10.7498/aps.66.125202

?國家自然科學基金(批準號:11675025)和國家自然科學基金重點項目(批準號:11135007)資助的課題.

?通信作者.E-mail:ning_cheng@iapcm.ac.cn

?2017中國物理學會Chinese Physical Society

http://wulixb.iphy.ac.cn

*Project supported by the National Natural Science Foundation of China(Grant No.11675025)and the Key Program of the National Natural Science Foundation of China(Grant No.11135007).

?Corresponding author.E-mail:ning_cheng@iapcm.ac.cn