李正吉 陳偉? 孫愛萍 于利明 王卓 陳佳樂 許健強李繼全 石中兵 蔣敏 李永高 何小雪 楊曾辰 李鑒
1) (核工業西南物理研究院,成都 610041)
2) (中國科學院等離子體物理研究所,合肥 230031)
3) (61287 部隊,成都 610000)
可控核聚變裝置是利用輕原子核聚合形成重原子核過程中產生的能量進行發電的裝置,其能量轉換效率是傳統化石能源的百萬倍,并且原料儲量豐富,運行條件安全、無污染,成為了許多國家新能源發展的重點目標.到目前為止,托卡馬克裝置被認為是最有希望實現可控核聚變要求的實驗裝置.對于磁約束聚變裝置托卡馬克和球形環等,等離子體的歸一化環向比壓βN是標志磁約束聚變裝置約束性能和經濟性的重要指標,βN=βt/(Ip/aBt),其中,βt=Pt/(Bt/2μ0) 為等離子體環向比壓;a,Bt,Ip,Pt和μ0分別是裝置小半徑、環向磁場強度、等離子體電流、等離子體壓強和真空磁導率.高βN有三大優勢: 首先,高βN有利于等離子體達到點火條件,即達到聚變三乘積nTτE的閾值,其中,n為等離子體密度,T為等離子體溫度,τE為能量約束時間,H98為高約束因子,q95是安全因子在磁通量95%位置的值),當三乘積到達閾值時才能發生聚變反應;其次,由于βN與聚變功率成正比因此高βN有利于提高聚變反應的燃燒效率;最后,βN還與等離子體產生的自舉電流份額(fBS)成正比(即,其中ε是反環徑比)[2],而等離子體自舉電流份額越高,越有利于托卡馬克裝置實現穩態運行[3].在DIII-D[4–7],C-Mod[8],JT-60U[9,10],JET[11]等裝置上,通過對中性束(NBI)加熱、離子回旋加熱(ICRF)、電子回旋加熱(ECRH)、低雜波電流驅動(LHCD)等裝置的加熱和電流驅動調控,并結合對等離子體的參數及剖面優化,得到了離子或電子溫度的內部輸運壘(ITB)和邊緣輸運壘(ETB)同時存在的雙輸運壘(DTB)模式,進而得到了高βN等離子體.在HL-2A 裝置大功率NBI,以及NBI 和LHCD 協同作用的情況下實現了離子溫度的ITB 和ETB 共存的DTB 高級運行模式[12,13],等離子體的芯部最高離子溫度(Ti,0)達到2.5 keV,等離子體儲能(WE)達到了42 kJ,βN可以長時間穩定在2.5 以上,并且實現了瞬態βN≈3.0 的高性能等離子體.
由于大型托卡馬克裝置放電跨越多個時間和空間尺度,涉及多種物理過程,而目前大部分計算模型通常只考慮了某個時空尺度下對應的物理過程,所以其結果很難反映托卡馬克等離子體總體的真實情況,為了更好地描述托卡馬克等離子體,集成模擬這一概念便應運而生,集成模擬能結合多個時空尺度以及對應物理過程的代碼,并使用循環迭代的方式使它們得出自洽的結果,能很好地反映托卡馬克等離子體放電過程.OMFIT 程序(one modeling framework for integrated tasks)是由美國通用原子能公司(General Atom)開發的集成模擬平臺,其具有可視化界面,并且對于python腳本有良好的支持,并能讀寫多種科學文件類型,已被用于DⅢ-D,CFETR,HL-3 和EAST 等托卡馬克裝置的集成模擬工作[14–18].國內集成模擬程序的發展相對美國和歐洲而言起步較晚,對于HL-2A 裝置的相關集成模擬工作剛剛展開.本文對HL-2A 裝置上單獨NBI 加熱條件下的高βN雙輸運壘等離子體進行OMFIT 集成模擬.使用OMFIT平臺的EFIT 代碼計算平衡,通過輸運代碼ONE TWO 耦合中性束代碼NUBEAM 計算加熱,并利用TGYRO 代碼耦合湍流輸運代碼TGLF 和新經典輸運代碼NEO 來計算輸運流和E×B剪切.本文第2 節為HL-2A 裝置高βN的物理實驗現象介紹,第3 節通過OMFIT 集成模擬平臺對高βN物理實驗現象進行模擬,并進一步分析雙輸運壘形成的原因,第4 節為總結.
在HL-2A 裝置[19]單獨NBI 加熱條件下實現了βN長時間(t=550—700 ms)維持在2.5 以上,且在ITB 和ETB 共存的DTB 基礎上實現了瞬態βN達到約3.05 的高性能等離子體.圖1 給出了HL-2A 裝置上一次典型高βN放電(炮號27055)的主要參數隨時間的演化,自上而下分別為Ip、電子線平均密度ne,l、NBI 注入功率PNBI、氘阿爾法Dα光譜信號和WE,以及對應的高約束因子H98和βN等,這些物理量也是下文中OMFIT 集成模擬分析的實驗樣本.本次實驗采用圓形等離子體截面下單零偏濾器氘等離子體位形放電,等離子體電流范圍Ip=149—160 kA,環向磁場強度Bt≈1.3 T.NBI加熱采用了同電流方向注入,其束粒子能量Eb和PNBI分別約為40 keV 和0.85 MW.在NBI 注入前(t≈ 502 ms),等離子體參數為Ip≈0 kA,ne,l≈2.0×1019m-3,WE≈26 kJ,對應的βN≈1.7 .隨著NBI 的注入,等離子體參數大幅上升,即等離子體約束性能迅速提高,通過觀察Dα信號的轉變,可以看到在t=545 ms 時刻附近等離子體從低約束模式(L 模)向高約束模式(H 模)轉換,等離子體進入無邊緣局域模(ELM)階段.無ELM 的H 模約維持了13 ms (t=547—560 ms),隨著ne,l,WE等參數的持續增加,在t=561 ms 時刻,Dα開始出現大幅度的ELM 特征,并在ne,l和WE等信號上發現與ELM 對應的垮塌.在t=602 ms 時刻,主要等離子體參數為Ip≈152 kA,ne,l≈2.1×1019m-3,等離子體的WE達到了最高~46 kJ,與此對應的βN=3.05,H98=1.65和τE≈40—45 ms .隨著時間的推移,大ELM 的爆發和ne,l的持續升高等使WE逐漸降 低,對應的H98和βN也持續 下降.在t=675 ms 附近,等離子體參數趨于穩定,該時刻等離子體參數為Ip≈154 kA,ne,l≈2.1×1019m-3,WE≈43 kJ,H98≈1.39,βN≈2.83 .通過分析與計算,發現560—700 ms 之間的ELM 頻率在 80—140 Hz 之間,WE的損失為 6%—9%,因此判斷其為I 型ELM.

圖1 HL-2A 裝置27055 次高 βN 放電實 驗主 要參數隨時間的演化 (a)等離子體電流 Ip ;(b)電子線平均密度 ne,l ;(c) NBI 注入和 低混雜 波電流 驅動功率PNBI和PLHCD ;(d)氘 阿爾法 輻射 Dα 信號;(e)等 離子體 儲能 WE ;(f)約束品質因子 H98 ;(g) βNFig.1.Evolution of discharge parameters in HL-2A high βN plasmas in shot 27055: (a) Plasma current;(b) average electron density;(c) power of NBI and LHCD;(d) deuterium α signal;(e) stored energy;(f) confinement improvement factor;(g) normalized beta.
在HL-2A 裝置上,通過電荷交換復合光譜儀(CXRS)診斷系統測量Ti和等離子體環向旋轉頻率ft[20],使用電子回旋輻射計(ECE)診斷系統測量電子溫度Te[21],通過HCOOH 激光干涉儀[22]和調頻連續波(FMCW)微波反射儀[23]診斷系統分別測量芯部和臺基區電子密度ne.圖2(a)—(d) 藍色和紅色三角形數據點和曲線分別展示了通過上述診斷 系統得 到的在t=602 ms (βN=3.05)和675 ms (βN=2.83)兩個時刻的Ti,Te,ne和ft等主要等離子體參數的徑向分布圖,其中ne的剖面通過將芯部和邊緣測量得到的密度剖面耦合而成.可以明顯地觀察到,除了圖2(c)的ne分布之外,圖2(a),圖2(b)和圖2(d)給出的Ti,Te和ft等剖面數據在602 ms 時刻的數值均明顯地高于675 ms 時刻.并且,在兩個時刻的Ti和ft剖面芯部和邊緣同時觀測到了ITB 和ETB,ITB 主要集中在歸一化小半徑ρ 為0.38 附近,而ETB 的位置則在ρ=0.9 附近.在602 ms 時刻,Ti和ft擁有更強的ITB,芯部離子溫度Ti,0=1.9 keV 和芯部環向旋轉頻率ft,0=21.5 kHz,而 在675 ms 時刻的Ti,0=1.7 keV,ft,0=18.2 kHz.

圖2 HL-2A 裝置高 βN 放電等離子體參數徑向分布圖 (a)離子溫度 Ti ;(b)電子溫度 Te ;(c)電子密度 ne ;(d)等離子體環向旋轉頻率 ft .Exp. 表示通過診斷測量得到的實驗數據點,Fit. 表示通過擬合得到的剖面數據.在t=602 ms 和675 ms 時刻的實驗和擬合數據分別用藍色和紅色三角形和實線表示Fig.2.Profiles of the high βN discharge profiles of (a) Ti,(b) electron temperature,(c) electron density and (d) toroidal rotation frequency of plasmas.Exp. and Fit. mean the experimental data from diagnostics systems and the fitting curves.The colors of blue and red present the data from t=602 ms and 675 ms,the triangulares and solid curves mean the data in experiment and fitting.
使用OMFIT 集成模擬代碼對27055 次放電中βN分別為3.05 和2.83 的t=602 ms 和675 ms時刻(圖1 中藍色和紅色虛線所指示的時刻)的實驗參數進行了演算.在模擬過程中,使用的等離子體Ti,Te,ne和ft等主要實驗剖面如圖2 所示.等離子體的q分布、大半徑R、小半徑a,以及等離子體平衡等主要通過EFIT 代碼得到.表1 為實驗數據與OMFIT 模擬得到的宏觀參數的對比(如ne,l,WE,H98和βN,以及歸一化電 子密度(ne,l/ne,G,其中,ne,G=Ip/(πa2) 是Greenwald 密度),進一步通過OMFIT 計算得到的fBS和等離子體內感li等實驗無法直接得到的數據.通過對比發現: 通過OMFIT 模擬得到的高βN條件下兩個時刻的參數與實驗值基本一致,證明了該模擬計算的準確性.此外,還通過該模擬得到了在βN=3.05 和2.83情況下,li分別為0.97 和0.91,fBS已分別達到了46%和45%,意味著HL-2A 裝置有實現混合運行[24]或者穩態運行[1]等托卡馬克先進運行模式的可能,fBS詳細的分布情況如圖3(c)和圖3(d) 所示.

表1 27055 次放電在t=602 ms 和675 ms 時刻的實驗(Exp.)與OMFIT 集成模擬結果(Mod.)的宏觀參數對比Table 1.Comparisons between the experimental data and OMFIT simulation results at t=602 ms and 675 ms in shot 27055.

圖3 通過OMFIT 集成模擬得到的(a)等離子體壓強(Press)和q 剖面,(b) NBI 加熱功率沉積密度 pNBI,(c) 602 ms 時刻的等離子體 電流密 度剖面(j),(d) 675 ms 時刻的 等離子 體電流 密度剖 面(j).其 中,圖(a) 中的實線為602 ms 時刻的剖面,虛線為675 ms 時刻的剖面;圖(b) 中的實線和虛線分別表示為NBI 的能量沉積在離子和電子的能量密度分布;圖(c)和圖(d)中的 jTotal,jNBI,jOhm和jBS 分別為總電流密度分布、中性束驅動電流密度分布、歐姆驅動電流密度分布和自舉電流密度分布Fig.3.Integrated simulation results from OMFIT about distributions of (a) plasma pressure and q profile,(b) NBI power deposition,(c) 602 ms plasma current density,(d) 675 ms plasma current density.The solid line represents data of 602 ms and the dashed line represents data of 675 ms in panel (a).The solid and dashed curves present the NBI power deposited on ions and electrons in panel (b).jTotal,jNBI,jOhm and jBS means the total current density,current density driven by NBI,Ohmic and bootstrap current densities in panels (c) and (d).
通過OMFIT 集成模擬得到的等離子體壓強,q,NBI 加熱功率沉積密度pNBI和電流密度j等參數的詳細分布情況分別在圖3(a)—(d)進行了展示.其中,圖3(a)中的實線和虛線分別代表了t=602 ms 和675 ms 時刻的壓強和q剖面分布圖,可以看到,在602 ms 時刻的等離子體芯部的壓強要明顯高于675 ms 時刻的數值.從圖3(c)和圖3(d)的總電流分布jTotal可發現,由于等離子體電流的離軸分布,使得圖3(a)中602 ms 時刻和675 ms時刻都形成了q0>2且qmin≥1.5 的反磁剪切q剖面,并且由于602 ms 時刻的總電流分布更加峰化,其qmin≈1.5 比675 ms 時刻的qmin(>1.5) 更低.兩個放電時刻fBS≈45%—46% 接近,但NBI 驅動電流的比例在602 ms 和675 ms 分別為 33%和28% .這可能與NBI 加熱過程中在兩個時刻的等離子體密度分布有關,如圖1(b)和圖2(c)所示,低密度等離子體有利于NBI 電流驅動.從圖3(b)可以發現,NBI的功率主要沉積于離子,沉積位置集中在歸一化小半徑ρ=0.38 以內,這與離子溫度ITB 的位置一致.對比圖3(a)—(d)四幅圖可以發現,NBI 的驅動電流、能量沉積剖面的最大值位置,與總電流分布最大值位置以及q剖面最小值的位置是重合的(在歸一化小半徑ρ=0.25 左右的位置),可見NBI 驅動電流誘導了反磁剪切的產生.因此,NBI主要加熱位置的離軸(主要加熱位置在ρ=0.25 左右)形成了離子溫度ITB 與q剖面的反磁剪切.
雖然602 ms 時刻和675 ms 時刻的等離子體擁有相同的NBI 加熱功率,以及相近的離子功率沉積剖面,但相較于675 ms 時刻,602 ms 時刻的離子溫度卻形成了更強的ITB,我們認為這是由于NBI 產生的快離子對等離子體輸運的影響造成了602 ms 時刻與675 ms 時刻等離子體剖面的差別[25,26].因此,本文使用了TGYRO 耦合湍流輸運代碼TGLF 和新經典輸運代碼NEO,對等離子體歸一化小半徑ρ=0—0.8 以內能量流的湍性輸運和新經典輸運模式進行了計算,并分析了沒有快離子情況下的輸運,計算模型是靜電模式,并包含了E×B剪切.圖4 從上到下依次展示了快離子比壓分布(βf)、離子湍性輸運能量流分布(Qi,Tur)、離子新經典輸運能量流分布(Qi,Neo),左邊的圖4(a)、圖4(c)、圖4(e)三幅圖和右邊的圖4(b)、圖4(d)、圖4(f)三幅圖分別為602 ms 時刻和675 ms 時刻對應的參數.圖4中βf=Pf/(Bt/2μ0),Pf為快離子壓強,G.B.為回旋波姆單位(GyroBohm units).從圖4(a)和圖4(b)可看到,602 ms 時刻的快離子比壓βf顯著強于675 ms 時刻.對比圖4(c)、圖4(d)和圖4(e)、圖4(f),上述兩個不同時刻對應的湍性輸運流均比新經典輸運流強度高出兩個量級,可見等離子體輸運以湍性輸運為主.此外,通過圖4(c)和圖4(d)可以看到,無論是602 ms 時刻,亦或是675 ms 時刻,有快離子時的離子熱輸運強度整體是弱于沒有快離子的,可見,NBI 產生的快離子對于湍流輸運有抑制作用.

圖4 快離子與能量輸運流示意圖 (a),(c),(e) 分別為602 ms 時刻的快離子分布、湍性輸運能量流、新經典輸運能量流;(b),(d),(f) 分別為675 ms 時刻的快離子分布、湍性輸運能量流、新經典輸運能量流.圖中紅線表示沒有快離子,藍線表示有快離子Fig.4.Fast ion and energy flux: (a),(c),(e) Fast ion,energy flux of turbulent transport,energy flux of neoclassic transport for 602 ms;(b),(d),(f) fast ion,energy flux of turbulent transport,energy flux of neoclassic transport for 675 ms;the red line means with fast ions,the blue line means with out fast ions.
對于快離子抑制湍流輸運的物理機制,這里討論線性與非線性兩種情況.首先對于快離子抑制湍流輸運的線性機制,目前認為有4 個可能的因素[27]:1)高能量粒子由于其軌道遠大于湍流的特征垂直波長,與湍流的相互作用較小.因此,高能量粒子可以稀釋熱離子的比重,從而減少了熱離子對湍流的驅動作用,進而表現為高能量粒子對離子溫度梯度不穩定性(ITG 模)起致穩作用.2)高能量粒子的存在可以增加局域的等離子體壓強梯度,從而改變局域的磁面距離,致使離子的磁場梯度漂移反向,從而抑制由于“壞曲率”驅動的ITG 不穩定性.3)與ITG 發生共振的具有相對較低能量的高能量粒子分布函數具有負梯度時,高能量粒子的朗道阻尼可能會顯著地抑制ITG 的不穩定性.4)高能量粒子經常伴隨的電磁效應,會引起磁力線彎曲而增加自由能,進而能夠顯著地抑制ITG 不穩定性.可見,快離子抑制湍流輸運的線性機制,主要就在于快離子對ITG 湍流的稀釋與抑制.對于快離子抑制湍流輸運的非線性機制,目前認為可能是微觀湍流與阿爾芬不穩定性共同驅動的帶狀流結構強化了剪切對湍流的抑制作用[28],但對這一非線性問題的理解還不夠充分,仍需進一步的研究.
離子熱輸運被快離子抑制之后,等離子約束得到改善,從而使得602 ms 時刻和675 ms 時刻的等離子體得到了更高的離子溫度剖面.但從圖4(c)和圖4(d)可以發現,無論是否存在快離子,在ITB以內(ρ<0.38) 602 ms 時刻的湍流輸運強度顯著弱于675 ms 時刻,可見除了NBI 產生的快離子,還有其他因素在影響輸運,考慮到圖2(d)中602 ms時刻比675 ms 時刻更峰化的旋轉頻率剖面以及輸運模型包含了E×B剪切效應,因此602 ms 時刻較弱的湍性能量流很可能是E×B剪切對湍流輸運產生的影響所導致的.為此,本文使用湍流輸運代碼TGLF 的回旋朗道流體輸運模塊[29]計算了E×B相關參數.圖5 展示了歸一化的E×B旋轉剪切率aγE/cs在歸一化小半徑ρ=0—0.8 的分布,其中γE≈(r/q)(dqVE×B·r)/dr是E×B旋轉剪切率,cs為離子聲波,VE×B為E×B剪切速度.由于E×B剪切流能夠有效地通過剪切作用降低湍流結構的相關長度從而抑制湍流輸運[30,31],從圖5 可以看到,無論是否有快離子的存在,602 ms時刻的 歸一化 旋轉剪切率aγE/cs整體都 是強于675 ms 時刻,尤其在歸一化小半徑ρ≤0.38 的ITB范圍以內.可見圖4(c)和圖4(d)中,ITB 以內(ρ<0.38) 602 ms 時刻的湍流輸運強度峰值顯著小于675 ms 時刻,是由于芯部更強的E×B剪切的存在,更好地抑制了芯部的湍流輸運,使602 ms 時刻擁有更強的ITB,進而得到了更高的等離子體比壓.

圖5 E×B 剪切率,圖中藍線為602 ms 時刻的剖面,紅線為675 ms 時刻的剖面,而虛線線表示沒有快離子,實線表示有快離子.Fig.5.E×B shear rate.The blue line represents data of 602 ms and the red line represents data of 675 ms,and the the dashed line means with out fast ions,the solid line means with fast ions.
本文討論的HL-2A 裝置27055 次放電實驗,在NBI 加熱條件下,獲得了穩定的高βN等離子體(βN>2.5,持續時間t≈150 ms),并且實現了瞬態βN=3.05,ne/ne,G≈0.6,H98≈1.65,fBS≈46%的高約束性能.使用OMFIT對βN=2.83和βN=3.05 兩種情況下的等離子體進行了集成模擬,得到了電流密度分布、NBI 能量沉積剖面、壓強和安全因子分布,并與實驗的宏觀參數進行了對比,模擬結果能夠較為準確地還原實驗參數.對輸運流進行了計算和分析,得到了離子溫度ITB 形成的原因:NBI 加熱使等離子體進入H 模,高βN的輸運以湍流輸運為主,而快離子作用和E×B剪切使得等離子體芯部的湍流輸運被抑制,從而形成了離子溫度ITB.并且,在H 模條件下,ITB 與ETB 相互協同形 成了高βN的DTB 離子溫 度剖面.對 比602 ms 時刻和675 ms時刻,可以看到602 ms 時刻在更低的密度條件下得到了更高的βN,所以我們進一步的工作是在同樣NBI 加熱下,計算和尋找在什么樣的等離子體密度區間內,能得到更強的溫度ITB、更高的穩態βN和更大的自舉電流份額,為HL-2M 裝置高比壓實驗提供依據.