PIC-MCC程序手冊全套

PIC-MCC程序手冊PIC_MCC的模擬方法和數(shù)值計(jì)算1.1PIC.MCC的模擬原理1.2PIC模型1.3MCC模型1.4模擬中所涉及的放電粒子1.5模擬中的碰撞數(shù)據(jù)分析程序正常運(yùn)行所需文件及文件意義主要輸岀文件的意義1PIC_MCC的模擬方法和數(shù)值計(jì)算的模擬原理目前在計(jì)算機(jī)模擬中大量采用低壓射頻放電模型來模擬材料的加工及改性。但在低壓情況下，粒子和中性氣體的碰撞不足，無法使其達(dá)到平衡，在這種情況下我們認(rèn)為放巨所產(chǎn)生的負(fù)離子及電子的速度已經(jīng)偏離了Maxwellian分布。因此，流體力學(xué)的模擬方法已經(jīng)無法準(zhǔn)確的解決此類問題。我們選用一種新的方法，即用運(yùn)動(dòng)分析的方法來解決低壓等離子反應(yīng)器中的物理和化學(xué)過程，用包含大量粒子的模型來解決Boltzmann方程。粒子模型入僅可以解決刻蝕在基板上粒子的能量問題，而且還可以很好的解決粒子刻蝕率和各向異性的問題。MonteCarlo算法與PIC模擬方法的有機(jī)結(jié)合就形成了PIC_MCC模擬方法。在PIC_MCC中，我們假設(shè)中性氣體是在時(shí)間與空間位置上的一種特定分布。PIC_MCC模擬中采用的是一維空間(Zn),速度方向?yàn)槿S(Vx,Vy,Vz)°PIC模型在Pic模型中(如圖a),粒子在電場力的作用下運(yùn)動(dòng)。粒子模擬只能解決少量粒子存在的模型，這個(gè)模型中的粒子數(shù)量遠(yuǎn)遠(yuǎn)小于真實(shí)情況下等離子體中的粒子數(shù)目，模擬中的每個(gè)粒子即超粒子表示106-109個(gè)粒子。在模擬中我們必須有足夠多的超粒子，以減少粒子離散及噪聲擾動(dòng)。超粒子與所劃分的網(wǎng)格點(diǎn)數(shù)之比必須大于1?在模擬中，我們主要解決MaxwelI方程以及F=ma=q(E+v*B)。電場可以通過MaxwalI方程解出。粒子在電場和磁場中所受的力可通過Newton_Lorentz方程解岀。

RFpowersource圖1.2.1模擬中所用的系統(tǒng)模型(a),系統(tǒng)模型z軸方向上網(wǎng)格點(diǎn)的劃分。0點(diǎn)為接射頻電極極板端，Zn為接地電極極板端。1.2.1PIC模型中一個(gè)時(shí)間步循環(huán)的數(shù)值計(jì)算(4)(5)(4)(5)圏1.2.2PIC-個(gè)時(shí)間步的循環(huán)運(yùn)算(1) 電荷密度Pj被分配到每個(gè)網(wǎng)格點(diǎn)j上，這個(gè)過程稱之為電荷分配。因此，先從連續(xù)的節(jié)點(diǎn)Zpi然后再到離散的節(jié)點(diǎn)Zj來計(jì)算電荷密度pjo電荷分配函數(shù)可用S(Zj-Zpi)來表示，包括零節(jié)點(diǎn)，第一節(jié)點(diǎn)和最后節(jié)點(diǎn)。圖1.2.3中所描述的是第一節(jié)點(diǎn)的電荷分配方式。這種分配方式把Zj節(jié)點(diǎn)上的j單胞和Zj+1節(jié)點(diǎn)上的(j+1)單胞這部分電荷看做帶電粒子或電荷云°這種帶電粒子看做是一種有限度的剛性電荷云，他們可以在通過彼此時(shí)不受束縛而自由運(yùn)動(dòng)，這種模型我們稱之為cloud_in_cell或者CIC.圖1.2.3圖1.2.3如果電荷粒子的密度是一定的，j和j+1之間的距離為AZ,那么電荷粒子qp,分配給節(jié)點(diǎn)j的電量為:AzZj+1IpiAzZj+1IpiAz分配給節(jié)點(diǎn)j+1的電量為:%+1二細(xì)%+1二細(xì)七0戶1 )Az因此，在Zpj上的電荷粒子q”分配給j節(jié)點(diǎn)的電荷密度pj為:pj=￡qpiS(Zj-Zp)Pi（2） 電荷密度可以用來計(jì)算網(wǎng)格點(diǎn)上的電場E。在靜電場模擬中，V*E=-eB/et20,所以E=-V0>,由一維條件下的Poisson方程可以得到:dz1電場可以由以下公式計(jì)算得出:

E=E==竺&"(3) E又按照函數(shù)S(Zj-Zpi)分配給網(wǎng)格點(diǎn)上的粒子。在一維的靜電場模型中,電場分配各網(wǎng)格點(diǎn)的電場為:Z）+i二Z）+i二piAzE.+帶電粒子所受的電場力為F=qE,-維靜電場模型中:&F"S（Z宀）.（4）運(yùn)動(dòng)方程可以計(jì)算出帶電粒子新的位置和速度。（4）運(yùn)動(dòng)方程可以計(jì)算出帶電粒子新的位置和速度。鼻Fdrdx—=v

dr在一維靜電場模型中，可用以下方程代替上述運(yùn)動(dòng)方程:小5—LW小5—LW尸ArJ+&—/J 二"+也/2At~因?yàn)閹щ娏Ｗ拥乃俣萔和位置X是不能同時(shí)確定的，所以leap-frog算法要采用不同模式原則。圖1.2.2leap_frog算法的網(wǎng)格點(diǎn)劃分示意圖。應(yīng)注意到初始條件下帶電粒子在時(shí)間t=0時(shí)的速度是需要改變的，把v（0）處的速度V退回到v（-△t/2）處，然后通過帶電粒子所受的電場力還計(jì)算t=0時(shí)的速度。（5） 檢查邊界條件，檢查粒子是否附著在極板上。初始時(shí)間t=0時(shí)粒子的位置和初始-△t/2時(shí)的速度已經(jīng)給出，帶電粒子的密度也可通過計(jì)算得出。圖1.2.2中所描述的（1）到（5）只是重復(fù)的循環(huán)，直到等離子體達(dá)到收斂。1.2.2邊界處等離子體粒子的模擬在射頻放電產(chǎn)生等離子體的模擬中，不僅要考慮中心等離子體處的粒子行為，也要模擬邊界處，即鞘層處的粒子行為。位勢方程的邊界條件可通過Gauss法則得岀：§sEdS二J/5+也土厶=0,%%S一離子體區(qū)域和上下兩極板的總面積Ao—下極板（接射頻電壓的極板）的表面區(qū)域面積A?—上極板（接地電極極板）的表面區(qū)域面積Oo—下極板（接射頻電壓的極板）的帶電粒子密度上極板（接地電極極板）的帶電粒子密度網(wǎng)格點(diǎn)的電勢可通過以下方程計(jì)算得岀：①丿T-2①/+①Pj（攵）2 %j=1,2,……，N-1,N為所劃分的網(wǎng)點(diǎn)。一維系統(tǒng)的邊界條件為：①N二°卩—%5。MCCmodeI1.3.1無碰撞的模擬方法pic模擬方法是一種碰撞模型。即使在低電壓情況下帶電粒子和中性氣體的碰撞也對維持放電起著非常重要的作用。碰撞可以將PIC和MC兩種方法結(jié)合起來進(jìn)行模擬運(yùn)算。PIC模擬的是所有粒子在同一時(shí)間步的運(yùn)動(dòng)。而MC方法模擬的是在碰撞中一些隨機(jī)粒子的行為，只對每個(gè)時(shí)間步中的部分粒子的行為進(jìn)行模擬，我們稱之為MCC模擬方法。1.41.5電子和中性粒子的碰撞1.5.1碰撞截面的數(shù)據(jù)我們認(rèn)為中性氣體（Ar,CL和奧）的分布是均勻的，它們的速度分配在室溫情況下（Tg=eV或者300K）是服從麥克斯韋分布的。因此，中性氣體粒子和皂子（平均Tg>2eV）相比所具有的能量很小，我們認(rèn)為它們是靜止的。碰撞截面我們用。

(g)來表示，g是粒子在碰撞之前的相對速度。在模型中所有中性粒子的碰撞截面數(shù)據(jù)和相應(yīng)的閾值在表中給出。TypeofReactionThreshold(eV)collisionArElasticscatteringe+Ar-*e+ArTota(electronice+ArTe+Ar*excitation1onizatione+Ar-*2e+Ar‘CF?Momentumtransfere+CRjTe+CF<Vibrations1excitatione+CFx->e+CF4Vibrationalexcitatione+CFLTe+CFaVibrationa1excitatione+CRTe+CF,Electronicexcitatione+CFqTe+CF<*Electronattachmente+CF—F+CF35Electronattachmente+CF—F+CF35Dissociatione+CFle+F+CF3'12Dissociativeionizatione+CF4-*2e+F+CF3'16Neutra1dissociatione+CFx-*e+F+CFa12Neutra1dissociatione+CFd-*e+2F+CF?17Neutra1dissociatione+CF〈Te+3F+CF18N2Momentumtransfere+N?-*e+N2excitatione+N?Te+板(Y)a…Ionizatione+卜2+N2'(Y)b

Ionizatione+N?T2+N/(Y)*(B?Ionizatione+N?T2+N/(Y)*(B?Z)18.1.5.2粒子碰撞后的速度計(jì)算方法一般情況下，我們所考慮的是兩個(gè)均勻球體粒子的之間的彈性碰撞。粒子在碰撞之前，我們設(shè)它們的質(zhì)量和速度分別為，m和M,v和V,它們之間的相對速度為g=v-V。在不考慮一般情況下的損失，我們認(rèn)為在碰撞之前的系統(tǒng)條件為：V=0,v=go碰撞之后的速度為V’，V’，v=go碰撞之后的速度為V’，V’，g'V'。因?yàn)榕鲎睬昂笙到y(tǒng)的總動(dòng)量守恒，在這種條件下，可以計(jì)算出質(zhì)心的速度Vg（圖）。v=v-Va/=—g=^-gT

m+Mniv=v-vCM=—-g=-^gni+M MfnM因?yàn)樵谫|(zhì)心坐標(biāo)系中兩個(gè)粒子的初始動(dòng)量大小相等，方向相反，所以,在系統(tǒng)中，兩個(gè)粒子所受的力大小相等，方向相反。所以，碰撞后的動(dòng)量也是大小相等，方向相反。碰撞之后粒子的速度方向仍然是平行的，但是卻偏離了一定角度x,如圖

圖粒子在質(zhì)心坐標(biāo)系下的運(yùn)動(dòng)軌跡。X為散射角圖粒子在質(zhì)心坐標(biāo)系下的運(yùn)動(dòng)軌跡。X為散射角要求出粒子碰撞之后的速度，就要先求出散射角X。如果可以求岀散射角X,那么可以通過以下公式來求出粒子碰撞之后的速度,v'=Vni+M[g(l-cosz)+hsinz]V=Vv'=Vni+M[g(l-cosz)+hsinz]V=Vrnni+M[g(l-cos/)+hsin/].叫=g±cos。_gxgycos(p+gg,sillq>"產(chǎn)一-~Tgx&cos低一gg,sin。n.= = ' g丄g=Jg：+g；+g：,g—g：+g；在質(zhì)心坐標(biāo)系中粒子運(yùn)動(dòng)軌跡所在的平面叫做碰撞平面。碰撞平面和參考平面之間的夾角為屮，參考平面是任意選取的，因此夾角屮為：(p=2^RR是在特定分布［0,1］之間的一個(gè)隨機(jī)數(shù)。兩個(gè)粒子碰撞之后的速度有兩個(gè)主要影響因素，第一個(gè)是角度屮，第二個(gè)是b。圖1.1.2碰撞參數(shù)屮和b系統(tǒng)要維持通量守恒，所以通過環(huán)面2ndb和立體角2nsinXdx的總粒子數(shù)相等。r2^bdb=-Fcr(g./)2庁sin/d/上式中負(fù)號表示如果增加b,x會(huì)相應(yīng)的減小。,、bdb

b(g./)=——.sin/d/db/dx在這里去絕對值，因?yàn)槠錇樨?fù)值。4了 x弓=jb(g?/)g=2才Jb(g,/)sin/d/.散射角x與兩粒子之間的電勢和速度有關(guān)系，不同的散射角，碰撞截面不同。dQ=sinzd/d(p

P（g?，）=團(tuán)邊sinTOC\o"1-5"\h\zJ 02再二Jb（g./）sin/d/=_R?

W （）不同的碰撞截面取決于粒子之間的相互作用力和電勢。例如，在電子和Ar的碰撞中，在計(jì)算電勢時(shí)可以忽略掉庫倫電勢的影響（彈性碰撞，激發(fā)，電離）。1+8￡不同情況下的碰撞截面由以下公式給出，1+8￡b(￡./)= 力,4了(1+4￡-4￡COS/)" ()￡=E/Eo是無量綱的能量，E是電子的相對能量，E。原子的單位能量（E°=eV）.2R散射角x與隨機(jī)數(shù)R及e之間的關(guān)系是:2Rcos/=1 .1+8￡(1-&)()這樣就可以求出散射角X,而粒子彈性碰撞之后的速度也可通過公式。和（）這樣就可以求出散射角X,而粒子彈性碰撞之后的速度也可通過公式。和（）求出。角度屮可通過公式（）求出。多數(shù)情況下的電子和中性粒子碰撞中，公式。和（）a)e_Arexcitatione+Ar—e+Ar*根據(jù)動(dòng)量和能量守恒公式,〃八'+MV=mv+MVIXV-Vy+E^l^v-V)2o以是非彈性碰撞的閾值。g^vV’為I，的標(biāo)量，MM()E=mv2/2是電子碰撞之前的能量。激發(fā)過程如果看做是一個(gè)彈性碰撞過程，碰撞前的速度為；何V。碰撞之后的速度可以由公式()和()求出。公式中所有的V’用；，來代替。b)electron_AionizationA是質(zhì)量為M的中性粒子。%(v)+4(V)tq(V)+勺(V)+A+(V)()A?為碰撞后的離子，&是初始電子，e,是發(fā)射出的電子。能量守恒公式，wv2MV2 m，,2MV^ + = +———+ +Et,2 2 2 2 2 ()孔為電離反應(yīng)的閾值。因?yàn)槿肷潆x子與電子的質(zhì)量比很大，我們可以認(rèn)為電子的動(dòng)量遠(yuǎn)小于中性粒子的動(dòng)量，電子與中性粒子碰撞之后，電子偏離，中性粒子變成離子，碰撞之后的軌跡仍然是未被干擾的。這種假設(shè)一方面使得在碰撞之前中性粒子的速度V'=V；另一方面，能量守恒公式表示為，nnf2_mv,2—~間-亍+^-()等式的左邊我們看做是電離能量的改變量△￡,需要找到一種方法把散射電子和發(fā)射電子區(qū)分開?！璲07WV*2R=2=(1-R)宜0E=ao+atantan-1—+tail-1—|-tan-1—'aa/aa=10.3,Elno=mv72為初始電子的能量。a：和a。單位是電子伏特。當(dāng)區(qū)分開后，我們可以從公式（）中來求得V’。言2㈤+烏）Vm（）類似于激發(fā)情況，碰撞之后散射電子的速度V’可有公式（）求出。式中V用I，來代替。碰撞之前，發(fā)射岀的電子的能量實(shí)際是不存在的，我們假設(shè)為，碰撞之后的速度可由公式?？伤愠觯藐?duì),來代替V,V’來代替V’。由于很多反應(yīng)的碰撞截面數(shù)據(jù)是無法在文獻(xiàn)中查詢到的。如果沒有碰撞截面的數(shù)據(jù)，我們就無法求出相應(yīng)的散射角X。總的碰撞截面可由公式（）得出aT=4na（g）,這樣可以得到sinXdXd^/4n。意味著碰撞之后的速度方向隨機(jī)，我們稱這種性質(zhì)的散射為各向同性。因?yàn)閑lectron_CF,和electron.Nz碰撞的碰撞截面數(shù)據(jù)是沒有的，因此我們認(rèn)為它們之間的碰撞是各向同性的。各向同性的散射角在區(qū)間[0,n],可表示為,cos/=1-2^.（）角度屮可從公式。求岀，而electron_CF.和electron.Nz彈性碰撞之后的速度可由公式（）和（）求出。與electron_Ar的彈性碰撞相似，M+m^M,g^Vo1.6lon_neutralcoIIisions1.6.1Cross_sectiondata（a）Co11isionsofAr+withneutrals因?yàn)樗蠥r*參與的反應(yīng)的碰撞截面數(shù)據(jù)均已給出，我們用null_collision方法來處理此類問題。主程序讀取文件，并計(jì)算gasdensitygasdensity=pertorr*pressure/gtemp主程序打開data/文件，讀取e和Ar的碰撞截面數(shù)據(jù)根據(jù)語句“來判斷主程序讀取哪些碰撞截面數(shù)據(jù)文件。.‘a(chǎn)rgon/CF4,)主程序讀取CF4ion與CF4反應(yīng)所需能量的文件(70,file="data/")(70,file="data/")(70,file="data/")Spec(i)所代表的粒子分別為：spec(1)-espec(2)-Arspec(3)-CF3+spec(4)-F-spec(5)-CF3-判斷語句if(nflag_col_spec(1).,由于此時(shí)只有e—種離子存在，因此調(diào)用子程序prob_e_ngas(),來模擬計(jì)算電子和所有中性氣體粒子的碰撞反應(yīng)。if(nfIag_coI_spec(1).thencalIprob_e_ngas()endif對模擬中氣體種類進(jìn)行判斷()，然后對不同的species以及spec(i)選擇調(diào)用不同的子程序來進(jìn)行碰撞模擬。if.'argon'.'argon/CF4'.'argon/CF4/*N2')thenif(nfIag_col_spec(2).thencalIprob_Ar_Ar0cc********************************************************\\\ccprobabiIityofvariousion-moIecuIecolIisions;betainf=3.poIar_cf4=polar_Ar=if.,argon/CF4,.,argon/CF4/N2,)thendois=3,nspeccc********************************************************\\\ccprobabiIityofCF4ions+CF4reactiveandeIIasticcoIIisionsif(nfIag_coI_spec(is),thenprob_col_spec(is)=e*sqrt(poIar_cf4*PI*(mass(is)+mass(0))/mass(is)/mass(0)/EPSO)*beta_inf*beta_inf** ratio_gas(2)*gden*dtis(is)cc********************************************************thenprob_col_ArCF=e*sqrt(poIar_cf4*PI*(mass(2)+mass(0))/mass(2)*/mass(0)/EPSO)*beta_inf*beta_inf*ratio_gas(2)*gden*dtis(2)endifcc********************************************************\\\ccprobabiIityofCF4ions+CF4reactiveandeIIasticcoIIisionsccincaseofpureCF4dischargeelseif.'CF4')thendois=2,nspecif(nfIag_coI_spec(is),thenprob_col_spec(is)=e*sqrt(poIar_cf4*PI*(mass(is)+*mass(0))/mass(is)/mass(0)/EPSO)*beta_inf*beta_inf** ratio_gas(2)*gden*dtis(is)endifenddoendifcc********************************************************argon/CF4/N2')thencalIprob_Ar_N2callprob_N2_N2endifcc********************************************************\\\ccwriteallprobabiIitieswrite(*,*)'prob',(prob_col_spec(i),i=1,nspec)write(*,*)if.'argon/CF4'.'argon/CF4/N2')thenwrite(*,*)'probCF4ions-Ar',(prob_col_CFAr(i),i=3,nspec)write(*,*)'probAr+一CF4',prob_col_ArCFwrite(*,*)'probAr+-N2chargeexchange',prob_col_ArN2write(*,*)'probN2+-ArchexandeIcol',prob_col_N2Arwrite(*,*)endifcc********************************************************thenif.'argon')thencalImc_e_ngas(-1,2,-1,-1,-1,-1)elseif.'argon/CF4')thenca11mc_e_ngas(0,2,3,4,5,-1)elseif.'argon/CF4/N2')thencaIImc_e_ngas(0,2,3,4,5,6)elseif.'CF4')thencaIImc_e_ngas(0,-1,2,3,4,-1)endifendifcc********************************************************argon'.'argon/CF4'.or.* .'argon/CF4/N2')thencc********************************************************thenca11mc_Ar_Ar0endifendifcc********************************************************argon/CF4/N2')thencallmc_N2_N2()endifcc********************************************************argon/CF4'.'argon/CF4/N2')thencc********************************************************thencalImc_ion_neut(3,0,dif_engy_cf3p,n_reac_cf3p,*num_col_ij_3,extra_col_ij_3,rate_col_ij_3)ca11mc_cf_Ar(3,2,num_col_ij_3,rate_col_ij_3)endifcc********************************************************thenca11tnc_ion_neut(4,0,dif_engy_fm,n_reac_fm,*num_col_ij_4,extra_col_ij_4,rate_col_ij_4)calImc_cf_Ar(4,2,num_col_ij_4,rate_col_ij_4)endifcc********************************************************thencalImc_ion_neut(5,0,dif_engy_cf3m,n_reac_cf3m,*num_col_ij_5,extra_col_ij_5,rate_col_ij_5)calImc_cf_Ar(5,2,num_col_ij_5,rate_col_ij_5)endifg********************************************************cF4')thenthencaIImc_ion_neut(2,0,dif_engy_cf3p,n_reac_cf3p,*num_col_ij_3,extra_col_ij_3,rate_col_ij_3)endifjtfLjKZsT>0^jfJfKJPL1J*C^*rZsCZfC>0^>*C>vr^*CZgCGG'^TTr'TTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTTP1thencaIImc_ion_neut(3,0,dif_engy_fm,n_reac_fm,*num_col_ij_4,extra_col_ij_4,rate_coI_ij_4)endifcc********************************************************thencalImc_ion_neut(4,0,dif_engy_cf3m,n_reac_cf3m,num_col_ij_5,extra_col_ij_5,rate_col_ij_5)g********************************************************mv'+MV'=mv+MV121-^v'-Vr+E[h=-^v-VYTakingintoconsiderationthatM+m^Mandg^vweobtainWhereE=mv2/2istheeIectronenergybeforecollision.TheexcitationprocessistreatedasifitwereaneIasticcolIisionwithpre_colIisionvelocitiesvandV.Thepost_colIisionvelocitiesaregivenbyEqs.-inwhichallvsarerepIacedbyva)lonization(electron_A)WhereAdenotestheneutraIparticlewithamass,whichcan

beassumedequaItothatoftheionM.Theprocessisrepresentedas烏（卩）+厶（卩）—馬頃）+勺頃'）+A+（V'）WhereA*denotesanion,eiistheincidenteIectron,e2istheejectedeIectron,andthesymbolsintheparenthesesdenotethevelocities.TheenergybaIanceequationismv2mv2MV2mv'2 + = 2 2 2mv%~T~2WhereEthisthethreshoIdenergyoftheionization.Becauseofthelargeion_to_eIectronmassratio,wecanassumethatthemomentumoftheincidenteIectronismuchlowerthanthemomentumoftheneutraIparticle,.theincidentelectronremovesaneIectronfromtheneutraI,andtheneutraIbecomesanion,continuingonitstrajectronundisturbed.Ononehand.thisassumptionmeansthatthecreatediontakestheveIocityanddirectionoftheneutraIparticIebeforecollision,.V'=V.Ontheotherhand,theenergybaIancecanberewrittenmv2廣mvl2mv'；2 " 2 2Theleft_handsideisknownastheexcessenergyAEafterionizationandweneedtofindanalgorithmhowtodivideitintothescatteredandejectedelectrons.IngeneraI,whenthereisnopubIishedworkonthedivisionoftheexcessenergy,itisdividedrandomlyintotwoEscat= -=R\Escat2mv.,2E/f-=（l_R）AEFortheeIectron_Arionzation,however,weusetheexpressionforsampIingtheEejbyarandomnumberR,Eej=%+qtan[/?（tan_l—+tan"—）-tan-1—]aa ac100%=—E機(jī)+10E-Eax= ——一％a=10.3,Where￡床二01撰/2istheenergyoftheincidenteIectron.Theunitsofaoanda1areelectronvolts.Oncetheexcessenergyisdivided,v'iscaIcuIatedfromEq.〃「一2（5與）Thepost_coIIisionveIocityofthescatteredeIectronv'isobtainedbyEqinwhichvisrepIacedbyv,Thepre_colIisionvelocityfortheejectedeIectron,whichdoesnotexistinreaIity,isassumedandthepost_coIIisionveIocityisobtainedagainfromEq.byreplacingvwithvejandv'withv'ej.FormanytypesofcolIisionnodataondifferentiaIcross_sectionisassumednottodependtothedefIectionangIex.ThetotaIcross_sectioncaIcuIatedfromEq.givesoT=4no(g).Therefore,theprobabiIityDefinedbybecomes8inxdxd。/4n.Thismeansthatthepost_coIIisionvelocitiestakerandomdirection.ThescatteringwiththispropertyiscalledtherearenosufficientdataforthedifferentiaIcross_sectionsofeIectron_CF4andeIectron_N2coIIisionsweassumethatthescatteringisisotropic.TheapproximationoftheinteractionpotentiaIwiththescreenedCouIombpotentiaImadeforcolIisionswithatomicArisnotvalidformolecularcaseofisotropicscatteringthedefIectionangIe,whichisintheintervaI[0,n],israndomlysampIedbycos/=l-2RTheangIe。srandomIysampIedfromandthepost_coIIisionvelocitiesforeIectron_CF4andeIectron_N2elasticcolIisionscanbefoundfromEqs.SimiIartoeIectron_AreIasticcoIIisions,theseequationscanbesomewhatsimpIifiedbyuseofM+m^Mandg^v.EIectronCF4andeIectronN2inelasticcolIisionsaretreatedsimiIartoeIectron_ArineIasticcolIisions,takingintoaccounttheisotropicscatteringaftercaseofionization,theexcessenergyisdividedbyarandomnumberusingEqs.andbetweenthescatteredandejectedelectron.IncaseofattachmenttheincidenteIectronisremovedfromthecalculationandthecreatednegativeiontakestheveIocityanddirectionoftheneutraIparticlebeforethecolIision.a、Isotropicscatteringi.Crosssectiondataii.Reactione+CFLe+CFqe+N2->e+N2iii.CaIcuIationofthepost_coIIisionvelocities(P=2^RI'M 。cos/=l-2A. 0昂=2風(fēng)汎(]_g冋iv.subroutinenewvel_e_isotrsubroutine newveI_e_isotr(sengy,sveI,svz,svy,svx,n__fIag,is1,is2)incIude"param_p1d"integern_flag,is1,is2reaI*8sengy,sveI,svx,svy,svz,rd,aIphareaI*8phi1,cosphi,sinphi,coschi,sinchi,up1,up2,up3,smag,r11,r12,r13,r21,r22,r23,r31,r32,r33TOC\o"1-5"\h\zcoschi= ()sinchi二dsqrt*coschi)phi1=*PI*rand_gen() ()cosphi=cos(phi1)sinphi二sin(phi1)8********************************************************mass(is1)*/mass(is2) ()sveI=svel*dsqrt ()sengy=sengy* ()endifr13=svzr23=svyr33=svxifthenr22=r22r22=r22/(smag+up1=up2=up3=elseup1=up2=up3=endifr12=r23*up3-r33*up2r22=r33*up1-r13*up3r32=r13*up2-r23*up1smag=dsqrt(r12*r12+r22*r22+r32*r32)r12=r12/(smag+r32=r32/(smag+r11=r22*r33-r32*r23r21=r32*r13-r12*r33r31=r12*r23-r22*r13svz=sveI*(r11*sinchi*cosphi+r12*sinchi*sinphi+r13*coschisvy=sveI*(r21*sinchi*cosphi+r22*sinchi*sinphi+r23*coschisvx=sveI*(r31*sinchi*cosphi+r32*sinchi*sinphi+r33*coschi2Rcc******************************************************2Rcos/=1 **1+85(1—R)**s(is2)v=v()彳廣f-〃“S2(1_cos同(m]+ mass(is1)*/mas2mm()sveI=sveI*dsqrt()sengy=sengy*r13=svzr23=svyr33=svxifthenup1=up2=up3=elseup1=up2=up3=endifr12=r23*up3-r33*up2r22=r33*up1-r13*up3r32=r13*up2-r23*up1smag=dsqrt(r12*r12+r22*r22+r32*r32)r12=r12/(smag+r22=r22/(smag+r32=r32/(smag+r11=r22*r33-r32*r23r21=r32*r13-r12*r33r31=r12*r23-r22*r13svz=sveI*(r11*sinchi*cosphi+r12*sinchi*sinphi+r13*coschi)svy=sveI*(r21*sinchi*cosphi+r22*sinchi*sinphi+r23*coschi)svx=sveI*(r31*sinchi*cosphi+r32*sinchi*sinphi*+r33*coschi)cc********************************************************vx(8**********************************************************csubroutinemaxwvel(svx,svy,svz,ip,np,nflag)incIude"param_p1d"reaI*8svx,svy,svz,svr,stemp,theta,randomreaI*8tvx,tvy,tvz,phiintegerip,np,nflagreaI*8sinphi,cosphi,sinthe,costhe,vel_trmcc********************************************************PI*rand_gen()random=rand_gen0dowhile random=rand_gen()enddovel_trm=dsqrt(e*stemp/mass(ip))svr=veI_trm*dsqrt*log(random))svx=svr*dcos(theta)svy=svr*dsin(theta)random=rand_gen0dowhile random=rand_gen()enddosvz=veI_trm*dsqrt*log(random))* *sin*PI*rand_gen())cc********************************************************mv'+MV1=mv+MV1 2i-^V'-vr+E[h=-^v-vyTakingintoconsiderationthatM+m^Mandg^vweobtainWhereE=mv2/2istheeIectronenergybeforecollision.TheexcitationprocessistreatedasifitwereaneIasticcolIisionwithpre_colIisionvelocitiesvandV.Thepost_colIisionvelocitiesaregivenbyEqs.-inwhichallv'sarerepIacedbyv.ii.b)lonization(eIectron_A)WhereAdenotestheneutraIparticlewithamass,whichcanbeassumedequaItothatoftheionM.Theprocessisrepresentedasq(y)+A(V)Tq(v')+e2(>')+A+(V')WhereA*denotesanion,eiistheincidenteIectron,e2istheejectedeIectron,andthesymboIsintheparenthesesdenotethevelocities.TheenergybaIanceequationis2 2 2 2 2WhereEthisthethreshoIdenergyoftheionization.BecauseoftheIargeion_to_eIectronmassratio,wecanassumethatthemomentumoftheincidenteIectronismuchIowerthanthemomentumoftheneutraIparticIe,.theincidenteIectronremovesaneIectronfromtheneutraI,andtheneutraIbecomesanion,continuingonitstrajectronundisturbed.Ononehand.thisassumptionmeansthatthecreatediontakestheveIocityanddirectionoftheneutraIparticIebeforecollision,.V'=V.Ontheotherhand,theenergybaIancecanberewrittenmv2_mvl2g'%Theleft_handsideisknownastheexcessenergyAEafterionizationandweneedtofindanalgorithmhowtodivideitintothescatteredandejectedeIectrons.IngeneraI,whenthereisnopubIishedworkonthedivisionoftheexcessenergy,itisdividedrandomlyintotwo

%=寫=戲mv.,2Eej==-=(1-R)-EFortheeIectron_Arionzation,however,weusetheexpressionforsampIingtheEejbyarandomnumberR,E.=%+atan[/?(tan_l—+tan-1—)-tan-1—]aaa2 1°°~知+10E-E

ax= ―。=10.3,WhereEi/mv?/?istheenergyoftheincidenteIectron.Theunitsofaoanda1areelectronvolts.Oncetheexcessenergyisdivided,v'iscaIcuIatedfromEq.2(Efh+Eej)mEth+￡可E時(shí)Thepost_coIIisionvelocityofthescatteredeIectronv'isobtainedbyEqinwhichvisrepIacedbyv.Thepre_colIisionvelocityfortheejectedeIectron,whichdoesnotexistin

reaIity,isassumedmMj~VVej=—mMj~VVej=—Vandthepost_coIIisionveIocityisobtainedagainfromEq.byandthepost_coIIisionveIocityisobtainedagainfromEq.byreplacingvwithandvwithvreplacingvwithandvwithvFormanytypesofcolIisionnodataondifferentiaIcross_sectionisassumednottodependtothedefIectionangIeEq.givesoTEq.givesoT=4no(g).Therefore,theprobabiIityDefinedbybecomessinxdxd。/4n.Thismeansthatthepost_coIIisionvelocitiestakerandomdirection.ThescatteringwiththispropertyiscalledDefinedbybecomessinxdxd。/4n.Thismeansthatthepost_coIIisionvelocitiestakerandomdirection.ThescatteringwiththispropertyiscalledtherearenoofthesufficientdataforthedifferentialcrosssectionseIectronCF4andeIectronN2colIisionsweassumethatscatteringisisotropic.TheapproximationoftheinteractionpotentiaIwiththescreenedCouIombpotentiaImadeforcolIisionswithatomicArisnotvalidformoIecuIarcaseofisotropicscatteringthedefIectionangIe,whichisintheintervaIofthe[0,n],israndomlysampIedbycos^f=1-27?Theangle。israndomlysampIedfromandthepost_coIIisionvelocitiesforeIectron_CF4andeIectron_N2eIasticcoIIisionscanbefoundfromEqs.SimiIartoeIectron_AreIasticcoIIisions,theseequationscanbesomewhatsimpIifiedbyuseofM+m^Mandg^v.EIectron_CF4andeIectron_N2ineIasticcoIIisionsaretreatedsimiIartoeIectron_ArineIasticcolIisions,takingintoaccounttheisotropicscatteringaftercaseofionization,theexcessenergyisdividedbyarandomnumberusingEqs.andbetweenthescatteredandejectedeIectron.IncaseofattachmenttheincidentelectronisremovedfromthecalculationandthecreatednegativeiontakesthevelocityanddirectionoftheneutraIparticlebeforethecollision.2.Ion_neutraIco11isions1)Cross_sectiondataa)ColIisionsofAr+withneutralsAr+_AreIasticisotropicscatteringandscatteringinbackwarddirection(tosimulatechargetransfer),Ar+_ArcollisionsAr+_N2colIisionsarepresentedin.TheAr+_CF4elasticcross_sectionisapproximatedtotheLangevincross_sectionforpoIarization.Sineethecross_sectiondataareavailableforalIthesecollisions,theAr+_neutraIcolIisionaretreatedbythenuIl_colIisionmethod.

Relativeenergy(eV)0.1 1 10Laboratoryionenergy(eV)FIGURE2.13Ar*Laboratoryionenergy(eV)FIGURE2.13Ar*+Archargetransfer(1)andAr*+Arelasticisotropicscattering(2)crosssectionsasafunctionofthelaboratoryionenergyineV[84],andAr*+N2chargetransfer(3)crosssectionasafunctionoftherelativeenergyineV[85].ReactionAr+A/TAr+ArIsotrAr+Ar*->Ar*+ArIsotrAr+A/TAr+ArIsotrAr+Ar*->Ar*+ArIsotrAr*+CF4TAL+CF4beta<=1isotrbeta>1anisotrAr*+N2->Ar+N2,maxwveIb)ColIisionsofN2+withneutralsi. N2*_N2andN2^_AreIasticandchargetransfercolIisionsareconsidered.FigureshowstheN2_N2eIasticisotropicscatteringandchargetransfercrosssectionsasafunctionoftheIaboratoryionenergyineV,andN2*_ArchargetransferandeIasticisotropicscatteringcrosssectionsasafunctionofthereIativeenergyineV.ThecrosssectionforN2*_N2elasticisotropicscatteringQiisapproximatedupto50eVtotheLangevincrosssectionforpolarizationscatteringandforenergieshigherthan50eVto(Qm-20CT)IikeithasbeendonefortheAr*_Arelasticisotropiccrosssections,respectively.ThedataforN2_AreIasticisotropicscatteringcrosssectionareestimatedfromthedatavaIiIabIeforAr*_ArcolIisionstakingintoaccountthatQ,1/如,where卩isthereducedmass.Sincethecross_sectiondataareavaiIabIeforaIIthesecolIisions,theN2_neutraIcolIisionsaretreatedbythenull_colIisionmethod.Relativeenergy(eV)0.1 1 10 100111111111111111111 1111111111018壹J0T9b1020N2++(N2,Ar) \1 11111111 1 11111111 1 111111110.1 1 10 100Laboratoryionenergy(eV)FIGURE2.14.N2++N2chargetransfer⑴andN2*+N2elasticisotropicscattering(2),crosssectionsasafunctionofthelaboratoryionenergyineV[84],andN2*+Archargetransfer(3)[85]andelasticisotropicscattering(4)crosssectionsasafunctionoftherelativeenergyineV.ii. N2*+N2TN2，+N2isotrN2+N2TN2+N2+maxwveIN2'+Ar-4N2*+ArisotrN；+ArHAr+N；cc**********************************************************cccSUBROUTINEMC_N2_N2-(N2+)+N2andArchargeexchangeCC******************************************************CC******************************************************cc**********************************************************ccc**********************************************************csubroutinemc_N2_N2()incIude,,param_p1d"integernum_col_extra,num_coI_extra1,num_coI_extra2integern_sion_extra,nnpreaI*8tempsigma1,tempsigma2,cs(10),engy,engyl,engy2,veI,* vneutx,vneuty,vneutz,dum,tempreaI*8beta_inf,beta,randomreaI*8gx,gy,gz,vel_reI,engy_reIreaI*8sigma6_1,sigma6_2,sigma6_3,sigma6_4reaI*8rand_genn_sion_extra=0beta_inf=cc*****************************************************mass(6)/e*(vx(6,j)*vx(6,j)+vy(6,j)*vy(6,j)*+vz(6,j)*vz(6,j))**thencs⑴二elsecs(1)=ratio_gas(3)*sigma6_1(engy)endifif(nfIag_coI ij_6(2).thencs⑵二elsecs(2)=ratio_gas(3)*sigma6_2(engy)endiftempsigma1=cs(1)/sigma_totaItempsigma2=(cs(1)+cs(2))/sigma_totaIkk=idint((zp(6(j)-zg(0))/dz)u=(zp(6,j)-zg(kk))/dzCC*********************************************************mass(6)/e*(vx(6,j)*vx(6,j)+vy(6,j)*vy(6,j)* +vz(6,j)*vz(6,j))energy_Iost_coI i=energy_Iost_col_i+(engy1-engy2)木super(6)cc********************************************************\\\cc(N2+)+Ar->eIasticandchargetransferthenj=nnp+num_coI_extra1+1dowhilesamplethevelocityoftheneutraIArvneutx=vneuty=vneutz=calImaxwveI(vneutx,vneuty,vneutz,2,0,1)engy1=.5d0*mass(6)*(vx(6,j)*vx(6,j)+vy(6,j)*vy(6,j)+vz(6,j)*vz(6,j))/ecc********************************************************d0*(mass(2)*mass(6)/(mass(2)+mass(6)))* *veI_reI*veI_reI/esigma_totaI=max_sigmav_spec(4)/veI_reIrandom=rand_gen()if(nflag_col_ij_6(3).thencs⑴二elsecs(1)=ratio_gas(1)*sigma6_3(engy_reI)endifif(nflag_col_ij_6(4).thencs⑵二elsecs(2)=ratio_gas(1)*sigma6_4(engy_reI)endiftempsigma1=cs(1)/sigma_totaItempsigma2=(cs(1)+cs(2))/sigma_totaIkk=idint((zp(6,j)-zg(0))/dz)u=(zp(6,j)-zg(kk))/dzcc*********************************************************mass(2)/e*(vx(2,k)*vx(2,k)+vy(2,k)*vy(2,k)+vz(2,k)*vz(2,k))zp(2,k)=zp(6,j)npart(2)=npart(2)+1ii=ii+1enddoendif8*********************************************************mass(6)/e*(vx(6,j)*vx(6,j)+vy(6,j)*vy(6,j)* +vz(6,j)*vz(6,j))energy_Iost_col_i=energy_Iost_col_i+(engyl~engy2)*super(6)elsecc********************************************************ThetotaIcross_sectionoftheion_moIecuIecolIisionsisderivedfromtheLangevin_HassemodeI.ThemethodcaleuIatestheprobabiIityforagivenreactionbyapplyingtheRice_Rampsperger_KasseI(RRK)thereactivecolIisionswithcorrespondingthresholdenergiesarepresentedintheAppendixcalculatedcross_sectionsfromthismethodareshownbelowinthesectiondiscussingCF3+/ CF4andF_CF4colIision,CF3"_ArandF_AreIasticcollisioncross_sectionsareapproximatedtotheLangevincross_sectionforpolarizationscattering.TheAr?_N?elasticcollisionsarenotconsideredinthepresentsimuIationssineethedensityofN2isIowcomparedtothatofArorCF4andconsequently,theprobabiIityforcollisionissmaII.SimilarIy,N2*_CF4andCF3T"_N2andF_N2co