負三角形變位型下剝離氣球模的非線性演化特征*

秦晨晨 牟茂淋? 陳少永

1) (四川大學物理學院,成都 610065)

2) (四川大學,高能量密度物理國家重點實驗室,成都 610065)

托卡馬克實驗中已經(jīng)實現(xiàn)了負三角形變位型下的高約束放電,其特點是具有較低的臺基,并伴隨幅值較小且頻率較高的邊界局域模.本文基于不同三角形變的托卡馬克平衡,研究了負三角形變位型條件下剝離氣球模的非線性演化特征.研究發(fā)現(xiàn),由于弱場側(cè)壞曲率區(qū)域增大,負三角形變位型會使剝離氣球模失穩(wěn);在非線性階段,負三角形變位型下的剝離氣球模壓強擾動分布在極向截面上擴展到了弱場側(cè)的頂部和底部區(qū)域,使得邊界局域模更早發(fā)生崩塌,同時,在負三角形變位型下,多種環(huán)向模數(shù)的擾動被激發(fā)并增長,故而具有更明顯的湍流輸運特性.

1 引言

高約束運行模式(high confinement regime,H 模)[1]因其較高的約束參數(shù)已成為未來聚變堆穩(wěn)態(tài)運行的基本模式之一,但高約束條件下伴隨的Ⅰ型邊界局域模(edge localized mode,ELM)[2]會引起等離子體邊緣區(qū)域臺基的周期性崩塌,崩塌過程中釋放的強粒子流和熱流會導(dǎo)致偏濾器熱負荷過載.ELM 現(xiàn)象通常被認為由剝離氣球模(peelingballooning mode,P-B 模)[3]激發(fā)產(chǎn)生,該不穩(wěn)定性由H 模放電時臺基區(qū)電流驅(qū)動的剝離模和壓強梯度驅(qū)動的氣球模耦合而成.

三角形變是等離子體形狀參數(shù)之一,在DIII-D[4],JET[5],ASDEX Upgrade[6]及其他托卡馬克裝置[7]上的H 模放電實驗發(fā)現(xiàn),正三角形變(positive triangularity,PT)位型放電時約束水平顯著提高.相關(guān)模擬表明,PT 可以解耦P-B 模,使臺基獲得更大壓強梯度和電流,從而有助于進入氣球模第二穩(wěn)定區(qū)[3,8].然而,在實驗和模擬中,負三角形變(negative triangularity,NT)位型具有截然不同的特性[9].一方面,NT 可以增加低約束模式(low confinement regime,L 模)向H 模轉(zhuǎn)換的閾值,進而可在L 模實現(xiàn)類似H 模的高比壓放電,同時不會出現(xiàn)ELM[9,10].另一方面,TCV 實驗發(fā)現(xiàn),NT位型下H 模放電的臺基壓強水平低于PT,但其ELM 頻率更高,幅值更小[11].此外,NT 也有其獨特的放電優(yōu)勢,如: 更大的外分界面、更靈活的偏濾器配置空間、更大的刮削層空間、更低的內(nèi)部極向線圈磁場需求,以及偏濾器室的更大泵送電導(dǎo)[12]等.相關(guān)模擬工作表明,NT 可使氣球模失穩(wěn),并且導(dǎo)致氣球模第二穩(wěn)定區(qū)通道完全關(guān)閉[12,13],但已有研究并未給出NT 位型下P-B 模的非線性演化特征,而P-B 模的非線性演化直接關(guān)系到ELM 爆發(fā)產(chǎn)生的熱流和粒子流大小,故本文對此進行了詳細研究.

本文通過Corsica 的TEQ 平衡模塊構(gòu)建了具有不同三角形變的托卡馬克平衡[14],并使用BOUT++代碼[15,16]對P-B 模不穩(wěn)定性進行數(shù)值模擬,模擬考慮了有限電阻效應(yīng)和逆磁效應(yīng),并分別討論了不同三角形變位型條件下P-B 模的線性和非線性特征.第2 節(jié)介紹了物理模型和托卡馬克平衡;第3 節(jié)詳細分析了線性和非線性模擬結(jié)果;最后在第4 節(jié)給出了結(jié)論.

2 物理模型和平衡參數(shù)

模擬使用了BOUT++三場模塊,在非圓截面平衡位型下求解簡化的三場雙流體方程組[15,16],這些方程包含了渦度U、壓強P和平行磁矢勢A//的演化.

方程中變量定義為

研究首先基于BOUT++模擬常用的cbm18_dens8 平衡構(gòu)造了一系列具有不同三角形變的平衡,該系列平衡基本參數(shù)與JET 裝置接近,其優(yōu)點是可以靈活設(shè)置托卡馬克實驗常見的形變參數(shù),而不用擔心數(shù)值誤差問題.具體參數(shù)如下: 小半徑a=1.2 m,大半徑R0=3.4 m,磁軸處環(huán)向磁場B0=2 T,等離子體電流IP=2 MA,拉長比κ=1.6,該拉長比參數(shù)屬于托卡馬克實驗常見參數(shù)范圍(如ASDEX Upgrade[6]和TCV[11]).通過CORSICA代碼TEQ 模塊[14]構(gòu)建一系列具有不同三角形變參數(shù)(- 0.3 ≤δ≤0)的托卡馬克平衡,如圖1 分別給出了無三角形變(δ=0)和NT 位型(δ=-0.3)的磁面形狀,圖中綠色區(qū)域為好曲率區(qū)域(?B2·?P ＞0),黃色區(qū)域為壞曲率區(qū)(?B2·?P ＜0)[9,18].模擬中所有平衡壓強分布相同,如圖2 所示,該壓強分布取自BOUT++模擬中常用的雙曲正切函數(shù).平衡通過Sauter 模型考慮了自舉電流的影響,自舉電流份額為0.5,同時保持總等離子體電流恒定[19-21],可以看出當臺基壓強不變時,自舉電流隨三角形變的變化并不明顯.平衡的芯部粒子密度取n0=3×1019m-3,密度分布由ne(ψn)=n0(P0(ψn)/P0(0))0.3給出,假設(shè)離子和電子的溫度和密度相等,即Te(ψn)=Ti(ψn)=P0(ψn)/2ne(ψn),其中ψn=(ψ-ψaxis)/(ψsep-ψaxis) 是歸一化磁通,ψaxis和ψsep分別是磁軸和最外閉合磁面處的磁通.

圖1 不同三角形變位型的磁面 (a) δ=0;(b) δ=-0.3,藍色實線代表 ψn=0.4 到1 的磁面,歸一化磁通間隔為ψn=0.1,其中 ψn=(ψ-ψaxis)/(ψsep-ψaxis),綠色區(qū) 域為壞曲率區(qū)域(? B2·?P ＞0),黃色區(qū)域為好曲率區(qū)(? B2·?P ＜0)[9,18],黑色虛線為中平面位置Fig.1.Comparison of magnetic surfaces with varied δ: (a)δ=0;(b) δ=-0.3.Blue lines represent the magnetic surfaces from ψn=0.4 to 1 with an interval of 0.1,ψn=(ψ-ψaxis)/(ψsep-ψaxis) is the normalized radial coordinate.The green areas show the unfavorable curvature regions where ? B2·?P ＞0 and yellow areas show the favorable curvature regions where ? B2·?P ＜0.Black dashed line shows the position of the midplane.

圖2 壓強剖面 P0 (黑色實線)和三角形變分別為 δ=-0.3(紅色實線),δ=0.0 (藍色點線)的平行電流剖面 J‖0Fig.2.The pressure P0 (black solid line) and parallel current J‖0 profiles for cases δ=0.0 (blue dotted line) and δ=0.3 (red solid line).

3 模擬結(jié)果分析

3.1 NT 對P-B 模線性增長率的影響

對于不同的三角形變(δ=-0.3—0)位型,通過P-B 模線性不穩(wěn)定性計算,圖3 給出了相應(yīng)的線性增長率模譜,由圖可知,隨著三角形變的降低,在中低n(n＜ 45)模區(qū)間,線性增長率增加;但在高n(n＞ 45)模區(qū)間,P-B 模不穩(wěn)定性受到抑制.由于磁面的形狀由極向電流磁線圈控制,因此不同的三角形變位型下平衡極向磁場有所不同.極向磁場的上下兩端受到反向電流磁線圈的作用而降低,磁面向弱場側(cè)移動.如圖1 所示,NT 位型的壞率區(qū)域(綠色區(qū)域)隨三角形變降低而增大,因為壞曲率區(qū)集中在弱場側(cè),具有更大的平均半徑,所以沿環(huán)向積分后具有更大的壞曲率空間,從而為PB 模提供更多自由能[9,22],故中低n模區(qū)間的PB 模線性增長率增大.與此同時,P-B 模不穩(wěn)定性的高n模部分還受到局域磁剪切的穩(wěn)定作用[23].定義外中平面局域磁剪切[21],其中,為局域安全因子,hθ為曲率半徑.圖4給出了不同三角形變位型的外中平面上局域磁剪切的徑向變化,由圖可知,隨三角形變降低,臺基區(qū)局域磁剪切的絕對值不斷增大,其對高n模的穩(wěn)定作用逐漸增強,故高n模區(qū)間的線性增長率受到抑制.

圖3 不同三角形變(δ=-0.3—0)位型的P-B 模線性增長率模譜Fig.3.Linear growth rates versus toroidal mode number for δ=-0.3 —0.

圖4 不同三角形變(δ=-0.3—0)位型外中平面上的局域磁剪切 sl 在徑向上的變化Fig.4.Profiles of local shear sl at the outer midplane for δ=-0.3 —0.

綜上所述,NT 位型引起的磁場結(jié)構(gòu)變化帶來更大的壞曲率區(qū)使P-B 模失穩(wěn),同時,外中平面上臺基區(qū)的局域磁剪切有助于穩(wěn)定高n模.該線性結(jié)果與已有模擬結(jié)果[9,13]一致.

3.2 P-B 模非線性演化特征

不同三角形變位型的托卡馬克磁場結(jié)構(gòu)不同,因此,由ELM 崩塌導(dǎo)致的湍流輸運過程特性也會有差異,所以,NT 位型下的P-B 模非線性模擬對未來的NT 位型實驗和對其中的ELM 物理機制的理解具有重要意義.在非線性模擬中,考慮了離子逆磁、有限電阻率和超電阻率效應(yīng),設(shè)置的初始擾動的環(huán)向模數(shù)為線性階段最不穩(wěn)定模式的模數(shù).ELM 能量損失通常用臺基能量損失(ΔWped)與臺基儲能(Wped)之比表示,定義為[16]

其中ψin是模擬的內(nèi)邊界;ψout是壓強梯度最大處.

接下來,通過δ=0 和 - 0.2 兩種位型的對比來闡明NT 位型的P-B 模非線性演化特征.圖5為不同三角形變位型下ELM 能量損失的對數(shù)值隨時間的演化,這里取ELM 能量損失的對數(shù)主要是為了更加清晰地展示ELM 崩塌的先后順序.由圖5 可知,隨三角形變降低,ELM 崩塌時間提前.圖6 為t=193τA時擾動壓強的環(huán)向平均在極向截面的分布,此時仍處于P-B 模的線性崩塌階段,可以發(fā)現(xiàn),隨著三角形變降低,極向截面上下兩端(圖中黑色虛線方框區(qū)域)的擾動壓強增強.由于NT 位型下壞曲率區(qū)域的擴展,P-B 模擾動在小截面上具有更大的增長區(qū)域,故其不穩(wěn)定性增長得更快,崩塌也發(fā)生得更早.

圖5 不同三角形變位型下ELM 能量損失的對數(shù)值隨時間的演化Fig.5.Time evolution of the logarithm of ELM size for different triangularity cases.

圖6 (a) δ=0 和(b δ=-0.2 位型下,t=193τA 時擾動壓強的環(huán)向平均在極向截面的分布.弱場側(cè)頂部和底部區(qū)域的黑色虛線框顯示了比較區(qū)域,黑色點線表示中平面位置Fig.6.Distribution of the toroidal-averaged pressure perturbation at the poloidal cross section at t=193τA for cases (a) δ=0 and (b) δ=-0.2.Black dashed frames at the top and bottom areas in the low field side show the regions for comparison.Black dotted line shows the position of the midplane.

圖7 分別給出了δ=0 和 - 0.2 兩種形變位型條件下擾動壓強在外中平面上的二維分布隨時間的演化.由于P-B 模的環(huán)向周期性,在BOUT++中只計算了圓環(huán)的四分之一區(qū)域.對比兩種位型下的擾動分布可以發(fā)現(xiàn),在t=100τA和t=200τA時,初始擾動始終占據(jù)主導(dǎo)地位,兩種位型下的擾動分布具有明顯的周期性;然而,在t=300τA時,無三角形變的壓強擾動仍基本保持初始擾動的周期,但NT 位型下的P-B 模擾動在環(huán)向的周期性則已被破壞,P-B 模具有更明顯的湍流輸運特性.

圖7 (a1)—(a3) δ=0 和(b1)—(b3) δ=-0.2 在 t=100, 200, 300τA 時在外中平面上的壓強擾動Fig.7.Pressure perturbation at t=100,200,300τA at the outer midplane for cases: (a1)—(a3) δ=0;(b1)-(b3) δ=-0.2.

通過對兩種位型下壓強梯度最大處的擾動做環(huán)向模式傅里葉分解可知,如圖8 所示,在ELM崩塌階段(t?200τA),初始擾動模式n=20(紅線)快速增長,臺基壓強擾動的負值部分向里移動(如圖7 所示),導(dǎo)致壓強臺基崩塌.隨著初始擾動模式的增長,其他環(huán)向模數(shù)的不穩(wěn)定模式也被激發(fā)起來,在湍流輸運階段(t?200τA),初始擾動模式幅值降低,n=0 模(藍線)快速增長,能量首先轉(zhuǎn)移到n=0 模,而n=0 模通常被認為是帶狀流,其與湍流相互作用發(fā)生能量耦合,但不會影響臺基儲能[24],隨后,部分低n模也相繼被激發(fā).不同之處在于,無三角形變位型下被激發(fā)的不穩(wěn)定模式增長并不明顯,n=20 的初始擾動模式仍然在較長時間里占主導(dǎo),故圖7(a3)中的壓強擾動仍然基本保持初始擾動的周期性;但在NT 位型下,如圖8(b),n=4和n=8 的模相繼增長,并依次占據(jù)主導(dǎo)地位,形成多種不穩(wěn)定模式共存的狀態(tài).在ELM 非線性模擬中[25],通常將擾動的這種演化狀態(tài)稱為湍流.相比于單一模式占主導(dǎo)的情況,由于模式間的相互競爭,初始擾動的周期性被破壞,故在NT位型下擾動表現(xiàn)出更加明顯的湍流輸運特性.

圖8 (a) δ=0 和 (b) δ=-0.2 位型下的環(huán)向模式演化Fig.8.Modes evolution for cases: (a) δ=0;(b) δ=-0.2.

4 總結(jié)與討論

本文考慮非理想效應(yīng)后,研究了NT 位型對P-B 模不穩(wěn)定性的影響,得到了P-B 模在NT 位型條件下的線性和非線性演化特征,可為進一步研究NT 與等離子體約束間的作用機制提供參考.

在線性階段,NT 位型壞曲率區(qū)域增加,導(dǎo)致氣球模驅(qū)動源增強,同時該位型中較大的局域磁剪切有利于穩(wěn)定高n模,故隨著NT 的減小,P-B 模的中低n模增長率增大,這使得NT 位型下的PB 模不穩(wěn)定性閾值降低,即邊緣臺基在較低水平時就有不穩(wěn)定性增長,引發(fā)臺基崩塌,從而限制了邊緣臺基的提升.相比之下,已有的對PT 的研究[3]發(fā)現(xiàn),PT 可以增加P-B 模不穩(wěn)定性閾值,進而提高等離子體約束.上述模擬結(jié)果與實驗[11]中觀察到的現(xiàn)象一致,即NT 位型下測得的邊緣臺基水平明顯低于PT 位型.

在非線性階段,隨著三角形變降低,較大的壞曲率區(qū)導(dǎo)致弱場側(cè)頂部和底部區(qū)域的壓強擾動增大,擾動分布區(qū)域的擴大有利于擾動的增長,減少了ELM 崩塌所需的時間,一定程度上有利于增加ELM 頻率;同時,NT 位型下多種環(huán)向不穩(wěn)定模式的增長,使得ELM 具有更明顯的湍流輸運特性,多模共存的湍流輸運狀態(tài)下的擾動傳播相對單一模式更弱[26],故其有利于ELM 能量損失的降低.值得一提的是,由于實驗中NT 與PT 位型導(dǎo)致的最明顯的差異之一在于邊緣臺基的不同,而臺基高度與ELM 頻率和幅值密切相關(guān),故下一步工作將基于真實托卡馬克參數(shù)(如: TCV,DⅢ-D),結(jié)合實驗測量的臺基剖面,分別模擬NT 和PT 位型對P-B 模的影響,以期為未來NT 位型下的高約束放電優(yōu)化提供理論依據(jù).