基于Matlab的二維PIC/MCC模型的實(shí)現(xiàn)

王俊杰 鄭錦華,2 魏新煦 吳 雙 許 璐

（1. 鄭州大學(xué)化工與能源學(xué)院熱能系統(tǒng)節(jié)能技術(shù)與裝備教育部工程研究中心，鄭州 450001；2. 鄭州瑞邦石油機(jī)械有限公司，鄭州 450001）

基于Matlab的二維PIC/MCC模型的實(shí)現(xiàn)

王俊杰1鄭錦華1,2魏新煦1吳 雙1許 璐1

（1. 鄭州大學(xué)化工與能源學(xué)院熱能系統(tǒng)節(jié)能技術(shù)與裝備教育部工程研究中心，鄭州 450001；2. 鄭州瑞邦石油機(jī)械有限公司，鄭州 450001）

為了實(shí)現(xiàn)限定條件下氣體放電的二維模擬，基于Matlab通過編程手段開發(fā)模擬程序,并用此模擬程序?qū)χ绷鞯蛪簵l件下氬氣的放電過程進(jìn)行驗(yàn)證。對(duì)電子速度進(jìn)行麥克斯韋初始化，并進(jìn)行結(jié)果驗(yàn)證。提取50個(gè)時(shí)間步長(zhǎng)內(nèi)的電子數(shù)變化趨勢(shì)，驗(yàn)證了氬氣放電過程中電子崩的發(fā)生；并且通過100個(gè)時(shí)間步長(zhǎng)結(jié)束時(shí)的電子位置分布，可以觀察到陽極鞘層的產(chǎn)生；根據(jù)解泊松方程至200步和300步結(jié)果的對(duì)比研究，驗(yàn)證在本文的模擬條件下200步解泊松方程是足夠的；對(duì)電子速度均方根處理，得到溫度分布示意圖，驗(yàn)證了直流氬氣放電時(shí)陽極板附近溫度高的實(shí)驗(yàn)現(xiàn)象。至此，基于Matlab開發(fā)的二維粒子/蒙特卡洛碰撞的耦合（PIC/MCC）程序，能夠進(jìn)行氣體放電模擬。

Matlab；PIC/MCC；直流放電；數(shù)學(xué)模型；麥克斯韋分布

高電壓絕緣氣體放電模擬一直以來是電工學(xué)領(lǐng)域的熱點(diǎn)，PIC/MCC模型是低壓條件下常用的放電模型，該模型自20世紀(jì)50年代被南部等人提出，至今幾十年被不斷完善[1-2]。很多國(guó)內(nèi)的研究者采用這種方法給電器技術(shù)領(lǐng)域帶來了發(fā)展，如：湖南大學(xué)的汪沨教授采用PIC模型對(duì)10%～90%SF6/N2混合氣體中絕緣子沿面放電特性進(jìn)行研究[3]；上海交通大學(xué)的研究團(tuán)隊(duì)采用改進(jìn)的MCC方法研究SF6和CO2混合氣體電子崩的參數(shù)[4]；沈陽工業(yè)大學(xué)研究團(tuán)隊(duì)采用PIC方法對(duì)斷路器短間隙氣體放電過程進(jìn)行模擬研究[5]。

氣體放電數(shù)學(xué)模型有：流體模型即N-S模型也稱宏觀模型；粒子（PIC）模型也稱微觀模型；二者混合模型三種[6]。本文針對(duì)低壓氣體放電進(jìn)行粒子（PIC）模擬，因?yàn)樵诘蛪簵l件下（＜10mtorr）流體模型并不適用[7]。本文首次采用Matlab編程的方式，實(shí)現(xiàn)二維PIC/MCC模型。并對(duì)氬氣直流條件下的放電過程進(jìn)行模擬驗(yàn)證。氬氣是微電子器件表面處理的常用氣體，并在各種放電條件下被研究。因此選擇氬氣作為驗(yàn)證氣體。

Matlab具有直接面向矩陣的特點(diǎn)，決定其在粒子模型實(shí)現(xiàn)方面有很大的優(yōu)勢(shì)。Matlab程序簡(jiǎn)潔且易于理解，并且可以與其他匯編語言進(jìn)行混合編程，解決編程問題手段多變[8]，故采用Matlab來實(shí)現(xiàn)PIC/MCC模型的程序開發(fā)。

二維放電模型與一維模擬[9]相比有以下不同：①二維模型需要考慮平板邊界條件對(duì)放電過程的影響，而一維模型無需考慮；②二維模擬能夠直觀的反映出在放電過程中，電子與離子的運(yùn)動(dòng)情況和各個(gè)參數(shù)的變化；③對(duì)于二維模型而言，在相同條件下的粒子放電模擬，要達(dá)到同樣的參數(shù)精度，網(wǎng)格數(shù)目是一維時(shí)的二次方倍，計(jì)算時(shí)間成倍的增加，工作效率降低[10]。

本文從離散玻爾茲曼方程、建立完善的二維數(shù)學(xué)模型，到程序流程建立以及編程實(shí)現(xiàn)，再到氬氣放電驗(yàn)證等角度，來討論直流氣體放電過程中的粒子模擬，并將模擬結(jié)果與實(shí)驗(yàn)現(xiàn)象以及其他研究者所取得的相關(guān)模擬結(jié)果進(jìn)行比較。

1 PIC/MCC模型的Matlab實(shí)現(xiàn)

對(duì)直流低壓條件下氬氣PIC/MCC放電模型的研究基于以下幾點(diǎn)假設(shè)：①忽略電子運(yùn)動(dòng)形成的磁場(chǎng)；②忽略帶電粒子間的庫倫作用，即忽略庫倫碰撞；③邊界為良導(dǎo)體，電子到達(dá)表面會(huì)被吸收；④為了達(dá)到驗(yàn)證目的并減少工作量，采用矩形放電區(qū)域，網(wǎng)格為正方形網(wǎng)格。

圖1為編程實(shí)現(xiàn)模型的計(jì)算流程圖，先要給各個(gè)粒子賦初始值。（其中包括粒子所帶電荷，粒子的位置，初始速度等）進(jìn)入第一個(gè)循環(huán)，按照面積權(quán)重將粒子所帶電荷分配到各自鄰近的網(wǎng)格節(jié)點(diǎn)上。接下來將網(wǎng)格節(jié)點(diǎn)的帶電荷與外電場(chǎng)進(jìn)行耦合，解泊松方程，得到各個(gè)節(jié)點(diǎn)的電場(chǎng)強(qiáng)度。帶電粒子在得到的新的電場(chǎng)的作用下運(yùn)動(dòng)，然后利用MCC方法處理離子之間的碰撞過程。最后判斷是否結(jié)束，沒有結(jié)束進(jìn)入下一個(gè)時(shí)間步長(zhǎng)，直到循環(huán)結(jié)束，并進(jìn)行結(jié)果分析。

圖1 Matlab實(shí)現(xiàn)模型的流程圖

1.1 滿足麥克斯韋初始速度分布的實(shí)現(xiàn)

按照前面提到的模擬流程，第一步要實(shí)現(xiàn)參數(shù)的初始化。在Ar氣放電的開始階段，粒子滿足麥克斯韋速度分布，需要對(duì)粒子速度進(jìn)行初始化即速度分布滿足方程為

式中，T為氣體溫度；m為電子質(zhì)量；v為粒子速度；k為玻爾茲曼常數(shù)。

圖2 統(tǒng)計(jì)結(jié)果與麥克斯韋分布的比較

Matlab具有產(chǎn)生滿足一定分布隨機(jī)數(shù)的內(nèi)部函數(shù)，例如：滿足正態(tài)分布隨機(jī)數(shù)的產(chǎn)生，平均分布隨機(jī)數(shù)的產(chǎn)生等，但是沒有產(chǎn)生滿足麥克斯韋分布的隨機(jī)數(shù)的函數(shù)?；谑剑?）麥克斯韋分布函數(shù)，采用舍選法[11]，編程實(shí)現(xiàn)滿足麥克斯韋分布函數(shù)的隨機(jī)速度的生成，并編程對(duì)隨機(jī)數(shù)進(jìn)行驗(yàn)證，統(tǒng)計(jì)在模擬壓力和溫度的條件下，隨機(jī)數(shù)的分布情況。對(duì)200個(gè)隨機(jī)數(shù)的統(tǒng)計(jì)結(jié)果進(jìn)行八次多項(xiàng)式擬合，圖2中的擬合結(jié)果與麥克斯韋分布曲線一致，驗(yàn)證了編程的正確性，說明利用此方法，可以完成粒子速度大小的初始化。

1.2 玻爾茲曼方程離散之PIC部分

玻爾茲曼方程是氣體動(dòng)力學(xué)基本方程，PIC/MCC方法是基于玻爾茲曼離散方程的一種模擬方法。根據(jù)文獻(xiàn)[1-2]玻爾茲曼方程可以離散為粒子部分和碰撞部分。其中粒子部分包括電荷分配模型、電場(chǎng)模型以及運(yùn)動(dòng)模型。碰撞部分用MCC模型進(jìn)行研究。

在粒子賦予初值之后，利用權(quán)重函數(shù)（面積權(quán)重）進(jìn)行二維電荷分配，實(shí)現(xiàn)電荷向網(wǎng)格節(jié)點(diǎn)的分配。圖3為權(quán)重向節(jié)點(diǎn)分配的示意圖，利用式（2）[1]計(jì)算節(jié)點(diǎn)A的面積權(quán)重。圖3中（Xp,Yp）為粒子點(diǎn)的位置坐標(biāo)，正方形單元網(wǎng)格的邊長(zhǎng)AB=Δh。與式（2）類似可以得到B、C、D點(diǎn)的權(quán)重式。

郜教授提出，乘的本質(zhì)是縮放。讀者可反復(fù)讀“乘的本質(zhì)是縮放”這一段，對(duì)乘法會(huì)有新的認(rèn)識(shí)。張景中院士的文章也提出，現(xiàn)代數(shù)學(xué)中，并不一定從加法出發(fā)引入乘法，“幾個(gè)幾”的說法，不是數(shù)學(xué)語言，而是直觀描述，從集合概念出發(fā)，建立乘法的概念。然而對(duì)于小學(xué)生來說，直觀可描述建立乘法概念也是很好的，作為教師當(dāng)然不應(yīng)該只滿足于這些，所以郜教授給我們提供的解讀也是非常必要的，對(duì)于解讀教材與教學(xué)設(shè)計(jì)應(yīng)是有幫助的。郜教授的解讀說明教師掌握本體性知識(shí)的重要性，提醒教師看教材要有一種居高臨下的能力。

圖3 面積權(quán)重分配示意圖

根據(jù)網(wǎng)格節(jié)點(diǎn)的電荷密度，利用二維泊松方程實(shí)現(xiàn)電場(chǎng)求解，具體方法是有限差分法。方程的離散形式為

式中，ρi,j是網(wǎng)格節(jié)點(diǎn)（Xi,Yj）（任意節(jié)點(diǎn)）處的電荷密度。在泊松方程的求解過程中，每個(gè)步長(zhǎng)內(nèi)用200步就可以實(shí)現(xiàn)網(wǎng)格節(jié)點(diǎn)電勢(shì)的求解（具體內(nèi)容在本文第三部分論述）。由泊松方程解得網(wǎng)格節(jié)點(diǎn)的電勢(shì)分布，電勢(shì)向網(wǎng)格節(jié)點(diǎn)線性差值，可以實(shí)現(xiàn)節(jié)點(diǎn)電場(chǎng)強(qiáng)度的求解。在編程的過程中，需要注意邊界線性差值公式與中間節(jié)點(diǎn)的區(qū)別。每一個(gè)時(shí)間步長(zhǎng)求解泊松方程，會(huì)降低程序的運(yùn)行效率，但根據(jù)相關(guān)文獻(xiàn)，采取每一步長(zhǎng)內(nèi)對(duì)方程求解，能提高方程的收斂性。

解泊松方程得到網(wǎng)格節(jié)點(diǎn)電場(chǎng)的分布，帶電粒子在極板電場(chǎng)和帶電粒子的耦合電場(chǎng)作用下運(yùn)動(dòng)。運(yùn)動(dòng)方程的離散為

式（4）、式（5）分別是位移和速度在x方向上的離散，而y方向以及vy,vz的求解，與x方向類似，其中p代表第p個(gè)粒子。這里是耦合電場(chǎng)，在計(jì)算vx時(shí)需要通過解得新一個(gè)步長(zhǎng)的電場(chǎng)分布，新一個(gè)步長(zhǎng)的電場(chǎng)的求解需通過式（4）計(jì)算出新的粒子位置，進(jìn)而算出t+Δt時(shí)刻電場(chǎng)的大小，再利用式（5），得到新時(shí)間步長(zhǎng)下的速度分布。這樣便實(shí)現(xiàn)了電場(chǎng)與運(yùn)動(dòng)方程的耦合。

1.3 玻爾茲曼方程離散之MCC模型

氣體在放電的過程中，離子會(huì)發(fā)生劇烈的碰撞。碰撞形式有：離子或電子與分子間的彈性碰撞、電荷交換碰撞、電離碰撞、激發(fā)碰撞等。在電場(chǎng)中，發(fā)生何種形式的碰撞可以用MCC方法進(jìn)行判斷。MCC模型假設(shè)碰撞發(fā)生在一個(gè)時(shí)間步長(zhǎng)的結(jié)束時(shí)刻，這樣對(duì)于編程是可以實(shí)現(xiàn)的。對(duì)粒子先后用PIC和MCC模型處理，便實(shí)現(xiàn)了粒子模型與碰撞模型的耦合[1-2,12-13]。在利用Matlab編程的過程中需注意碰撞過程中電子的消失與出現(xiàn)的問題，可以利用Matlab提供的空函數(shù)和賦值函數(shù)來實(shí)現(xiàn)。

MCC模型在一個(gè)時(shí)間步長(zhǎng)內(nèi)，一個(gè)粒子最多只能處理一次碰撞，這就需要考慮時(shí)間步長(zhǎng)的選取不能太大。參考birdsall提出的方法[1,12]，選取時(shí)間步長(zhǎng)Dt=10-8s。

在放電過程中發(fā)生碰撞的概率，采用以下式確定：

式中，n為粒子密度；v為速度；ε為電子能量；σi為總碰撞截面，是各種類型的碰撞截面的總和。碰撞截面是與能量有關(guān)的量，隨著粒子能量的變化而變化。意味著，粒子在不同時(shí)刻，所具有的總能量不同，決定發(fā)生碰撞的概率不同。

根據(jù)動(dòng)量和能量守恒定律，中子和電子碰撞前的速度分別為V, v，碰撞后的速度分別為V', v'。當(dāng)處理電子與中性粒子發(fā)生碰撞時(shí)，假設(shè)M+m≈M，V'≈V，即二者質(zhì)量不在同一數(shù)量級(jí)，忽略中性粒子碰撞前后的速度差。電子碰撞后的速度可由下式計(jì)算得到[12-13]：

式中，h為笛卡爾分量；χ為散射角。散射角的計(jì)算用各項(xiàng)同性散射假設(shè)，參考式cosχ=1-2R2，R2為[0, 1]之間的隨機(jī)數(shù)。關(guān)于散射角的處理可以參考其他論文的方法[1-2]。

2 模擬條件與結(jié)果分析

模擬條件：上表面加100V的直流電壓，極間距為1mm，極板寬度為1mm；初始粒子權(quán)重為1010；氣體壓力：1Pa；初始時(shí)刻溫度：300K。其余3個(gè)壁面接地，進(jìn)行模擬。

模擬氬氣放電過程，只考慮電子與中性粒子之間的彈性碰撞、電離碰撞和激發(fā)碰撞。經(jīng)過100個(gè)時(shí)間步長(zhǎng)之后，得到如圖4所示的電子圖4（a）和離子圖4（b）位置分布，與實(shí)拍實(shí)驗(yàn)照片圖4（c）對(duì)比。由圖可以看到，在陽極附近區(qū)域有很多電子積累，電子的密度也呈現(xiàn)有一定規(guī)律的不均勻分布。

圖4 氬氣放電過程模擬

圖4（a）為電子的位置分布，可以看到明顯的陽極鞘層的形成。在放電初始階段，由于電子的質(zhì)量比離子質(zhì)量小的多，在電場(chǎng)的作用下，高速向陽極運(yùn)動(dòng)被陽極大量吸收。離子的滯后性，會(huì)在電場(chǎng)停留積累。隨著時(shí)間的推移，電子被吸收的越多，留下的正離子也越多，這時(shí)離子與陰極板形成的電場(chǎng)會(huì)中和掉一部分的極板電場(chǎng)，在放電區(qū)域與陰極板間會(huì)形成電壓梯度很大的鞘層區(qū)域，稱為陰極鞘層。由于離子運(yùn)動(dòng)速度比電子慢的多，甚至?xí)陉枠O附近出現(xiàn)離子正電勢(shì)大于基板正電勢(shì)的情況，這樣便在陽極附近形成與陰極類似的區(qū)域，稱為陽極鞘層，與圖4（c）比較，可以看出電子分布與實(shí)驗(yàn)圖片吻合。

圖4（b）為離子的位置分布。經(jīng)過一段時(shí)間的運(yùn)動(dòng)，離子被陰極吸附。由于陰極鞘層的作用，陽極會(huì)發(fā)生二次電子發(fā)射，來維持等離子體區(qū)域的電子平衡，這是等離子體形成的重要過程。本文模擬沒有考慮二次電子發(fā)射的影響，這是模型需要改進(jìn)的地方。實(shí)驗(yàn)照片圖4（c）中，可以明顯看出陰極板上方有很明亮的放電區(qū)域，即為陰極鞘層。此現(xiàn)象與圖4（b）中離子分布的結(jié)果一致。

在二維模擬電子和離子分布的情況下，受邊界0電荷條件的限制，電場(chǎng)分布并不是理想的直線分布。分布在不同位置的電子所受到的電場(chǎng)力是不同的。速度和位置的不同，決定粒子能量不同，發(fā)生碰撞的概率也不相同。因此二維模擬相比與一維更能反映真實(shí)電場(chǎng)中的情況。

取模擬步長(zhǎng)的前50步，得到模擬電子數(shù)在前50個(gè)步長(zhǎng)中的變化趨勢(shì)，如圖5所示。模擬區(qū)域中電子數(shù)目是反映放電激烈程度的主要參數(shù)。

在模擬開始階段電子數(shù)目有一些波動(dòng)，這是電子運(yùn)動(dòng)到陽極板和粒子電離出電子的雙重作用的結(jié)果。在模擬進(jìn)行到一段時(shí)間后，由于陰極附近離子的增加，被陰極吸收的電子數(shù)將顯著減少，這時(shí)原子電離占主要作用，電子數(shù)目迅速增加。此結(jié)果與相關(guān)論文[9]中的模擬結(jié)果是一致的。

圖5 電子數(shù)隨時(shí)間步長(zhǎng)的變化曲線

根據(jù)步長(zhǎng)結(jié)束時(shí)電勢(shì)分布圖（如圖6所示），可以看出，電勢(shì)分布在二維區(qū)域中的分布是光滑的曲面。與300個(gè)循環(huán)的結(jié)果對(duì)比發(fā)現(xiàn)，200個(gè)循環(huán)解本模型中的泊松方程是足夠的。并且能夠明顯地看到邊界電勢(shì)與中間部位的電勢(shì)分布趨勢(shì)是不一樣的。中間區(qū)域變化較平滑，即電場(chǎng)強(qiáng)度比較接近。而在放電區(qū)域邊緣，電勢(shì)變化較劇烈，即場(chǎng)強(qiáng)變化較大。再一次驗(yàn)證了二維模擬與一維模擬是有顯著區(qū)別的，故不能忽略二維邊界條件對(duì)放電過程的影響。

圖6 步長(zhǎng)結(jié)束時(shí)電勢(shì)分布圖

對(duì)速度均方根處理，得到50個(gè)時(shí)間步長(zhǎng)結(jié)束時(shí)，速度的對(duì)數(shù)分布示意圖如圖7所示。由于速度與溫度的正相關(guān)性，速度分布可以反映出溫度分布。由圖可以看出在距離陽極板0.4mm左右，橫向中間的位置，半徑為0.13mm的區(qū)域內(nèi)溫度最高，電子

圖7 速度的對(duì)數(shù)（溫度）分布

運(yùn)動(dòng)速度最快，放電最激烈。所以陽極板容易受到周圍放電反應(yīng)的影響。此結(jié)果與氬氣放電反應(yīng)的實(shí)驗(yàn)相吻合，也驗(yàn)證了模型的正確性。

3 結(jié)論

本文基于Matlab編程實(shí)現(xiàn)了二維氣體放電的PIC/MCC模型，并對(duì)直流低壓條件下氬氣放電過程進(jìn)行了模擬驗(yàn)證：①編程實(shí)現(xiàn)麥克斯韋隨機(jī)速度分布，實(shí)現(xiàn)粒子速度的初始化。其次，通過觀察100個(gè)時(shí)間步長(zhǎng)后的二維電子與離子位置分布圖，可以看到陽極鞘層的形成，此現(xiàn)象與實(shí)際氬氣放電過程相符。模擬過程考慮了邊界條件對(duì)模型的影響，相較與一維模擬，該方法更接近真是的實(shí)驗(yàn)條件；②比較200個(gè)循環(huán)結(jié)束時(shí)的電勢(shì)分布與300個(gè)循環(huán)時(shí)的泊松方程的解，證明泊松方程200個(gè)循環(huán)的求解是足夠的。從此，根據(jù)模擬結(jié)果50個(gè)時(shí)間步長(zhǎng)結(jié)束時(shí)超電子數(shù)達(dá)到110，并有直線增加的趨勢(shì)，可以說明電子崩的發(fā)生。最后，根據(jù)溫度分布示意圖可以得到溫度最高區(qū)域在陽極附近，此結(jié)果與氣體放電實(shí)驗(yàn)相吻合。至此，可以證明通過Matlab編程基本實(shí)現(xiàn)了PIC/MCC模型。

[1] Birdsall C K. Particle-in-Cell Charged-Particle simulations, plus Monte Carlo collision with neutral atoms, PIC-MCC[J]. IEEE Trans.Plasma Sci, 1991, 19(2): 65-85.

[2] Nanbu K. Probability theory of electron-molecule, ion-molecule, molecule-molecule, and coulomb collisions for particle modeling of materials processing plasma and gases[J]. IEEE Trans. on Plasma Science, 2000, 28(3): 971-989.

[3] 汪沨, 肖曉林, 張憲標(biāo), 等. 基于PIC法SF6/N2混合氣體中絕緣子沿面放電特性研究[J]. 電工技術(shù)學(xué)報(bào), 2011, 26(8): 220-226.

[4] 吳變桃, 肖登明. 用改進(jìn)的蒙特卡羅法模擬SF6和CO2混合氣體電子崩參數(shù)[J]. 電工技術(shù)學(xué)報(bào), 2007, 22(1): 13-16.

[5] 李靜, 曹云東, 王爾智, 等. 斷路器短間隙氣體擊穿過程的粒子模擬[J]. 中國(guó)電機(jī)工程學(xué)報(bào), 2010(16): 125-130.

[6] 汪沨, 李錳, 潘雄峰, 等. 基于FEM-FCT算法的SF6/N2混合氣體中棒-板間隙電暈放電特性的仿真研究[J]. 電工技術(shù)學(xué)報(bào), 2013, 28(9): 261-267.

[7] 馮明, 張鑒, 黃慶安. ICP刻蝕中等離子體分布的模擬[J]. 電子器件, 2006, 28(3): 529-531.

[8] 張志涌. MATLAB教程[M]. 5版. 北京: 北京航空航天大學(xué)出版社, 2003.

[9] 李靜, 曹云東, 鄒積巖, 等. 直流氣體放電過程中的離子分子碰撞模型[J]. 高電壓技術(shù), 2009(7): 1677-1682.

[10] 邵福球. 等離子體離子模擬[M]. 北京: 科學(xué)出版社, 2002.

[11] 劉沛華, 魯華祥, 龔國(guó)良, 等. 基于FPGA的高速任意分布偽隨機(jī)數(shù)發(fā)生器[J]. 應(yīng)用科學(xué)學(xué)報(bào), 2012, 30(3): 306-310.

[12] Nanbu K. Particle modeling of nonequilibrium plasmas and gases for materials processing[J]. Vacuum Science Technology, 2004(5): 1-57.

[13] Birdsall C K. Plasma physics via computer simulation[M]. New York: McGraw-Hill, 1985.

Realization of Two-Dimensional PIC/MCC Model based on Matlab

Wang Junjie1Zheng Jinhua1,2Wei Xinxu1Wu Shuang1Xu Lu1
(1. Engineering Research Center of Energy-saving Technology & Equipment of Thermal Energy System, Ministry of Education, School of Chemical Engineering and Energy, Zhengzhou University, Zhengzhou 450001;
2. Zhengzhou Reborn Petroleum Machinery Co., Ltd, Zhengzhou 450001)

In order to achieve the two-dimensional simulation of gas discharge under limited conditions, a simulation program based on Matlab has been developed in the paper, and utilized to verify the discharge process of argon under direct current and low pressure. Electronic speed was initialized through Maxwell equations and the results were verified. The occurrence of electron avalanche in the process of argon discharge was verified by extracting the change trend of the number of electrons in the 50 time steps, and the production of anode sheath was observed through the distribution of electron position at the end of the 100 time step. According to the solution of Poisson’s equation to the 200 step and 300 step, the solution of Poisson's equation in the 200 step was sufficient under the simulated conditions in this paper. The diagram of temperature distribution given by the root mean square treatments of the electron velocity successfully verified the experimental phenomena that the temperature near the anode plate was higher in the DC argon discharge. Thus, the two-dimensional PIC/MCC program based on Matlab can satisfactorily simulate the gas discharge.

Matlab; PIC/MCC; DC discharge; mathematical model; maxwell distribution

王俊杰（1992-），男，碩士研究生，研究方向?yàn)榈入x子體放電的模擬研究。