高能X射線光譜分析工具—SasalPy

郝明月,梁貴云,王 凱,毛俊捷

(1.河北大學,保定 071002; 2.中國科學院國家天文臺,北京 100101; 3.清華大學,北京 100084)

1 引 言

高能X射線研究是天文學家探究高能天體物理過程的主要研究手段，常見的X射線輻射源如太陽或者恒星的星冕和耀斑、超新星遺跡、激變變星、X射線雙星和黑洞周圍的吸積盤、星系、活躍星系核、星系團等.隨著現代天文的蓬勃發展，尤其是XMM-Newton、Chandra和Hitomi等天文望遠鏡的觀測為天文學家研究高能X射線、探究和認知宇宙高能物理過程提供了重要的觀測數據.

建立在各種物理條件下的等離子體模型是天文學家解析反演上述天文望遠鏡觀測到的數據的重要工具，是現代天文學研究工作的基礎.目前，被廣大天文學家熟知的等離子體模型有CHIANTI[1]、APEC[2]、Raymond-Smith[3]、MEKAL[4]、SPEX[5]、XSTAR[6]、Cloudy[7].其中，由Landi等[8]人提出的CHIANTI主要關注的研究對象是太陽的極紫外輻射光譜；由Smith等[2]人提出的APEC[2]最早用于計算輻射率和冷卻曲線，是當前美日學者分析X射線能譜數據的首選方法；由Raymond[3]提出的Raymond-Smith等離子體模型適用于星際介質；由Mewe等[4]人提出的MEKAL是高能X射線研究早期被廣泛應用的等離子模型；目前為天文學家廣知的XSPEC[9]、SPEX[5]都包涵MEKAL的基本模型，XSPEC是由Arnaud等[9]提出的光譜擬合工具，它除了包含MEKAL模型外，同時還集成了APEC[2]和XSTAR[6]模型，由Kallman和McCray提出的XSTAR主要用于研究由光電離等離子體產生的X射線能譜；SPEX是由Kaastra等人[5]提出的光譜擬合工具包，它包含的碰撞電離平衡模型延續了Mewe等人[4]的工作，同時它還包含非電離平衡模型和有限的光電離平衡模型；對于光電離等離子體來講，由Ferland等人[7]提出的Cloudy模型也被天文學家應用，近期該模型也向X射線波段擴展.

本工作中描述了一種全新的等離子體光譜分析模型，即SasalPy，它是為了滿足新一代X射線望遠鏡海量數據和極高能譜分辨需求提出的新的等離子體光譜分析工具，該模型主要是基于Python語言并在SASAL[10](Spectral Analysis Sys-tem for Astrophysical and Laboratory)模型的IDL版本基礎上發展的新的X射線光譜分析工具包，以激發/退激發、復合、電離三類原子物理過程為基礎構建熱碰撞電離平衡(Collisional Ionisation Equilibrium，簡稱CIE)模型.對于熱等離子體來講，其微觀物理過程構成具體包含電子和質子的激發和退激發、光子激發和受激輻射、輻射復合和雙電子復合、三體復合和碰撞電離等原子物理過程，熱等離子體廣泛存在于恒星星冕、星系、星系團、超新星遺跡等天體中.SASALDB是SasalPy數據交互的基礎數據庫，涵蓋所有天體高豐度元素，即從氫(Z=1)到鋅(Z=30)30號元素的相關原子參數.SASALDB數據庫的持續更新和發展有助于等離子體模型SasalPy精確性提升.SasalPy仍在開發當中，下文主要從SasalPy當前已經完成的碰撞電離平衡模型的理論、主要框架和功能展開介紹.第二部分簡要介紹SasalPy的碰撞電離平衡模型的理論基礎，第三部分主要介紹碰撞電離平衡的總體框架和主要功能以及與其他等離子體模型的對比研究，第四部分是對SasalPy的碰撞電離平衡模型研究的結果與討論.

2 碰撞電離平衡模型

假設等離子體光學薄并且無外部輻射場的影響，碰撞電離平衡狀態下等離子體的電離過程和復合過程處于動態平衡狀態，即由復合捕獲的自由電子數等于由電離產生的自由電子數.該狀態下主導的電離過程為直接碰撞電離和激發自電離，主導的復合過程為輻射復合和雙電子復合.當等離子體處于低溫時，電荷交換過程的影響不可以被忽略.電荷交換過程不在本文研究范圍，后續的研究中會加入到模型中.

2.1 電離平衡

電離平衡主要是描述不同電荷態離子相對于離子對應元素總量在一定溫度和一定密度下的相對含量，或者稱其為電荷態分布.電離平衡過程可用梁貴云等[10]工作中第二部分(1)-(6)式描述的速率方程表述，通過求解方程即可得到特定溫度和密度條件下該電荷態的分布情況.不考慮光致電離過程與電荷交換過程的速率方程可表述為：

(1)

(2)

(3)

(4)

2.2 線輻射

線輻射對于熱等離子體來說起著至關重要的作用.在特定情況下，線輻射在總的輻射量中占據主導地位.由碰撞激發產生的線輻射可以分為離子外層電子的激發和輻射衰變兩個階段，而由復合過程引起的線輻射是在電子與離子發生復合后接續發生輻射衰變產生的.在發生輻射之前，離子需要先通過碰撞激發、吸收光子、內殼層電子電離以及輻射復合或雙電子復合過程達到激發態，隨后處于激發態的離子通過輻射回到基態或者亞穩態，產生線輻射.假設原子或離子先激發至激發態j，隨后經過輻射躍遷回到基態或能量較低的能級i，單位時間處于能級j的電子由于自發輻射躍遷至較低能級i的速率為Aji(單位：s-1)，對于不同能級的輻射躍遷，輻射躍遷速率不同.線輻射率可以表述為：

εji=ΔEnjAji

(5)

ΔE是激發能量，即能級i與能級j之間的能級能量差，nj是經過歸一化的處于能級j的相對離子數密度，即nj=Nj/∑Nk，Nj是能級j的能級布居，∑Nk表示所有求解目標能級的布居總和.

2.3 連續輻射

連續輻射主要包括軔致輻射、復合輻射和雙光子輻射.軔致輻射可以概述為等離子體中的自由電子在離子的庫倫場作用下做加速運動產生的輻射.Nozawa等[11]工作中(23)式很好的描述了相對論條件下的軔致輻射，本文(6)式在此式基礎上應用了Itoh等[12]計算的相對論條件下的岡特因子：

(6)

(7)

(8)

其中，h是普朗克常數，c是光速，k是玻爾茲曼常數，λ是輻射光子對應的波長，Ael是離子對應的元素豐度，FZj是離子對應溫度T下的電荷態分布數據，Zj為原子序數，此處j特指對應的離子，即離子電荷量加1，gzj，itoh為岡特因子，me為電子質量.本文(9)式描述的是采用了Sutherland等[13]計算的岡特因子的軔致輻射：

(9)

(10)

由Itoh等[12]計算的相對論岡特因子適用的溫度范圍為106K≤T≤108.5K，適用的波長范圍為14.39 ?≤λ≤4549.81 ?.本文中對于超出以上溫度和波長范圍的輻射采用由Sutherland計算[13]的岡特因子.

復合輻射是由具有一定動能的自由電子經過高電荷態離子近旁時，被離子捕獲復合成較低電荷態離子，同時產生的輻射.復合輻射(自由-束縛)與光電離(束縛-自由)互為逆過程，借助復合輻射截面(自由-束縛碰撞截面)與光電離截面的Milne關系[14]：

(11)

其中，σbf為光電離截面，σfb為復合輻射截面，ve為電子的速度，Eγ是輻射光子的能量，ω0是復合前離子的基態的統計權重，ωi是復合后離子的能級i的統計權重.

復合輻射過程產生的輻射可由光電離截面求得.本文中的復合輻射計算包含了Verner等[15]計算的基態光電離截面和Karzas等[16]計算的激發態光電離截面.復合輻射可由下式表述：

(12)

(13)

雙光子輻射是指電子由激發態輻射兩個光子回到基態，該過程對于較低電子密度的類氫和類氦等離子體尤為重要.光學薄條件下的雙光子輻射過程產生的輻射采用Young等人[17]工作中(11)式表述，雙光子輻射由下式表述：

(14)

其中，Aji為輻射躍遷速率，Ael為元素的相對豐度，FZj為電荷態為j的離子在溫度T下對應的電荷態分布數據，Nj(X+q)為處于能級j的帶+q電荷的離子的數密度，λ為波長，ne為電子密度，λ0為激發態與基態間能量差對應的波長，λ0/λ表示波長為λ的光子的能量與對應的激發態與基態間能量差的比值，φ(λ0/λ)為譜分布函數[18，19]，即雙光子躍遷速率隨光子能量的分布.

3 碰撞電離平衡模型結構及其應用

圖1為SasalPy的碰撞電離平衡模型的理論概述圖，碰撞電離平衡模型主要包括電離平衡計算或者電荷態分布計算(Charge stage distribution，簡稱CSD)、線輻射率計算(Line emissivity)、連續輻射率計算(Continuum emissivity)、合成光譜(Spectrum)以及輻射損失(Radiative loss)計算，對應的程序模塊依次為getCsdCie、getLineEmiss、getContEmiss、getSpectrum和getRadLoss.圖中E-E是電子碰撞激發過程，Ph-E是光子激發過程，P-E是質子激發過程，DI是直接碰撞電離，EAI是激發自電離過程，RR是輻射復合過程，DR是雙電子復合過程，TBR是三體復合過程.

圖1 SasalPy碰撞電離平衡模型理論概述圖.數字序號表示不同程序接口與模型最上層模塊的關聯性

Abundance是天體元素豐度，微分發射度(Differential emission measure，簡稱DEM)和發射量度(Emission measure，簡稱EM)是理論光譜與觀測光譜建立聯系的關鍵參數，用于描述天體X射線輻射源溫度結構與發射量度之間的關系[20，21].Response由觀測望遠鏡的響應幾率和有效面積組成，取決于望遠鏡的性能.Broaden表示光譜展寬的方法，一般采用高斯展寬.光譜合成計算基于線輻射率和連續輻射率計算結果，輻射率經由電離平衡(CSD)、元素豐度(Abundance)與DEM或EM加權求和，而后經過展寬(Broaden)或與儀器的響應(Response)卷積便可得到模擬的合成光譜.圖2演示SasalPy獲取光譜的流程(以Fe XIV為例)：

圖2 SasaPy獲取Fe XIV離子光譜的演示

1.第一步(stage 1)豐度數據設定，對應圖1中的序號5；

2.第二步(stage 2)電離平衡數據設定或計算，對應圖1中的序號1；

3.第三步(stage 3)獲取特定溫度、密度網格有限波長范圍內的線輻射率數據，通過參數read選擇計算(read=0)或者讀取(read=1)相應數據，對應圖1中的2；

4.第四步(stage 4)調用光譜合成模塊(getSpectrum)計算，計算前需設定波長網格(lambgrid)和譜線展寬(instrBroaden)，通過參數doContinuum設定是否包含連續輻射計算，photons參數用于設定輸出單位形式(photons=1時單位為光子形式)；

5.第五步(stage 5)顯示光譜結果.

以上從程序層面簡要描述了getSpectrum獲取光譜的過程，下面將進一步描述各個模塊的實現情況.

3.1 電離平衡計算

電離平衡計算模塊(getCsdCie)核心計算過程主要在于求解由微觀物理過程中的速率系數構成的速率方程，由碰撞等離子體主導的電離平衡計算主要考慮的微觀物理過程主要有電子、質子與粒子的碰撞激發及其逆過程、直接碰撞電離、激發自電離、雙電子復合、輻射復合和電荷交換.電荷交換過程暫未包含到模型中來，它將會在下一步工作中包含到模型中.SasalPy使用Scipy中的bicg求解速率方程AN=0，以一維單位向量作為bicg求解時的初始猜測，減少求解所需的迭代次數.電離平衡條件下的電荷態分布結果表明，元素從中性至裸離子的離化產生限于一定的溫度范圍，定義該溫度范圍為該元素的電荷態分布的有效溫度范圍.圖3給出了高能天體物理研究中所關注的30號元素的電荷態分布的有效溫度范圍，溫度低于有效范圍，元素將以中性粒子形式存在，元素電離度低；溫度高于有效范圍，元素將以裸離子形式存在，元素電離度高.

圖3 原子序數30以內的元素的電荷態分布的有效溫度范圍

除了了解元素的電荷態分布的有效溫度范圍外，了解元素對應的所有電荷態離子在電離平衡條件下存在分布的溫度范圍也很重要，該溫度范圍有助于模型判斷特定溫度下在計算中需要考慮的離子，參考該溫度范圍，程序在運算過程中自動判斷選擇需要加載的離子參數，可以有效節省運算時間.SasalPy取0.001作為電荷態分布的最小值，保留離子的電荷態分布大于0.001的溫度區間內的數據.

圖4是SasaPy計算的氧和鐵的電荷態分布與Chianti數據庫[22]、AtomDB數據庫[23]以及梁貴云等人[24]計算的電荷態分布的對比.其中，氧元素在低溫區域(<105K)隨著溫度的降低，本文計算的電荷態分布與Chianti[22]、AtomDB[23]和梁貴云等[24]的結果隨著溫度的降低表現出的差異性逐漸增大，此溫度區間主要介于紅外波段與紫外波段之間，氧元素低電荷態離子(O I、O II、O III)對X射線波段影響很小，此處的差異可以忽略；SasalPy計算的位于高溫區域(>105K)的氧離子與Chianti[22]、AtomDB[23]的結果符合的很好，與梁貴云等[24]表現出的差異主要是因為其計算結果對應的溫度網格較為稀疏，對比中由插值計算帶來的影響，但梁貴云等[24]的實際計算值與SasalPy能夠很好的符合.由于Chianti[22]、AtomDB[23]和梁貴云等[24]計算的電荷態分布結果對應的溫度范圍僅介于104—109K之間，超出此溫度范圍的結果由外推法插值計算獲得，SasalPy計算的鐵元素的電荷態分布與前述三種結果在低溫區小于104K處表現出的較大差異主要由插值計算引起，溫度越低則插值精度越低，差異越大.如上述氧的討論中所述，處于更低溫度的低電荷態鐵離子對X射線波段的影響可以忽略.Urdampilleta等人[25]認為電荷態分布數據之間差異小于10-20%是好的符合.按照相同的評定標準，在低溫區(104—105K)：SasalPy計算的Fe II、Fe III和Fe IV的電荷態分布與AtomDB、Chianti和梁貴云等人的結果符合的較好，本文計算的Fe V的結果與Chianti和梁貴云等的結果符合的很好，與AtomDB結果的離子比例峰值相差24.8%；在高溫區(107.5—109K)：本文計算的Fe XXV、Fe XXVI和Fe XXVII與AtomDB、Chianti和梁貴云等符合的很好；在溫度區域105—107.5內：SasalPy計算的結果與AtomDB、Chianti和梁貴云等人的結果差異明顯，其中Fe XV的離子比例峰值達到了梁貴云等人結果的2.2倍.

圖4 元素氧和鐵的電荷態分布.Sa.表示由SasalPy計算的電荷態分布結果，Ch.表示Chianti數據庫中的電荷態分布結果(chianti.ioneq)，At.表示AtomDB數據庫中3.0.7版本電荷態分布結果，L14.表示梁貴云等計算的電荷態分布結果，圖中底部描述的是Sa.與其余三者的差異，即不同模型結果的差值.不同的顏色表示不同的電荷態的離子

電離平衡數據即電荷態分布作為線輻射、連續輻射計算的重要參數，該參數中存在的誤差會影響光譜的計算結果.電離平衡取決于電離過程和復合過程，電離速率和復合速率基于理論計算和實驗測量產生，理論計算時能級的選取和能級數據的限制以及實驗測量誤差均會對電離和復合速率產生影響，同時不同的理論模型也會引入一定的系統誤差，這就導致了各個模型收錄或者計算的電離平衡數據存在差異[26].為了進一步對電離平衡數據進行差異性探討，本文選取氧和鐵元素的類氦離子的電離平衡數據進行對比，O VII和Fe XXV對X射線波段的貢獻高于附近其他電荷態的離子.具體做法是，以SasalPy計算的電離平衡數據為基礎，選取各離子電離平衡數據大于0.001的溫度范圍，對不同的電離平衡數據進行對比，見圖5.為了方便對比，本文在溫度范圍內以插值的方法統一溫度網格.

圖5 氧和鐵類氦離子不同數據來源的電離平衡數據對比.Sa.、Ch.、At.和L14.表示的含義與圖4相同

SasalPy計算的類氦氧離子(O VII)的電荷態分布與Chianti、AtomDB以及梁貴云等計算的結果符合的很好，對比中表現出的振蕩式的差異曲線主要由于插值計算引起，由于梁貴云等計算的結果采用的溫度網格較其余三者稀疏，差異曲線表現出的振蕩程度最強，表明更加精細的溫度網格有助于提高電荷態分布結果的精確性.SasalPy計算的類氦鐵離子電荷態分布在曲線上升部分與Chianti、AtomDB以及梁貴云等的結果差異明顯，本文計算的結果與Chianti和AtomDB最大相差64.7%，與梁貴云等的結果最大相差65.2%；對于曲線下降部分本文結果與Chianti和梁貴云等的結果最大相差13%，與AtomDB最大相差21.6%.Heuer等[26]以PyAtomDB為基礎通過數值模擬的方法指出系統誤差與溫度偏移存在著聯系.類氦鐵離子的電荷態分布曲線表明SasalPy計算的電荷態分布曲線峰值溫度相比AtomDB、Chianti、梁貴云等的結果向低溫區產生了偏移.SasalPy采用的電離平衡計算方法與梁貴云等[10]人的方法相似，主要不同的地方在解AN=0線性方程時采用的方法不同，SasalPy主要應用Python第三方庫Scipy中的bicg函數求解方程，而梁貴云等人采用的是IDL語言中的LINBCG函數求解方程.此外，SasalPy與梁貴云等人采用了相同的電離速率數據和復合速率數據，本文計算的電荷態分布與梁貴云等的結果的差異可能由求解線性方程組的方法不同引起.

SASALDB數據庫除了本文計算的電離平衡數據外，還收錄了多種電離平衡數據供用戶選擇，詳見表1.

表1 SASALDB收錄的其他電荷態分布數據

3.2 線輻射率計算

線輻射率計算模塊(getLineEmiss)旨在計算線輻射率隨波長或者能量的分布，計算核心在于求解某一電荷態離子某一能級j的歸一化的能級布居，求解能級布居時需考慮電子和質子與離子的碰撞激發及去激發、輻射復合和雙電子復合、直接碰撞電離以及激發自電離的影響.解得的能級布居與輻射躍遷速率、光子能量的乘積即是線輻射率，如(5)式所示.對于熱等離子體來講，線輻射的貢獻是不容忽略的.原子參數的不同影響線輻射率計算，例如能級能量、輻射躍遷速率、碰撞強度的不同會對線強產生顯著影響.圖6是類氦碳離子C V在10—60 ?波段內對線輻射貢獻較大的四種躍遷在單一溫度(106K)低密度條件下的線輻射結果.圖中標注的四種躍遷根據波長由短到長對應的上能級分別是1s3p1P1、1s2p1P1、1s2p3P1、1s2s3S1，躍遷對應的下能級均為1s21S0.表2列出了四種躍遷在SASALDB與Chianti數據庫[22]中對應的能級、波長、輻射躍遷速率(A Value)和差異以及四種躍遷對應的有效碰撞強度.

表2 C V四種躍遷對應的輻射躍遷速率和有效碰撞強度

圖6 C V四種躍遷的線輻射率.Ch.C V是使用ChiantiPy基于Chianti數據庫計算的結果；Sa.1 C V是SasalPy基于SASALDB數據庫計算的結果；Sa.2 C V是SasalPy使用Chianti數據庫中C V的能級數據與輻射躍遷速率計算的結果；Sa.3 C V表示的結果中SasalPy使用了Chianti數據庫C V的碰撞強度；Sa.4 C V表示的結果中SasalPy同時使用了Chianti數據庫中C V的能級能量、輻射躍遷速率以及碰撞強度數據

從圖6中可以看到，SasalPy與ChiantiPy對于躍遷1s21S0-1s3p1P1的線輻射率計算結果符合的很好.躍遷1s21S0-1s2p1P1對應的Ch.C V、Sa.1 C V、Sa.2 C V線輻射符合的較好，Sa.3 C V由于使用Chianti數據庫的碰撞強度，線輻射率相比Sa.1 C V提升了4.74%；同理，Sa.4 C V由于碰撞強度不同的影響相比Sa.2 C V提升了30.95%；Sa.2 C V由于使用Chianti數據庫的能級能量和輻射躍遷速率數據，線輻射率相比Sa.1 C V降低了0.078%，變化可以忽略；Sa.4 C V由于能級能量和輻射躍遷速率不同的影響相比Sa.3 C V提升了24.9%.對于躍遷1s21S0-1s2s3S1，Sa.3 C V相比Sa.1 C V提升了0.78%，Sa.4 C V相比Sa.2 C V降低了1.18%；Sa.2 C V相比Sa.1 C V提高了4.03%，Sa.4 C V相比Sa.3 C V降低了16.25%.對于躍遷1s21S0-1s2p3P1，碰撞強度數據的不同導致Sa.3 C V相比Sa.1 C V提升了18.5%，Sa.2 C V和Sa.4 C V與ChiantiPy符合的很好.對比結果表明，原子參數如輻射躍遷速率、有效碰撞強度對于較強躍遷的線輻射率計算會有顯著影響.

3.3 連續譜輻射率計算

連續譜輻射率計算模塊(getContEmiss)同以上描述的電離平衡計算模塊和線輻射率計算模塊是碰撞電離模型的核心計算模塊，進一步的光譜合成計算模塊和輻射損失計算模塊將依賴于這三個模塊的運算結果.連續輻射主要包括軔致輻射、復合輻射以及雙光子輻射.

圖7為SasalPy與ChiantiPy在溫度為106K和低密度(1 cm-3)條件下的連續譜結果對比，SasalPy與ChiantiPy使用了相同的豐度數據[27]和電荷態分布數據，電荷態分布數據使用Chianti數據庫[22]中的chianti.ioneq文件.可以看到，SasalPy與Chianti在軔致輻射和復合輻射計算上符合的很好，SasalPy的軔致輻射在應用Nozawa等[11]計算的相對論岡特因子的基礎上，用Sutherland等[13]提出的岡特因子進行了修正，方法與Chianti選用的方法相同.對于復合輻射，SasalPy和ChiantiPy應用相同的光電離截面，即Karzas等計算的光電離截面[16].軔致輻射及復合輻射存在微量的差異，一方面是因為SasalPy和ChiantiPy在處理相關系數時約化取值存在不同，例如SasalPy中光速取值為2.997825×1010cm·s-1，而ChiantiPy中取值為2.99792458×1010cm·s-1，再如根據hc/λ對波長進行單位轉換時(KeV與?相互轉換)，SasalPy對hc的取值為12.398，ChiantiPy對hc的取值為12.39841875；另一方面是因為SasalPy與ChiantiPy在電荷態分布數據插值方法上的不同，為了避免插值計算時負值的出現，SasalPy在電荷態分布數據插值取值時將其轉換為以10為底的對數進行插值計算，ChianitPy則是將電荷態分布數據轉換為以e為底的對數進行插值計算，最終的結果在合理的誤差范圍.對連續譜差異貢獻最大的主要是雙光子輻射.雙光子輻射主要包括類氫和類氦離子的貢獻.引起雙光子輻射差異的主要是類氦離子，類氫碳離子和類氦碳離子在雙光子輻射中貢獻最多，以下以C V和C VI為例分析雙光子輻射的差異.雙光子輻射計算中，SasalPy與ChiantiPy使用相同的譜分布數據[18，19]，即相同的譜分布函數φ(λ0/λ)，見(14)式.

圖7 SasalPy與ChiantiPy連續譜對比.圖中Ch.total和Sa.total分別為ChaintiPy和SasalPy總的連續譜結果；Ch.ff和Sa.ff對應的是軔致輻射；Ch.fb和Sa.fb對應的是復合輻射；Ch.tp和Sa.tp對應的是雙光子輻射

圖8左圖是類氫碳離子C VI的雙光子輻射對比，Ch.C VI tp是ChiantiPy的計算結果，ChiantiPy在106K時對應的經插值后C VI的電荷態分布數據為0.5855，Sa.C VI tp a是SasalPy的計算結果，SasalPy在106K時對應的經插值后的C VI電荷態分布數據為0.5481，僅由插值方法不同引起的差異為6.4%.Sa.C VI tp b是SasalPy使用與ChiantiPy相同的電荷態分布數據0.5855計算的結果，SasalPy與ChiantiPy在使用相同的差值方法后結果符合的很好，說明C VI的雙光子輻射差異僅由插值方法的不同引起.圖8右圖是類氦碳離子C V的雙光子輻射對比，Ch.C V tp是ChiantiPy的計算結果，ChiantiPy在106K時對應的經插值后C V的電荷態分布數據為0.2677，Sa.C V tp a是SasalPy的計算結果，SasalPy在106K時對應的經插值后的C V電荷態分布數據為0.2744，此時僅由差值方法引起的差異為2.5%，但是與由能級數據和輻射躍遷速率引起的差異相比，可以忽略不計.Sa.C V tp b是SasalPy使用與ChiantiPy相同的輻射躍遷速率和能級數據計算的結果.可以看到，SasalPy采用與ChiantiPy相同的輻射躍遷速率與能級數據得到的結果與ChiantiPy符合的很好，但Sa.C V tp b與Ch.C V tp依然存在18.4%的差異，這是由于SasalPy與ChiantiPy的數據庫收錄的C V的電離、復合以及碰撞強度數據不同導致的.相比由輻射躍遷速率與能級數據不同引起的差異，由C V的電離、復合以及碰撞強度數據不同導致的差異可以忽略.對于C V，SasalPy的數據庫中目前僅收錄了主量子數n=2-4的49個能級的數據[28]，Chianti數據庫收錄了主量子數n=2-8的577個能級[29-31].圖9中SasalPy與ChiantiPy的數據庫收錄的C V的輻射躍遷速率在10-16 ?波段的對比表明SasalPy的輻射躍遷速率多數小于Chia-ntiPy，SasalPy在求解能級布居時自發輻射對能級布居的影響遠小于ChiantiPy，故由SasalPy計算的C V的能級布居結果要遠大于ChiantiPy，進而使得由SasalPy計算的C V的雙光子輻射Sa.C V tp a高于ChiantiPy的結果Ch.C V tp.

圖8 C VI與C V由不同等離子體模型計算的雙光子連續譜對比.上圖是類氫碳離子C VI的雙光子輻射對比，Ch.C VI tp是ChiantiPy的計算結果，Sa.C VI tp a是SasalPy的計算結果，Sa.C VI tp b是SasalPy使用與ChiantiPy相同的電離平衡數據計算的結果；下圖是類氦碳離子C V的雙光子輻射對比，Ch.C V tp是ChiantiPy的計算結果，Sa.C V tp a是SasalPy的計算結果，Sa.C V tp b是SasalPy使用與ChiantiPy相同的輻射躍遷速率和能級數據計算的結果

圖9 C V在10-60 ?波段來自不同數據庫的輻射躍遷速率對比

3.4 合成光譜

圖10 6.199-123.985 ?波段內的合成光譜.光譜b的溫度為106.125 K，DEM數值為0.23，電子密度為2×108 cm-3；光譜c的溫度為106.875 K，DEM數值為6.99；光譜d的溫度為107.525 K，DEM數值為0.53，光譜c和d的電子密度同為1012 cm-3；光譜a為光譜b、c、d的合成譜，DEM的單位見右圖.右圖中三個顏色區域由左至右依次對應低溫區、中溫區和高溫區

3.5 輻射損失

輻射損失計算模塊(getRadLoss)主要用于計算由微觀物理過程引起的等離子體能量損失隨溫度的分布.輻射損失計算的核心計算主要是線輻射率的計算和連續輻射率的計算，在考慮連續輻射的影響時，僅考慮軔致輻射和復合輻射的貢獻.隨著X射線光譜儀性能的不斷提升，使得我們可以觀測到更為精細的天體結構信息.與此同時，新的觀測數據向等離子體模型提出了挑戰.新的觀測數據捕捉到的精細天體結構信息需要合理的等離子體模型和更加精確的原子數據進行解析.不同的等離子體模型由于其本身采用的理論和近似方法以及原子數據的不同，在輻射損失計算時會表現出明顯的差異性[33].圖11是SasalPy與不同的等離子體模型的輻射損失或冷卻曲線的對比.

圖11 不同等離子體模型之間的輻射損失對比.Sa.表示SasalPy的計算結果，Ch.表示Chianti的計算結果，Sp.表示SPEX的計算結果，Cl.表示Cloudy的計算結果，At.表示PyAtomDB的計算結果，Sc.表示Schure等人的計算結果

圖11中SPEX、PyAtomDB、Cloudy、Schure四種等離子體模型的輻射損失或冷卻曲線源自于Stofanova等人[33]的工作.表3匯總了不同等離子體模型在輻射損失計算中使用的豐度數據和電荷態分布數據.

表3 不同等離子體模型計算輻射損失使用的豐度數據與電離平衡數據

SASALDB數據庫目前更新了類碳[34]、類氮[35]、類氧[36]等電子系以及類氫離子和類氦離子[37]的理論碰撞激發數據，即能級數據、輻射躍遷速率以及碰撞強度.通過NIST數據庫[38]的實驗能級數據對SASALDB中的離子的波長進行了修正.SasalPy的輻射損失結果是基于更新后的數據庫得到的，輻射損失計算用到的具體離子見圖12.將輻射損失以0.1 keV為界限將輻射損失結果分為兩部分進行討論.可以看到，當溫度大于0.1 keV時，SasalPy的輻射損失結果與PyAtomDB在這一部分符合的很好.SasalPy與Cloudy在0.1-10 keV范圍符合的很好，大于10 keV時差異增大，差異最大達到29.4%.在高溫區(大于0.1 keV)，連續輻射中的軔致輻射主導地位加深，Stofanova等人[33]認為，由于不同的等離子體模型的軔致輻射采用不同修正方法計算的岡特因子，例如非相對論、半相對論和相對論修正等，會導致不同程度的差異.SasalPy的軔致輻射采用Nozawa等[11]計算的經過相對論修正的岡特因子，這一點與PyAtomDB和Chianti的輻射損失結果采用的岡特因子相同.SasalPy與PyAtomDB的輻射損失結果在高溫區符合的很好，差異最大為8.25%，出現在溫度最高處.Chianti的輻射損失在高溫區略低于SasalPy，溫度最高處與SasalPy的差異為7.24%，此處差異最大.在0.1-10 keV范圍內主導輻射損失的為線輻射，SasalPy與SPEX、Schure等[39]差異變化明顯是因為不同等離子體模型采用的線輻射計算的基本參數不同，這其中輻射躍遷速率和碰撞強度數據的不同導致的影響不可忽略.SasalPy在高溫區與SPEX、PyAtomDB、Chianti以及Schure等人[39]的輻射損失表明，SasalPy在X射線波段的數據是足夠精確的.

圖12 SasalPy輻射損失計算用到的離子.藍色的五角星表示SASALDB原始離子數據，紫色方塊表示更新后的離子數據；q是離子所帶電荷量

當溫度小于0.1 keV時，輻射損失主要由線輻射里的碰撞激發過程主導.在0.002-0.1 keV范圍內，SasalPy的輻射損失相比于其他等離子體模型最小.在0.001-0.015 keV范圍氫原子的貢獻占據主導地位，Stofanova等人[33]認為氫原子的碰撞強度的不同對氫原子對輻射損失的影響較大，他們換用Chianti數據庫[1]新的氫原子碰撞強度，氫原子對輻射損失的貢獻在溫度為2 eV處的減少了3倍，得到了與PyAtomDB與Cloudy符合較好的結果.SASALDB中氫原子的碰撞強度采用的是未修正過的Anderson等[40]計算的數據.因此，SasalPy在0.015 keV左右相比其他等離子體模型要低很多.隨著SASALDB數據庫的碰撞激發數據的更新，SasalPy在小于0.1 keV的溫度范圍將會有所改善，相關工作將會在后續工作中開展.經過上述分析，從不同等離子體模型的輻射損失結果表現出的差異性來看，等離子體模型和原子數據的持續精確地更新研究是尤為重要的.

4 結 論

提出了一種新的基于Python語言的適用于高能X射線的等離子體模型SasalPy，簡要闡述了SasalPy的碰撞電離平衡理論，主要就碰撞電離平衡模型總體框架和電離平衡計算、線輻射率計算、連續輻射率計算、合成光譜及輻射損失計算五個重要模塊作了詳細描述，從功能性或差異性方面說明了SasalPy的可用性和精確性.此外，電離平衡計算對氧元素和鐵元素對比討論表明，由SasalPy的電離平衡計算的氧元素和鐵元素的部分溫度區域的結果與其他等離子模型高度相符.另外，通過對比來自不同數據源的類氦離子電荷態分布數據，可以為用戶選取電離平衡數據提供參考.線輻射計算中對C V四種躍遷產生的線輻射率的討論表明，原子參數如躍遷速率、碰撞強度等的不同會對線輻射率和線強產生影響.對于這種現象，SasalPy的未來版本將支持用戶選擇并測試不同的原子物理數據，填補目前其他等離子體模型不具備的功能應用空缺.SasalPy與ChiantiPy的連續譜的對比表明，由輻射躍遷速率數據的不同產生的影響對類氦碳離子的雙光子輻射不可忽略.采用Capella星冕的低、中、高三個不同溫度區域的體微分發射度，實現了光譜合成應用.SasalPy與PyAtomDB、SPEX、Chianti、Cloudy以及Schure等人的輻射損失對比表明，Sasalpy在X射線能段的數據與其他等離子體模型極度相近.SasalPy作為自主研發的新一代等離子體模型將會為由我國主導的下一代X射線高分辨能譜空間望遠鏡[41](Hot Universe Baryon Surveyor，簡稱HUBS)的高分辨率X射線光譜分析提供有力支持，該模型在高能X射線研究領域具有較高的應用價值.