高冬平 丁耀根 劉濮鯤 杜朝海 趙 鼎

①(中國(guó)科學(xué)院電子學(xué)研究所中國(guó)科學(xué)院高功率微波源與技術(shù)重點(diǎn)實(shí)驗(yàn)室 北京 100190)②(中國(guó)科學(xué)院研究生院 北京 100039)

1 引言

速調(diào)管自1939年問(wèn)世以來(lái)，作為一種高功率、高增益和高效率的微波放大器件，已經(jīng)在科學(xué)研究、國(guó)防建設(shè)和工業(yè)生產(chǎn)等諸多領(lǐng)域得到了廣泛的應(yīng)用。速調(diào)管工作穩(wěn)定可靠，使用壽命長(zhǎng)，在大功率微波應(yīng)用領(lǐng)域有著不可替代的作用。但是，速調(diào)管的工作機(jī)理很復(fù)雜，實(shí)際研究過(guò)程中通常依靠多次制作實(shí)驗(yàn)樣管來(lái)修正設(shè)計(jì)，不僅研制周期長(zhǎng)，而且耗資巨大。為更加準(zhǔn)確地設(shè)計(jì)速調(diào)管，國(guó)內(nèi)外的許多研究機(jī)構(gòu)開發(fā)出一系列模擬軟件，比如基于1維(1D)模型的KLY6，JPNDSK；基于2.5維(2.5D)模型的FCI[1]，Arsenal. MSU，中國(guó)科學(xué)院電子所的LSP[2]，TESLA[3]，和Magic2D[4]；以及基于3維(3D)模型的Magic3D[4,5]。1維模型通常采用圓盤模型。該模型忽略了電子的徑向運(yùn)動(dòng)和電子間的徑向作用力，認(rèn)為每個(gè)圓盤中包含的電子具有全同的運(yùn)動(dòng)狀態(tài)。1D模型能較好地模擬聚焦良好且電子注波動(dòng)小的速調(diào)管。但是對(duì)于導(dǎo)流系數(shù)較高，電子注波動(dòng)較大的速調(diào)管，電子的徑向運(yùn)動(dòng)已經(jīng)不可忽略，故此時(shí)得到的模擬結(jié)果與實(shí)測(cè)值偏離較大。3維粒子模擬程序(PIC：Particle-In-Cell simulation)可以完全考慮粒子在3個(gè)方向的運(yùn)動(dòng)信息，通過(guò)求解全3維的Maxwell方程來(lái)推進(jìn)粒子的運(yùn)動(dòng)。理論上講，該方法最接近速調(diào)管的實(shí)際物理情形。然而，3維粒子模型由于是全時(shí)域粒子模擬代碼，需要巨大的計(jì)算資源，這嚴(yán)重影響該模型的實(shí)用性。

常規(guī)的速調(diào)管具有軸對(duì)稱結(jié)構(gòu)，可以忽略電子的角向運(yùn)動(dòng)對(duì)互作用系統(tǒng)的影響，從而將全3維模型簡(jiǎn)化為2.5D(2維位移，3維速度)模型。此外對(duì)于導(dǎo)流系數(shù)較大、電子注徑向波動(dòng)較劇烈的情形，2.5D粒子模型能夠準(zhǔn)確地考慮電子的徑向運(yùn)動(dòng)，很好地處理漂移管壁對(duì)電子的截獲，全面準(zhǔn)確地分析輸出腔入口和出口位置處電子注的群聚狀態(tài)和電子注的運(yùn)動(dòng)電流。這些問(wèn)題都是1維模型所不能處理的?，F(xiàn)有的大多數(shù)2.5D模型都是基于端口近似的頻域穩(wěn)態(tài)模型，不能研究系統(tǒng)的瞬態(tài)變化過(guò)程；此外Magic2D采用全數(shù)值的時(shí)域模型，雖然能夠較全面地模擬軸對(duì)稱結(jié)構(gòu)的互作用系統(tǒng)，但是消耗資源較大，運(yùn)算時(shí)間太長(zhǎng)，缺乏實(shí)用性。本文改進(jìn)的2.5D模型是一個(gè)將端口近似和PIC相結(jié)合的時(shí)域粒子模型。該模型可以模擬互作用系統(tǒng)隨時(shí)間的演化過(guò)程，從而有助于分析諧振腔間隙電壓從開始建立至達(dá)到穩(wěn)定的復(fù)雜過(guò)程；同時(shí)，由于采用端口近似模擬高頻諧振腔，減少了高頻場(chǎng)的計(jì)算時(shí)間，比全數(shù)值的Magic2D高效，更具有實(shí)用性。

本文將著重介紹用于研究速調(diào)管中復(fù)雜注波互作用過(guò)程的2.5D粒子模型。該模型是首次將端口近似方法和PIC結(jié)合起來(lái)的時(shí)域粒子模型。大量的研究表明端口近似方法能夠比較準(zhǔn)確地模擬諧振腔對(duì)互作用系統(tǒng)的影響[1?4]，所以該模型具有高效、準(zhǔn)確和節(jié)省計(jì)算資源等優(yōu)點(diǎn)。在每一個(gè)時(shí)間步中，該模型可以選擇采用面積權(quán)重或體積權(quán)重將粒子電荷恰當(dāng)?shù)胤峙涞娇臻g網(wǎng)格各點(diǎn)上，再根據(jù)電流連續(xù)性方程將電流分配到網(wǎng)格各邊上；采用時(shí)域有限差分方法求解麥克斯韋方程組以精確解出電子注在漂移管中產(chǎn)生的空間電荷場(chǎng)；采用端口近似來(lái)模擬高頻諧振腔場(chǎng)，諧振腔場(chǎng)分布可以采用解析場(chǎng)(漂移頭間隙電場(chǎng)分布采用雙曲余弦函數(shù)近似表示)或數(shù)值場(chǎng)(通過(guò)Superfish或CST等軟件仿真求得)兩種形式；最后，采用Boris旋動(dòng)格式求解洛倫茲力方程來(lái)更新粒子的速度和位置，并在柱坐標(biāo)系下對(duì)粒子的速度推進(jìn)作了修正。C++語(yǔ)言特有的類設(shè)計(jì)也使得KLY2D具有很強(qiáng)的擴(kuò)展性，它提供了多種模擬邊界，有金屬邊界、對(duì)稱軸邊界和漂移管右端的吸收邊界。本文將在最后給出對(duì)具體速調(diào)管實(shí)例的計(jì)算結(jié)果與相關(guān)討論。

2 物理模型

2.1 運(yùn)動(dòng)粒子電荷與電流的PIC分配

電荷PIC分配采用面積加權(quán)[6]的算法來(lái)處理，具體過(guò)程如圖1所示。設(shè)該宏粒子所處的位置為xn，格點(diǎn)1所占的面積為a，格點(diǎn)2所占的面積為b，格點(diǎn)3所占的面積為c，格點(diǎn)4所占的面積為d，軸、徑向的網(wǎng)格步長(zhǎng)分別為hz，hr，該網(wǎng)格總面積為S=h?h，另外再定義= (xn?)?z0，zr=(xn?xj,k)?r0，其中z0、r0分別為z、r方向的單位向量。

圖1 PIC電荷分配到網(wǎng)格點(diǎn)

根據(jù)面積權(quán)重分配原則，4個(gè)網(wǎng)格點(diǎn)的電荷密度為

在經(jīng)過(guò)一個(gè)時(shí)間步長(zhǎng)dt以后，宏粒子由初始位置A(xn)運(yùn)動(dòng)到了終點(diǎn)位置B(xn*)，位移為AB(如圖2所示)，下面僅討論該宏粒子位移處于一個(gè)網(wǎng)格之內(nèi)的情形(該宏粒子位移穿過(guò)多個(gè)網(wǎng)格的情形可以分開成多段采用疊加原理進(jìn)行處理)。

另Δw=wn*?wn，=(wn*+wn)/2，根據(jù)式(1)的電荷分配原則和電流連續(xù)性方程?Q/?t+∫○J?dS=0，對(duì)于圖2中節(jié)點(diǎn)1有ΔQj,k/Δt ++Ir,j,k+1/2=0，同理對(duì)于節(jié)點(diǎn)2、3和4可以列出類似的方程，聯(lián)立求解這組方程可以獲得相應(yīng)網(wǎng)格邊電流密度，

圖2 PIC電流分配到網(wǎng)格各邊

這樣在每個(gè)時(shí)間步推進(jìn)以后，都可以根據(jù)以上分析來(lái)確定空間各個(gè)網(wǎng)格點(diǎn)的電荷密度以及各網(wǎng)格邊上的電流密度，從而為下一步利用FDTD方法來(lái)求解麥克斯韋方程組作好準(zhǔn)備。

2.2 粒子注的空間電荷場(chǎng)

本粒子模型在一個(gè)時(shí)間步中的推進(jìn)流程如圖3所示：首先從右下角的方框開始，將電荷、電流分別加權(quán)分配于各網(wǎng)格點(diǎn)和各網(wǎng)格邊，然后結(jié)合一定的邊界條件采用時(shí)域有限差分算法來(lái)求解完整的麥克斯韋方程組，從而得到各個(gè)網(wǎng)格點(diǎn)的空間電荷場(chǎng)，即圖3中的(Esc,Bsc)，再結(jié)合高頻諧振腔場(chǎng)(Ecav,Bcav)得到施加于宏粒子的總電磁場(chǎng)，從而得到作用于宏粒子上的力，求解洛倫茲力方程，得到宏粒子新的速度與位置，接下來(lái)應(yīng)用粒子在邊界上的發(fā)射、吸收條件(某個(gè)邊界作為粒子發(fā)射面使用，某個(gè)邊界作為吸收邊界使用，這些都可以在輸入文件中定義)。隨著時(shí)間步的推進(jìn)，當(dāng)一個(gè)高頻周期T到來(lái)時(shí)，再根據(jù)Ramo定理來(lái)計(jì)算諧振腔新的間隙電壓，一直迭代到所有的諧振腔間隙電壓都達(dá)到穩(wěn)定為止。

粒子注在漂移管空間所產(chǎn)生的電磁場(chǎng)即空間電荷場(chǎng)，程序中采用時(shí)域有限差分方法[7,8]求解Maxwell旋度方程來(lái)獲得空間電荷場(chǎng)隨時(shí)間的演化。整個(gè)求解過(guò)程是在圓柱形漂移管中進(jìn)行的，忽略了具體諧振腔的形狀，僅僅是一個(gè)漂移管模型，高頻諧振腔對(duì)粒子注的影響將在后面討論。

圖3 粒子的推進(jìn)流程

2.3 高頻諧振腔場(chǎng)

速調(diào)管的高頻諧振腔采用端口近似(portapproximation)來(lái)處理，如圖4所示。

圖4 端口近似的示意圖

若已知諧振腔某模式(本程序中暫時(shí)先考慮TM010模)下的電磁場(chǎng)分布形式，如電場(chǎng)分布函數(shù)為fz(r, z)和fr(r, z)，磁場(chǎng)分布函數(shù)為fθ(r, z)，實(shí)際的諧振腔間隙電磁場(chǎng)值可以由諧振腔的間隙電壓求得。本程序中的場(chǎng)形分布函數(shù)可以采用解析場(chǎng)和數(shù)值場(chǎng)兩種形式來(lái)表征。解析場(chǎng)是根據(jù)諧振腔的特征結(jié)構(gòu)參數(shù)由雙曲余弦級(jí)數(shù)求和表示，這樣得出來(lái)的諧振腔間隙場(chǎng)形與實(shí)際情形比較接近。數(shù)值場(chǎng)是先根據(jù)諧振腔的幾何參數(shù)利用第3方軟件，如CST微波工作室或Superfish來(lái)計(jì)算出諧振腔的諧振頻率，R/Q，Q值，以及場(chǎng)形分布函數(shù)，而具體網(wǎng)格點(diǎn)的電、磁場(chǎng)值再通過(guò)插值方法來(lái)得到。

運(yùn)動(dòng)粒子經(jīng)過(guò)諧振腔時(shí)將會(huì)在諧振腔中激起感應(yīng)電流，如果間隙之間的距離遠(yuǎn)小于電子波長(zhǎng)，此時(shí)可以略去滯后效應(yīng)，從能量守恒的觀點(diǎn)出發(fā)可以導(dǎo)出感應(yīng)電流定律(Ramo’s theorem)[9]如下：

其中iind(t)為諧振腔間隙的感應(yīng)電流，Vgap(t)為間隙電壓，這里的Vgap(t)為縱向電場(chǎng)Ez沿諧振腔軸線的積分，即在式(3)右邊的積分符號(hào)內(nèi)，J(r,t)為電子注的運(yùn)動(dòng)電流，Ecav(r,t)為諧振腔間隙電場(chǎng)，二者都是與位置r有關(guān)的矢量。

根據(jù)等效電路模型有

其中in為感應(yīng)電流的第n次諧波分量，Zn為第n次諧波腔的間隙阻抗QL為諧振腔的有載品質(zhì)因數(shù)，fr為諧振腔的諧振頻率)，若為基波腔則取n=1，若為2次諧波腔則取n=2。

綜上，先在每個(gè)時(shí)間步中計(jì)算各諧振腔的感應(yīng)電流，當(dāng)一個(gè)高頻周期結(jié)束后再將該周期內(nèi)的感應(yīng)電流作傅里葉分解以得到其基波分量和高次諧波分量，再通過(guò)式(4)求得新的間隙電壓，直至各腔間隙電壓達(dá)到穩(wěn)定值為止。

2.4 粒子的推進(jìn)

粒子的推進(jìn)采用Boris[10]旋動(dòng)格式實(shí)現(xiàn)，粒子的推進(jìn)包括粒子的速度和位置推進(jìn)，在Particle-In-Cell 算法中我們采用中心時(shí)間差分，粒子位置的更新是從(n?1)Δt到nΔt，粒子速度的更新是從(n?1/2)Δt 到(n+1/2)Δt ，粒子的位置和速度的更新是相差半個(gè)時(shí)間步長(zhǎng)的，這即是傳統(tǒng)的蛙跳算法。

考慮相對(duì)論效應(yīng)，采用粒子的歸一化動(dòng)量u=γv，粒子的運(yùn)動(dòng)方程可寫為

其中m，q，x，v分別為粒子質(zhì)量、電量、位置和速度；γ為相對(duì)論因子；E，B對(duì)應(yīng)于當(dāng)前粒子位置處的電場(chǎng)強(qiáng)度與磁感應(yīng)強(qiáng)度，通過(guò)插值就可以由網(wǎng)格點(diǎn)上記錄的場(chǎng)值得到其它任意位置處的電磁場(chǎng)分量。

經(jīng)過(guò)推導(dǎo)，粒子速度由時(shí)間步n?1/2的n?1/2u更新為n+1/2的n+1/2u推進(jìn)流程為

式中2t=t? t。粒子的位置推進(jìn)形式為

其中γn由u?根據(jù)來(lái)近似確定。

3 計(jì)算結(jié)果及相關(guān)討論

下面給出以中國(guó)科學(xué)院電子學(xué)研究所研制的S波段50 MW高峰值功率速調(diào)管為算例的研究結(jié)果。該速調(diào)管應(yīng)用于電子直線加速器，工作頻率為2856 MHz，電子注直流電壓為309 kV，電子注電流為376 A，采用均勻磁場(chǎng)聚焦。該速調(diào)管由6個(gè)諧振腔組成，采用單間隙輸出腔，其余腔也均為單間隙諧振腔，輸入腔與輸出腔間隙中心間的距離為0.49 m。本程序計(jì)算時(shí)采用的輸入功率為275 W(與實(shí)驗(yàn)時(shí)一致)，下面給出計(jì)算結(jié)果并對(duì)其進(jìn)行分析討論。

圖5為最后一個(gè)時(shí)間步漂移管中的宏粒子分布，由于計(jì)算的需要本文給模擬區(qū)域加了一段必要的過(guò)渡段，即輸入腔的間隙中心的z軸位置為0.08 m，各間隙中心間的距離均保持不變，縱軸r為宏粒子的徑向位置對(duì)漂移管半徑R0的歸一化值，宏粒子是從發(fā)射面隨機(jī)發(fā)射出來(lái)的。從圖5中可以看到，在末前腔和輸出腔出現(xiàn)了明顯的電子注群聚塊，有很強(qiáng)的密度調(diào)制。

圖6表示在最后一個(gè)時(shí)間步漂移管空間中宏粒子歸一化軸向動(dòng)量γβz=γvz/c 隨軸向距離z的分布，可知在輸出腔(其間隙中心對(duì)應(yīng)圖6中的0.57 m處)中大部分電子由于受到間隙電場(chǎng)的減速作用，電子速度銳減，同時(shí)實(shí)現(xiàn)了從電子動(dòng)能到微波能量的轉(zhuǎn)移。

圖7為在最后一個(gè)時(shí)間步漂移管空間中宏粒子歸一化徑向動(dòng)量γβr=γvr/c 隨軸向距離z的分布。隨著運(yùn)動(dòng)距離的增加，宏粒子徑向動(dòng)量的零散也就越大，特別是到了輸出腔的位置，由于有很強(qiáng)的徑向電場(chǎng)作用，宏粒子的徑向速度波動(dòng)很大。

圖8為宏粒子的相軌跡，它反映了宏粒子在依次受到各個(gè)諧振腔的調(diào)制而逐漸出現(xiàn)群聚的整個(gè)變化過(guò)程，它可以作為速調(diào)管結(jié)構(gòu)及各項(xiàng)參數(shù)選擇是否合理的重要參考之一。從圖8中可以很明顯地看出，在該參數(shù)結(jié)構(gòu)下宏粒子出現(xiàn)了比較理想的群聚。

圖9為電子注歸一化運(yùn)動(dòng)電流的前4次諧波分量隨軸向距離z的變化曲線，隨著軸向距離的增大，電子注的群聚就越明顯，運(yùn)動(dòng)電流的幅值也就越大，其中基波分量幅值最大，其大小也決定了速調(diào)管效率的高低，二階及以上的諧波分量相對(duì)較小。

圖10為同一工作電壓(309.0 kV)下由KLY2D，LSP以及KLY6程序計(jì)算出的傳輸特性曲線，可見(jiàn)KLY2D與KLY6的差異比較大，由后者求出的最大輸出功率高于前者的計(jì)算結(jié)果，而且可以注意到隨著輸入功率的增加KLY6程序預(yù)測(cè)的輸出功率很快就達(dá)到了飽和，其原因是1維模型中電子動(dòng)量都集中在軸向，高估了速調(diào)管中的互作用效率，而由KLY2D計(jì)算得到的輸出功率變化要相對(duì)平緩很多。從圖10中還可以看出，相同輸入功率下KLY2D的計(jì)算結(jié)果要略低于LSP的結(jié)果，其原因?yàn)檩敵銮恢须娮幼⒃谳S向和徑向都會(huì)出現(xiàn)嚴(yán)重的超越和交叉，LSP圓環(huán)模型所依賴的層流性假設(shè)不再成立，而KLY2D則不存在這一問(wèn)題。

圖5 r-z平面宏粒子的分布

圖6 宏粒子的軸向動(dòng)量隨運(yùn)動(dòng)距離z的變化

圖7 宏粒子的徑向動(dòng)量隨運(yùn)動(dòng)距離z的變化

圖8 宏粒子的相軌跡

圖9 電子注歸一化運(yùn)動(dòng)電流的諧波分量隨運(yùn)動(dòng)距離z的變化

圖10 由KLY2D、LSP和KLY6計(jì)算的傳輸特性曲線

當(dāng)輸入功率為275 W時(shí)，由程序計(jì)算得到的輸出功率為51.5 MW，效率為44.3%，而對(duì)應(yīng)樣管實(shí)驗(yàn)測(cè)得的輸出功率為50.9 MW，由此可見(jiàn)，誤差是很小的。

4 結(jié)束語(yǔ)

本文所介紹的粒子模型首次將粒子模擬思想與等效電路端口近似相結(jié)合。與1維模型相比，它考慮了宏粒子的徑向運(yùn)動(dòng)，模擬了大量的宏粒子的運(yùn)動(dòng)，這樣更接近真實(shí)的物理情形。采用端口近似方法，能夠準(zhǔn)確地模擬諧振腔高頻場(chǎng)對(duì)注波互作用系統(tǒng)的影響，同時(shí)比全數(shù)值的3D模型節(jié)約運(yùn)算資源。基于該模型我們開發(fā)出了國(guó)內(nèi)首個(gè)速調(diào)管專用粒子模擬程序KLY2D。其運(yùn)算時(shí)間與通常的粒子模擬程序相比大為縮短，從而極大地提高了數(shù)值模擬的效率。此外，本文的2.5D模型又是一個(gè)時(shí)域的粒子模型，能夠記錄下所有粒子的瞬態(tài)信息，可以模擬互作用系統(tǒng)隨時(shí)間的演化過(guò)程，從而有助于分析諧振腔間隙電壓從開始建立至達(dá)到穩(wěn)定的復(fù)雜過(guò)程，這是傳統(tǒng)的頻域代碼無(wú)法與之相比的。KLY2D運(yùn)行穩(wěn)定，計(jì)算結(jié)果可信度高，與本單位研發(fā)的實(shí)際管型實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比驗(yàn)證，結(jié)果表明程序與實(shí)測(cè)值吻合得很好。

此模型還有很多需要完善的地方，接下來(lái)我們將考慮注波互作用區(qū)域入口處電子注的速度零散，使其更接近真實(shí)的物理情形。同時(shí)，我們將引入群聚腔為雙間隙耦合腔的模型，使其能夠用于擴(kuò)展互作用速調(diào)管的模擬。

[1] Shintake T. Klystron simulation and design using the Field Charge Interaction (FCI) code[J]. Nuclear Instruments and Methods in Physics Research, 1995, Section A(363): 83-89.

[2] 趙鼎. 速調(diào)管非線性注波互作用程序的研究[D]. [博士論文],中國(guó)科學(xué)院電子學(xué)研究所，2007.Zhao D. Research on nonlinear beam-wave interaction program for klystrons[D]. [Ph.D.dissertation], Institute of Electronics, Chinese Academic of Sciences, 2007.

[3] Chernyavskiy I A, Cooke S J, and Vlasov A N, et al.. Parallel simulation of independent beam-tunnels in multiple-beam klystrons using TESLA[J]. IEEE Transactions on Plasma Science, 2008, 36(3): 670-681.

[4] Goplen B, Ludeking L, and Smithe D, et al..User-configurable MAGIC for electromagnetic PIC calculations[J]. Computer Physics Communications, 1995,87(1-2): 54-86.

[5] 鄒峰, 薛謙忠, 劉濮鯤. 大回旋電子注雙磁會(huì)切電子槍的數(shù)值模擬[J]. 電子與信息學(xué)報(bào), 2008, 30(9): 2276-2278.Zou F, Xue Q Z, and Liu P K. Numerical simulation of large orbit gyrotron electron beam double CUSP gun[J]. Journal of Electronics & Information Technology, 2008, 30(9):2276-2278.

[6] Verboncoeur J P, Landon A B, and Glad N T. An object-oriented electromagnetic PIC code[J]. Computer Physics Communications, 1995, (87): 199-211.

[7] Huang Z Y and Pan G W. Universally applicable uniaxial perfect matched layer formulation for explicit and implicit finite difference time domain algorithms[J]. IEEE Transaction on Microwaves, Antennas & Propagation, 2008,2(7): 668-676.

[8] Wu D G, Chen J, and Liu C R. An efficient FDTD method for axially symmetric LWD environments[J]. IEEE Transactions on Geoscience and Remote Sensing, 2008, 46(6):1652-1656.

[9] 丁耀根. 大功率速調(diào)管的理論與計(jì)算模擬[M]. 北京：國(guó)防工業(yè)出版社, 2008: 33-43.Ding Y G. Theory and Computer Simulation of High Power Klystron[M]. Beijing: National Defense Industry Press, 2008:33-43.

[10] Birdsall C K and Langdon A B. Plasma Physics Via Computer Simulation[M]. Bristol: Adam Hilger, 1991: 58-63.