熱離子能量轉換器銫電弧工況下電子勢能的分布特性

鐘武燁，呂 征，鄭劍平，楊啟法

(中國原子能科學研究院 反應堆工程技術研究部，北京 102413)

熱離子反應堆電源是采用熱離子能量轉換方式將反應堆裂變釋放的熱能直接轉換成電能的能量轉換裝置，其具有比功率高、無機械轉動、所需比散熱面積小等優點，是國際上最具優勢的高效空間核電源之一[1]。目前，具有實際應用價值的熱離子能量轉換器的結構形式是充銫的工作在電弧工況的熱離子二極管[2]，其通常采用難熔金屬作電極對，發射極被加熱到高溫(～1 500 ℃)產生熱電子發射，充入電極間隙的銫蒸氣在其中電離形成等離子體，電子輸運通過等離子體區域后達到接收極，經外電路的負載做功后重新回到發射極，形成回路。熱離子轉換器可看作以電子作為介質的熱機[3]，電極間隙的勢能分布則決定該熱機的能量來源。因此，電子在電極間隙的勢能分布是研究熱電轉換特性的關鍵。

目前，對關于電極間隙的電子勢能分布已有較深入的了解，但一些基礎問題仍有待進一步深入研究，如熱離子轉換器運行在電弧工況不同區域(包括飽和區、阻塞區)時，對于電極鞘層的電子勢能分布的認識仍不夠深入，忽略對阻塞區虛擬陰極的討論[4]；在討論阻塞區虛擬陰極時又無法給出數值上可解的方法[5]；缺乏對接收極電勢分布特征的研究等。導致該問題的原因主要是現有的模型并沒有很好地處理等離子體輸運及其邊界條件的匹配問題。本文基于動理學方法建立電弧工況下的等離子體輸運模型，對等離子體鞘層的電子勢能分布進行分析，并將等離子體鞘層處理為輸運的邊界條件，自編程實現等離子體輸運方程與邊界條件的自適應解耦求解。

1 理論建模

電子勢能分布的研究需建立等離子體的輸運模型。當熱離子轉換器工作在電弧工況(圖1)時，在電極間隙中，等離子體與電極表面相接觸的區域形成電勢突然躍變的、非準中性的等離子體鞘層，其尺度與德拜長度相當，遠小于電極間隙。因此，電極間隙的等離子體可分為電極附近的鞘層和中間的準中性等離子體區[6]。在具體的建模過程中，在準中性等離子體區建立帶電粒子的輸運模型，而等離子體鞘層則處理為輸運的邊界條件。

1.1 電弧工況等離子體輸運模型

熱離子能量轉換器電極間隙的等離子體屬于低溫弱電離的等離子體，此時必須考慮電子、離子以及中性原子3組分同時存在的多粒子碰撞輸運問題。等離子體中帶電粒子體系間存在復雜的相互作用，要精確描述等離子體的行為是極其困難的。目前，只能根據不同的條件和研究的問題，采用不同的近似方法[7]。文獻給出的建模方法可歸為唯象模型和基于介觀尺度的動理學模型兩大類[8]。唯象模型作為定性分析是有效的，作為定量分析計算仍有差距。本文的輸運方程采用動理學的模型。動理學模型是輸運過程嚴格的定量理論。動理學方法研究其輸運過程的基本思路是從粒子碰撞遷移所滿足的玻爾茲曼方程出發，然后根據轉換器等離子體輸運的具體特點對分布函數和碰撞輸運方程作合理簡化，建立輸運控制方程組，并給出輸運系數。

由于電極間隙遠小于電極的軸向尺寸，故只需建立電極間隙準中性等離子體輸運的一維模型，在建模過程中以發射極作為原點建立坐標系。中性銫原子按理想氣體處理：靠近電極的中性原子溫度等于電極溫度，電極間隙各點的溫度隨電極溫度呈線性變化，不必考慮中性原子的輸運問題，但必須考慮由于中性原子的存在導致的散射效應。等離子體輸運的描述的目的是獲得電極間隙等離子體濃度n、電子溫度Te、負電勢V、電子電流密度je及其能量流密度Se、離子電流密度ji及其能量流密度Si等分布狀態。

圖1 熱離子能量轉換器的工況分類及工作原理Fig.1 Operation mode and working principle of thermionic energy converter

首先，建立電子、離子、中性原子3組分共存條件下的玻爾茲曼方程：

式中：β=e,i,a(e、i、a分別表示電子、離子、中性原子)；α=e,i；v為速度；r為位置；e為電量；下標in表示非彈性碰撞；Jαβ為碰撞積分項；fα為分布函數并假定其偏離平衡態不遠。通過對玻爾茲曼方程的簡化，可建立等離子體參量的通量及連續性方程[9-10]。

1) 通量方程。帶電粒子的通量包括引起帶電粒子輸運的濃度梯度、溫度梯度以及弱電場力的貢獻，電子、離子的電流密度分別為：

為考慮離子電流密度的貢獻，總電流密度jtot包括電子電流密度和離子電流密度的部分。

jtot=je-ji

對于帶電粒子的能量通量，考慮了宏觀運動的動能、勢能以及熱擴散的貢獻。

式中：k為玻爾茲曼常數；De,i、μe,i、ke,i、βe,i分別為濃度擴散系數、遷移率、熱導率、熱擴散比；Rie為離子受電子的作用力。

2) 連續性方程。帶電粒子通量的連續性是由體電離效應決定的。

式中：Γ為電離函數；σ0(Te)為與電子溫度相關的等效電離截面；Na為原子濃度；n(Te)為熱平衡態下的等離子體濃度。

對于能量通量的連續性，重粒子在電離/復合過程中吸收/放出的能量由電子提供/帶走，忽略電子的輻射ΔSrad以及e-i、e-a的碰撞損失ΔSei、ΔSea，則有：

ΔSei-ΔSea=-Γ(Eion-eV)

式中，Eion為銫離子的第一電離能。

上述的通量及連續性方程組成了等離子體輸運的控制方程組，其中包含了je、ji、Se、Si、n、V、Te7個變量。對于其中的輸運系數，當分布函數偏離平衡態不遠時，基于動理學理論可對分布函數作Sonine多項式展開，從而將輸運系數分級展開，取其中的前3項即可使截斷誤差小于10%[11]。

1.2 電弧工況等離子體輸運的邊界條件

圖2 電弧工況飽和區可能的電子勢能分布Fig.2 Possible distribution of electron potential energy at saturation region of ignited mode

為獲得等離子體輸運微分方程的邊界條件，需建立帶電粒子通量、能量通量的邊界值與等離子體邊界參量的關系，其中的關鍵點是等離子體鞘層的電子勢能分布。在電弧工況的不同區域，等離子體鞘層電子勢能的分布狀態是有差異的，所有可能的勢能分布狀態如圖2、3所示，其中VE,C分別為發射極、接收極鞘層電勢的躍變，φE,C分別為發射極、接收極功函數，Ud、U分別為電弧壓降及輸出電壓。對于發射極鞘層，其分布在飽和區是單調下降的(圖2)；而在阻塞區則是先增后降(圖3)，形成了電子發射的勢壘ΔVE。另一方面，對于接收極鞘層，其電勢分布并不依賴于具體的工作區域，其正負性取決于邊界條件本身的匹配。以下依據電極鞘層內電子勢能的分類對邊界條件進行討論，并將接收極表面的電勢定義為零電勢參考點。

圖3 電弧工況阻塞區可能的電子勢能分布Fig.3 Possible distribution of electron potential energy at obstructed region of ignited mode

1) 發射極鞘層的邊界條件

電弧工況發射極鞘層的邊界條件需區分為飽和區與阻塞區兩種情況進行討論。

在飽和區，發射極附近的電子勢能是單調下降且躍變較大，由此產生顯著的附加電場促進發射，即產生肖脫基效應，從而獲得較大的輸出電流密度[12]。邊界通量需考慮帶電粒子的熱發射以及來自于等離子體的反射。

jeE=(1-r1)jSE(EE)-

式中：jeE、jiE、SeE為發射極鞘層與等離子體交界面上的通量；TeE為交界面上的電子溫度；nE為交界面上的等離子體濃度；jSE(EE)為熱電子的肖脫基發射項；r1、r2為考慮了發射的動理學反射特性而引入的參數[13]。

在阻塞區，發射極鞘層的電子勢能是先增后降的，即熱電子發射的勢壘增大了ΔVE，形成虛擬陰極，由于未促進發射附加電場，肖脫基效應也就不存在[14]，此時的發射電流密度是飽和發射電流密度jS0。

2kTeE+je(Ud+ΔVE)

2) 接收極鞘層的邊界條件

當接收極的電勢躍變VC>0時(對應于圖2a與圖3a)，有：

式中：jeC、jiC、SeC為接收極鞘層與等離子體交界面上的通量；f0與f1為描述接收極鞘層內電子運動的各向異性而引入的函數[15]。

當接收極的電勢躍變VC<0時(對應于圖2b與圖3b)，有：

2 模型求解與算法流程

根據上述建立的輸運方程和采用的邊界條件，輸運系數依賴于等離子體參量本身，因此輸運方程本身是非線性的，故本文建立的模型只能數值求解。另一方面，輸運的邊界條件依賴于工作參數所處的工況區域(飽和區或阻塞區)，但實際計算中并不能根據輸入條件(工作參數)預知其工況區域。因此，在數值求解的過程中，應使得邊界條件與輸運方程自適應匹配求解。本文采用的算法流程如圖4所示。

1) 輸入計算所必要的參數：基本物理參量、熱離子轉換器的運行參數，如電極溫度(TE、TC)、銫溫TCs、電極的吸銫功函數(φE、φC)、電極間隙d、總電流密度(輸出電流密度)jtot等。

圖4 熱離子轉換器電弧工況等離子體輸運算法流程Fig.4 Algorithm flow of plasma transportation of TEC at ignited mode

2) 給定接收極鞘層與準中性等離子體區交界處的電子溫度TeC和接收極鞘層的電勢躍變VC作為迭代初始值，根據VC的正負性選擇相應的接收極邊界條件模型，計算此邊界的jeC、jiC、nC及SeC、SiC。

3) 以流程2計算得到的參量為邊界值，采用龍格-庫塔方法對輸運微分方程進行積分，得到準中性等離子體區Te、n、je、ji、V、Se、Si的分布，特別地，記發射極邊界上的電子電流密度積分為je_inte(TeC,VC)、離子電流密度積分為ji_inte(TeC,VC)。

4) 先假定工作點處于阻塞區，根據阻塞區發射極邊界條件計算非單調的電勢躍變情況VE及ΔVE。若解存在，則認定虛擬陰極的存在(即處于阻塞區)；若解不存在，則認定處于飽和區(即虛擬陰極不存在，發射極電勢躍變是單調的)，并據此計算相應的電勢躍變VE。

5) 根據邊界條件的自適應判斷，計算相應工況下發射極邊界上的電子電流密度je_cal(TeC,VC)、離子電流密度ji_cal(TeC,VC)，并與積分所得邊界值聯立成方程組je_cal(TeC,VC)-je_inte(TeC,VC)=0與ji_cal(TeC,VC)-ji_inte(TeC,VC)=0。采用牛頓迭代法求解方程組，直至收斂到設定的精度。若發散，則重新選擇初始條件進行計算。

6) 輸出計算結果。

3 結果討論

由于程序的輸入包含電極溫度、電極功函數、銫溫、電極間隙、總電流密度(輸出電流密度)等7個參數，在結果討論中不適宜采用控制變量的方式。Bullis和Rasor[16]在其綜述中指出，盡管熱離子轉換器的工作參數有差異，其輸出特征可采用離子裕度β及參量pCsd來區分。其中，離子裕度β由發射極的發射特性決定，用來表征熱發射離子中和熱發射電子的程度，其臨界值為1；參量pCsd由輸運過程決定，用來表征電子輸運過程的散射強度，其臨界值為20 mile·Torr(1 Torr≈133.322 Pa，1 mil=25.4 μm。

為在較寬的參數范圍內驗證計算模型的準確性，受限于電極材料，本文研究的熱離子轉換器是電極對在給定的電極溫度下運行，僅在可變的銫溫下進行輸出，此時離子裕度β被限定而參量pCsd可在較寬的范圍取值。為便于與文獻的數據進行比較[17]，程序計算中輸入的工作參數列于表1，選定的電極對材料為Re-Mo，限定電極的工作溫度TE=1 700 K、TC=925 K，電極間隙d=0.25 mm。根據表1，由于電極特性的制約，離子裕度β限制在小于1的區域，而參量pCsd則覆蓋了較寬的范圍。將上述工作參數作為程序的輸入，獲得本文所有的計算結果。

表1 程序計算中輸入的工作參數Table 1 Input working parameter for procedure

1) 計算結果的有效性判定

首先，為判定計算結果的合理性，需將等離子體的德拜長度rD(表2)與特征尺寸(電極間隙)進行比較，德拜長度均遠小于電極間隙，說明電極間隙形成了等離子體，所計算的工況都處于電弧工況；并且隨著總電流密度的增大，工況狀態逐漸從電弧工況的阻塞區進入飽和區，在飽和區發射極鞘層電場值滿足EE遠大于0。因此，計算結果驗證了邊界條件分區設置的合理性。

表2 計算結果有效性判定所需的數據Table 2 Datum needed to judge validity of computing result

2) 伏安特性

圖5示出本文計算得到的伏安特性及其與文獻值的對比，計算值與實驗值符合較好。圖6比較了TCs=551 K與591 K時的伏安特性，計算表明，當參量pCsd<20 mile·Torr時，飽和區的輸出電流密度表現出較明顯的飽和特征；而當參量pCsd>20 mile·Torr時，飽和電流密度的特征不明顯。

圖5 本文計算的伏安特性曲線與文獻值的對比Fig.5 Comparison of U-I characteristic between this paper and lecture results

圖6 本文計算的不同銫溫時伏安特性的比較Fig.6 Computed U-I characteristic in this paper at different temperatures of cesium

3) 電子勢能分布

當熱離子轉換器的銫溫為551 K時，電極間隙電子勢能的變化如圖7所示，涵蓋了電子從發射極費米面熱發射輸運到接收極費米面的全過程。整體而言，電子勢能在發射極鞘層的躍變大于接收極鞘層的躍變，而在準中性等離子體區則相對平緩(其中躍變發生的尺度僅有德拜長度的量級，為便于展示，圖片對躍變的尺度進行放大)。

對于特定的工作點(例如jtot=6 A/cm2)，電子始末位置的能級差，即兩電極費米面的能級差，構成了對外輸出電壓。兩電極表面勢能的差用于維持電弧狀態，屬于固有的內部壓降。

圖7 TCs=551 K時電子勢能隨電極間位置的變化Fig.7 Dependence of electron potential energy on location of interelectrode at TCs=551 K

對于發射極鞘層的電勢躍變，當電流密度較小時(jtot=1 A/cm2)，發射極鞘層的電子勢能是先增后減的非單調躍變，在發射極的鞘層內形成了附加的發射勢壘(即所謂的虛擬陰極)，相當于增加了發射極功函數，從而抑制了電子發射，但發射極功函數的提升有利于獲得更大的電壓輸出，此時轉換器工作在阻塞區。隨著電流密度的增大，發射極鞘層的電子勢能變為單調向下躍變，其量值隨電流密度而增大，此勢能躍變使電子熱發射的過程中得到電場的激勵(即肖脫基效應)，同時該電場較好地阻礙等離子體中的電子返回發射極；為維持大電流密度，就必然需要更高的電弧壓降(內壓增大)，從而輸出電壓也就降低了，此時轉換器工作在飽和區。

圖8 TCs=551 K時發射極鞘層電子電流密度的構成Fig.8 Consistency of electron current density in sheath of emitter at TCs=551 K

整體而言，在輸出電流密度遞增的過程中，發射極鞘層的電子勢能分布經歷了由非單調變成單調的過程。對發射極鞘層電子電流密度的構成分析可進一步說明產生此分布特征的機理，如圖8所示，在小電流密度的阻塞區，從等離子體返回的電子電流密度超過凈電子電流密度，相當于發射極鞘層“累積”了電子，這將增加電子發射的勢壘，形成了虛擬陰極；在大電流密度的飽和區，熱發射的電子電流密度遠大于從等離子體返回的電子電流密度，此時在發射極附近產生強烈的電離(圖9a)，發射極附近產生大量離子，形成了促進電子發射的單調向下的電子勢能分布。

對于接收極鞘層的電勢躍變，在輸出電流密度遞增的過程中，其電子勢能分布經歷了從向下躍變(VC<0)到向上躍變(VC>0)的過程。產生這種分布的原因是在輸出電流密度從0開始增大的過程中，實驗觀察到最先發生電離的區域是接收極的鞘層[18]，由此從擴散工況進入飽和工況的阻塞區，因此需要向下躍變的電子勢能分布促使接收極附近發生電離(圖9b)。而隨著輸出電流密度的增大，電離的峰值逐漸向發射極移動，接收極一側趨向于電離平衡(圖9a)，此時接收極鞘層表現為準中性等離子體鞘層的性質，即熱速率較大的電子運動到接收極表面逸走而留下富余的正離子，從而使電子勢能向上躍變。

以上關于電極間隙電子勢能的分布特性，在轉換器的銫溫為591 K時也同樣存在(圖10)。

圖9 TCs=551 K時電離速率隨電極間位置的變化Fig.9 Dependence of rate of ionization on location of interelectrode at TCs=551 K

圖10 TCs=591 K時電弧工況不同工作點下對應的電子勢能分布Fig.10 Dependence of electron potential energy at different working points of ignited mode at TCs=591 K

4) 根據輸出功率特性對電子勢能分布特性進行確認

為確認上述的熱離子轉換器電極間隙電子勢能分布特性，進行最大輸出功率的討論。一般而言，當電池的內外壓降相等時有最大的電功率輸出，這是電池的一基本屬性。根據電子勢能分布，可獲得輸出電壓、電弧壓降與電流密度的關系，如圖11所示，電流密度越大，輸出電壓越小，而電弧壓降則越大，兩者的交點為內外壓相等的工作點。同時，根據計算的伏安特性，可獲得輸出電功率密度與電流密度的關系(圖12)，當轉換器銫溫為551、591 K時，在對應內外壓相等的工作點jtot=5、14 A/cm2下獲得最大的輸出電功率。這說明了本文關于熱離子轉換器的計算結果是自洽的，并與一般電池的功率特性相一致。

圖11 輸出電壓、電弧壓降與輸出電流密度的關系Fig.11 Dependence of output voltage and plasma voltage drop on output current density

圖12 輸出電功率密度與電流密度的關系Fig.12 Dependence of output power density on output current density

4 結論

本文基于動理學方法建立了描述熱離子轉換器電極間隙電弧工況等離子體輸運的方程及相應的邊界條件，并設計了一套基于牛頓迭代法的計算程序實現其自適應解耦求解，詳細討論了等離子體的電子勢能分布特性。

1) 電子勢能在發射極鞘層的躍變大于接收極鞘層的躍變，而在準中性等離子體區則相對平緩。隨著輸出電流密度的增大，發射極鞘層的電子勢能躍變特征將從先增后減的非單調變為單調下降，這是由發射極鞘層內電子電流密度構成的差異導致的。

2) 同時接收極鞘層的電子勢能則從向下躍變變為向上躍變，這是由接收極鞘層附近的電離狀態決定的。

3) 基于等離子體鞘層電子勢能分布的分析，改進了等離子體輸運的邊界條件模型，所計算的伏安特性曲線與文獻值較為吻合，說明本文的改進模型是合理的。