• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    離子推力器空心陰極放電模型研究進(jìn)展

    2016-02-08 19:24:41谷增杰賈艷輝
    真空與低溫 2016年6期
    關(guān)鍵詞:方法模型

    谷增杰,郭 寧,賈艷輝

    (蘭州空間技術(shù)物理研究所 真空技術(shù)與物理重點(diǎn)實(shí)驗(yàn)室,蘭州 730000)

    離子推力器空心陰極放電模型研究進(jìn)展

    谷增杰,郭 寧,賈艷輝

    (蘭州空間技術(shù)物理研究所 真空技術(shù)與物理重點(diǎn)實(shí)驗(yàn)室,蘭州 730000)

    為了解空心陰極放電物理過(guò)程的各種建模方法,調(diào)研了流體模型、PIC模型、混合模型和降級(jí)模型等多種建模方法的研究現(xiàn)狀,分析了各種模型的特點(diǎn)、應(yīng)用范圍和實(shí)施要點(diǎn),為促進(jìn)建模仿真作為解決空心陰極研制、試驗(yàn)過(guò)程中遇到的具體工程問(wèn)題的技術(shù)途徑,選擇有效、合理的建模方法提供依據(jù)和參考。

    離子推力器;空心陰極;物理模型;建模方法;

    0 引言

    空心陰極是一種依靠氣體放電而工作的真空電子器件,具有發(fā)射電子束流密度大、功耗低、可靠性高、壽命長(zhǎng)等優(yōu)點(diǎn)。作為一種高效率的宇航級(jí)電子源,空心陰極既是離子推力器和霍爾推力器的核心器件[1-2],又在空間站主動(dòng)電位控制系統(tǒng)[3-4]、空心陰極推力器[5]、電動(dòng)繩系等離子體接觸器[6]等領(lǐng)域得到廣泛應(yīng)用。并且作為一種工業(yè)等離子體源,空心陰極在相關(guān)領(lǐng)域有廣闊的應(yīng)用前景。

    空間用空心陰極的結(jié)構(gòu)主要由陰極管、發(fā)射體、加熱器、觸持極和絕緣器等幾部分組成,受空心陰極結(jié)構(gòu)尺寸小、工作溫度高以及試驗(yàn)真空潔凈度要求高等因素限制,采用相關(guān)等離子體診斷技術(shù)研究空心陰極工作物理過(guò)程存在困難,有必要通過(guò)建模仿真,結(jié)合試驗(yàn)數(shù)據(jù)研究空心陰極放電機(jī)理。

    通過(guò)調(diào)研國(guó)內(nèi)外空心陰極放電模型,對(duì)不同模型的特點(diǎn)、應(yīng)用范圍和實(shí)施要點(diǎn)進(jìn)行分析,有針對(duì)性的解決空心陰極研制過(guò)程中的特定問(wèn)題,選擇合適的建模方法以及開展建模仿真具有鋪墊作用。

    1 國(guó)外研究現(xiàn)狀

    1.1 流體方法

    基于流體理論的等離子體建模方法,是通過(guò)等離子體密度、流速、溫度等局部平均量的演化過(guò)程進(jìn)行數(shù)值研究的,即通過(guò)數(shù)值求解每種粒子的概率密度函數(shù)隨時(shí)間演化的玻爾茲曼方程,得到等離子體參數(shù)空間分布特性。

    美國(guó)NASA/JPL針對(duì)空心陰極放電產(chǎn)生的低溫等離子體開發(fā)了兩種流體模型:IROrCa2D模型為發(fā)射體區(qū)二維軸對(duì)稱時(shí)變流體模型;OrCa2D模型為全局的二維軸對(duì)稱非時(shí)變流體模型,計(jì)算區(qū)域?yàn)榘l(fā)射體區(qū)、孔腔區(qū)和羽流區(qū)。2005年之后經(jīng)過(guò)多次模型修正,仿真結(jié)果與實(shí)測(cè)結(jié)果吻合度逐步提升。

    1.1.1 IROrCa2D模型

    2004年,Mikellides等[7-9]開發(fā)IROrCa2D模型的最初目的是為了研究轟擊發(fā)射體的離子電流密度,并分析這種轟擊效應(yīng)是否會(huì)引起發(fā)射體性能降級(jí)。模型利用質(zhì)量守恒、動(dòng)量守恒和能量守恒關(guān)系建立控制方程,邊界條件包括發(fā)射體內(nèi)表面、陰極頂內(nèi)表面、發(fā)射體上游陰極管內(nèi)表面、節(jié)流孔入口以及工質(zhì)注入邊界,采用有限差分法離散守恒方程,迭代求解直至等離子體密度、電子溫度等參數(shù)達(dá)到穩(wěn)態(tài)。該模型能夠得到等離子體參數(shù)空間分布特性以及離子和電子的流場(chǎng)。與試驗(yàn)結(jié)果對(duì)比顯示,模型所得節(jié)流孔邊界處電子溫度1.8 eV低于實(shí)測(cè)值3.4 eV,與節(jié)流孔附近區(qū)域由雙流不穩(wěn)定性引起的異常加熱機(jī)制以及電荷碰撞交換過(guò)程有關(guān)。

    1.1.2 OrCa2D模型

    2005年,Mikellides等[10]首次公布了OrCa2D模型的細(xì)節(jié),主要目標(biāo)是評(píng)估觸持極濺射腐蝕和壽命,因此其計(jì)算區(qū)域擴(kuò)展到了節(jié)流孔和羽流區(qū)。該模型是在IROrCa2D模型基礎(chǔ)之上發(fā)展而來(lái),控制方程考慮了中性氣體動(dòng)力學(xué),增加了除離子慣性以外的所有與離子聲波和電離碰撞有關(guān)的時(shí)變項(xiàng),邊界條件增加了節(jié)流孔壁面、陽(yáng)極壁面和自由擴(kuò)散邊界,模型算法多采用隱式格式,提升了代碼運(yùn)行效率,減少了有顯示格式帶來(lái)的數(shù)值振蕩。與試驗(yàn)結(jié)果對(duì)比表明,仿真所得羽流區(qū)電阻率對(duì)電場(chǎng)強(qiáng)度的貢獻(xiàn)小于電子壓力梯度作用力,導(dǎo)致羽流區(qū)等離子體電勢(shì)分布與試驗(yàn)所測(cè)軸線方向單調(diào)遞增的結(jié)果不一致,這種差異也與雙流不穩(wěn)定性引起的異常電阻率增大有關(guān)。

    1.1.3 流體模型的修正

    2005年,Mikellides等[11]針對(duì)模型中由雙流不穩(wěn)定性引起的等離子體電阻率增大,進(jìn)而導(dǎo)致空心陰極羽流區(qū)等離子體處于高電勢(shì)狀態(tài)的現(xiàn)象,對(duì)OrCa2D模型進(jìn)行修正,研究反常電阻、中性氣體黏性、觸持極壁面邊界及壁面溫度等因素對(duì)羽流區(qū)等離子體電勢(shì)分布的影響。模型修正前后的仿真結(jié)果以及NEXIS空心陰極等離子體診斷數(shù)據(jù)對(duì)比可知,用離子聲波反常電阻率代替經(jīng)典等離子體電阻率,節(jié)流孔和觸持極附近等離子體電勢(shì)提升約40%,仿真所得電子溫度和等離子體電勢(shì)與試驗(yàn)結(jié)果一致性均得到提高。

    2006年,Mikellides等[12]用IROrCa2D模型仿真NSTAR和NEXIS兩種規(guī)格空心陰極發(fā)射體區(qū)等離子體過(guò)程時(shí)發(fā)現(xiàn),NSTAR空心陰極(陰極管直徑0.635 cm,節(jié)流孔直徑0.1 cm,放電電流12 A,流率0.42 mg/s)節(jié)流孔電流密度是NEXIS空心陰極(陰極管直徑1.5 cm,節(jié)流孔直徑0.3 cm,放電電流25 A,流率0.54 mg/s)的4.5倍,中性氣體密度前者是后者的6倍,小規(guī)格NSTAR空心陰極內(nèi)等離子體密度峰值(1.7 e 21 m-3)比NEXIS(2.5 e 20 m-3)高一個(gè)量級(jí)。較高的等離子體密度導(dǎo)致發(fā)射體表面德拜長(zhǎng)度足夠小,使得發(fā)射體表面鞘層變形凸入發(fā)射體表面的微孔,增大發(fā)射體有效發(fā)射面積,引起發(fā)射體發(fā)射電子增強(qiáng)現(xiàn)象。模型修正前后與試驗(yàn)結(jié)果對(duì)比表明,對(duì)于發(fā)射體區(qū)等離子體密度較大的空心陰極,建模時(shí)須考慮發(fā)射體發(fā)射增強(qiáng)修正系數(shù)。

    2006~2007年,Mikellides等[13-14]對(duì)OrCa2D模型升級(jí),主要在計(jì)算網(wǎng)格、觸持極邊界、中性氣體模型等方面進(jìn)行改進(jìn),從原來(lái)的直線劃分網(wǎng)格改進(jìn)為直線-曲線混合網(wǎng)格,網(wǎng)格數(shù)量減少了60%,計(jì)算精度和速度都有了很大提升。模型增加了觸持極外表面邊界,但是沒有考慮觸持極和陰極頂之間的間隙。仿真對(duì)比分析了全流體的中性氣體模型和混合無(wú)碰撞流體中性氣體模型下中性氣體密度和流場(chǎng)的差別,后者比前者所得節(jié)流孔區(qū)域中性氣體密度高出約4倍。模型還研究了等離子體密度梯度作用下的高電子遷移速度引起的碰撞截面積增大對(duì)等離子體密度的影響,考慮該因素可以更好的擬合等離子體密度參數(shù),但是所得近羽流區(qū)電子溫度比實(shí)測(cè)值低2~5倍。

    2007~2008年,Mikellides等[15-18]針對(duì)NSTAR離子推力器30 352 h長(zhǎng)壽命試驗(yàn)中放電陰極和中和陰極濺射腐蝕問(wèn)題,用建模仿真的方法開展空心陰極濺射腐蝕機(jī)理研究。中和陰極濺射腐蝕主要發(fā)生在節(jié)流孔內(nèi)部,放電陰極濺射腐蝕主要發(fā)生在觸持極頂外表面。采用流體方法建立發(fā)射體區(qū)和孔腔區(qū)的仿真模型,通過(guò)模擬轟擊壁面的離子電流密度和離子能量計(jì)算節(jié)流孔孔壁濺射產(chǎn)額,即可預(yù)測(cè)節(jié)流孔孔壁濺射腐蝕深度。把計(jì)算區(qū)域擴(kuò)展到羽流區(qū),仿真計(jì)算觸持極附近區(qū)域的穩(wěn)態(tài)等離子體電勢(shì),與發(fā)射探針診斷試驗(yàn)所得該區(qū)域的時(shí)變電勢(shì)脈動(dòng)疊加,即為觸持極附近區(qū)域高能離子能量來(lái)源,進(jìn)而計(jì)算由高能離子轟擊引起的觸持極濺射腐蝕速率,預(yù)測(cè)觸持極壽命。

    2009年,Mikellides[19]在上述放電陰極觸持極濺射腐蝕機(jī)理研究的流體模型基礎(chǔ)上,在中性氣體擴(kuò)散方程中加入氣體黏性流動(dòng)項(xiàng),研究氣體黏性對(duì)中性氣體和離子流動(dòng)的影響,結(jié)果表明,中性氣體黏性對(duì)流場(chǎng)影響很大,陰極內(nèi)壓力提高約40%,氣體壓力梯度作用力與離子拖拽作用相當(dāng)甚至超過(guò)后者。2009~2010年,Mikellides等[20-21]用考慮中性氣體黏性的流體模型仿真NEXT離子推力器中和陰極在不同工況、節(jié)流孔直徑和觸持極頂厚度條件下的放電性能,并通過(guò)改變節(jié)流孔邊界的形狀,研究節(jié)流孔孔壁被濺射腐蝕以后空心陰極放電特性的變化。

    2010~2013年,Mikellides等[22]在模型初始化條件和中性氣體流動(dòng)連續(xù)性等方面改進(jìn)了OrCa2D模型,啟動(dòng)仿真程序不再需要預(yù)處理參數(shù),中性氣體從發(fā)射體區(qū)黏滯流到羽流區(qū)的分子流平滑過(guò)渡,解決了先前中性氣體流動(dòng)過(guò)渡區(qū)域的不連續(xù)問(wèn)題。NASA/JPL用改進(jìn)的模型模擬了100 A大電流空心陰極放電等離子體特性,模型預(yù)測(cè)離子聲波形成的反常電阻率和雙流不穩(wěn)定性產(chǎn)生的高能離子,仍然會(huì)是引起空心陰極觸持極濺射腐蝕的重要原因。

    1.2 PIC方法和混合方法

    等離子體粒子建模方法,是通過(guò)直接在相空間中研究等離子體分布函數(shù)的演化,在高速計(jì)算機(jī)上通過(guò)跟蹤大量帶電粒子在其自洽場(chǎng)和外加電磁場(chǎng)中的運(yùn)動(dòng)模擬等離子體動(dòng)力學(xué)特性。PIC方法則在粒子方法基礎(chǔ)上引入權(quán)重粒子,用一個(gè)模擬粒子來(lái)代替分布在相空間小范圍內(nèi)的粒子團(tuán)。

    混合方法嘗試把流體方法和PIC方法的優(yōu)點(diǎn)結(jié)合起來(lái),既能獲得流體模型快的運(yùn)算速度,又能得到PIC模型的計(jì)算準(zhǔn)確度。通?;旌夏P头指畛蓭讉€(gè)子模型,PIC模型用于求解如碰撞頻率和遷移系數(shù)等微觀參數(shù),這些參數(shù)是流體模型所需要的輸入?yún)?shù),流體模型求解粒子密度和電磁場(chǎng),這些參數(shù)又是PIC模型追蹤粒子運(yùn)動(dòng)所需要的輸入條件。

    2001~2004年,Crofton等[24-25]用混合方法建模仿真T5和T6空心陰極羽流,用PIC方法模擬重粒子(氙原子和氙離子)的行為,用流體方法建立連續(xù)方程、動(dòng)量方程和能量方程模擬電子的行為,直接模擬蒙特卡洛方法模擬粒子之間的碰撞。該模型能夠模擬從節(jié)流孔到羽流區(qū)的等離子體電離、碰撞電荷交換、束流中帶電粒子復(fù)合等過(guò)程,研究羽流等離子體參數(shù)對(duì)節(jié)流孔邊界、輸運(yùn)系數(shù)等初始化條件的敏感性[26]。

    2015年,Kubota等[28]探索混合方法在空心陰極發(fā)射體區(qū)、孔腔區(qū)和羽流區(qū)全局化放電模型研究領(lǐng)域的適用性,離子的行為用PIC方法模擬,碰撞電荷交換過(guò)程和彈性碰撞過(guò)程用直接蒙特卡洛模擬,電子的行為用遷移—擴(kuò)散模型模擬,由于中性原子密度顯著高于離子,代表中性原子的模擬粒子的權(quán)重取得比離子大,用中性原子模擬粒子質(zhì)量的減少表示電離過(guò)程。初步研究了發(fā)射體溫度、碰撞電荷交換對(duì)電子溫度、等離子體電勢(shì)、離子和電子流動(dòng)的影響,與美國(guó)NASA/JPL的NSTAR空心陰極羽流區(qū)等離子體診斷試驗(yàn)所得數(shù)據(jù)對(duì)比表明,當(dāng)發(fā)射體峰值溫度為1 900 K時(shí),等離子體密度分布特性基本吻合。

    1.3 降級(jí)方法

    空心陰極降級(jí)模型通常是基于0維的唯像分析方法或1維的流體方法建立的放電模型,多以解析格式的能量平衡、電流平衡、粒子平衡等作為控制方程。一般不需要依賴試驗(yàn)數(shù)據(jù)作為輸入條件,較為容易實(shí)現(xiàn)自洽運(yùn)算。

    1981~1982年,Siegfried等[29-31]基于唯像分析方法建模研究節(jié)流型空心陰極發(fā)射電子電流組分,考慮陰極管內(nèi)表面、發(fā)射體發(fā)射表面、陰極頂內(nèi)表面、節(jié)流孔內(nèi)表面和陰極頂外表面的電流平衡關(guān)系,得到各表面對(duì)空心陰極發(fā)射電子電流的貢獻(xiàn)比重。該模型可以預(yù)測(cè)發(fā)射體區(qū)壓力、放電電流對(duì)有效發(fā)射長(zhǎng)度、發(fā)射體溫度、等離子體密度和電勢(shì)的影響,與試驗(yàn)測(cè)量數(shù)據(jù)吻合度能夠滿足工程需求。

    1992年,Salhi等[32]為了研究空心陰極發(fā)射體區(qū)等離子體過(guò)程,以解析表達(dá)形式的電流、壓力、能量平衡方程為基礎(chǔ),建立了一種0維模型。該模型可以預(yù)測(cè)不同尺寸特性的空心陰極使用氙氣和汞作為工質(zhì)氣體時(shí),發(fā)射體溫度、電子溫度、等離子體密度和電勢(shì)隨放電電流和發(fā)射體區(qū)壓力的變化趨勢(shì)。

    1994年,Mandell等[33]和Katz等[34]以孔腔區(qū)電離平衡和能量平衡關(guān)系為基礎(chǔ),建模研究空心陰極節(jié)流孔內(nèi)等離子體過(guò)程。節(jié)流孔內(nèi)部等離子體過(guò)程很大程度上決定了空心陰極的放電性能,該模型能夠預(yù)測(cè)孔腔區(qū)電子溫度、離子電流密度、等離子體密度和電勢(shì)以及孔腔區(qū)電壓降等參數(shù)隨工質(zhì)流率、放電電流、節(jié)流孔尺寸的變化趨勢(shì)。此外,還把仿真結(jié)果與試驗(yàn)結(jié)果結(jié)合,通過(guò)比較觸持極電壓與等離子體電勢(shì)的關(guān)系,研究了空心陰極的羽狀模式和斑狀模式的工作模式轉(zhuǎn)換過(guò)程。1996年,Katz等[35]用該模型模擬了空心陰極等離子體接觸器孔腔區(qū)等離子體特性,研究工況參數(shù)對(duì)節(jié)流孔下游邊界發(fā)射電子電流的影響規(guī)律,并結(jié)合試驗(yàn)研究真空艙內(nèi)壁狀態(tài)、陽(yáng)極位形對(duì)接觸器伏安特性的影響。

    Capacci等[36]總結(jié)了Siegfried等[30]的唯像分析建模思路,首次用解析格式的薩哈方程、電流平衡、擴(kuò)散平衡、能量平衡和壓力平衡方程建立空心陰極發(fā)射體區(qū)、孔腔區(qū)和觸持極區(qū)比較完整的降級(jí)模型,把微觀等離子體參數(shù)和宏觀放電參數(shù)結(jié)合起來(lái),給定空心陰極特征尺寸、放電電流、工質(zhì)流率,即可仿真計(jì)算放電電壓、發(fā)射體溫度以及三個(gè)區(qū)域平均等離子體參數(shù),但是限于試驗(yàn)測(cè)量結(jié)果有限,仿真所得觸持極電壓能夠與試驗(yàn)結(jié)果大致吻合,未測(cè)量的等離子體參數(shù)無(wú)法與仿真結(jié)果對(duì)比,無(wú)法判定模型中的假設(shè)條件合理性[37-39]。

    2002年,Rossetti等[40]開發(fā)了一種空心陰極發(fā)射體區(qū)的1維模型,基于傳熱學(xué)的發(fā)射體熱平衡方程、等離子體能量平衡方程、中性氣體連續(xù)性方程、電流平衡方程為控制方程,結(jié)合微觀等離子體物理參數(shù),求解發(fā)射體溫度沿軸向的分布特性隨工質(zhì)流率的變化規(guī)律。

    乞求皮特發(fā)慈悲的想法讓我怒從膽邊生,沖動(dòng)之下,我抬腳朝皮特的側(cè)身踢去。事與愿違,他抓住我的腳,向前一拽,我一下失去了平衡,背部著地,仰面朝天狠狠地摔在地上,只好把腳抽回來(lái),掙扎著站起身。

    2002年,Domonkos[41]開發(fā)了一種0維的粒子和能量平衡模型仿真計(jì)算空心陰極發(fā)射體區(qū)和孔腔區(qū)等離子體參數(shù)。兩個(gè)仿真計(jì)算區(qū)域分別以離子密度守恒關(guān)系、電流平衡關(guān)系和能量平衡關(guān)系作為控制方程,除離子擴(kuò)散方程為微分格式外,其他關(guān)系式均為解析格式。模型中把電子與中性原子之間的激發(fā)碰撞系數(shù)、激發(fā)能量作為自由參數(shù),通過(guò)仿真預(yù)測(cè)自由參數(shù)在一定范圍內(nèi)變化時(shí)對(duì)等離子體密度、電子溫度等的影響,分析了空心陰極功耗對(duì)孔腔區(qū)的原子激發(fā)能量、激發(fā)碰撞系數(shù)以及發(fā)射體區(qū)等離子體電勢(shì)的敏感性。

    2002~2003年,Katz等[42-44]為了預(yù)測(cè)空心陰極發(fā)射體壽命,研究節(jié)流孔孔壁濺射腐蝕機(jī)理,建立了一種1維模型模擬孔腔區(qū)和發(fā)射體區(qū)放電過(guò)程,每個(gè)計(jì)算區(qū)域分別用質(zhì)量連續(xù)方程、離子和電子的動(dòng)量方程、擴(kuò)散方程、能量方程和電流連續(xù)方程作為控制方程。模型中,中性氣體流動(dòng)為黏滯流,考慮了離子與中性原子的電荷碰撞交換過(guò)程,電離過(guò)程用碰撞電離方程描述,不再需要求解薩哈方程。模型能夠預(yù)測(cè)發(fā)射體材料的蒸發(fā)、擴(kuò)散和沉積及由孔腔區(qū)離子徑向擴(kuò)散轟擊節(jié)流孔壁面造成的材料濺射。這種模型適用于離子推力器中和陰極和空間等離子體接觸器研究。

    2004年,Katz等[45]為研究空心陰極觸持極濺射腐蝕機(jī)理,探索高能離子的能量來(lái)源,在上述孔腔區(qū)和發(fā)射體區(qū)模型基礎(chǔ)之上,建立了一種1維的羽流區(qū)模型。模型忽略了離子和電子的慣性,用雙極擴(kuò)散方程描述粒子的擴(kuò)散運(yùn)動(dòng),電場(chǎng)強(qiáng)度用歐姆定律求解,控制方程還包括電子能量平衡方程和中性氣體擴(kuò)散方程,用朗繆爾條件確定了羽流區(qū)雙鞘層特性。仿真所得羽流區(qū)等離子體電勢(shì)和電子溫度沿軸向分布特性與試驗(yàn)結(jié)果能夠吻合,該降級(jí)模型預(yù)測(cè)了羽流區(qū)粒子密度梯度對(duì)等離子體電勢(shì)的影響,以及羽流區(qū)雙鞘層的位置和結(jié)構(gòu)特性,但是等離子體的二維運(yùn)動(dòng)和轟擊觸持極的離子的密度和能量還需要更復(fù)雜的二維模型深入研究[46]。

    2013~2015年,Albertoniy等[47]和Pedrini等[48-49]為預(yù)測(cè)空心陰極放電特性開發(fā)了一種涉及發(fā)射體區(qū)和孔腔區(qū)的0維模型,每個(gè)計(jì)算區(qū)域分別以離子平衡、能量平衡、電流平衡和壓力平衡關(guān)系為主要控制方程,并加入空心陰極關(guān)鍵結(jié)構(gòu)的傳熱學(xué)模型,迭代求解等離子體參數(shù)、放電電壓、發(fā)射體溫度等參數(shù)。該模型能夠分析空心陰極功率耗散分布情況,關(guān)鍵結(jié)構(gòu)尺寸和工況參數(shù)對(duì)放電性能的影響,預(yù)測(cè)由材料蒸發(fā)損耗決定的發(fā)射體壽命。與試驗(yàn)結(jié)果對(duì)比表明,空心陰極功耗以及發(fā)射體溫度、放電電壓隨放電電流、工質(zhì)流率的變化趨勢(shì)與試驗(yàn)結(jié)果吻合較好,該簡(jiǎn)化模型可以為空心陰極設(shè)計(jì)參數(shù)選擇提供依據(jù)。2015年,Korkmaz等[50]簡(jiǎn)化了上述Pedrini的模型,孔腔區(qū)離子平衡方程改為徑向擴(kuò)散方程,傳熱學(xué)模型簡(jiǎn)化為節(jié)流孔孔壁能量平衡方程,仿真所得參數(shù)變化趨勢(shì)仍能與試驗(yàn)結(jié)果吻合。

    2 國(guó)內(nèi)研究現(xiàn)狀

    國(guó)內(nèi)開展空心陰極工程研制的單位主要有蘭州空間技術(shù)物理研究所和上??臻g推進(jìn)研究所,采用建模仿真研究空心陰極放電等離子體特性和放電參數(shù)相對(duì)較少,試驗(yàn)研究仍然是空心陰極放電性能研究的主要手段。

    郭寧等[51]和賈艷輝等[52]根據(jù)空心陰極發(fā)射體材料(六硼化鑭)的蒸發(fā)損耗特點(diǎn),采用數(shù)學(xué)建模方法研究發(fā)射體發(fā)射熱電子的電流密度與發(fā)射體材料蒸發(fā)速率的關(guān)系,最終得到空心陰極發(fā)射體壽命與放電電流的關(guān)系,對(duì)空心陰極的壽命進(jìn)行了預(yù)測(cè)[53]。將模型預(yù)測(cè)的LHC-5L空心陰極發(fā)射體損耗量與LIPS-200離子推力器空心陰極壽命試驗(yàn)3 500 h試驗(yàn)測(cè)量值對(duì)比,證明該模型適用于空心陰極發(fā)射體蒸發(fā)損耗壽命預(yù)測(cè)。

    陳梅[54]總結(jié)了Mandell和Katz等的空心陰極降級(jí)模型,分別用離子連續(xù)性方程、電流連續(xù)性方程和能量守恒方程作為發(fā)射體區(qū)和孔腔區(qū)的控制方程,采用牛頓迭代法進(jìn)行計(jì)算,得到電子溫度和密度隨空心陰極節(jié)流孔直徑和流率的變化規(guī)律。馮亮[55]建立了空心陰極發(fā)射體表面雙鞘層的1維數(shù)學(xué)模型,利用連續(xù)性方程、粒子守恒方程、能量守恒方程和泊松方程,結(jié)合初始條件和邊界條件,研究有無(wú)虛陰極兩種情況,發(fā)射體表面附近區(qū)域等離子體和發(fā)射表面之間的相互作用,兩個(gè)模型主要改變發(fā)射電子電流密度和壁面電勢(shì),其他參數(shù)保持不變。王穎[56]根據(jù)空心陰極的自洽熱平衡條件,建立發(fā)射體表面雙鞘層模型以及自持工作和加熱時(shí)的結(jié)構(gòu)模型,利用ANSYS分析了空心陰極的溫度場(chǎng)。安秉健[57]使用Comsol軟件建立全局二維空心陰極放電模型,模型涉及了等離子體放電、流體流動(dòng)、傳熱等多因素耦合物理過(guò)程,并根據(jù)實(shí)驗(yàn)結(jié)果進(jìn)行了修正,最終模擬結(jié)果與試驗(yàn)測(cè)量結(jié)果的吻合較好。

    3 總結(jié)與展望

    國(guó)內(nèi)外已有多種空心陰極建模方法,包括流體方法、PIC方法、混合方法和降級(jí)方法,通過(guò)分析各種模型建模特點(diǎn),結(jié)合空心陰極試驗(yàn)反饋修正模型時(shí)需要考慮的因素,得出結(jié)論:

    (1)二維流體方法是應(yīng)用比較成熟的一種空心陰極建模方法,能夠模擬空心陰極發(fā)射體區(qū)、孔腔區(qū)和羽流區(qū)等離子體參數(shù),結(jié)合等離子體診斷試驗(yàn),還可研究與空心陰極發(fā)射體壽命、觸持極濺射腐蝕等多種工程問(wèn)題。但是,流體模型建模和算法設(shè)計(jì)比較復(fù)雜,程序代碼運(yùn)行時(shí)間較長(zhǎng),由于計(jì)算網(wǎng)格尺寸往往遠(yuǎn)大于德拜長(zhǎng)度,使得仿真結(jié)果不能反應(yīng)金屬壁面鞘層、發(fā)射體區(qū)與孔腔區(qū)之間雙鞘層、羽流區(qū)雙鞘層的特性,而這些鞘層結(jié)構(gòu)對(duì)研究空心陰極放電機(jī)理至關(guān)重要。此外,流體模型需要以發(fā)射體溫度等作為邊界條件,無(wú)法脫離試驗(yàn)結(jié)果獨(dú)立開展仿真研究;

    (2)PIC方法和混合方法在空心陰極放電模型領(lǐng)域應(yīng)用較少。但是,粒子方法建模思路清晰明了,用計(jì)算機(jī)模擬跟蹤單個(gè)粒子的運(yùn)動(dòng),再對(duì)微觀粒子進(jìn)行統(tǒng)計(jì)平均,由此得到宏觀物體的物質(zhì)特性和運(yùn)動(dòng)規(guī)律。從原則上講,粒子方法既考慮了粒子微觀運(yùn)動(dòng),又考慮了等離子體固有的集體效應(yīng),這種方法最能反映等離子體粒子運(yùn)動(dòng)規(guī)律,尤其對(duì)于物理規(guī)律尚不清晰的問(wèn)題,可以用粒子模擬方法初步探索,幫助設(shè)計(jì)實(shí)驗(yàn)并預(yù)示結(jié)果。粒子方法在航天器帶電、離子推力器和霍爾推力器工作性能以及電推進(jìn)系統(tǒng)羽流特性等研究領(lǐng)域得到廣泛應(yīng)用,未來(lái)粒子建模方法也將會(huì)逐步在空心陰極放電等離子體研究中發(fā)揮作用;

    (3)基于唯像分析的空心陰極降級(jí)模型建模方法較為簡(jiǎn)便,運(yùn)算量小,多以分區(qū)域等離子體參數(shù)平均分布假設(shè)為前提,因此不能得到等離子體空間分布特性。但是,這種基于唯像分析的建模方法,能夠自洽完成運(yùn)算,不需要先驗(yàn)數(shù)據(jù)作為輸入條件,以平均等離子體參數(shù)為中介,研究空心陰極關(guān)鍵結(jié)構(gòu)特性、工況參數(shù)對(duì)空心陰極放電性能的影響,有效支撐空心陰極設(shè)計(jì)工作。

    以建模仿真作為研究空心陰極放電特性的技術(shù)途徑,有助于縮短空心陰極研制周期、降低研發(fā)成本,提升仿真設(shè)計(jì)能力和對(duì)空心陰極設(shè)計(jì)結(jié)果的預(yù)示能力。但是,只借助一種建模方法很難解決與空心陰極有關(guān)的所有問(wèn)題。在空心陰極研制階段,關(guān)注的主要問(wèn)題是結(jié)構(gòu)參數(shù)和工況參數(shù)對(duì)空心陰極放電性能、功耗、效率等的影響,基于降級(jí)方法的簡(jiǎn)單模型即可預(yù)示上述參數(shù)的關(guān)系;在空心陰極試驗(yàn)階段,關(guān)注的主要問(wèn)題是空心陰極的零部件壽命、放電穩(wěn)定性、濺射腐蝕速率等問(wèn)題,基于流體方法的二維流體模型能夠滿足上述問(wèn)題的研究需求;在空心陰極應(yīng)用階段,關(guān)注的主要問(wèn)題是空心陰極與工作環(huán)境的匹配性和耦合工作特性等問(wèn)題,基于PIC方法或混合方法的放電模型在空心陰極與離子推力器放電室、空心陰極與霍爾推力器耦合工作研究中是比較理想的方法。

    鑒于目前國(guó)內(nèi)離子推力器和霍爾推力器的技術(shù)發(fā)展對(duì)新型和新規(guī)格空心陰極工程應(yīng)用的迫切性,可從降級(jí)方法的空心陰極放電模型入手,與空心陰極性能試驗(yàn)相結(jié)合,首先滿足空心陰極快速研制的需要;再通過(guò)流體方法或粒子方法詳細(xì)研究空心陰極和工作環(huán)境耦合條件下的放電特性,使得建模仿真逐步達(dá)到滿足空心陰極工程研制所需水平。

    [1]郭寧,顧佐,邱家穩(wěn),等.空心陰極在空間技術(shù)中的應(yīng)用[J].真空,2005,42(5):32-35.

    [2]郭寧,邱家穩(wěn),顧佐,等.20 cm氙離子發(fā)動(dòng)機(jī)空心陰極3 000h壽命試驗(yàn)[J].真空與低溫,2006,12(4):204-206.

    [3]周昊澄,張?zhí)炱?LEO大型載人航天器主動(dòng)電位控制技術(shù)進(jìn)展[J].真空與低溫,2014,20(4):243-247.

    [4]田立成,石紅,李娟,等.航天器表面充電仿真計(jì)算和電位主動(dòng)控制技術(shù)[J].航天器環(huán)境工程,2012,29(2):144-149.

    [5]Angelo Niko Grubi?i?.Microthrusters Based on the T5 and T6 Hollow Cathode[D].Southampton:University of Southamp?ton,2009.

    [6]李安杰.等離子體接觸器電子收集雙鞘層模型數(shù)值解和空心陰極研制[D].北京:中國(guó)科學(xué)院空間科學(xué)與應(yīng)用研究中心,1996.

    [7]MikellidesI,KatzI,GoebelD,etal.Theoreticalmodelofahol?lowcathodeinsertplasma[J].AIAA,2004,3817:11-14.

    [8]GoebelDM,JamesonKK,WatkinsRM,etal.Hollowcathode theory and experiment I Plasma characterization using fast miniature scanning probes[J].Journal of Applied Physics,2005,98(11):113302-113309.

    [9]MikellidesIG,KatzI,GoebelDM,etal.Hollowcathodetheory and experiment II A two-dimensional theoretical model of the emitter region[J].Journal of Applied Physics,2005,98(11):113303.

    [10]Mikellides I G,Katz I,Goebel D M,et al.Theoretical model of a hollow cathode plasma for the assessment of insert and keeper lifetimes[M].Pasadena,CA:Jet Propulsion Laborato?ry,NationalAeronauticsandSpaceAdministration,2005.

    [11]Mikellides I G,Katz I,Goebel D M.Numerical simulation of the hollow cathode discharge plasma dynamics[C]//IEPC Pa?per IEPC-2005-200,29th International Electric Propulsion Conference,Princeton,2005.

    [12]Mikellides I G,Katz I,Goebel D M,et al.Plasma processes inside dispenser hollow cathodes[J].Physics of Plasmas(1994-present),2006,13(6):504-510.

    [13]Mikellides I,Katz I,Jameson K,et al.Driving processes in the orifice and near-plume regions of a hollow cathode[J]. AIAAPaper,2006,5151:9-12.

    [14]MikellidesIG,KatzI,GoebelDM,etal.Evidenceofnonclas?sical plasma transport in hollow cathodes for electric propul?sion[J].Journalofappliedphysics,2007,101(6):301-311.

    [15]Mikellides I G,Katz I,Goebel D,et al.Partially-Ionized Gas and Associated Wear in Electron Sources for Ion Propulsion,II:DischargeHollowCathode[J].AIAAPaper,2007,7:5168.

    [16]Mikellides I G,Katz I,Goebel D,et al.The Partially-Ionized Gas and Associated Wear in Electron Sources for Ion Propul?sion,II:Discharge Hollow Cathode[J].AIAA Paper,2007,7:5192.

    [17]Mikellides I G,Katz I.Wear mechanisms in electron sources for ion propulsion,I:Neutralizer hollow cathode[J].Journal ofPropulsionandPower,2008,24(4):855-865.

    [18]Mikellides I G,Katz I,Goebel D M,et al.Wear mechanisms in electron sources for ion propulsion,II:Discharge hollow cathode[J].Journal of Propulsion and Power,2008,24(4):866-879.

    [19]Mikellides I G.Effects of viscosity in a partially ionized chan?nel flow with thermionic emission[J].Physics of Plasmas(1994-present),2009,16(1):501-511.

    [20]MikellidesIG,SnyderJS,GoebelDM,etal.NeutralizerHol?low Cathode Simulations and Comparisons with Ground Test Data[C]//31st International Electric Propulsion Conference. 2009:20-24.

    [21]Mikellides I G,Goebel D M,Snyder J S,et al.The discharge plasma in ion engine neutralizers:Numerical simulations and comparisons with laboratory data[J].Journal of Applied Physics,2010,108(11):113308.

    [22]MikellidesIG,GoebelDM,JornsBA,etal.NumericalSimu?lations of the Partially Ionized Gas in a 100 A LaB 6 Hollow Cathode[J].IEEE Transactions on Plasma Science,2015,43(1):173-184.

    [23]JornsBA,MikellidesIG,GoebelDM.Temporalfluctuations in a 100 A LaB6 hollow cathode[C]//Proceedings of the 33rd International Electric Propulsion Conference,Washington,DC,2013:385.

    [24]CroftonMW,BoydID.Plumemeasurementandmodelingre?sults for a hollow cathode micro-thruster[R].AIAA Paper,2001.

    [25]Crofton M,Boyd I.The T6 Hollow Cathode:Measurements and Modeling[C]//34 th AIAA Plasmadynamics and Lasers Conference,2003.

    [26]Boyd I D,Crofton M W.Modeling the plasma plume of a hol?low cathode[J].Journal of Applied Physics,2004,95(7):3285-3296.

    [27]Levko D,Krasik Y E,Vekselman V,et al.Two-dimensional model of orificed micro-hollow cathode discharge for space application[J].PhysicsofPlasmas(1994-present),2013,20(8):083512.

    [28]Kubota K,Oshio Y,Watanabe H,et al.Hybrid-PIC Simula?tiononPlasmaFlowofHollowCathode[C]//IEPC,2015.

    [29]Siegfried D E,Wilbur P J.Phenomenological model describ?ing orificed,hollow cathode operation[J].AIAA Journal,1983,21(1):5-6.

    [30]Siegfried D E,Wilbur P J.A model for mercury orificed hol?low cathodes-Theory and experiment[J].AIAA journal,1984,22(10):1405-1412.

    [31]Siegfried D E.Phenomenological Model for Orificed Hollow Cathodes[D].USA:ColoradoStateUniversity,1983.

    [32]Salhi A,Turchi P J.A first-principles model for orificed hol?low cathode operation[C]//Joint Propulsion Conference and Exhibit,1992.

    [33]Mandell M J,Katz I.Theory of hollow operation in spot andplume modes[C]//AIAA,ASME,SAE,and ASEE,Joint Pro?pulsionConferenceandExhibit,30th,Indianapolis,IN,1994.

    [34]Katz I,Mandell M J,Patterson M,et al.Sensitivity of IHollow Cathode Performan,ce to Design and Operating Parameters [J].AnnArbor,1999,1001:48109.

    [35]KatzI,GardnerB,JongewardG,etal.Amodelofplasmacon?tactor behavior in the laboratory[C]//AIAA,Aerospace Sci?encesMeetingandExhibit,34th,Reno,1996.

    [36]CapacciM,MinucciM,SeveriA.Simplenumericalmodelde?scribing discharge parameters in orificed hollow cathode de?vices[J].AIAAPaper,1997,2791:6-9.

    [37]FearnDG,PhilipCM.Aninvestigationofphysicalprocesses inahollowcathodedischarge[J].AiaaJournal,1973,11(2):131-132.

    [38]DelcroixJL,TrindadeAR.Hollowcathodearcs[J].Advances inElectronicsandElectronPhysics,1974,35:87-190.

    [39]Siegfried D E,Wilbur P J.An investigation of mercury hollow cathodephenomena[R].AIAA,1978.

    [40]Rossetti P,Paganucci F,Andrenucci M.A Hollow Cathode Model for Application to the Electric Propulsion[C]//38th JointPropulsionConference,Indianapolis,2002:7-10.

    [41]DomonkosMT.Aparticleandenergybalancemodeloftheor?ificed hollow cathode[M].National Aeronautics and Space Administration,GlennResearchCenter,2002.

    [42]Katz I,Anderson J R,Polk J E,et al.A model of hollow cath?odeplasmachemistry[R].AIAAPaper,2002.

    [43]Katz I,Anderson J R,Polk J E,et al.Model of hollow cathode operationandlifelimitingmechanisms[C]//JPropulsionLabo?ratory,InternationalElectricPropulsionConference,2003.

    [44]Katz I,Anderson J R,Polk J E,et al.One-dimensional hol?low cathode model[J].Journal of Propulsion and Power,2003,19(4):595-600.

    [45]Katz I,Mikellides I,Goebel D.Model of the plasma potential distribution in the plume of a hollow cathode[M].Pasadena,CA:Jet Propulsion Laboratory,National Aeronautics and SpaceAdministration,2004.

    [46]Goebel D M,Katz I.Fundamentals of electric propulsion:ion andHallthrusters[M].JohnWiley&Sons,2008.

    [47]Albertoni R,Pedrini D,Paganucci F,et al.A Reduced-Or?derModelforThermionicHollowCathodes[J].IEEETransac?tionsonPlasmaScience,2013,41(7):1731-1745.

    [48]PedriniD,AlbertoniR,PaganucciF,etal.Theoreticalmodel of a lanthanum hexaboride hollow cathode[J].IEEE Transac?tionsonPlasmaScience,2015,43(1):209-217.

    [49]PedriniD,AlbertoniR,PaganucciF,etal.ModelingofLaB 6 hollow cathode performance and lifetime[J].Acta Astronauti?ca,2015,106:170-178.

    [50]KorkmazO,CelikM.GlobalNumericalModelfortheAssess?ment of the Effect of Geometry and Operation Conditions on Insert and Orifice Region Plasmas of a Thermionic Hollow Cathode Electron Source[J].Contributions to Plasma Phys?ics,2014,54(10):838-850.

    [51]郭寧,江豪成,顧佐.離子發(fā)動(dòng)機(jī)空心陰極壽命預(yù)測(cè)[J].真空,2009,46(4):114-116.

    [52]賈艷輝,郭寧,李娟,等.LIPS-200離子推力器關(guān)鍵部組件壽命分析[J].空間控制技術(shù)與應(yīng)用,2013,39(4):54-58

    [53]郭寧,江豪成,高軍,等.離子發(fā)動(dòng)機(jī)空心陰極失效形式分析[J].真空與低溫,2005,2005,11(4):239-242

    [54]陳梅.空心陰極的數(shù)值模擬及實(shí)驗(yàn)研究[D].哈爾濱:哈爾濱工業(yè)大學(xué),2007.

    [55]馮亮.電推力器空心陰極雙鞘層數(shù)值模擬[D].哈爾濱:哈爾濱工業(yè)大學(xué),2009.

    [56]王穎.自持空心陰極雙鞘層數(shù)值模擬及熱分析[D].哈爾濱:哈爾濱工業(yè)大學(xué),2010.

    [57]安秉健.自持空心陰極關(guān)鍵尺寸影響與設(shè)計(jì)準(zhǔn)則研究[D].哈爾濱:哈爾濱工業(yè)大學(xué),2014.

    REVIEW ON DISCHARGE MODELS OF ION THRUSTER HOLLOW CATHODE

    GU Zeng-jie,GUO Ning,JIAYan-hui
    (Science and Technology on Vacuum Technology and Physics Laboratory,Lanzhou Insitute of Physics,Lanzhou 730000,China)

    Aimed at understanding the modeling methods about the discharge physical process of ion thruster hollow cathode,various modeling methods were reviewed,including fluid models,Particle-in-Cell(PIC)models,hybrid models and reduced-order models.The physical models’characteristics,application areas and implementation outlines were analyzed,which would serve the modeling methods selection of hollow cathode modeling and simulation with theoretical basis.It would promote modeling solutions to the engineering problems that encountered in the development and experiment process.

    ion thruster;hollow cathode;physical model;modeling method

    V439+.1

    A

    1006-7086(2016)06-0324-07

    10.3969/j.issn.1006-7086.2016.06.003

    2016-05-25

    谷增杰(1989-),男,山東菏澤人,碩士研究生,主要從事空心陰極技術(shù)研究。E-mail:guzengjie@foxmail.com。

    猜你喜歡
    方法模型
    一半模型
    重要模型『一線三等角』
    重尾非線性自回歸模型自加權(quán)M-估計(jì)的漸近分布
    學(xué)習(xí)方法
    可能是方法不對(duì)
    3D打印中的模型分割與打包
    用對(duì)方法才能瘦
    Coco薇(2016年2期)2016-03-22 02:42:52
    FLUKA幾何模型到CAD幾何模型轉(zhuǎn)換方法初步研究
    四大方法 教你不再“坐以待病”!
    Coco薇(2015年1期)2015-08-13 02:47:34
    賺錢方法
    一级毛片 在线播放| 日日啪夜夜爽| 亚洲综合精品二区| 丝袜美足系列| 99久久综合免费| 日日摸夜夜添夜夜爱| 国产精品一区二区在线不卡| 免费观看a级毛片全部| 熟女电影av网| 国产av一区二区精品久久| 国产日韩欧美亚洲二区| 99久久综合免费| 九色成人免费人妻av| 一本色道久久久久久精品综合| 亚洲欧美精品自产自拍| 大香蕉久久网| 国产成人午夜福利电影在线观看| 在线天堂中文资源库| 在线观看免费日韩欧美大片| 久久ye,这里只有精品| a级毛片黄视频| 丝袜喷水一区| 肉色欧美久久久久久久蜜桃| 18禁裸乳无遮挡动漫免费视频| 午夜福利视频精品| 老熟女久久久| 亚洲经典国产精华液单| 色吧在线观看| 日本欧美国产在线视频| 久久久a久久爽久久v久久| 国产精品.久久久| 国产精品一国产av| av国产精品久久久久影院| 精品视频人人做人人爽| 久久免费观看电影| 国产精品欧美亚洲77777| 老司机影院成人| av在线老鸭窝| 高清黄色对白视频在线免费看| 七月丁香在线播放| 高清欧美精品videossex| 伊人久久国产一区二区| 久久国产精品大桥未久av| 天天躁夜夜躁狠狠久久av| 亚洲欧美一区二区三区国产| 女性生殖器流出的白浆| 18禁动态无遮挡网站| 亚洲综合精品二区| av免费在线看不卡| 桃花免费在线播放| 少妇的丰满在线观看| 成人无遮挡网站| 久久久久久人妻| 国产在线免费精品| 国产精品秋霞免费鲁丝片| 少妇猛男粗大的猛烈进出视频| 亚洲国产欧美日韩在线播放| 免费看光身美女| 国产精品一国产av| 日韩成人av中文字幕在线观看| 麻豆精品久久久久久蜜桃| 国产成人午夜福利电影在线观看| 蜜臀久久99精品久久宅男| 91久久精品国产一区二区三区| 国产精品99久久99久久久不卡 | 校园人妻丝袜中文字幕| 男人舔女人的私密视频| 亚洲av电影在线进入| 日韩不卡一区二区三区视频在线| 精品久久久精品久久久| 一级爰片在线观看| 一级,二级,三级黄色视频| 有码 亚洲区| 一本久久精品| 亚洲少妇的诱惑av| 免费观看在线日韩| 欧美少妇被猛烈插入视频| 97人妻天天添夜夜摸| 免费人成在线观看视频色| 成人毛片60女人毛片免费| 777米奇影视久久| 国产精品嫩草影院av在线观看| 久久久久久久亚洲中文字幕| 国产成人精品久久久久久| av又黄又爽大尺度在线免费看| 欧美 亚洲 国产 日韩一| 国产成人精品一,二区| 中文字幕另类日韩欧美亚洲嫩草| 七月丁香在线播放| 精品少妇内射三级| 乱人伦中国视频| 成年动漫av网址| 制服丝袜香蕉在线| av天堂久久9| 99热这里只有是精品在线观看| 久久av网站| 国产 一区精品| 国产免费视频播放在线视频| 三上悠亚av全集在线观看| 伦理电影大哥的女人| 日韩不卡一区二区三区视频在线| 亚洲av免费高清在线观看| 国产又爽黄色视频| 亚洲精品,欧美精品| 日本av手机在线免费观看| 亚洲精品456在线播放app| 亚洲精品美女久久av网站| 最新中文字幕久久久久| 又黄又粗又硬又大视频| 日韩精品免费视频一区二区三区 | 9191精品国产免费久久| 久久精品熟女亚洲av麻豆精品| 美女大奶头黄色视频| www.熟女人妻精品国产 | 又大又黄又爽视频免费| 成年人午夜在线观看视频| 日韩欧美一区视频在线观看| 日本色播在线视频| 在线观看免费高清a一片| 老熟女久久久| 91午夜精品亚洲一区二区三区| 久久ye,这里只有精品| 成人二区视频| 一级毛片黄色毛片免费观看视频| www.熟女人妻精品国产 | av免费在线看不卡| 亚洲精品成人av观看孕妇| 人人妻人人爽人人添夜夜欢视频| 国产精品无大码| 卡戴珊不雅视频在线播放| 97在线人人人人妻| 亚洲一区二区三区欧美精品| 色视频在线一区二区三区| 韩国精品一区二区三区 | 韩国av在线不卡| 欧美人与善性xxx| 亚洲欧洲国产日韩| 热99国产精品久久久久久7| 亚洲精品久久成人aⅴ小说| 在线观看免费日韩欧美大片| 久久精品国产综合久久久 | 国产高清国产精品国产三级| 亚洲美女搞黄在线观看| 午夜福利在线观看免费完整高清在| 高清av免费在线| 一级黄片播放器| 99国产精品免费福利视频| 日韩精品有码人妻一区| 国产片特级美女逼逼视频| 久久毛片免费看一区二区三区| 男人操女人黄网站| 一区二区三区乱码不卡18| 男人添女人高潮全过程视频| av免费观看日本| 最后的刺客免费高清国语| 成人综合一区亚洲| 国产黄频视频在线观看| 咕卡用的链子| 亚洲精品一二三| 美女主播在线视频| 国产精品熟女久久久久浪| 两个人免费观看高清视频| 国产精品无大码| av片东京热男人的天堂| 欧美日韩av久久| 人妻一区二区av| 大香蕉久久成人网| 国语对白做爰xxxⅹ性视频网站| 欧美97在线视频| 2022亚洲国产成人精品| 午夜老司机福利剧场| 精品视频人人做人人爽| 精品福利永久在线观看| 亚洲成国产人片在线观看| 免费av中文字幕在线| 中国国产av一级| 久久午夜福利片| 成年人午夜在线观看视频| 青春草国产在线视频| av在线app专区| 激情视频va一区二区三区| 超碰97精品在线观看| 国产日韩欧美亚洲二区| 99re6热这里在线精品视频| 国产成人精品无人区| 日韩精品有码人妻一区| 人人妻人人添人人爽欧美一区卜| a 毛片基地| 久久久国产欧美日韩av| 免费人妻精品一区二区三区视频| av电影中文网址| 国产不卡av网站在线观看| 欧美bdsm另类| 亚洲欧美中文字幕日韩二区| 成人亚洲精品一区在线观看| 久久久久人妻精品一区果冻| 国产一区有黄有色的免费视频| av在线观看视频网站免费| 精品少妇内射三级| 欧美人与善性xxx| 亚洲一级一片aⅴ在线观看| 亚洲欧美日韩另类电影网站| 老司机影院毛片| √禁漫天堂资源中文www| 欧美成人精品欧美一级黄| 成年美女黄网站色视频大全免费| 啦啦啦中文免费视频观看日本| 免费黄频网站在线观看国产| 婷婷色综合大香蕉| 成年人午夜在线观看视频| 日韩,欧美,国产一区二区三区| 日韩中文字幕视频在线看片| 天天躁夜夜躁狠狠躁躁| 免费人妻精品一区二区三区视频| 成年人午夜在线观看视频| 在线天堂中文资源库| 啦啦啦在线观看免费高清www| 蜜桃国产av成人99| 91午夜精品亚洲一区二区三区| 在线观看免费高清a一片| 日韩欧美一区视频在线观看| 最近最新中文字幕免费大全7| 一边摸一边做爽爽视频免费| 国产男女超爽视频在线观看| 国产综合精华液| 亚洲欧美色中文字幕在线| 91aial.com中文字幕在线观看| 国内精品宾馆在线| 国产白丝娇喘喷水9色精品| 97在线视频观看| 老熟女久久久| 最黄视频免费看| 曰老女人黄片| 免费高清在线观看日韩| 在线观看美女被高潮喷水网站| 免费少妇av软件| 一区在线观看完整版| 亚洲四区av| 不卡视频在线观看欧美| 国国产精品蜜臀av免费| 99久久人妻综合| 亚洲一区二区三区欧美精品| 国产精品久久久久成人av| 久久精品国产亚洲av涩爱| av又黄又爽大尺度在线免费看| 欧美日韩亚洲高清精品| 亚洲精品自拍成人| 国产精品熟女久久久久浪| 九草在线视频观看| 国产片内射在线| 老女人水多毛片| 欧美人与性动交α欧美精品济南到 | 欧美xxⅹ黑人| 亚洲高清免费不卡视频| 亚洲国产日韩一区二区| 亚洲精品国产色婷婷电影| 涩涩av久久男人的天堂| 制服诱惑二区| 国产精品一区www在线观看| 久久久久久久久久人人人人人人| 亚洲欧美日韩卡通动漫| 一本一本久久a久久精品综合妖精 国产伦在线观看视频一区 | 在线观看三级黄色| 国产亚洲av片在线观看秒播厂| 如日韩欧美国产精品一区二区三区| 亚洲精品日本国产第一区| 伦理电影大哥的女人| 在线观看国产h片| 在线观看一区二区三区激情| 精品久久久久久电影网| 美女视频免费永久观看网站| 汤姆久久久久久久影院中文字幕| 久久久久精品性色| 成年美女黄网站色视频大全免费| 国产成人免费无遮挡视频| kizo精华| 多毛熟女@视频| 男女午夜视频在线观看 | 日本91视频免费播放| av免费观看日本| 侵犯人妻中文字幕一二三四区| 中文字幕另类日韩欧美亚洲嫩草| 中文字幕av电影在线播放| 久久久欧美国产精品| 亚洲国产欧美在线一区| 乱码一卡2卡4卡精品| 日本色播在线视频| 亚洲精品中文字幕在线视频| 亚洲精品视频女| 婷婷色综合www| 国产白丝娇喘喷水9色精品| 视频中文字幕在线观看| 午夜免费男女啪啪视频观看| 多毛熟女@视频| tube8黄色片| 人人妻人人澡人人看| 超碰97精品在线观看| 男女啪啪激烈高潮av片| 九色亚洲精品在线播放| 国产成人精品无人区| 精品久久久精品久久久| 精品国产一区二区三区久久久樱花| 校园人妻丝袜中文字幕| 免费在线观看完整版高清| 男女边吃奶边做爰视频| 午夜免费观看性视频| 亚洲一级一片aⅴ在线观看| 韩国精品一区二区三区 | 又黄又爽又刺激的免费视频.| 国产欧美日韩综合在线一区二区| 国产成人欧美| 香蕉国产在线看| 少妇 在线观看| 男女下面插进去视频免费观看 | 这个男人来自地球电影免费观看 | 嫩草影院入口| 久久精品人人爽人人爽视色| 高清黄色对白视频在线免费看| 一区二区三区乱码不卡18| 国产精品.久久久| 日韩成人av中文字幕在线观看| 蜜桃在线观看..| 日韩av免费高清视频| 日本91视频免费播放| 国产麻豆69| 女人被躁到高潮嗷嗷叫费观| 69精品国产乱码久久久| 天堂俺去俺来也www色官网| 男的添女的下面高潮视频| 中文字幕精品免费在线观看视频 | 亚洲精品国产av成人精品| 少妇被粗大的猛进出69影院 | 国产精品人妻久久久影院| av电影中文网址| 亚洲图色成人| 大香蕉久久网| 久久久久精品久久久久真实原创| 99久久综合免费| 九九爱精品视频在线观看| 成人亚洲精品一区在线观看| 亚洲av男天堂| 国产精品一区www在线观看| 国产免费又黄又爽又色| 视频在线观看一区二区三区| 久久久久人妻精品一区果冻| 飞空精品影院首页| 九草在线视频观看| 日本vs欧美在线观看视频| 国产毛片在线视频| 男女边吃奶边做爰视频| 另类亚洲欧美激情| 欧美另类一区| 国产成人午夜福利电影在线观看| 亚洲av电影在线观看一区二区三区| 两个人免费观看高清视频| 国产精品免费大片| 日韩欧美一区视频在线观看| 国产激情久久老熟女| 国产色婷婷99| 国产男女内射视频| 热99久久久久精品小说推荐| 精品第一国产精品| 欧美人与善性xxx| 久久综合国产亚洲精品| 日韩,欧美,国产一区二区三区| 久久久久久久精品精品| 成年人免费黄色播放视频| 久久99热这里只频精品6学生| 欧美激情国产日韩精品一区| xxxhd国产人妻xxx| 免费大片黄手机在线观看| 久久精品久久精品一区二区三区| 免费大片黄手机在线观看| 国产av国产精品国产| 免费看光身美女| 亚洲,一卡二卡三卡| 日本黄色日本黄色录像| 黑人欧美特级aaaaaa片| 国产在线视频一区二区| 丝袜美足系列| 你懂的网址亚洲精品在线观看| 亚洲欧美清纯卡通| 寂寞人妻少妇视频99o| 欧美 亚洲 国产 日韩一| 亚洲欧美一区二区三区黑人 | 丁香六月天网| 国产成人精品婷婷| 考比视频在线观看| 在线观看美女被高潮喷水网站| 99国产精品免费福利视频| 久久青草综合色| 精品人妻熟女毛片av久久网站| 日韩成人av中文字幕在线观看| 欧美日本中文国产一区发布| av卡一久久| 欧美最新免费一区二区三区| av不卡在线播放| 在线看a的网站| 在线 av 中文字幕| 亚洲激情五月婷婷啪啪| 亚洲av综合色区一区| 国产 一区精品| 纯流量卡能插随身wifi吗| 精品国产乱码久久久久久小说| 欧美亚洲日本最大视频资源| 日韩伦理黄色片| 免费黄网站久久成人精品| 日本vs欧美在线观看视频| 久久久a久久爽久久v久久| xxx大片免费视频| 在线观看免费视频网站a站| 精品一区在线观看国产| 人人妻人人爽人人添夜夜欢视频| 国产1区2区3区精品| 亚洲av电影在线进入| 夜夜骑夜夜射夜夜干| 国产福利在线免费观看视频| 亚洲人与动物交配视频| 国产精品成人在线| 天天躁夜夜躁狠狠躁躁| 亚洲成国产人片在线观看| 久久狼人影院| 亚洲内射少妇av| 精品国产一区二区三区四区第35| 有码 亚洲区| 赤兔流量卡办理| 制服人妻中文乱码| 伦理电影大哥的女人| 一级毛片我不卡| 亚洲一码二码三码区别大吗| 少妇人妻精品综合一区二区| 香蕉国产在线看| 久久久久久人妻| 中国三级夫妇交换| h视频一区二区三区| 制服人妻中文乱码| 日韩熟女老妇一区二区性免费视频| 咕卡用的链子| 母亲3免费完整高清在线观看 | 纯流量卡能插随身wifi吗| 日本vs欧美在线观看视频| 巨乳人妻的诱惑在线观看| 欧美少妇被猛烈插入视频| 色5月婷婷丁香| 91在线精品国自产拍蜜月| 成人亚洲精品一区在线观看| www日本在线高清视频| 日韩av在线免费看完整版不卡| 青春草国产在线视频| 精品国产乱码久久久久久小说| av黄色大香蕉| 国产男女内射视频| 亚洲欧美成人综合另类久久久| 成人手机av| 国产老妇伦熟女老妇高清| 国产免费福利视频在线观看| 欧美国产精品va在线观看不卡| 永久网站在线| 国内精品宾馆在线| 久久久久久久久久人人人人人人| 在线观看美女被高潮喷水网站| 久久久久久久久久成人| 人人妻人人添人人爽欧美一区卜| 欧美激情国产日韩精品一区| 狂野欧美激情性xxxx在线观看| 亚洲图色成人| 免费在线观看黄色视频的| 国内精品宾馆在线| 夫妻性生交免费视频一级片| 国产激情久久老熟女| 欧美变态另类bdsm刘玥| 精品久久国产蜜桃| 男女无遮挡免费网站观看| 成人国语在线视频| 欧美成人午夜精品| 老司机影院毛片| 亚洲成人手机| 欧美激情国产日韩精品一区| 丰满乱子伦码专区| 亚洲国产看品久久| 国产欧美另类精品又又久久亚洲欧美| 自线自在国产av| 久久鲁丝午夜福利片| 国产综合精华液| 午夜福利,免费看| 久久久久人妻精品一区果冻| 午夜免费男女啪啪视频观看| 一级毛片电影观看| 人妻 亚洲 视频| 夫妻午夜视频| 精品一品国产午夜福利视频| 在线观看三级黄色| 国产色婷婷99| 日韩av免费高清视频| 亚洲国产看品久久| 日本猛色少妇xxxxx猛交久久| 国产淫语在线视频| 国产一区亚洲一区在线观看| 亚洲伊人色综图| 日韩不卡一区二区三区视频在线| 夜夜骑夜夜射夜夜干| 久久精品夜色国产| 亚洲国产精品一区二区三区在线| 精品卡一卡二卡四卡免费| 少妇高潮的动态图| 久久午夜福利片| 满18在线观看网站| 丝袜在线中文字幕| 亚洲精品色激情综合| 熟女av电影| 精品少妇内射三级| 日韩av在线免费看完整版不卡| 中文字幕免费在线视频6| 啦啦啦在线观看免费高清www| 高清视频免费观看一区二区| 婷婷色综合www| av在线观看视频网站免费| 亚洲国产精品国产精品| 久久ye,这里只有精品| 亚洲经典国产精华液单| 最近的中文字幕免费完整| 1024视频免费在线观看| 秋霞在线观看毛片| 只有这里有精品99| 高清在线视频一区二区三区| 色哟哟·www| 大香蕉久久网| 国产精品99久久99久久久不卡 | 韩国精品一区二区三区 | 色5月婷婷丁香| 欧美日韩视频高清一区二区三区二| 免费看不卡的av| 久久精品久久精品一区二区三区| 最近最新中文字幕大全免费视频 | 一级毛片 在线播放| 18+在线观看网站| videossex国产| 国产精品国产三级专区第一集| 高清在线视频一区二区三区| 在线观看国产h片| 99视频精品全部免费 在线| 麻豆乱淫一区二区| 七月丁香在线播放| 国产毛片在线视频| 日本爱情动作片www.在线观看| 国产精品人妻久久久久久| 少妇的逼好多水| 午夜福利乱码中文字幕| 精品一区在线观看国产| 午夜免费男女啪啪视频观看| 成人亚洲欧美一区二区av| 国产免费一区二区三区四区乱码| 国产免费一级a男人的天堂| 久久精品夜色国产| 亚洲av国产av综合av卡| 久久久久国产精品人妻一区二区| 亚洲精品乱码久久久久久按摩| 咕卡用的链子| 人体艺术视频欧美日本| 嫩草影院入口| 久久人妻熟女aⅴ| 曰老女人黄片| av在线观看视频网站免费| 国产精品久久久久久久久免| 亚洲av成人精品一二三区| 国产亚洲精品第一综合不卡 | 在线观看www视频免费| 最黄视频免费看| 少妇人妻精品综合一区二区| 大香蕉久久成人网| 国产成人欧美| 久久精品国产自在天天线| 男女午夜视频在线观看 | 亚洲三级黄色毛片| 自拍欧美九色日韩亚洲蝌蚪91| 少妇被粗大的猛进出69影院 | 18禁动态无遮挡网站| 国产成人91sexporn| 巨乳人妻的诱惑在线观看| 免费av中文字幕在线| 人妻一区二区av| 国精品久久久久久国模美| 日韩一本色道免费dvd| av在线app专区| 精品一区在线观看国产| 中文字幕人妻丝袜制服| 丝袜在线中文字幕| 97在线视频观看| 国产有黄有色有爽视频| 如何舔出高潮| 精品一区在线观看国产| 免费观看av网站的网址| 亚洲精品国产av蜜桃| 伦理电影免费视频| 久久人妻熟女aⅴ| 狠狠精品人妻久久久久久综合| 丰满乱子伦码专区| 在线天堂最新版资源| 国产成人欧美| 一级毛片电影观看| 日本免费在线观看一区| 纯流量卡能插随身wifi吗| 搡老乐熟女国产| 在线看a的网站| 男人舔女人的私密视频| 一级毛片电影观看| 免费观看av网站的网址| 男人舔女人的私密视频| 欧美性感艳星| 全区人妻精品视频| 女人被躁到高潮嗷嗷叫费观| 精品一区在线观看国产| 永久网站在线| 国产深夜福利视频在线观看| 国产日韩一区二区三区精品不卡|