秦晨晨 牟茂淋? 陳少永

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

2) (四川大學(xué),高能量密度物理國(guó)家重點(diǎn)實(shí)驗(yàn)室,成都 610065)

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

1 引言

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

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

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

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

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

方程中變量定義為

研究首先基于BOUT++模擬常用的cbm18_dens8 平衡構(gòu)造了一系列具有不同三角形變的平衡,該系列平衡基本參數(shù)與JET 裝置接近,其優(yōu)點(diǎn)是可以靈活設(shè)置托卡馬克實(shí)驗(yàn)常見的形變參數(shù),而不用擔(dān)心數(shù)值誤差問題.具體參數(shù)如下: 小半徑a=1.2 m,大半徑R0=3.4 m,磁軸處環(huán)向磁場(chǎng)B0=2 T,等離子體電流IP=2 MA,拉長(zhǎng)比κ=1.6,該拉長(zhǎng)比參數(shù)屬于托卡馬克實(shí)驗(yàn)常見參數(shù)范圍(如ASDEX Upgrade[6]和TCV[11]).通過CORSICA代碼TEQ 模塊[14]構(gòu)建一系列具有不同三角形變參數(shù)(- 0.3 ≤δ≤0)的托卡馬克平衡,如圖1 分別給出了無(wú)三角形變(δ=0)和NT 位型(δ=-0.3)的磁面形狀,圖中綠色區(qū)域?yàn)楹们蕝^(qū)域(?B2·?P ＞0),黃色區(qū)域?yàn)閴那蕝^(qū)(?B2·?P ＜0)[9,18].模擬中所有平衡壓強(qiáng)分布相同,如圖2 所示,該壓強(qiáng)分布取自BOUT++模擬中常用的雙曲正切函數(shù).平衡通過Sauter 模型考慮了自舉電流的影響,自舉電流份額為0.5,同時(shí)保持總等離子體電流恒定[19-21],可以看出當(dāng)臺(tái)基壓強(qiáng)不變時(shí),自舉電流隨三角形變的變化并不明顯.平衡的芯部粒子密度取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,藍(lán)色實(shí)線代表 ψn=0.4 到1 的磁面,歸一化磁通間隔為ψn=0.1,其中 ψn=(ψ-ψaxis)/(ψsep-ψaxis),綠色區(qū) 域?yàn)閴那蕝^(qū)域(? B2·?P ＞0),黃色區(qū)域?yàn)楹们蕝^(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 壓強(qiáng)剖面 P0 (黑色實(shí)線)和三角形變分別為 δ=-0.3(紅色實(shí)線),δ=0.0 (藍(lán)色點(diǎn)線)的平行電流剖面 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 對(duì)P-B 模線性增長(zhǎng)率的影響

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

圖3 不同三角形變(δ=-0.3—0)位型的P-B 模線性增長(zhǎng)率模譜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 位型引起的磁場(chǎng)結(jié)構(gòu)變化帶來更大的壞曲率區(qū)使P-B 模失穩(wěn),同時(shí),外中平面上臺(tái)基區(qū)的局域磁剪切有助于穩(wěn)定高n模.該線性結(jié)果與已有模擬結(jié)果[9,13]一致.

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

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

其中ψin是模擬的內(nèi)邊界;ψo(hù)ut是壓強(qiáng)梯度最大處.

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

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

圖6 (a) δ=0 和(b δ=-0.2 位型下,t=193τA 時(shí)擾動(dòng)壓強(qiáng)的環(huán)向平均在極向截面的分布.弱場(chǎng)側(cè)頂部和底部區(qū)域的黑色虛線框顯示了比較區(qū)域,黑色點(diǎn)線表示中平面位置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 兩種形變位型條件下擾動(dòng)壓強(qiáng)在外中平面上的二維分布隨時(shí)間的演化.由于P-B 模的環(huán)向周期性,在BOUT++中只計(jì)算了圓環(huán)的四分之一區(qū)域.對(duì)比兩種位型下的擾動(dòng)分布可以發(fā)現(xiàn),在t=100τA和t=200τA時(shí),初始擾動(dòng)始終占據(jù)主導(dǎo)地位,兩種位型下的擾動(dòng)分布具有明顯的周期性;然而,在t=300τA時(shí),無(wú)三角形變的壓強(qiáng)擾動(dòng)仍基本保持初始擾動(dòng)的周期,但NT 位型下的P-B 模擾動(dòng)在環(huán)向的周期性則已被破壞,P-B 模具有更明顯的湍流輸運(yùn)特性.

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

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

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

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