HL-2A 装置高βN 双输运垒实验的集成分析*

李正吉 陈伟† 孙爱萍 于利明 王卓 陈佳乐 许健强李继全 石中兵 蒋敏 李永高 何小雪 杨曾辰 李鉴

1) (核工业西南物理研究院,成都 610041)

2) (中国科学院等离子体物理研究所,合肥 230031)

3) (61287 部队,成都 610000)

1 引言

可控核聚变装置是利用轻原子核聚合形成重原子核过程中产生的能量进行发电的装置,其能量转换效率是传统化石能源的百万倍,并且原料储量丰富,运行条件安全、无污染,成为了许多国家新能源发展的重点目标.到目前为止,托卡马克装置被认为是最有希望实现可控核聚变要求的实验装置.对于磁约束聚变装置托卡马克和球形环等,等离子体的归一化环向比压β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 节为总结.

2 HL-2A 装置高βN 物理实验

在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.

3 高βN 等离子体放电的OMFIT 集成模拟

使用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.

4 总结

本文讨论的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 装置高比压实验提供依据.