高頻感應(yīng)等離子體流動(dòng)模型與微型電推進(jìn)器流場特性分析

李亦非, 付宸聰, 劉 偉, 蔡國飆, 王偉宗

(北京航空航天大學(xué), 北京 100191)

引 言

“航天探索, 動(dòng)力先行”. 與利用燃料氧化燃燒產(chǎn)生推力的化學(xué)推進(jìn)相比, 電推進(jìn)利用電離推進(jìn)劑加速工質(zhì)形成射流, 具有比沖高、 推力小、 控制精度高、 使用壽命長等優(yōu)點(diǎn), 尤其適用于航天器位置保持、 軌道轉(zhuǎn)移、 軌道維持、 姿態(tài)控制和深空探測等太空任務(wù)[1]. 根據(jù)推進(jìn)劑加速機(jī)制的不同, 可將電推進(jìn)分為3類: 電熱推進(jìn)[2-4]、 電磁推進(jìn)[5-7]和靜電推進(jìn)[8]. 射頻離子推力器基于感性耦合放電技術(shù), 無放電電極, 不存在電子轟擊式推力器的陰極壽命問題, 同時(shí)沒有電子回旋共振式離子推力器需要外部靜磁場和固定推進(jìn)劑流量以達(dá)到電子回旋共振的限制, 具有結(jié)構(gòu)簡單、 工作壽命長、 推力動(dòng)態(tài)范圍大、 性能調(diào)節(jié)響應(yīng)靈敏[9]等特點(diǎn), 是國際電推進(jìn)領(lǐng)域的研究熱點(diǎn)之一.

無論是高推力精度要求的太空任務(wù)(如歐洲航天局的Darwin， GAIA， LISA等), 還是協(xié)同運(yùn)作的立方星星座任務(wù)(如SpaceX的StarLink), 都迫切需要微推進(jìn)系統(tǒng)[10]. 微型射頻離子推力器滿足體積小、 質(zhì)量小、 推力小、 響應(yīng)快速、 控制精確的要求, 除此之外結(jié)構(gòu)簡單、 工作壽命長, 是最具有發(fā)展?jié)摿Φ耐屏ζ髦?

目前, 對微型射頻離子推力器的研究主要集中在歐洲和美國, 典型樣機(jī)有德國吉森大學(xué)的RIT-10[11-12]， RIT-4[10,13], 美國賓夕法尼亞州立大學(xué)的MRIT[14], 美國BUSEK公司的BIT-3. 其中RIT-10于2001年在ARTEMIS通信衛(wèi)星上得到正式應(yīng)用, 用于執(zhí)行南北位置保持任務(wù)[15], BIT-3預(yù)計(jì)2020年之后將搭載太空發(fā)射系統(tǒng)升空[16]. 國內(nèi)研制微型射頻離子推力器起步較晚, 相關(guān)研究較少, 最早是由中國科學(xué)院空間科學(xué)技術(shù)中心對RIT-15采用氙作為工質(zhì)進(jìn)行性能實(shí)驗(yàn)測試[17]. 中科院力學(xué)所的μRIT-2.5已進(jìn)行飛行驗(yàn)證[18].

射頻離子推力器電離室內(nèi)的感性耦合放電等離子體特性和推力器的性能密切相關(guān). 感性耦合放電等離子體數(shù)值仿真可以提供比實(shí)驗(yàn)更加豐富的結(jié)果, 甚至可以揭示在實(shí)驗(yàn)中無法觀察到的現(xiàn)象, 因此得到國內(nèi)外的廣泛關(guān)注. 針對射頻離子推力器內(nèi)部感性耦合等離子體的數(shù)值仿真模型主要有3種: 全局模型[19-20]、 流體模型[21-22]、 粒子模型[23-24]. 流體模型與全局模型相比可以得到等離子體特性參數(shù)的空間分布, 與粒子模型相比節(jié)省計(jì)算時(shí)間與空間消耗, 具有顯著優(yōu)勢.

在射頻離子推力器的工作過程中, 放電室內(nèi)的物理場包括流場、 等離子體濃度場、 射頻線圈施加的電磁場. 各個(gè)物理場之間呈現(xiàn)強(qiáng)耦合關(guān)系: 流場影響中性粒子及放電室氣壓的分布; 電磁場影響等離子體中的電離過程, 改變帶電粒子的濃度分布; 等離子體通過Joule和非Joule加熱效應(yīng), 影響流場分布. 微型射頻離子推力器流體模型考慮了上述多物理場的耦合以獲得等離子體特征參數(shù)分布, 并揭示放電過程中質(zhì)量、 動(dòng)量和能量傳遞的機(jī)制.

文獻(xiàn)中對大尺寸感性耦合放電有較多的流體模型研究[25-26], 本文對依據(jù)縮尺準(zhǔn)則設(shè)計(jì)的放電室內(nèi)徑40 mm的微型射頻離子推力器進(jìn)行分析, 建立了低氣壓、 小尺度下的感性耦合等離子體流體模型, 探究了等離子體放電特征參數(shù)隨推進(jìn)劑工質(zhì)氣壓、 放電吸收功率、 射頻頻率以及線圈匝數(shù)等因素的變化規(guī)律. 結(jié)果表明, 推進(jìn)劑工質(zhì)氣壓、 放電吸收功率是調(diào)節(jié)微型射頻離子推力器性能的主要因素.

1 等離子體流體模型

1.1 流體條件驗(yàn)證

流體模型將等離子體中各種粒子看作連續(xù)流體, 用質(zhì)量、 動(dòng)量、 能量方程描述等離子體行為. 在使用流體模型分析推力器放電過程之前, 首先須明確在所設(shè)計(jì)推力器工況下流體模型是否適用.

當(dāng)工質(zhì)氣壓在在13.3 Pa(100 mtorr)以上時(shí), 等離子體內(nèi)碰撞占據(jù)主導(dǎo), 可以將帶電粒子看作連續(xù)流體. 但微型射頻離子推力器的工質(zhì)氣壓一般為0.1～1 Pa, 故須預(yù)先分析流體模型在此情況下是否適用. Mikellides等[27]的研究表明, 當(dāng)粒子的平均自由程遠(yuǎn)小于放電室的特征長度時(shí), 碰撞會(huì)占據(jù)主導(dǎo), 具體表述見下式

式中,Kn是粒子的Knudsen數(shù),λmfp是粒子的平均自由程,L是放電室的特征長度. 射頻離子推力器內(nèi)部電子Knudsen數(shù)較大, 但電子由于會(huì)被鞘層反射, 碰撞機(jī)會(huì)增大, 故電子Knudsen數(shù)不具代表性. 考慮離子與中性粒子的碰撞, 其平均自由程如下式所示[22].

式中,nn是放電室內(nèi)部中性粒子數(shù)密度,Qin是離子-中性粒子散射碰撞截面, 單位為m2[28].

式中,cr是由約化質(zhì)量定義的單電荷離子與中性粒子的相對熱運(yùn)動(dòng)速度.

式中,kB是Boltzman常數(shù),Ti是離子溫度,mi是離子質(zhì)量. 中性粒子數(shù)密度通過中性流的Maxwell平均通量方程計(jì)算.

由雙柵系統(tǒng)近似公式可得到有效透明度為

由上述公式, 可得所設(shè)計(jì)微型射頻離子推力器放電室內(nèi)離子Knudsen數(shù)隨工質(zhì)流量的變化規(guī)律. 同時(shí), 根據(jù)理想氣體狀態(tài)方程可得放電室內(nèi)部壓強(qiáng)隨工質(zhì)流量的變化規(guī)律, 如圖1所示.

圖1 流體條件驗(yàn)證Fig. 1 Fluid condition verification

由圖1可知, 在推力器工質(zhì)氣壓為0.133 Pa時(shí), 工質(zhì)流量約為0.6 mg/s, 對應(yīng)Kn為0.37, 小于0.4. 因此可以推斷, 當(dāng)工質(zhì)氣壓在0.133 Pa以上時(shí)， 采用流體模型是合理的.

1.2 幾何構(gòu)建

圖2為按照所設(shè)計(jì)微型射頻離子推力器建立的二維軸對稱簡化模型的幾何構(gòu)造, 考慮了放電室入口與出口間工質(zhì)流動(dòng). 軸對稱處理會(huì)導(dǎo)致現(xiàn)實(shí)中的螺旋型線圈變?yōu)檩S對稱的環(huán)形, Tsay[22]用實(shí)驗(yàn)驗(yàn)證了采用兩種不同構(gòu)型線圈得到結(jié)果的最大誤差不超過5%, 證明了二維軸對稱模型的合理性.

由圖2可知, 模型包括等離子體、 石英放電室、 銅線圈和真空4個(gè)部分, 物理性質(zhì)見表1, 其中等離子體的物理性質(zhì)由等離子體模塊計(jì)算得出.

圖2 二維軸對稱模型Fig. 2 Two-dimensional axisymmetric model

表1 各部分材料的物理性質(zhì)

1.3 氙氣流場

中性粒子的數(shù)密度以及壓強(qiáng)、 速度分布對電離影響很大, 故須對氙氣流場進(jìn)行計(jì)算. 對于射頻離子推力器, 電離度在10-3的量級, 離子對流場的影響可以忽略. 采用動(dòng)量方程和連續(xù)性方程對氙氣流場進(jìn)行描述.

式中,ρ是流體密度,t是時(shí)間,u是流體速度矢量,p是流體壓強(qiáng),I是單位矩陣,K是黏性應(yīng)力張量,F是力矢量.

壁面邊界條件為無滑移, 入口設(shè)置為流量邊界條件, 出口設(shè)置為壓力邊界條件.

研究表明, 背景氣體溫度對等離子體特征參數(shù)影響不大[29]. 故微型射頻離子推力器流體模型不考慮等離子體的加熱效應(yīng), 背景氣體溫度設(shè)置為室溫300 K.

1.4 等離子體

等離子體模塊通過導(dǎo)入碰撞截面數(shù)據(jù), 在時(shí)域計(jì)算電子碰撞反應(yīng)以及重物質(zhì)反應(yīng), 最終得到等離子體特性參數(shù)分布. 選取化學(xué)反應(yīng)體系見表2, 其中，ΔΣ表示反應(yīng)過程中的能量損耗，kf表示反應(yīng)速率， 1～5的反應(yīng)速率由其碰撞截面數(shù)據(jù)計(jì)算. 碰撞截面數(shù)據(jù)來自文獻(xiàn)[31].

表2 模型中考慮的反應(yīng)

流體模型將粒子看作連續(xù)流體.電子濃度方程與能量方程為

式中,μe是電子遷移率,De是電子擴(kuò)散系數(shù),με是電子能遷移率,Den是電子能擴(kuò)散系數(shù), 在電子能量分布函數(shù)滿足Maxwell分布的情況下, 通過Einstein關(guān)系計(jì)算[31].

式中,νm是彈性碰撞頻率. 用漂移-擴(kuò)散方程計(jì)算重物質(zhì)輸運(yùn), 適用于除電子之外的所有粒子.

式中,ωk是第k種物質(zhì)的質(zhì)量分?jǐn)?shù),jk是第k種物質(zhì)的擴(kuò)散通量矢量,Rk是第k種物質(zhì)的產(chǎn)率系數(shù). 采用Poisson方程計(jì)算等離子體的電勢V.

式中,ρq是等離子體的體電荷密度. 邊界條件設(shè)置參考對感性耦合等離子體放電的研究[29].初始電子密度設(shè)置為為1017m-3, 初始平均電子能量為 5 eV. 在壁面處不考慮二次電子發(fā)射, 壁面邊界條件如下[29]

式中,n表示表面的法向量,veth是電子的熱運(yùn)動(dòng)速度. 對靜電邊界條件, 工質(zhì)入口與出口處設(shè)置為零電荷與絕緣, 壁面處設(shè)置為接地.

1.5 電磁場

電磁場模塊通過在頻域求解Maxwell方程組, 得到磁場以及感應(yīng)電場的變化規(guī)律.

J=σE+jωD+σv×B+Je

將上述4式化簡合并后, 可以得到

式中, j是虛數(shù)單位,ω是射頻源的角頻率,σ是等離子體電導(dǎo)率,ε0是真空介電常數(shù),A是磁矢勢,μ0是真空磁導(dǎo)率,Je是應(yīng)用外部電流. 由上式求解磁矢勢, 根據(jù)磁矢勢計(jì)算出磁場與感應(yīng)電場.

1.6 多物理場耦合

等離子體模塊與磁場模塊通過等離子體電導(dǎo)與變壓器模型形成耦合關(guān)系. 首先根據(jù)反應(yīng)進(jìn)度實(shí)時(shí)計(jì)算出等離子體的電導(dǎo)率

式中,q是元電荷,me是電子質(zhì)量,νe是電子碰撞頻率. 其次引入變壓器模型, 該模型將線圈看作變壓器的初級, 等離子體看作次級.[32]

圖3 變壓器模型Fig. 3 Transformer model

通過等離子體內(nèi)的感應(yīng)電流密度以及感應(yīng)電場來計(jì)算等離子體吸收的功率密度如下

式中,J是等離子體內(nèi)的感應(yīng)電流密度,E是等離子體內(nèi)的感應(yīng)電場強(qiáng)度.

2 結(jié)果與討論

2.1 特定工況計(jì)算結(jié)果

利用構(gòu)建的微型射頻離子推力器流體模型, 對工質(zhì)流量為0.33 mg/s、 放電室壓強(qiáng)為0.8 Pa、 吸收功率為20 W、 線圈射頻頻率為2.5 MHz工況下的射頻放電等離子體進(jìn)行了數(shù)值仿真.

圖4為達(dá)到穩(wěn)態(tài)時(shí)推力器內(nèi)部流場. 由圖可知, 速度與壓強(qiáng)的最大值均發(fā)生在工質(zhì)入口處, 這是因?yàn)槿肟谔幍牧魍娣e很小, 氙氣工質(zhì)進(jìn)入時(shí)產(chǎn)生堆積形成了較大的速度、 壓強(qiáng)峰值. 同理, 在工質(zhì)出口可以看到類似的速度、 壓強(qiáng)極值. 由圖4(b)可知, 放電室內(nèi)的工質(zhì)氣壓均勻保持在0.8 Pa, 為等離子體的均勻分布提供了條件.

(a) Velocity field

(b) Pressure field圖4 流場計(jì)算結(jié)果Fig. 4 Flow field calculation results

圖5為等離子體電子密度與電子溫度分布. 由圖5(a)可知, 電子密度在軸線處有最大值, 在壁面處有最小值, 這是因?yàn)楸诿嫣庪x子與電子通過表面反應(yīng)再結(jié)合產(chǎn)生中性粒子, 造成離子電子損失. 由圖5(b)可知, 電子溫度在線圈附近有最大值, 在軸線處有非零最小值, 這是因?yàn)橼吥w效應(yīng)導(dǎo)致功率大多沉積在線圈附近, 使電子溫度較高.

電勢與感應(yīng)電場的分布如圖6所示. 由圖可知, 電勢的峰值出現(xiàn)在軸線處; 而電場強(qiáng)度模在壁面處達(dá)到最大值, 向軸向逐漸減小, 在軸線處趨近于零.

圖5與圖6中的等離子體特性參數(shù)以及電場分布與文獻(xiàn)[21]結(jié)果有相同趨勢, 證明結(jié)果合理性.

(a) Electron density

(b) Electron temperature圖5 等離子體特性參數(shù)分布Fig. 5 Distribution of plasma characteristic parameters

(a)Potential

(b)Electric field圖6 電場分布Fig. 6 Distribution of electric field

2.2 工質(zhì)氣壓對等離子體特性參數(shù)影響

在吸收功率為20 W, 射頻頻率為2.5 MHz的工況下, 工質(zhì)氣壓對徑向等離子體特征參數(shù)的影響如圖7所示. 由圖7(a)可知, 電子密度隨工質(zhì)氣壓的增大而增加. 這是因?yàn)殡S著工質(zhì)氣壓增大, 中性粒子數(shù)密度會(huì)增加, 電子與中性粒子碰撞頻率增加, 加快了電離反應(yīng)的速率, 從而導(dǎo)致離子密度的上升, 由于等離子體的準(zhǔn)中性, 電子密度也隨之增加. 由圖7(b)可知, 電子溫度隨工質(zhì)氣壓的減小而升高, 這是因?yàn)楣べ|(zhì)氣壓減小后, 高能電子與中性粒子以及離子之間的碰撞頻率降低, 損失能量減少, 表現(xiàn)為電子溫度上升.

(a) Electron density

(b) Electron temperature圖7 工質(zhì)氣壓對等離子體特性參數(shù)影響Fig. 7 Influence of propellent pressure on plasma characteristic parameters

2.3 吸收功率對等離子體特性參數(shù)影響

在工質(zhì)氣壓為1.33 Pa, 射頻頻率為2.5 MHz的工況下, 吸收功率對徑向等離子體特征參數(shù)的影響如圖8所示. 由圖8(a)可知, 電子密度隨吸收功率的增大而增加, 這是因?yàn)橛懈嗟碾姽β视捎陔婋x生成等離子體. 相比之下, 電子溫度則表現(xiàn)出功率無關(guān)性, 幾乎不隨吸收功率的增大而改變, 如圖8(b)所示.

2.4 射頻頻率對等離子體特性參數(shù)影響

在吸收功率為20 W, 工質(zhì)氣壓為1.33 Pa時(shí), 選取2.5， 13.56， 50， 100， 200 MHz等5個(gè)常見頻率, 探究射頻頻率對等離子體特征參數(shù)的影響, 結(jié)果如圖9所示. 可以看出電子密度隨射頻頻率的增加略有下降, 與文獻(xiàn)[23]中粒子模型計(jì)算結(jié)果趨勢一致. 電子溫度隨頻率增大幾乎不變. 柵極附近的等離子體密度直接決定了柵極系統(tǒng)引出的束流大小, 進(jìn)而影響推力器的推力, 為使柵極附近鞘層等離子體密度盡可能大, 設(shè)計(jì)時(shí)選擇較小的射頻頻率. 由于射頻頻率對等離子體特性參數(shù)的影響很小, 且在推力器工作過程中改變其射頻頻率較為困難, 射頻頻率一般不作為調(diào)節(jié)推力性能的參數(shù).

2.5 線圈匝數(shù)對等離子體特性參數(shù)影響

線圈匝數(shù)通過影響電磁場進(jìn)而影響等離子體特征參數(shù)分布. 在吸收功率為20 W, 射頻頻率為2.5 MHz, 工質(zhì)氣壓為1.33 Pa時(shí), 選取線圈匝數(shù)為3, 5, 7, 9匝, 探究線圈匝數(shù)對等離子體特征參數(shù)的影響. 圖10為不同線圈匝數(shù)下磁場分布. 由圖可知, 隨著線圈匝數(shù)的增加, 磁通密度更加緊密地聚集在放電室內(nèi)部靠近壁面的位置, 但放電室內(nèi)磁場分布幾乎不變. 表3為不同匝數(shù)下等離子體特征參數(shù), 其中，ne代表電子密度，Te代表電子溫度. 結(jié)果表明, 線圈匝數(shù)對電離室內(nèi)平均電子密度以及平均電子溫度幾乎沒有影響.

(a) Electron density

(b) Electron temperature圖8 吸收功率對等離子體特性參數(shù)影響Fig. 8 Influence of absorption power on plasma characteristic parameters

圖9 射頻頻率對等離子體特性參數(shù)影響Fig. 9 Influence of RF frequency on plasma characteristic parameters

(a) Number of coils n=3

(b) Number of coils n=5

(c) Number of coils n=7

(d) Number of coils n=9

表3 線圈匝數(shù)對等離子體特性參數(shù)影響

3 結(jié)論

本文建立了微型射頻離子推力器感性耦合放電等離子體流體模型, 對特定工況下的推力器等離子體特性參數(shù)分布進(jìn)行了仿真研究, 探究了等離子體放電特征參數(shù)隨推進(jìn)劑工質(zhì)氣壓、 放電吸收功率、 射頻頻率以及線圈匝數(shù)等因素的變化規(guī)律. 結(jié)論如下:

(1)對于所設(shè)計(jì)內(nèi)徑為40 mm的微型射頻離子推力器, 在工質(zhì)氣壓為0.133 Pa以上時(shí), 可以采用流體模型進(jìn)行仿真分析.

(2)在吸收功率、 射頻頻率和線圈匝數(shù)保持不變時(shí), 工質(zhì)氣壓增大會(huì)使電子密度增加, 電子溫度降低.

(3)在工質(zhì)氣壓、 射頻頻率和線圈匝數(shù)保持不變時(shí), 吸收功率增大會(huì)使電子密度增加, 電子溫度幾乎不變.

(4)在其他參數(shù)相同時(shí), 射頻頻率和線圈匝數(shù)對等離子體特征參數(shù)的影響可以忽略.

(5)由上述結(jié)論可以看出, 推進(jìn)劑工質(zhì)氣壓、 放電吸收功率是調(diào)節(jié)微型射頻離子推力器性能的主要因素, 該研究為綜合調(diào)控微型射頻離子推力器的工作性能奠定了良好的基礎(chǔ).

致謝本研究獲得國家自然科學(xué)基金(51977110)、 中央高校基本科研業(yè)務(wù)費(fèi)(NT2020007)資助項(xiàng)目資助.