張洪銘 吳靜 李佳鮮 姚列明 徐江城吳巖占 劉旗艷 郭鵬程

1) (電子科技大學(xué)物理學(xué)院,成都 611731)

2) (核工業(yè)西南物理研究院聚變研究所,成都 610225)

磁約束托卡馬克裝置等離子體物理參數(shù),如等離子體電流、約束磁場(chǎng)、安全因子的實(shí)時(shí)監(jiān)測(cè)對(duì)高約束模式(H 模)的穩(wěn)態(tài)運(yùn)行至關(guān)重要.本文介紹了基于HL-2A 實(shí)驗(yàn)數(shù)據(jù)和集成建模手段,重建先進(jìn)高約束運(yùn)行模式(H 模)下等離子體位形,電流密度和安全因子參數(shù)剖面.通過(guò)搭建動(dòng)力學(xué)平衡位形集成模擬平臺(tái),結(jié)合工作流快速處理手段,采用真實(shí)實(shí)驗(yàn)數(shù)據(jù)和集成模擬模型,得到了HL-2A 高約束模式放電條件的離子和電子溫度、密度,電流密度和安全因子q 剖面,重建了高約束模式放電物理圖像和托卡馬克內(nèi)部磁面位形和等離子體邊界參數(shù)連續(xù)分布,計(jì)算了等離子體歐姆電流、自舉電流和射頻電流等剖面及各自份額,給出了先進(jìn)運(yùn)行模式下安全因子q 徑向分布.通過(guò)集成模擬,發(fā)現(xiàn)HL-2A 該炮號(hào)參數(shù)條件下臺(tái)基寬度約7.52 cm,壓強(qiáng)梯度的方向在ρ(r/a)=0.1 改變,在ρ=0.7 附近數(shù)值達(dá)到最大,可能是負(fù)剪切產(chǎn)生的內(nèi)部輸運(yùn)壘位形(ITB)引起的.通過(guò)等離子體磁場(chǎng)和電流等剖面重建和實(shí)時(shí)監(jiān)測(cè),可以評(píng)估HL-2A 高約束放電實(shí)驗(yàn)的品質(zhì),為HL-2M 高比壓（βn）先進(jìn)模式的穩(wěn)態(tài)運(yùn)行提供參考.

1 引言

在磁約束聚變裝置中,高約束模式[1](H 模)相比于低約束模式(L 模),能大幅度提高等離子體性能[2].在高約束模式下,等離子體達(dá)到同等水平的運(yùn)行參數(shù)所需的裝置規(guī)模要小,可大大提高未來(lái)聚變能的經(jīng)濟(jì)性.磁約束托卡馬克出現(xiàn)H 模首先在德國(guó)的磁約束托卡馬克裝置ASDEX 發(fā)現(xiàn),一直是各個(gè)聚變裝置的重要研究?jī)?nèi)容,后來(lái)通過(guò)不斷開(kāi)展放電物理實(shí)驗(yàn)及理論模擬研究,逐漸掌握低約束模式到高約束模式(L-H 模)的轉(zhuǎn)換機(jī)理及對(duì)等離子體性能的影響.H 模放電的典型特征之一是壓強(qiáng)剖面(溫度、密度相關(guān))在等離子體邊界附近具有一定臺(tái)基,即邊界附近壓強(qiáng)剖面突然變陡峭,壓強(qiáng)梯度迅速增大.臺(tái)基的高度和寬度與等離子體性能緊密相關(guān),臺(tái)基高度越高,也更容易得到芯部參數(shù)更高的溫度密度剖面.臺(tái)基高度和寬度受到動(dòng)力氣球模和剝離氣球模不穩(wěn)定性的限制[3],可由EPED程序進(jìn)行預(yù)測(cè)[4,5].邊界附近壓強(qiáng)梯度的變化較大產(chǎn)生較大自舉電流,影響等離子體電流的剖面分布,進(jìn)而影響安全因子q的分布.電流剖面和安全因子q徑向分布是研究H 模以及比H 模先進(jìn)運(yùn)行模式(參數(shù)更高)的重要分析指標(biāo).但是電流分布,尤其是自舉電流分布及安全因子q分布、芯部數(shù)據(jù),難以從實(shí)驗(yàn)診斷實(shí)時(shí)獲取,一般需要通過(guò)輸運(yùn)模擬程序模擬分析得到.研究人員需要選取物理模型(粒子和流體模型),進(jìn)行集成建模,將診斷獲取的參數(shù)進(jìn)行整合,重建電流分布剖面和安全因子剖面,對(duì)托卡馬克穩(wěn)態(tài)運(yùn)行進(jìn)行實(shí)時(shí)監(jiān)控.

磁約束聚變研究涉及很多前沿科學(xué)和工程技術(shù)難題,在磁約束托卡馬克裝置中,放電實(shí)驗(yàn)主要依靠診斷設(shè)備獲取相關(guān)信號(hào)(等離子體溫度和密度,電子溫度和密度,旋轉(zhuǎn)速度等),并轉(zhuǎn)換為便于理解和分析的實(shí)驗(yàn)數(shù)據(jù).由于診斷技術(shù)、診斷空間、經(jīng)費(fèi)投入等因素,使得實(shí)驗(yàn)獲取的診斷數(shù)據(jù)難以滿(mǎn)足放電實(shí)驗(yàn)實(shí)時(shí)控制需求.比如,溫度密度剖面數(shù)據(jù),長(zhǎng)期難以得到完整的剖面數(shù)據(jù),由于分辨率低,獲取的剖面無(wú)法反描述H 模放電的臺(tái)基參數(shù)等.因此,需要利用有限的實(shí)驗(yàn)數(shù)據(jù)和理論模型分析,建立綜合的集成模擬平臺(tái),獲取更豐富的放電信息,更加全面評(píng)估放電品質(zhì).

美國(guó)通用原子公司(General Atom,GA)開(kāi)發(fā)了基于DIII-D 實(shí)驗(yàn)數(shù)據(jù)的磁面重建程序EFIT[6,7],可以獲得磁場(chǎng)位形,利用外部的電磁診斷設(shè)備(磁通環(huán)和磁探針)獲取磁信號(hào),對(duì)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行全面深入分析,準(zhǔn)確重建等離子體平衡位形和相關(guān)的物理參數(shù).由于缺乏芯部診斷數(shù)據(jù),EFIT 平衡反演程序難以給出準(zhǔn)確的電流及壓強(qiáng)剖面分布.可以通過(guò)增加內(nèi)部診斷手段,利用可見(jiàn)光譜的動(dòng)態(tài)斯塔克效應(yīng)(motional Stark effect,MSE)得到等離子體電流剖面分布[8]、等離子體壓強(qiáng)(可見(jiàn)光譜診斷等手段得到的等離子體離子密度和溫度)的診斷數(shù)據(jù),利用EFIT 完成了更為準(zhǔn)確的磁面重建和剖面參數(shù)分析.EFIT 應(yīng)用于磁約束托卡馬克裝置,成為分析實(shí)驗(yàn)數(shù)據(jù)的主要工具,可以重建等離子體邊界和磁場(chǎng)位形參數(shù),也能給出等離子體壓強(qiáng)、電流密度及安全因子分布剖面.隨著聚變研究的深入,基于EFIT 程序無(wú)法提供更準(zhǔn)確的剖面分布,尤其是不同電流份額及分布(比如自舉電流)條件下的磁面重建結(jié)果,難以滿(mǎn)足高約束模式H 模,尤其是先進(jìn)運(yùn)行模式的深入研究需求.

托卡馬克裝置DIII-D 發(fā)展了集成模擬框架程序OMFIT (one modeling framework for integrated tasks)[9-11].通過(guò)開(kāi)發(fā)工作流,根據(jù)目的將不同程序進(jìn)行集成運(yùn)行,可以開(kāi)展動(dòng)力學(xué)平衡位形重建,磁流體穩(wěn)定性研究,高次諧波快波電流驅(qū)動(dòng)優(yōu)化,磁流體和微觀湍流的相互作用以及自洽的平衡輸運(yùn)模擬等物理內(nèi)容研究.DIII-D 托卡馬克裝置研究人員提出設(shè)想并使用該程序進(jìn)行了實(shí)踐工作,充分顯示了OMFIT 框架在托卡馬克聚變研究的眾多領(lǐng)域的超強(qiáng)能力[12].

HL-2A 裝置作為國(guó)內(nèi)實(shí)現(xiàn)偏濾器位形高約束模式(H 模)放電實(shí)驗(yàn)的托卡馬克裝置,使中國(guó)成為歐美日之后成功實(shí)現(xiàn)高約束模式運(yùn)行的國(guó)家[13,14].HL-2A 具有相對(duì)較高的加熱/驅(qū)動(dòng)功率,開(kāi)展了與ITER(國(guó)際熱核聚變實(shí)驗(yàn)堆)物理相關(guān)研究.該研究基于HL-2A 高約束模式相關(guān)物理實(shí)驗(yàn),利用OMFIT 集成模擬平臺(tái),搭建了包含EFIT,ONETWO[15,16],TGYRO[17,18]等核心程序的具備開(kāi)展動(dòng)力學(xué)平衡位形重建的工作流,通過(guò)實(shí)驗(yàn)診斷數(shù)據(jù)給出的離子/電子的溫度/密度剖面,開(kāi)展了高約束模式(H 模)放電實(shí)驗(yàn)的集成模擬分析,重建了HL-2A高約束模式下等離子體的磁面位形,給出了準(zhǔn)確的等離子體邊界,提供了完整的溫度演化、密度剖面以及安全因子剖面.基于以上物理參數(shù),本文分析了等離子體電流成分,分別給出歐姆電流、自舉電流、驅(qū)動(dòng)電流占總電流的份額和剖面分布.基于OMFIT 模擬程序和HL-2A H 模實(shí)驗(yàn)數(shù)據(jù),考慮物理模型自洽性,基于物理量之間的數(shù)據(jù)智能傳遞耦合,得到了托卡馬克等離子體電流剖面和磁面的重建.本研究結(jié)合HL-2A 的實(shí)驗(yàn)結(jié)果、整合多個(gè)實(shí)驗(yàn)反演程序和模擬程序進(jìn)行集成建模分析(OMFIT),整合托卡馬克裝置大量的物理數(shù)據(jù),得到實(shí)時(shí)托卡馬克裝置物理整體剖面參數(shù)并進(jìn)行分析.

本文首先介紹理論分析模型和建立集成模擬程序OMFIT 框架及文中研究所需主要模塊,描述HL-2A 高約束放電集成模擬工作流設(shè)計(jì),然后選取HL-2A 實(shí)驗(yàn)高約束模式診斷物理結(jié)構(gòu)參數(shù),介紹了基于OMFIT 框架的整合,對(duì)實(shí)驗(yàn)數(shù)據(jù)和動(dòng)力學(xué)平衡反演程序計(jì)算,介紹了重建后的HL-2A磁面位形,電流(包含自舉電流份額)和安全因子q剖面,最后通過(guò)等離子體磁場(chǎng)和電流等剖面重建和實(shí)時(shí)監(jiān)測(cè),評(píng)估了HL-2A 高約束放電實(shí)驗(yàn)的品質(zhì).

2 集成模擬模型

2.1 理論模型

作為綜合模擬平臺(tái),集成模擬程序OMFIT 能夠用于各種物理程序之間的耦合組成復(fù)雜的流程模塊以分析物理問(wèn)題,其特點(diǎn)在于其對(duì)托卡馬克實(shí)驗(yàn)數(shù)據(jù)的實(shí)時(shí)整合.目前該平臺(tái)已經(jīng)整合了物理模擬程序,如EFIT,EPED,ONETWO和TGYRO 等.在自洽的平衡與輸運(yùn)預(yù)測(cè)模型基礎(chǔ)上,利用OMFIT集成模擬平臺(tái),搭建了包含EFIT,ONETWO,TGYRO 等核心程序的具備開(kāi)展動(dòng)力學(xué)平衡位形重建的建模過(guò)程,采用智能處理手段工作流針對(duì)HL-2A 高約束模式H 模放電進(jìn)行了研究.能量源和能量損失用平衡反演程序ONETWO 可以計(jì)算,包括輔助加熱和輻射損耗以及離子和電子之間的能量交換.湍流輸運(yùn)和新經(jīng)典輸運(yùn)特性用模擬代碼TGYRO,TGLF[19]和NEO[20]分別計(jì)算,得到能量、粒子和環(huán)向角動(dòng)量的輸運(yùn)通量.TGYRO 通過(guò)匹配能量通量和目標(biāo)通量來(lái)計(jì)算穩(wěn)態(tài)溫度分布,目標(biāo)通量由相關(guān)的能量源和損失項(xiàng)的體積積分給出,磁場(chǎng)分布由二維(2-D)EFIT 平衡反演程序求解.

本文采用固定邊界求解理想磁流體力學(xué)平衡Grad-Shafranov(G-S)方程來(lái)獲得平衡剖面演化特性.磁流體方程組為(1)式—(5)式,其中(1)式為粒子守恒方程,(2)式為動(dòng)量守恒方程,(3)式為能量守恒方程,(4)式為等離子體在電磁場(chǎng)中的安培定律,(5)式為法拉第定律.(6)式為受力平衡條件.以上為理想條件下的磁流體MHD 方程組,用來(lái)描述磁約束等離子體中范圍較大的粒子運(yùn)動(dòng)狀態(tài).等離子體的每個(gè)組分的粒子處于平衡態(tài),粒子的速度分布函數(shù)基本符合麥克斯韋分布,可以描述等離子體在電磁場(chǎng)中的平衡分布:

模擬條件假設(shè)為:磁流體重力相對(duì)于電磁力的約束很小,等離子體在磁場(chǎng)中達(dá)到平衡條件:?t v=0,(2) 式可以寫(xiě)成(7)式,假設(shè)等離子體的宏觀運(yùn)動(dòng)速度為0 m/s,即v=0 m/s,寫(xiě)成(8)式,分別用J和B來(lái)乘等式兩邊,可以得到(10)式:

磁場(chǎng)的洛倫茲力和熱壓力之間的平衡可以從(8) 式得到,在輸運(yùn)的尺度上作為力的平衡方程,包含電流密度、壓力和電磁力之間的作用.在時(shí)間尺度上作為與時(shí)間無(wú)關(guān)的平均速度為零的平衡理想磁流體方程的解.通過(guò)求解該方程可以得到等離子體的位形(包含壓強(qiáng)、磁場(chǎng)、粒子密度、溫度等)以及約束等離子體的外場(chǎng)(磁場(chǎng),電場(chǎng))隨空間坐標(biāo)變化.由 (9) 式可以得到壓強(qiáng)相等的磁面輪廓,由 (10) 式可得到磁面上的電流密度及徑向電流導(dǎo)致的電磁場(chǎng)力平衡方程.從動(dòng)量守恒方程出發(fā),假設(shè)等離子體的宏觀運(yùn)動(dòng)速度趨近于平衡,推導(dǎo)出理想條件下的MHD 平衡方程.考慮壓力張量是非負(fù)的,得到的平衡物理量能描述線性理想磁流體不穩(wěn)定性.作為壓力、電磁作用和電流密度之間的平衡的解,從中能判斷出約束等離子體的外場(chǎng)條件和絕對(duì)等離子體的形狀.在柱坐標(biāo)系f(R,φ,Z)中,用函數(shù)F(F=RBf)和磁通函數(shù)ψ來(lái)表達(dá)磁場(chǎng)B.φ作為大環(huán)方向的角度,ψ作為等離子體的極向磁通被準(zhǔn)確的用來(lái)描述磁面,其作為封閉環(huán)形磁面內(nèi)的極向磁通量的總和.在磁約束托卡馬克等離子體,圓截面的等離子體方程通常采用柱坐標(biāo)系,方向之間有如(11)式 的關(guān)系.R為徑向參數(shù),φ為角向參數(shù),Z為軸向參數(shù).由磁場(chǎng)B無(wú)旋性,如(12)式所示.考慮磁流體的對(duì)稱(chēng)性,那么在軸對(duì)稱(chēng)平衡中,將 (11) 式兩邊同時(shí)求散度.由于軸對(duì)稱(chēng)平衡中 0=?·(?ψ×?φ),(13) 式變?yōu)?14)式:

考慮到對(duì)稱(chēng)性,標(biāo)量函數(shù)F和ψ都與環(huán)向角度φ無(wú)關(guān),作為Z和R的函數(shù).將求解磁場(chǎng)B的(11) 式兩邊同時(shí)求旋度,考慮?×(?T)=0和?×(?×t)=?(?·t)-?2t,(15) 式可以推導(dǎo)出(16)式,代入安培定律,得到等離子體電流密度J.從(8)式—(10)式可知,由于磁面是平衡的等壓強(qiáng)磁面,等離子體壓強(qiáng)P只和極向磁通函數(shù)ψ有關(guān),即P=P(ψ).由于P的函數(shù)只考慮ψ,由(10)式得到(18)式,將(17)式中的等離子體電流密度代入(18)式,于是得到(19)式.通過(guò)前面的推導(dǎo),已經(jīng)得到J和B的式,代入MHD 式中.通過(guò)求解磁流體平衡方程,可以得到等離子體壓強(qiáng)與等離子體電流分布,得到平衡位形參大小數(shù),可在托卡馬克裝置中描述等離子體平衡電流位形和磁場(chǎng)位形分布:

將(19)式和(20)式代入(8)式中可得:

其中(?φ·?φ=R-2).比較等式(19),得到MHD磁流體平衡方程:

2.2 HL-2A 高約束放電集成模擬

基于磁流體平衡方程(23),將HL-2A 初始條件和邊界條件下,電子溫度、離子溫度和電子密度等進(jìn)行擬合處理,代入模擬程序,具體的工作流如圖1所示.采用EFIT 程序計(jì)算得到HL-2A 平衡位形重建的平衡物理參數(shù),代入集成模擬程序,將實(shí)驗(yàn)診斷數(shù)據(jù)(電子溫度Te,離子溫度Ti,電子密度ne和旋轉(zhuǎn)信息Vr)導(dǎo)入動(dòng)力學(xué)剖面計(jì)算程序,計(jì)算外部驅(qū)動(dòng)電流、歐姆電流和自舉電流,給出演化電流剖面信息,計(jì)算源和損失項(xiàng).源項(xiàng)包括粒子源、能量源、外部電流驅(qū)動(dòng)電源和動(dòng)量源,損失項(xiàng)為輻射損失.粒子源主要來(lái)自中性束注入氫及同位素、充氣和壁材料濺射.能量源主要包括輔助加熱,例如中性束和波加熱.動(dòng)量源主要來(lái)自非對(duì)稱(chēng)中性束注入引發(fā)的等離子體旋轉(zhuǎn).最終得到源和損失項(xiàng)的參數(shù)以及新的平衡參數(shù)、壓強(qiáng)信息與所有電流份額組成的平衡電流剖面信息.

圖1 HL-2A 高約束模式(H 模)放電模擬循環(huán)示意圖Fig.1.Design and analysis of HL-2A high confinement model (H model) discharge integrated simulation workflow.

當(dāng)輸出端存在電流演化時(shí),利用EFIT 反演平衡參數(shù)剖面來(lái)更新的平衡參數(shù)剖面,即壓力梯度和極向電流梯度,得到新的平衡參數(shù),包含壓強(qiáng)信息與所有電流份額組成的平衡電流剖面信息,代入到EFIT 程序,求解G-S 方程,得到自洽的二維平衡方程.將自洽的平衡和源項(xiàng)的參數(shù)輸入到計(jì)算程序,利用源和損失項(xiàng)以及更新的平衡參數(shù),進(jìn)行動(dòng)理學(xué)剖面的演化,采取臺(tái)基區(qū)固定,演化芯部的動(dòng)力學(xué)剖面,得到新的芯部的動(dòng)理學(xué)剖面,繼而將新平衡動(dòng)力學(xué)剖面(電子溫度Te,離子溫度Ti,電子密度ne和旋轉(zhuǎn)信息Vr)輸入到程序中,進(jìn)行迭代,直到相鄰迭代之間的誤差最小,表明得到了收斂自洽的解.通過(guò)比較實(shí)驗(yàn)數(shù)據(jù)與模擬的結(jié)果,完成HL-2A 磁面和安全因子與電流剖面重建.

2.3 HL-2A 裝置參數(shù)

本模擬采用HL-2A 幾何位形和物理參數(shù):大半徑R=1.65 m、小半徑a=0.4 m、環(huán)向磁場(chǎng)BT=2.8 T、環(huán)向等離子體電流IP=480 kA,具有偏濾器位形,圓截面,H 模放電實(shí)驗(yàn).HL-2A 具有相對(duì)較高的加熱功率,其電流加熱系統(tǒng)和輔助加熱包括歐姆加熱、離子回旋共振加熱、中性束加熱和電子回旋共振加熱以及低雜波系統(tǒng).HL-2A 裝置的電子回旋系統(tǒng)總規(guī)模為5 MW,包括68 GHz/3 MW/1 s 的一次諧波系統(tǒng)和140 GHz/2 MW/3 s 的二次諧波系統(tǒng).HL-2A 中性束注入系統(tǒng)注入功率為1.1 MW,束能量為44—55 keV、脈沖寬度為2 s,包括4 套大功率離子源、中性束線部件、中性束電源和診斷測(cè)量系統(tǒng)[22].

基于HL-2A 實(shí)驗(yàn)參數(shù)進(jìn)行先進(jìn)運(yùn)行模式(高約束模式,H 模)動(dòng)力學(xué)平衡位形重建,將實(shí)驗(yàn)測(cè)量的離子和電子溫度、離子密度、旋轉(zhuǎn)角速度剖面參數(shù)作為平衡物理量導(dǎo)入到EFIT 平衡反演程序中.電子密度和溫度,離子溫度和旋轉(zhuǎn)角速度的參數(shù)剖面輸入將是集成模擬分析重點(diǎn).在對(duì)放電參數(shù)選取的時(shí)候,選擇了高約束模式(H 模)的電子密度和溫度、離子溫度,通過(guò)數(shù)據(jù)庫(kù)選擇臺(tái)基明顯的H 模放電參數(shù),以及輸運(yùn)壘比較明顯的物理參數(shù),即#37012 炮放電參數(shù):等離子體電流為153799 A,磁場(chǎng)為1.3 T,等離子體位置坐標(biāo)小半徑a=0.376 m,大半徑R=1.53 m,采用兩束中性束加熱,加熱功率分別為600 keV和725 keV.該次放電炮號(hào)的放電情況如圖2所示,在882 s 后出現(xiàn)了高約束H模.加熱參數(shù)(加熱功率、電子溫度、離子溫度等)通過(guò)相關(guān)的計(jì)算得到:比壓βn=2.3,離子溫度Ti=1.58 keV,電子溫度Te=0.51 keV.設(shè)置好集成模擬程序初始等離子體位形重建物理量和邊界等離子體參數(shù)后,就可以進(jìn)入到反演電流剖面的計(jì)算過(guò)程.

圖2 HL-2A 炮號(hào)#37012 放電參數(shù)Fig.2.Shot #37012 parameters of HL-2A discharge.

3 分析和討論

3.1 HL-2A 實(shí)驗(yàn)數(shù)據(jù)

本文選取托卡馬克裝置HL-2A 參數(shù)#37012放電參數(shù)如下:模擬采用初始放電脈沖的平頂段處的電子溫度、離子溫度以及電子密度.從HL-2A 導(dǎo)入的放電數(shù)據(jù)繪制的曲線如圖3所示,通過(guò)將高能中性原子注入來(lái)將電子和離子進(jìn)行加熱,其中離子加熱為主,離子溫度高于電子溫度,選取芯部電子和離子溫度(大半徑R=1.65 m),對(duì)噪聲進(jìn)行平滑處理.電子密度的分布是等離子體芯部區(qū)域(ρ=0)到邊緣區(qū)域(ρ=±1,+1 代表右側(cè)弱場(chǎng)側(cè),—1 代表左側(cè)強(qiáng)場(chǎng)側(cè))的數(shù)據(jù),在處理時(shí)采用右半部分弱場(chǎng)側(cè)(0 ＜ρ＜ 1)的參數(shù).通過(guò)建立絕對(duì)坐標(biāo)的位置和磁軸的位置相對(duì)應(yīng)變換關(guān)系,得到磁面的位置來(lái)確定邊緣位置,分別對(duì)應(yīng)歸一化后的0 ＜ρ＜ 1位置.選取離散的數(shù)據(jù)集來(lái)構(gòu)造一個(gè)解析函數(shù),讓原來(lái)的離散函數(shù)盡可能接近擬合曲線來(lái)對(duì)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行擬合,通過(guò)迭代返回在x(磁零點(diǎn))處的結(jié)果,同時(shí)在0＜ρ＜1 平均生成離散點(diǎn),將數(shù)據(jù)進(jìn)行3 次樣條數(shù)據(jù)插值,流程圖如圖4所示.能量源和能量損失采用ONETWO 程序計(jì)算,將電子溫度、離子溫度與電子密度剖面導(dǎo)入ONETWO模塊.TGYRO 由ONETWO和EFIT 計(jì)算所得到的已歸一化的參數(shù)剖面進(jìn)行迭代,計(jì)算離子溫度和電子溫度、電子密度和等離子體環(huán)向旋轉(zhuǎn)速度及安全因子q的穩(wěn)態(tài)物理參數(shù)的剖面,從而計(jì)算芯部的輸運(yùn)通量.計(jì)算選取的HL-2A 的設(shè)計(jì)參數(shù)為:βn=2.4,2.3;ne=2.95×1019m—3;Ti=1.53,1.58 keV;Te=0.53,0.51 keV;PNBI1=600 kW;PNBI2=725 kW;fbs/fohm=34.6/53.3;fbeam=12.1%;a=0.376 m;R=1.667 m;Ip=153799 A;B=1.3 T.

圖3 HL-2A 托卡馬克平衡物理參數(shù)實(shí)驗(yàn)數(shù)據(jù) (a)等離子體中離子溫度Ti徑向分布;(b)電子溫度Te徑向分布;(c)電子密度ne徑向分布Fig.3.(a) Radial distribution profile of ion temperature Tiin HL-2A Tokamak plasma;(b) radial distribution profile of electron temperature Tein HL-2A Tokamak plasma;(c) radial distribution profile of electron density nein HL-2A tokamak plasma.

圖4 HL-2A 實(shí)驗(yàn)數(shù)據(jù)處理，物理參數(shù)擬合徑向分布和OMFIT 耦合流程Fig.4.Flow chart for getting fitting curve processing by using discrete and interpolation of HL-2A experimental data.

3.2 結(jié)果和討論

基于HL-2A 實(shí)驗(yàn)數(shù)據(jù)放電數(shù)據(jù),進(jìn)行實(shí)驗(yàn)擬合并將等離子體參數(shù)耦合到OMFIT 集成模擬平臺(tái),求解MHD 磁流體平衡方程,計(jì)算電流密度(帶電粒子)和能量的輸運(yùn)信息,建立邊緣和芯部溫度徑向分布函數(shù),給出邊緣溫度和計(jì)算出的臺(tái)基溫度剖面和電子密度徑向分布.芯部的新經(jīng)典輸運(yùn)和湍流輸運(yùn)由動(dòng)理學(xué)剖面的演化計(jì)算,通過(guò)得到芯部的動(dòng)理學(xué)剖面,固定臺(tái)基區(qū)的參數(shù)進(jìn)行迭代獲得.

圖5給出了重建后的電子溫度、離子溫度和電子密度,通過(guò)集成模擬重建電流剖面曲線與安全因子參數(shù)的剖面.通過(guò)整合等離子體離子溫度,電子溫度和密度,等離子體旋轉(zhuǎn)角速度,移動(dòng)斯塔克效應(yīng)測(cè)量的法拉第旋轉(zhuǎn)角等參數(shù),重建先進(jìn)的高約束模式下等離子體整體輪廓參數(shù)(壓強(qiáng)、電流密度、安全因子)演化物理圖像.

圖5 采用集成模擬程序OMFIT 重建的 (a) HL-2A 的電子溫度Te和離子溫度Ti剖面;(b)電子密度ne剖面Fig.5.Reconstruction of HL-2A parameters by OMFIT integrated simulation code:(a) Electron temperature Teand ion temperature Tiprofiles;(b) electron density ne profile.

圖6所示為穩(wěn)態(tài)等離子體參數(shù)剖面的實(shí)驗(yàn)擬合與集成模擬結(jié)果比較.當(dāng)達(dá)到穩(wěn)態(tài)時(shí),模擬得到的電子溫度、離子溫度、電子密度和旋轉(zhuǎn)因子,與#37012 次放電采用的集成模擬程序OMFIT 重建的該參數(shù)剖面進(jìn)行對(duì)比.在芯部(ρ ＜0.8)區(qū)域,湍流輸運(yùn)和新經(jīng)典輸運(yùn)采用TGYRO 計(jì)算得到,邊界附近處(0.8＜ρ＜1.0)的溫度采取實(shí)驗(yàn)數(shù)據(jù)作為初始值進(jìn)行迭代計(jì)算,芯部的電子溫度初始值為0.6 keV,芯部的離子溫度初始值為1.2 keV.從圖6可以看到,通過(guò)模擬計(jì)算得到的電子和離子溫度空間分布、電子密度和旋轉(zhuǎn)速度等模擬結(jié)果與實(shí)驗(yàn)曲線在演化過(guò)程中有少許差異,實(shí)驗(yàn)數(shù)據(jù)擬合誤差及模擬模型的近似有關(guān),但整體趨勢(shì)基本一致.

圖6 HL-2A 托卡馬克等離子體電子密度(a)、電子溫度(b)、離子溫度(c)和旋轉(zhuǎn)角速度(d)隨徑向剖面的集成模擬和實(shí)驗(yàn)擬合結(jié)果對(duì)比Fig.6.HL-2A Tokamak electron density (a),electron temperature(b),ion temperature (c),rotational angular velocity,and(d) versus radial profile by OMFIT and experimental fitting in shot #37012.

圖7所示為HL-2A 托卡馬克磁面位形、壓強(qiáng)分布、安全因子q、壓強(qiáng)梯度和極向電流梯度等物理參數(shù)剖面.將剖面參數(shù)導(dǎo)入集成建模工作流中,通過(guò)輸運(yùn)程序ONETWO和TGYRO 與平衡程反演序EFIT，迭代得到重建平衡參數(shù)剖面。通過(guò)反演平衡結(jié)果，可以重建電流密度與安全因子剖面.采用EFIT 平衡反演程序,結(jié)合偏振干涉儀測(cè)量的法拉第旋轉(zhuǎn)角（pitch angele）,反演獲得安全因子q分布,最終計(jì)算電流密度分布（Jdistributtion）.從圖7可以看到,壓強(qiáng)梯度在芯部(ρ=0)下降,后又上升,在ρ=0.7 附近有一個(gè)最大值,對(duì)應(yīng)著較強(qiáng)的內(nèi)部輸運(yùn)壘(internal transport barrier,ITB),內(nèi)部輸運(yùn)壘的形成有利于高約束模式的維持,起到凈空雜質(zhì)的作用,相當(dāng)于在芯部和邊緣之間建立了一個(gè)密度和溫度梯度勢(shì)壘.

圖7 集成模擬程序OMFIT 得到的HL-2A 托卡馬克等離子體磁場(chǎng)位形結(jié)構(gòu)圖 (a);等離子體壓強(qiáng)（單位Pa）徑向剖面分布(b);壓強(qiáng)梯度P’（無(wú)量綱）徑向分布(c);安全因子q（無(wú)量綱）隨著徑向變化分布剖面(d);極向電流梯度FF′（A/m3）徑向分布剖面(e)Fig.7.HL-2A Tokamak plasma magnetic field configuration obtained by integrated simulation code OMFIT (a);pressure profile(unit Pa) (b);profile distribution of the pressure gradient(Dimensionless quantity) (c);profile distribution of safety factor q(Dimensionless quantity) obtained (d);profile distribution of the polar current gradient FF′(unit A/m3)(e).

HL-2A 較強(qiáng)的內(nèi)部輸運(yùn)壘所形成的臺(tái)基區(qū)(pedestal region)[21-22],可以從離子溫度0.7 ＜ρ＜ 0.9之間形成的密度快速下降區(qū)域得到(如圖5所示).磁場(chǎng)的x點(diǎn)在z=—0.4 m 附近,該結(jié)果有利于設(shè)計(jì)偏濾器的位形,引導(dǎo)最后一個(gè)閉合磁面的帶電粒子流動(dòng)來(lái)轟擊偏濾器靶板.磁面的位形接近圓形,壓力圖在邊緣附近有梯度的變化.由實(shí)驗(yàn)結(jié)果圖5(a)(b)可以看到,與低約束模式（L 模）對(duì)比,邊緣密度、溫度剖面變陡相當(dāng)于壓強(qiáng)剖面中增加了一個(gè)臺(tái)階,因而在同樣的加熱功率下,總儲(chǔ)能增大,能量約束時(shí)間隨之延長(zhǎng).約束改善源于局域性的反常輸運(yùn)的減弱,即輸運(yùn)壘的建立,或邊緣附近湍流輸運(yùn)的減弱,使邊緣區(qū)的密度、電子和離子溫度剖面都變陡.

高約束模式(H 模)首次在德國(guó)托卡馬克裝置ASDEX 上發(fā)現(xiàn),具有重要的物理意義[23].高約束模式(H 模式)相對(duì)于低約束模式(L 模),粒子的約束密度極限和能量約束時(shí)間都有非常大的提升.當(dāng)能量連續(xù)注入到等離子體中時(shí),邊緣剪切流會(huì)在注入功率達(dá)到某一閾值時(shí)增強(qiáng),然后在邊界區(qū)域形成輸運(yùn)壘,等離子體就從L 模轉(zhuǎn)換成H 模.國(guó)際上其他的有較大加熱功率的聚變裝置都發(fā)現(xiàn)了H 模約束,在球馬克和仿星器上也都有發(fā)現(xiàn),說(shuō)明H 模的出現(xiàn)與偏濾器或限制器的位形和加熱方式無(wú)關(guān)[24],并且與約束裝置類(lèi)型無(wú)關(guān).

當(dāng)?shù)入x子體達(dá)到H 模以后,邊界的溫度和密度不斷上升,從而形成了稱(chēng)為臺(tái)基的結(jié)構(gòu),即邊界輸運(yùn)壘(ETB),此刻的等離子體的能量約束水平已大大提高,約束性能也大幅提高.文獻(xiàn)給出了ITER臺(tái)基區(qū)的預(yù)測(cè)圖[25],可得到臺(tái)基區(qū)的壓強(qiáng)和密度溫度預(yù)測(cè),臺(tái)基結(jié)構(gòu)處的計(jì)算可以用修正tanh 函數(shù)[26,27]來(lái)擬合:

其中,XSYM為臺(tái)基中心的徑向坐標(biāo),X為臺(tái)基區(qū)密度的徑向坐標(biāo),w是臺(tái)基的寬度,a是小半徑,為0.376 m,電子密度ne=2.95×1019m-3,A+B為臺(tái)基的高度[28].如圖3(b)所示,HL-2A 的該炮號(hào)產(chǎn)生的臺(tái)基區(qū)的寬度在0.8 ＜ρ＜ 1.0,臺(tái)基寬度約w=0.376 m× 0.2=7.52 cm.

托卡馬克聚變裝置中臺(tái)基的寬度一般僅僅為幾毫米到幾厘米,但對(duì)等離子體的約束有及其大的作用.在等離子體的芯部隨著壓力梯度的持續(xù)增長(zhǎng)達(dá)到某一臨界值的時(shí)候,湍流運(yùn)輸和各種不穩(wěn)定模式(例如捕獲電子模TEM)都會(huì)出現(xiàn)閾值效應(yīng)的爆發(fā)導(dǎo)致梯度變平,然后芯部等離子體的相關(guān)參數(shù)不能無(wú)限的增長(zhǎng).等離子體約束性能隨著臺(tái)基高度增大而加強(qiáng),使得等離子體的聚變功率提高[29].

H 模具有較長(zhǎng)的能量約束時(shí)間,較高的自舉電流份額和較大的粒子約束密度,被視為未來(lái)ITER實(shí)驗(yàn)聚變堆的基本運(yùn)行模式.研究人員對(duì)等離子體的能量約束時(shí)間和相關(guān)的參數(shù)進(jìn)行了許多實(shí)驗(yàn)[30],在L 模下的約束時(shí)間可以用(26)式定標(biāo)率來(lái)估計(jì):

在H 模具有高臺(tái)基時(shí),其H 模約束下約束時(shí)間可根據(jù)(27)式來(lái)估計(jì):

式中,Ip為等離子體電流,BT為環(huán)向磁場(chǎng)強(qiáng)度,M為原子質(zhì)量,R為大半徑,HL-2A 的大半徑為1.65 m,a為小半徑,HL-2A 的小半徑為0.376 m,κ為等離子體拉長(zhǎng)比,Pt是等離子體得到的凈的加熱功率.有了ITER98 的定標(biāo)率之后[31],可以推導(dǎo)出能量約束增強(qiáng)因子,為托卡馬克裝置獲得能量的約束時(shí)間和定標(biāo)率所預(yù)測(cè)的能量時(shí)間之間的比例,體現(xiàn)約束狀態(tài)的好壞.

圖8所示為采用集成模擬得到的電流密度剖面(電流份額)與安全因子q剖面.從圖8可以看出,安全因子q最小值為1.8,該值較低可能引發(fā)危險(xiǎn)性的低nMHD 模式,從而導(dǎo)致因?yàn)镸HD 不穩(wěn)定性的放電破裂.圖8(a)給出了中性束(beam)電流,自舉電流(BS),歐姆電流(Ohmic),射頻電流(Rf)和全電流(Tot)剖面分布,中性束加熱效果良好,可以看出在中心區(qū)(0 ＜ρ＜ 0.4)所占的份額最多,該區(qū)域中性束電流密度增大,中性束電流密度在中心區(qū)強(qiáng)所占份額最多,隨著徑向分布演化,在邊緣區(qū)域(0.7 ＜ρ＜ 1.0)逐漸減弱到趨近0,同時(shí)自舉電流密度(BS)份額大幅增大,而托卡馬克穩(wěn)態(tài)運(yùn)行需要大幅提高自舉電流的份額.通過(guò)計(jì)算得到的電流份額和空間分布可以對(duì)HL-2A 實(shí)驗(yàn)控制和等離子體位形設(shè)計(jì)提供參考.

圖8 (a)磁約束托卡馬克裝置shot #37012 由集成模擬程序OMFIT 計(jì)算的整體電流密度、包含中性束電流密度、自舉電流密度、歐姆電流密度和射頻電流密度及份額分布曲線;(b)OMFIT 程序和平衡反演程序EFIT 計(jì)算的安全因子q 剖面對(duì)比圖Fig.8.(a) Magnetic confinement Tokamak device shot #37012 total current densities calculated by the integrated simulation program OMFIT,including neutral beam current density,bootstrap current density,Ohmic current density,and RF current density with percentages of respective current densities in total current densities;(b) comparison of the safety factor q profiles calculated by the OMFIT program and the balance inversion code EFIT.

通過(guò)自舉電流分布可知,在0.2 ＜ρ＜ 1.0 區(qū)域,自舉電流在整體電流的份額在上升,而歐姆電流和中性束的電流份額在下降,對(duì)比電流驅(qū)動(dòng)剖面與安全因子q剖面發(fā)現(xiàn),自舉電流的改變與q的剖面的梯度的變化有聯(lián)系.總電流由歐姆電流Iohm、自舉電流IBS、中性束驅(qū)動(dòng)電流INBI和射頻感應(yīng)電流IRF組成.歐姆電流、自舉電流、中性束驅(qū)動(dòng)電流和射頻感應(yīng)電流份額通過(guò)計(jì)算分別為53.2%,34.6%,12.0%和0.2%.通過(guò)與圖3對(duì)比,發(fā)現(xiàn)中性束NBI注入以后,在1.6 m＜R＜1.7 m 處離子溫度上升明顯,證明加熱芯部效果良好,中性束粒子在注入過(guò)程中慢化的過(guò)程更加接近芯部.電子密度受到離子加熱的影響,在該區(qū)域溫度也快速上升,形成臺(tái)基區(qū)的內(nèi)部輸運(yùn)壘位形(ITB).

從圖7(c)壓強(qiáng)梯度隨徑向變化規(guī)律可以看到,中性束注入初期,電子溫度很快升高,形成一個(gè)中空的電流分布,后期的高功率NBI 使等離子體壓強(qiáng)梯度迅速加大,中心區(qū)的粒子密度很快提高,這兩種效應(yīng)都使自舉電流比例加大,負(fù)磁剪切位形得以長(zhǎng)時(shí)間維持,這種運(yùn)行模式能同時(shí)具有中心負(fù)磁剪切和非常峰化的密度分布,具有較長(zhǎng)的約束時(shí)間和較高的約束因子H98.在0.1 ＜ρ＜ 0.2 區(qū)間壓強(qiáng)梯度方向發(fā)生變化(由正變負(fù)),邊緣溫度的增加使等離子體壓強(qiáng)梯度在邊緣區(qū)不斷增加,如圖7(c)所示,溫度分布在邊緣區(qū)變陡,極向旋轉(zhuǎn)速度在此區(qū)域明顯加大,于是在此區(qū)域形成邊緣輸運(yùn)壘ETB(edge transport barrier).在穩(wěn)態(tài)托卡馬克JT-60U,TFTR和DIII-D 都觀測(cè)到約束改善模式H 模式和各種負(fù)剪切位形放電.內(nèi)部輸運(yùn)壘位形(ITB)均是以負(fù)磁剪切位形為基礎(chǔ)的改善約束模式,負(fù)磁剪切位形是在等離子體中心形成負(fù)的磁剪切區(qū),安全因子q的分布函數(shù)不再是單調(diào)下降分布,具有高的自舉電流份額,如圖8(a)所示,自舉電流大幅提升.在電流份額中,自舉電流所占份額對(duì)未來(lái)自持穩(wěn)態(tài)放電至關(guān)重要.一般可通過(guò)中性束加熱和低雜波驅(qū)動(dòng)等手段,大幅增加自舉電流份額.

穩(wěn)定狀態(tài)的托卡馬克等離子體中的可能產(chǎn)生扭曲模不穩(wěn)定性,需有環(huán)向磁場(chǎng)BT和具有上限的等離子體電流.在環(huán)形等離子體中,這個(gè)極限電流可用等離子體邊界的安全因子q＞ 1 滿(mǎn)足扭曲不穩(wěn)定性的穩(wěn)定條件,q稱(chēng)作安全因子,Bp是等離子體電流所產(chǎn)生的極向磁場(chǎng).安全因子可由兩種磁場(chǎng)來(lái)計(jì)算:

其中,r為位置坐標(biāo),可由r/a=ρ來(lái)計(jì)算,BT是環(huán)向磁場(chǎng),Bp是極向磁場(chǎng).

通過(guò)圖8(b)中q剖面的變化,模擬得到的安全因子q剖面并非單調(diào)分布曲線.在0.2 ＜ρ＜ 0.9區(qū)域,通過(guò)實(shí)驗(yàn)數(shù)據(jù)集成模擬得到的q剖面出現(xiàn)了波動(dòng),集成模擬得到的q最小值略小于2.0,在ρ=0處q≈ 2.0,在芯部略有波動(dòng).實(shí)驗(yàn)擬合得到安全因子q分布單調(diào)遞增,中間無(wú)明顯漲落趨勢(shì)變化.集成模擬得到的安全因子q的剖面比單獨(dú)采用EFIT由磁診斷數(shù)據(jù)擬合反演得到的q的剖面更加趨于合理.EFIT 是利用外部的電磁診斷設(shè)備(磁通環(huán)和磁探針)獲取的磁信號(hào),對(duì)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行了全面、深入的分析,準(zhǔn)確重建等離子體平衡位形和相關(guān)的形狀參數(shù).由于缺乏芯部診斷數(shù)據(jù),難以給出較為豐富和更加準(zhǔn)確的電流剖面及壓強(qiáng)剖面分布.

集成模擬OMFIT 是對(duì)芯部的數(shù)據(jù)進(jìn)行計(jì)算,能較為豐富地給出電流剖面及壓強(qiáng)剖面分布,即集成模擬得到的q剖面比單獨(dú)采用EFIT 由磁診斷數(shù)據(jù)反演得到的安全因子q分布更加合理.由圖8(b)可知,EFIT 得到的q變化趨勢(shì)平緩,而結(jié)合等離子體參數(shù),輸運(yùn)參數(shù)集成模擬得到的安全因子剖面q分布更加豐富.在0.3 ＜ρ＜ 0.7 區(qū)域,安全因子增加明顯.

4 結(jié)論

以EFIT,ONETWO,TGYRO 等核心模擬程序?yàn)榛A(chǔ),搭建了具備開(kāi)展動(dòng)力學(xué)平衡位形重建的集成模擬平臺(tái)OMFIT.針對(duì)HL-2A 先進(jìn)高約束模式(H 模)放電實(shí)驗(yàn)的數(shù)據(jù),采用集成建模手段,重建了HL-2A 放電實(shí)驗(yàn)的磁面位形,給出準(zhǔn)確的等離子體邊界物理量剖面,反演了完整的離子和電子溫度徑向分布、電子密度徑向分布、等離子體旋轉(zhuǎn)角速度徑向分布以及安全因子q徑向分布.

基于實(shí)驗(yàn)數(shù)據(jù)和模擬重建的結(jié)果,給出等離子體電流成分構(gòu)成,包括歐姆電流、自舉電流、驅(qū)動(dòng)電流占總電流的份額,并給出了這些電流詳細(xì)剖面分布.基于集成模擬計(jì)算的HL-2A 高約束模式放電集成模型動(dòng)力學(xué)剖面(電子溫度、離子溫度、電子密度、旋轉(zhuǎn))與實(shí)驗(yàn)診斷擬合圖像基本一致,驗(yàn)證了該集成模型.采用OMFIT 集成模擬平臺(tái),整合了物理模擬和數(shù)據(jù)重建程序EFIT,ONETWO,TGYRO 等,得到豐富的q剖面以及電流剖面參數(shù).在托卡馬克等離子體運(yùn)行過(guò)程中,對(duì)于安全因子q測(cè)量,研究人員通常采用運(yùn)動(dòng)斯塔克效應(yīng)光譜診斷來(lái)獲取法拉第旋轉(zhuǎn)角和徑向電場(chǎng)分布.但是對(duì)于可見(jiàn)光譜診斷,面向芯部的光譜測(cè)量遇到了芯部等離子體溫度高,原子被完全電離,光譜信號(hào)弱,邊緣等離子體受到來(lái)自偏濾器和外壁內(nèi)壁的反射雜散光干擾,因此獲取的安全因子圖像誤差較大,這對(duì)穩(wěn)態(tài)托卡馬克反饋控制安全有潛在風(fēng)險(xiǎn).從集成模擬思路出發(fā),基于HL-2A 實(shí)驗(yàn)數(shù)據(jù)和物理模擬模型整合,重建了非單調(diào)變化安全因子q的徑向分布,比單獨(dú)通過(guò)EFIT 平衡反演程序得到的單調(diào)變化安全因子徑向分布更為豐富.下一步計(jì)劃對(duì)集成模擬程序部分模塊進(jìn)行實(shí)時(shí)耦合和自洽性驗(yàn)證,得到HL-2M 相應(yīng)的等離子體參數(shù)的實(shí)時(shí)剖面信息演化,實(shí)現(xiàn)對(duì)穩(wěn)態(tài)運(yùn)行實(shí)時(shí)反饋控制.