FAR3d程序?qū)AST上阿爾芬不穩(wěn)定性識別過程中的應(yīng)用

王曦輝,黃 娟,趙海林,J. Varela,付 靜,孫延旭,史 唱,王書松

(1.中國科學(xué)院 合肥物質(zhì)科學(xué)研究院,安徽 合肥 230031;2.中國科學(xué)技術(shù)大學(xué),安徽 合肥 230026;3.Universidad Carlos Ⅲ de Madrid, 28911 Leganes, Madrid, Spain)

在未來的燃燒等離子體實(shí)驗(yàn)中,抑制或控制高能粒子引起的不穩(wěn)定性是一個(gè)至關(guān)重要的問題。在托卡馬克裝置中,高能粒子主要來源于中性束注入(neutral beam injection, NBI)、離子回旋共振加熱(ion cyclotron resonance heating, ICRH)等輔助加熱手段,以及聚變產(chǎn)物之一的α粒子(3.5 MeV)[1]。由于這些粒子攜帶著極高的能量,且特征速度接近阿爾芬速度,因此易與等離子體背景中的阿爾芬波共振傳遞能量,從而激發(fā)阿爾芬不穩(wěn)定性,導(dǎo)致快離子再分布或丟失,對裝置的加熱效率和安全性造成損害。因此,為有效控制高能粒子的行為,保障裝置的安全穩(wěn)定運(yùn)行,提高聚變能源的經(jīng)濟(jì)性和可持續(xù)性,對于阿爾芬不穩(wěn)定性的研究顯得尤為重要。

盡管阿爾芬不穩(wěn)定性在理論上已有較完整的物理模型,但在實(shí)際實(shí)驗(yàn)中,其激發(fā)和穩(wěn)定性仍需進(jìn)一步研究。通常情況下,阿爾芬不穩(wěn)定性是根據(jù)其激發(fā)方式和不穩(wěn)定特征來區(qū)分的,例如由于環(huán)效應(yīng)引起的TAE[2],在安全因子極值處被激發(fā)的RSAE(reversed shear Alfven eigenmode)[3],以及由于有限比壓引起的BAE(Beta-induced Alfven eigenmode)[4]等。實(shí)驗(yàn)中,磁探針、電子回旋輻射等診斷均可用來對不穩(wěn)定性進(jìn)行觀測[5]。但對于不穩(wěn)定性的激發(fā)機(jī)制的分析還需將診斷信息應(yīng)用到具體物理模型中,所以,各種用于研究不穩(wěn)定性的程序被開發(fā)出來,其中有采用完整粒子模型的GTC[6]、采用磁流體和動理學(xué)混雜模型的M3D-K[7]以及本征值程序NOVA[8]等。FAR3d是Spong、Garcia、Varela等開發(fā)的用來學(xué)習(xí)計(jì)算磁流體不穩(wěn)定性和阿爾芬本征模(AE)不穩(wěn)定性以及各種動理學(xué)效應(yīng)對它們的影響的模擬程序[9],該程序已在國外各裝置上有了廣泛的應(yīng)用。如在日本的仿星器裝置LHD上,其被用來研究由快離子激發(fā)的TAE不穩(wěn)定性的特征,模擬了快離子比壓對AE不穩(wěn)定性的影響[10]。在西班牙的TJ-Ⅱ仿星器裝置上,研究了多種快離子不穩(wěn)定性如TAE、EPM等共同存在時(shí)的特征,并通過掃描磁剪切的數(shù)值模擬,解釋了實(shí)驗(yàn)中的TAE頻率掃頻現(xiàn)象[11]。之后,針對DⅢ-D裝置中的高極向比壓放電,FAR3d同樣用于識別不穩(wěn)定性模式,在實(shí)驗(yàn)可達(dá)到的參數(shù)區(qū)間內(nèi)進(jìn)行參數(shù)學(xué)習(xí),提出可優(yōu)化不穩(wěn)定性的運(yùn)行參數(shù)區(qū)間[12]。以外,FAR3d還被用來研究高能粒子對EIC(energetic-ion-driven resistive interchange mode)事件的影響[13]。

目前在全超導(dǎo)托卡馬克裝置EAST實(shí)驗(yàn)中,已觀測到了大量的AE不穩(wěn)定性現(xiàn)象。如中性束注入條件下激發(fā)的TAE[14],由雜質(zhì)引起的撕裂模與BAE[15],磁島和測地聲模耦合引起的BAE[16]等。本文針對EAST高βN運(yùn)行模式條件下觀察到的AE不穩(wěn)定性,首次將FAR3d程序應(yīng)用到EAST上,并結(jié)合相關(guān)擾動測量進(jìn)行驗(yàn)證以及細(xì)致分析,為進(jìn)一步不穩(wěn)定性產(chǎn)生機(jī)制的研究和有效控制提供物理依據(jù)。

1 實(shí)驗(yàn)研究

1.1 EAST及快離子不穩(wěn)定性診斷

EAST裝置為全超導(dǎo)托卡馬克裝置,大半徑R=1.85 m,小半徑a=0.45 m。EAST裝置的位型與ITER類似,均為非對稱的偏濾器位型。在EAST實(shí)驗(yàn)中,快離子主要由ICRH和NBI加熱產(chǎn)生。目前,EAST上的NBI加熱系統(tǒng)擁有兩條束線,可穩(wěn)定地在65 keV束壓條件下提供共計(jì)5 MW以上的功率[17]。離子回旋系統(tǒng)目前擁有兩條天線,可穩(wěn)定提供2.5 MW以上的輸出功率[18]。

在快離子不穩(wěn)定性的診斷方面,EAST裝置上配備了最高可探測到500 kHz頻率的高頻磁探針,來測量不穩(wěn)定性造成的高頻磁擾動,并通過磁探針在環(huán)向空間上的分布可以計(jì)算得到不穩(wěn)定性的環(huán)向模數(shù)。EAST上的電子回旋輻射診斷則可以通過測量不穩(wěn)定性造成的電子回旋輻射擾動,借由其多道外差系統(tǒng),來提供厘米量級的空間分辨和微秒量級的時(shí)間分辨,給出不穩(wěn)定性的位置和頻率信息。除此以外,EAST中還配有軟X射線診斷、束發(fā)射光譜診斷等對快離子不穩(wěn)定性進(jìn)行測量。在上述的不穩(wěn)定性的直接測量手段以外,EAST中還配有多種診斷來提供快離子的信息。其中,快離子損失探針可以測量不穩(wěn)定性引起的快離子損失,以及FIDA(fast ion D-Alpha)診斷可以反演得到快離子的速度空間分布函數(shù)。本文中主要使用了高頻磁探針診斷和電子回旋輻射診斷,基于測到的不穩(wěn)定性的模數(shù)、頻率、位置等信息在實(shí)驗(yàn)上對不穩(wěn)定性做初步分析。

1.2 實(shí)驗(yàn)條件

93910炮放電的縱場Bt為1.6 T,等離子體電流Ip=400 kA,電子密度平臺大約在2.5×1019m-3,在t=6.5 s時(shí)刻,q95=4.0,βN=1.8,有效電荷數(shù)Zeff=2.5。該放電使用了總功率為1 MW的低雜波加熱和總功率為5 MW的NBI注入加熱。具體放電參數(shù)演化如圖1所示。圖1e為所選時(shí)間范圍內(nèi)的高頻磁探針診斷的頻譜圖,可看到在t=6.0 s后,磁探針信號出現(xiàn)了較大的擾動,多支不穩(wěn)定性已被激發(fā)。為避免多個(gè)模式的干擾以及診斷測量的時(shí)間分辨限制,本文選擇了該次放電t=6.5 s的時(shí)刻點(diǎn)來進(jìn)行分析,此時(shí)具有掃頻特性的不穩(wěn)定性已消失,僅剩一支90 kHz附近的頻率較穩(wěn)定的模式。

a——等離子體電流;b——輔助加熱,藍(lán)色線為4.6 GHz低雜波,橙色為兩條NBI束線;c——電子密度;d——中子產(chǎn)率;e——6～7 s時(shí)的高頻磁探針頻譜圖1 93910炮放電基本參數(shù)Fig.1 Basic parameters of shot 93910

2 FAR3d模擬程序簡介

2.1 程序概述和物理模型

FAR3d程序是由Spong、Garcia、Varela等開發(fā)的一用來學(xué)習(xí)研究線性MHD和AE不穩(wěn)定性的模擬工具[19-20]。FAR3d采用gyro-fluid模型,通過求解一系列的簡化的線性電阻MHD方程組[21],并在其中添加了快離子密度和平行方向的動量的變化,把與快離子相關(guān)的效應(yīng)加入到了程序的物理模型中。如通過加入線性的波-粒共振引入朗道阻尼,加入熱離子平行動量的響應(yīng)耦合測地聲模對不穩(wěn)定性的影響[22]。

(1)

(2)

(3)

(4)

FAR3d程序中將等離子體的熱成分和快成分分為兩部分。其中熱成分可寫為:

(5)

(6)

(7)

(8)

以渦量方程(式(6))為例,介紹模型中所考慮的效應(yīng)。式(6)中第1項(xiàng)表示的是平衡的環(huán)向旋轉(zhuǎn)帶來的影響;第2、3項(xiàng)為電流在磁場B方向的變化的影響;第4項(xiàng)為等離子體壓強(qiáng)梯度項(xiàng);第5項(xiàng)為快離子的密度梯度驅(qū)動項(xiàng);第6項(xiàng)為熱離子的有限拉莫半徑效應(yīng)項(xiàng);第7項(xiàng)為朗道阻尼項(xiàng);第8項(xiàng)表示雙流體效應(yīng)影響。此外,式(7)中第3、4項(xiàng)表示聲波帶來的壓縮效應(yīng)。FAR3d模型中,通過引入快離子密度和平行方向的動量的變化,把與快離子相關(guān)的效應(yīng)加入到了程序的物理模型中。描述快離子成分的關(guān)系式如下。

(9)

(10)

目前FAR3d模型中對于快離子部分的處理,僅引入了高能粒子的平行速度的變化,對平行能量較高的粒子引起的不穩(wěn)定性有較好的模擬效果,而對于快離子垂直方向的運(yùn)動引起的不穩(wěn)定性的模擬,還需后續(xù)版本的繼續(xù)升級來獲得更合理的結(jié)果。

2.2 程序輸入與輸出

FAR3d程序的輸入由3部分組成:平衡文件;輸入?yún)?shù)表單,其中主要包含各種效應(yīng)的開關(guān)選項(xiàng);還有就是等離子體的宏觀參數(shù)剖面,如熱離子的密度、溫度、快離子的密度剖面等。FAR3d現(xiàn)在僅支持VMEC格式的平衡文件輸入,目前EAST上使用的是EFIT的平衡[25],在使用FAR3d程序的過程中,需將EFIT平衡轉(zhuǎn)換為VMEC平衡再作輸入[26]。針對93910炮放電,本文中使用了EFIT、ONETWO/NUBEAM程序來對t=6.5 s時(shí)刻進(jìn)行動理學(xué)反演,得到EFIT格式的平衡輸入文件。

進(jìn)行動理學(xué)平衡反演的過程以及EAST裝置中的診斷可提供絕大多數(shù)FAR3d程序所需要的輸入?yún)?shù)剖面。圖2示出了診斷測量參數(shù)擬合過后的結(jié)果,即FAR3d中的輸入剖面。其中包括湯姆遜散射(TS, Thomson scattering)和偏振干涉儀診斷(POINT, polarimeter-interferometer)提供的電子溫度Te和密度ne剖面,電荷交換復(fù)合光譜(CXRS, charge exchange recombination spectroscopy)提供的離子溫度Ti剖面。在ONETWO/NUBEAM輸出中可得到快離子部分的能量Tf以及密度nf剖面。

a——安全因子剖面;b——快離子溫度與密度剖面;c——電子與離子密度剖面;d——電子與離子溫度剖面圖2 FAR3d模擬中的輸入?yún)?shù)剖面Fig.2 FAR3d input profiles

目前,FAR3d程序有兩種運(yùn)行模式,一種是initial value solver,盡管在FAR3d程序中該模式運(yùn)算速度很快,但對每步計(jì)算只能給出一個(gè)最不穩(wěn)定模式的信息。另一種運(yùn)行模式是eigensolver,該模式可計(jì)算出給定頻率/增長率范圍內(nèi)的所有的可能為不穩(wěn)定的模式,且包括其計(jì)算出的所有模式的本征函數(shù),而不僅只有增長率最高的主導(dǎo)模式。由于本文中所選擇的放電存在多種不穩(wěn)定性,所以本文所有結(jié)果均由eigensolver模式產(chǎn)出。在eigensolver模式下,程序可給出其發(fā)現(xiàn)的所有不穩(wěn)定性的頻率、增長率以及模結(jié)構(gòu),從中可得到各模式的環(huán)向和極向模數(shù)。

2.3 AE不穩(wěn)定性分析

1) 實(shí)驗(yàn)特征

如圖1e所示,在t=6.0～6.5 s之間,磁探針探測到強(qiáng)烈的不穩(wěn)定性現(xiàn)象,且具有隨時(shí)間向上掃頻的特性,頻率范圍為70～100 kHz,且還有一支約90 kHz的模式幾乎保持頻率不變。t=6.6 s開始,在50～60 kHz左右的較低頻率區(qū)域也出現(xiàn)了頻率范圍較寬的不穩(wěn)定現(xiàn)象,隨時(shí)間頻率略微降低。同時(shí)伴隨著一頻率較高但強(qiáng)度較弱的向上掃頻的模式。

通過EAST中的高頻磁探針在環(huán)向上的陣列,計(jì)算出各點(diǎn)信號的相位差即可得到不穩(wěn)定性的環(huán)向模數(shù)信息,如圖3a所示,兩支模式的環(huán)向模數(shù)均為n=1～2。結(jié)合圖3b所示的ECE外差系統(tǒng)的觀測結(jié)果,可得t=6.3 s時(shí)刻開始出現(xiàn)的頻率從70 kHz上升至100 kHz的模式,其環(huán)向模數(shù)為n=2,徑向位置在ρ=0.46。另一個(gè)在t=6.6 s時(shí)刻出現(xiàn)的60 kHz附近的模式同樣在徑向上位于ρ=0.46附近。需注意的是,該次放電為低縱場放電,電子回旋輻射基頻較低,ECE診斷測到的電子回旋輻射信號已經(jīng)處于3次諧頻范圍,所以其信號較弱。且該次放電的環(huán)向磁探針陣列僅有兩個(gè)高頻磁探針可用,所以對于環(huán)向模數(shù)的計(jì)算可能會有一定的誤差。在下一小節(jié)中,通過FAR3d程序?qū)=6.5 s時(shí)刻的不穩(wěn)定性做線性模擬分析。

a——環(huán)向磁探針陣列;b——ECE診斷測量圖3 不穩(wěn)定性的環(huán)向模數(shù)與位置Fig.3 Toroidal mode number and radial location of instability

2) FAR3d程序的模擬結(jié)果

FAR3d模擬顯示只有n=2的模式是不穩(wěn)定的。圖4為FAR3d在eigensolver運(yùn)行模式下的模擬輸出結(jié)果,增長率按阿爾芬時(shí)間進(jìn)行了歸一。其中紅色方框處為增長率最高的主導(dǎo)模式,頻率為87 kHz,藍(lán)色方框處為次級主導(dǎo)模式,頻率為62 kHz。圖5為主導(dǎo)模式和次級主導(dǎo)模式靜電勢的本征函數(shù)。從圖5a可看出,主導(dǎo)模式本征函數(shù)較寬,由n/m=2/3和n/m=2/4的模式耦合而成,位置在歸一化半徑ρ=0.45附近,這與實(shí)驗(yàn)觀測結(jié)果是一致的。圖5b展示了頻率為62 kHz次級主導(dǎo)模式,其本征函數(shù)非常窄,是由單一的極向模式所組成,其模數(shù)為n/m=2/4,位置在ρ=0.55附近。

圖4 開啟/關(guān)閉FLR效應(yīng)時(shí)FAR3d的模擬結(jié)果的增長率與頻率信息Fig.4 Modes’ growth rate and frequency of FAR3d results with/without FLR effect activated

a——主導(dǎo)模式(圖4中紅色方框);b——次級主導(dǎo)模式(圖4中藍(lán)色方框)圖5 主導(dǎo)模式和次級主導(dǎo)模式的本征函數(shù)Fig.5 Eigen functions of dominant and subdominant modes

為進(jìn)一步確定各模式的信息,本文使用StellGap程序計(jì)算出圖6中的阿爾芬連續(xù)譜[27]。圖中紅色橫線表示主導(dǎo)模式的位置,藍(lán)色橫線表示次級主導(dǎo)模式所處的位置。從中可看出,主導(dǎo)模式位于譜間隙的底部,而次級主導(dǎo)模式則位于連續(xù)譜內(nèi)。因此可判斷,主導(dǎo)模式為n/m=2/3,2/4的TAE,頻率為87 kHz,次級主導(dǎo)模式則為n/m=2/4的EPM。需注意,連續(xù)譜中的主導(dǎo)模式和連續(xù)譜線有少量的重疊,這是由于StellGap程序與FAR3d中對于聲波的耦合方法不同所導(dǎo)致的。在圖5a中2/4尖銳的模結(jié)構(gòu)則是由于模擬中未開啟FLR效應(yīng),因此FAR3d程序會計(jì)算出此處存在一較窄的模結(jié)構(gòu)并且與當(dāng)前模式混合。如圖7所示,在FLR效應(yīng)帶來的阻尼下,不穩(wěn)定性中較窄的模結(jié)構(gòu)已消失。同時(shí)也可看出,該模式并不是位于連續(xù)譜中的模式。至此,FAR3d模擬成功得到了與實(shí)驗(yàn)觀測一致的不穩(wěn)定性且識別了主導(dǎo)模式為TAE和次級主導(dǎo)模式為EPM的不穩(wěn)定性。模擬使用的宏觀參數(shù)剖面以及平衡均為t=6.5 s時(shí)刻,但在實(shí)驗(yàn)觀測中60 kHz附近的低頻模式直到t=6.6 s才出現(xiàn),而且本文的模擬中并未進(jìn)行非線性模擬分析,所以對于模擬中得到的次級模式EPM是如何演化并成為主導(dǎo)模式的,還需后續(xù)進(jìn)一步分析。

圖6 t=6.5 s時(shí)刻的阿爾芬連續(xù)譜以及不穩(wěn)定性所處的位置Fig.6 Alfven continuum and location of instability at t=6.5 s

圖7 開啟FLR效應(yīng)后的主導(dǎo)模式的本征函數(shù)Fig.7 Eigen function of dominant mode with FLR effect included

需要提到的是,即使FAR3d模擬中會計(jì)算出多種不穩(wěn)定模式,但由于驅(qū)動源所能提供的自由能量的限制,若多種模式均由同一種源來獲得能量,如某處的快離子壓強(qiáng)梯度,那么在自由能獲取的競爭中,通常增長率最高的主導(dǎo)模式會在競爭中勝出從而有最大的可能表現(xiàn)在實(shí)驗(yàn)中。對于不穩(wěn)定性從出現(xiàn)到各自的模式演化過程,想要得到模式更準(zhǔn)確的信息,還需要后續(xù)再進(jìn)行進(jìn)一步的非線性模擬分析。除此以外,由于FAR3d模擬程序僅從數(shù)值上得到收斂的結(jié)果,因此程序可能會計(jì)算出一些僅從計(jì)算上合理但無實(shí)際物理意義的“人造”解。所以,對于程序的計(jì)算結(jié)果,需要加入合理的物理理解,來排除干擾,解釋正確的物理現(xiàn)象。

在識別EAST等離子體中的不穩(wěn)定性模式時(shí),并不需將所有的動理學(xué)效應(yīng)項(xiàng)目均加入到計(jì)算中,可根據(jù)放電參數(shù)的不同,先使用FAR3d程序中各個(gè)動理學(xué)項(xiàng)對不穩(wěn)定性模式的影響進(jìn)行簡單的評估,再考慮是否需要將這些項(xiàng)加入到計(jì)算中。以有限拉莫半徑效應(yīng)為例,在圖4可看到,在93910炮放電中,開啟有限拉莫半徑效應(yīng)后,等離子體的不穩(wěn)定性在增長率上僅有微小的差別,對于模式的識別并沒有影響。這是因?yàn)镕LR效應(yīng)對于低n模式影響較小,而對于高n模數(shù)的模式影響較大。因此,在93910炮的計(jì)算中,選擇關(guān)閉有限拉莫半徑效應(yīng),以使計(jì)算更加快速。而對于高n的模式,FLR效應(yīng)則可能會帶來強(qiáng)烈的阻尼效應(yīng)。

3 總結(jié)與展望

針對EAST高βN運(yùn)行區(qū)間觀測到的快離子相關(guān)不穩(wěn)定性,在EAST上應(yīng)用FAR3d程序,并結(jié)合快離子不穩(wěn)定性診斷測量進(jìn)行相關(guān)驗(yàn)證和分析。基于93910炮放電,t=6.5 s的時(shí)刻點(diǎn)進(jìn)行線性模擬分析,成功重復(fù)出了實(shí)驗(yàn)中出現(xiàn)的不穩(wěn)定性模式,且判斷為m/n=3/2、4/2的TAE以及m/n=4/2的EPM,頻率與位置信息均與實(shí)驗(yàn)測量結(jié)果吻合。成功將FAR3d程序應(yīng)用到了EAST裝置中。并且評估了低n模數(shù)的AE不穩(wěn)定性,有限拉莫軌道半徑效應(yīng)對不穩(wěn)定性的影響較弱。在實(shí)驗(yàn)中,與磁探針提供的不穩(wěn)定性環(huán)向模數(shù)結(jié)合,適時(shí)的關(guān)閉該效應(yīng)有利于提高計(jì)算速度。本文僅討論了t=6.5 s時(shí)刻的AE不穩(wěn)定性。對于不穩(wěn)定性的掃頻現(xiàn)象,以及其激發(fā)機(jī)理的研究,后續(xù)會進(jìn)行進(jìn)一步的模擬分析。在未來實(shí)驗(yàn)中,將結(jié)合ICRF和NBI協(xié)同加熱產(chǎn)生的快離子激發(fā)的相關(guān)不穩(wěn)定性進(jìn)行深入分析,同時(shí)進(jìn)一步發(fā)展FAR3d程序相應(yīng)的物理模型,對快離子不穩(wěn)定性的有效控制提供實(shí)驗(yàn)依據(jù)和指導(dǎo)。