H-1NF仿星器標準磁場位形分析與高能量離子運動軌道模擬

向 東，夏彥博，宋澤宸，曹錦佳，杜 丹，管滬楠，龔學余

(1.南華大學 核科學技術學院，湖南 衡陽 421001；2.南華大學 數理學院，湖南 衡陽 421001)

仿星器和托卡馬克裝置是目前國際上研究最多的兩大磁約束可控熱核聚變裝置類型。由于建造仿星器的成本高、技術難度大等諸多原因，從20世紀70年代開始，大多數核聚變研究都聚焦于托卡馬克。目前中國擁有EAST、HL-2M、J-TEXT等托卡馬克裝置，還沒有仿星器裝置，因而在仿星器方面的理論研究相對較少。然而，由于仿星器內等離子體全部依靠外部磁場約束，易于實現自然穩態運行，激發大破裂的可能性非常低，這些優點在聚變界一直備受青睞[1-3]。特別是在2015年12月，目前世界上最大、設計最復雜、版本最新的德國W7-X仿星器成功實現高參數等離子體放電而轟動全球[4]。W7-X仿星器設計的遠景目標是實現1 800 s超長脈沖穩態運行[5]。目前，建造中國聚變工程實驗堆(CFETR)已被列為國家聚變能源發展的路線圖[6]。仿星器作為未來聚變反應堆備選方案之一，其相關研究同樣也受到了國內科學界的關注。西南交通大學和日本國立聚變研究所合作在中國設計籌建一個新型準軸對稱仿星器CFQS[7-8]，浙江大學、中國科技大學等高校聚變研究團隊也開始了關于仿星器的理論研究[9-10]。為了在國內廣泛開展H-1NF仿星器的實驗和理論研究，本文對H-1NF仿星器的磁場線圈組成和分布特征進行分析，研究其磁場位形特點及單粒子運動軌道特征，并模擬計算高能量離子在一種標準磁場位形中的典型運動軌道。

1 H-1NF仿星器磁場線圈分布

H-1NF仿星器的磁場沿大環方向三周期120°對稱分布[11-13]，如圖1所示，其磁場線圈系統共包括36餅環向場線圈(TFC)，1個中心環導體，即極向場線圈(PFC)，1條螺旋線圈(HCW)，1對內垂直場線圈(IVFC)，1對外垂直場線圈(OVFC)。其中，環向場線圈、極向場線圈、螺旋線圈和內垂直場線圈均安裝在仿星器真空腔體內，外垂直場線圈安裝在真空腔體外。

a——三維視圖；b——俯視圖圖1 H-1NF仿星器磁場線圈分布Fig.1 Distribution of magnetic field coil of H-1NF stellarator

36餅環向場線圈的大小規格相等，每餅環向場線圈內部繞有10圈銅導線，其平均等效半徑RTFC=0.38 m。這些環向場線圈的圓心不在同一平面上，若采用大柱坐標(R，φ，Z)來描述各線圈的分布位置，則每個環向場線圈的圓心坐標(Rc，φc，Zc)滿足以下關系式：

Rc=R0+ρscos(3φc)

(1)

Zc=ρssin(3φc)

(2)

φc=φ′-0.009 7sin(3φ′)

(3)

(4)

其中：R0為大半徑，R0=1.00 m；ρs為擺動半徑，ρs=0.22 m。

中心極向場線圈位于Z=0.00 m的赤道面上，其內部繞有36圈銅導線，其線圈等效半徑RPFC=1.00 m。1對內垂直場線圈分別對稱分布于赤道平面Z=±1.07 m的平面上，每個內垂直場線圈內部繞有16圈銅導線，其線圈等效半徑RIVFC=0.72 m。1對外垂直場線圈分別對稱分布于赤道平面Z=±0.70 m的平面上，其內部繞有12圈(通過外部引線可設置成4、8、12圈3種，其標準配置為8圈)銅導線，其線圈等效半徑ROVFC=2.13 m。螺旋線圈內部并排4條銅導線，它們沿環向呈3周期纏繞中心極向場線圈，纏繞的等效小半徑rHCW=0.10 m。

H-1NF仿星器線圈標準接線方式是將環向場線圈、極向場線圈、內垂直場線圈、外垂直場線圈中的銅導線全部串聯，由ABB電源的主回路供電，而螺旋線圈由ABB電源的從回路供電。通過設置ABB電源主、從兩路輸出的電流強度，可獲得不同的磁場位形。H-1NF仿星器標準運行模式是外垂直場線圈設置為8圈銅導線，且螺旋線圈不供電。在標準運行模式下，當主路輸入的銅導線電流強度為13 889 A時，若將中心極向場線圈等效看作是一圈圓電流，則其電流強度約為5.0×105A，此時磁軸上的磁感應強度約為1 T。本文以此種H-1NF仿星器標準運行模式下磁場位形為例，分析其位形特點，并模擬計算研究高能量離子在該標準磁場位形中的典型運動軌道特征。

2 標準配置下的磁場位形

H-1NF仿星器運行時，磁面將圍繞中心極向場線圈沿環向呈三周期麻花狀分布，令最外層磁面為等離子體邊界，繪制出其等離子體位形如圖2所示。基于H-1NF仿星器標準運行模式參數，采用HELIAC程序追蹤計算各磁面的磁場線[13]，分別繪制出不同環向角φ=0°、15°、30°、45°、60°、75°、90°、105°、120°等真空磁場的龐加萊圖，如圖3所示。由圖3可看出，仿星器真空磁面沿環向螺旋狀分布，這與托卡馬克磁面為環向對稱分布有很大不同，且仿星器的磁軸不在一平面上。另外，圖3a和圖3i顯示出φ=0°和120°磁面的極向截面完全相同，這說明H-1NF仿星器的磁場分布是沿大環方向3周期120°對稱。

圖2 H-1NF仿星器內等離子體位形Fig.2 Plasma shape in H-1NF stellarator

雖然仿星器磁面沿環向螺旋狀分布不具有環向對稱性，但是不同環向角的磁面極向截面除了沿磁軸旋轉一定角度外，其形貌特征相近。因此，基于磁軸的位置和磁面沿環向角的變化規律，以磁軸為旋轉軸，將圖3所示的不同環向角的極向截面的磁場龐加萊圖在原有截面坐標(R，Z)內，圍繞磁軸旋轉一定的角度，并將坐標原點平移至磁軸的位置，定義旋轉后的極向截面坐標為(R′，Z′)，此截面稱為旋轉坐標系下的等效標定極向截面。因此得到了一種旋轉坐標系下等效標定極向截面上的磁場龐加萊圖，圖3所示的各種極向截面的磁場龐加萊圖經過旋轉變換后對應的等效標定極向截面的磁場龐加萊圖如圖4所示。

在H-1NF仿星器標準運行模式下，螺旋線圈不通電，而環向場線圈、極向場線圈、內垂直場線圈和外垂直場線圈均可看作是一個等效的單匝圓形線圈導體。圓形導體產生的磁場可用第二類橢圓積分來表示，而這種積分可采用Hastings型多項式來計算[14]。因此，基于圓形導體產生磁場的計算方法和場強疊加原理，可實現H-1NF仿星器內任一點的磁感應強度計算。根據這一原理，HELIAC程序計算磁場以實現任一磁場線的追蹤。本文采用此方法，利用Matlab編程計算出環向角φ=0°時極向截面的磁場分布特征(圖5)。圖5a為旋轉坐標系下環向角φ=0°時等效標定極向截面的磁場龐加萊圖；圖5b為圖5a中示意磁面(極向截面為一個環)隨極向角變化的磁感應強度分布；圖5c為旋轉坐標系下Z′=0直線上沿R′方向的磁感應強度分布。

a——φ=0°；b——φ=15°；c——φ=30°；d——φ=45°；e——φ=60°； f——φ=75°；g——φ=90°；h——φ=105°；i——φ=120°圖3 不同環向角的極向截面的磁場龐加萊圖Fig.3 Poincaré map of magnetic field in poloidal cross-section with different toroidal angles

3 高能量離子的典型運動軌道模擬

忽略重力場等次要因素，離子在真空磁場B中的運動遵循動力學方程：

(5)

其中：m、q和v分別為離子的質量、電荷和運動速度；t為運動時間。

計算運動軌道的數值算法采用四階龍格庫塔算法[15]。另外，數值精度還與每個離散時間步Δt所對應的空間尺寸有關，因此計算前必須估算離子在計算區域內磁場中的拉莫爾半徑(rL=mv⊥/qB)，以保證離子在磁場中運動1個回旋周期的數值模擬經過了足夠多的時間步數演算。在磁場中運動離子的漂移速度與離子的能量、荷質比、投擲角和背景磁場均有關系。若離子的初始位置不同，其所在位置的磁場就可能不同，且在運動中由于離子位置處于不停地變化中，因此其軌道上每一點的磁場會不停地變化，投擲角也會不斷地改變。對于低能量離子，其漂移速度很小，回旋運動的導心基本會沿磁力線走，因此，以下主要研究高能量離子的運動軌道特征。

a——φ=0°；b——φ=15°； c——φ=30°；d——φ=45°；e——φ=60°； f——φ=75°；g——φ=90°；h——φ=105°；i——φ=120°圖4 旋轉坐標系下不同環向角的等效標定極向截面的磁場龐加萊圖Fig.4 Poincaré map of magnetic field in equivalent calibrated poloidal cross-section with different toroidal angles in rotating coordinate system

圖5 環向角φ=0°時極向截面的磁場分布特征圖5 Distribution characteristic of magnetic field in poloidal cross-section with toroidal angle φ=0°

不失一般性，本文令測試離子為質子，能量為10 keV，假設其初始位置的導心在環向角為φ=0°時圖5a所示磁面上的點A，其極向截面坐標為(R=1.30 m，Z=0.00 m)，相應旋轉坐標系中等效標定極向截面坐標為(R′=0.06 m，Z′=0.00 m)。令離子的初始投擲角為η，則cosη=v∥/v，v∥為離子速度平行磁場方向的分量。通行粒子軌道和捕獲粒子軌道間的臨界投擲角ηc滿足：

(6)

其中：Bs為離子初始位置的磁感應強度；Bm為離子運動軌道上的最大磁感應強度。在0°～90°范圍內，當離子的初始投擲角η<ηc時為通行粒子，反之為捕獲粒子。根據圖5b的磁場數據可初步理論估計ηc，但因為仿星器磁場分布不具有環向對稱性，在不同環向角對應的旋轉坐標系中標定極向截面上的磁感應強度沿極向角分布會有一定的差異，所以實際的ηc與根據圖5b簡單理論的估計值不同。經過精細計算，得到ηc=57°。若令離子的初始投擲角η=30°，則模擬計算出其運動全軌道如圖6所示，顯然這是一種典型的通行粒子軌道。若將離子的初始投擲角增加到η=60°，則可得到一種典型的捕獲粒子軌道，其模擬計算出的全軌道如圖7所示。

由于仿星器磁面不具有環向對稱性，因而無法像研究托卡馬克中的粒子軌道一樣將三維粒子軌道直接投影到二維極向截面上觀察。從圖6、7中無法直接觀察離子運動軌道相對于磁面的位置，也不能直接判斷離子是否飛出最外層閉合磁面。然而，H-1NF仿星器內不同環向角的極向截面沿磁軸旋轉一定角度后形貌特征相近(圖4)。因此，基于磁軸位置和和磁面沿環向角的旋轉變化規律，可將高能量離子的三維運動軌道投影到類似圖4所示的旋轉坐標系中不同環向角的等效標定極向截面上，從而可更加直觀清晰地顯示出高能量離子在仿星器磁場位形中的運動軌道特征。雖然旋轉坐標系中不同環向角的等效標定極向磁面不完全相同(圖4a～i)，但是差別不大。為了簡化，統一取φ=0°的標定極向磁面作為背景，分別繪制出圖6、7所示的三維粒子軌道對應在旋轉坐標系中等效標定極向截面的二維軌道投影(圖8、9)。從圖8可看出，仿星器中的通行粒子軌道與托卡馬克中的通行粒子軌道特征相似，但仿星器中的通行粒子在旋轉等效標定極向截面上繞幾圈后才閉合。從圖9可看出，仿星器中的捕獲粒子與托卡馬克中的捕獲粒子類似地形成了香蕉軌道，但在旋轉等效標定極向截面上仿星器中香蕉軌道并未閉合，且軌道逐漸向外側漂移，最終有可能飛出最外層閉合磁面導致粒子損失。

a——三維視圖；b——俯視圖圖6 高能量通行粒子全軌道模擬圖Fig.6 Full-orbit simulation of energetic passing particle

a——三維視圖；b——俯視圖圖7 高能量捕獲粒子全軌道模擬圖Fig.7 Full-orbit simulation of energetic trapped particle

圖8 高能量通行粒子軌道在旋轉坐標系下的 等效標定極向截面上的軌道投影Fig.8 Orbit projection of energetic passing particle on equivalent calibrated poloidal cross-section in rotating coordinate system

采用上述方法對不同初始位置的測試離子進行大量試算，以進一步證明此方法的有效性。例如，設置測試離子的初始位置的導心是環向角為φ=0°時圖5a所示磁面上的點B，其旋轉坐標系中的標定極向截面坐標為(R′=0.00 m，Z′=-0.12 m)，初始投擲角為η=60°，則該測試離子在旋轉標定極向截面的軌道投影如圖10所示。雖然該測試離子的初始投擲角仍為60°，但其所在初始位置的磁感應強度變大，則其臨界投擲角ηc相應變大，從而滿足條件η<ηc，此時該測試離子就變為通行粒子。

圖9 高能量捕獲粒子軌道在旋轉坐標系下的 等效標定極向截面上的軌道投影Fig.9 Orbit projection of energetic trapped particle on equivalent calibrated poloidal cross-section in rotating coordinate system

圖10 高能量粒子軌道在旋轉坐標系下的 等效標定極向截面上的軌道投影Fig.10 Orbit projection of energetic particle on equivalent calibrated poloidal cross-section in rotating coordinate system

4 結論

本文采用HELIAC程序計算了H-1NF仿星器在標準運行模式下不同環向角的真空磁場龐加萊圖，研究其磁場位形特點，并模擬計算出高能量離子在該標準磁場位形中典型的通行粒子和捕獲粒子三維運動全軌道。基于H-1NF仿星器標準磁場位形的磁軸位置和磁面沿環向角的旋轉變化規律，以磁軸為旋轉軸，按照旋轉規律將不同環向角的極向截面旋轉后得到旋轉坐標系下的等效標定極向截面，并將高能量離子的三維運動軌道投影到此等效標定極向截面上，從而能更加清晰地顯示出高能量離子在該磁場位形中的運動軌道特征。結果表明，H-1NF仿星器中的通行粒子軌道特征與一般托卡馬克中的通行粒子軌道特征相似，不同之處是仿星器中的通行粒子軌道在等效標定極向截面上繞幾圈后才閉合；H-1NF仿星器中的捕獲粒子類似于一般托卡馬克中的捕獲粒子形成的香蕉軌道特征，但H-1NF仿星器中香蕉軌道并未閉合，且軌道逐漸向外側漂移，最終可能導致粒子損失。

感謝澳大利亞國立大學Boyd Blackwell博士和Clive Michael博士為本工作提供的HELIAC程序代碼及建議。