王樹青 廖 鵬 彭 麗 楊 超 印茂偉

1 MeV介質(zhì)壁質(zhì)子直線加速器數(shù)值模擬優(yōu)化

王樹青 廖 鵬 彭 麗 楊 超 印茂偉

(西南科技大學(xué)國防科技學(xué)院 四川綿陽 621010)

采用自主開發(fā)的電磁粒子模擬軟件對中國工程物理研究院研制的1 MeV介質(zhì)壁質(zhì)子加速器進(jìn)行數(shù)值仿真，并對影響介質(zhì)壁加速器的主要參數(shù)進(jìn)行了系統(tǒng)分析，得到了可靠的優(yōu)化方案。模擬結(jié)果顯示:保持加速梯度20 MV/m的情況下，HGI越厚，得到的束流能量峰值越高，束流值越大，需綜合考慮加速器尺寸和所需的束流和束能，以選取合適的HGI厚度;峰值電流隨加速梯度的減小而增大，即在特定要求，可以通過損失部分加速梯度得到更大的輸出電流;輸出束流峰值隨輸入電流大小線性增長，且輸入束流每增加3 mA，對應(yīng)輸出峰值增加0.12 A，增長比率為1:40;輸出束流峰值隨輸入束流半徑的增加而減少，隨法蘭內(nèi)半徑的增大而增加;隨著延時間距的增大，束流先增大后減小，最佳延時為21，21.5，21.7 ns，且電流值為不考慮延時的1.73倍。

介質(zhì)壁加速器 電磁粒子模擬 質(zhì)子治療

介質(zhì)壁質(zhì)子加速器進(jìn)行放射治療的優(yōu)點表現(xiàn)如下:首先，加速器體積大大縮小，采用固態(tài)元器件所建立的介質(zhì)壁加速器的加速梯度理論上可以超過100 MV/m，因此輸出質(zhì)子束的加速器的長度在2 m左右，考慮到束線等因素，加速器的長度可以控制在5 m左右，這對于輸出相同能量的傳統(tǒng)加速器而言(約70 m)，長度大大縮短。其次，目前質(zhì)子治療系統(tǒng)由于加速器系統(tǒng)龐大不能實現(xiàn)移動調(diào)節(jié)，因此采用結(jié)構(gòu)相對簡單的治療機(jī)架實現(xiàn)病灶的定位調(diào)節(jié)，定位機(jī)架的成本超過加速器主體所需費用。而采用介質(zhì)壁質(zhì)子加速器進(jìn)行放射治療，由于加速器主體較小，可以采用通過調(diào)節(jié)加速器主體來實現(xiàn)定位調(diào)節(jié)的方案，因此成本大為降低，初步測算不到目前傳統(tǒng)質(zhì)子加速器治療系統(tǒng)成本的一半。上述優(yōu)點將有利于提高治療效率，降低治療成本，使質(zhì)子放射治療技術(shù)為一般患者所接受，有利于質(zhì)子放射治療技術(shù)的應(yīng)用和推廣。鑒于醫(yī)用加速器誘人的市場前景，國外的相關(guān)企業(yè)，如:LLNL and UC Davis，Compact Particle Acceleration Corporation等，已開展了相關(guān)的技術(shù)和應(yīng)用研究。在國內(nèi)，介質(zhì)壁加速器的研究工作主要集中在中國工程物理研究院流體物理研究所，其主要工作包括:微波介質(zhì)陶瓷材料的制備與加工、SiC光導(dǎo)開關(guān)的研究工作。然而，對介質(zhì)壁加速器系統(tǒng)化數(shù)值仿真的研究報道甚少，筆者對此進(jìn)行了相關(guān)研究。

1 1 MeV介質(zhì)壁質(zhì)子加速器整體布局

1 MeV介質(zhì)壁加速器由注入器(離子源、匹配段)、踢束器(Kicker)、介質(zhì)壁加速段(DWA)、漂移段和束流診斷系統(tǒng)組成，如圖1所示［2］。離子源引出的質(zhì)子束經(jīng)匹配段(LEBT)聚焦后，由Kicker將一個脈沖寬度150 ns左右的束團(tuán)導(dǎo)入介質(zhì)壁加速段，在5 cm的長度內(nèi)加速到約1 MeV。

圖1 介質(zhì)壁加速器總體布局Fig.1 The overall layout of DWA

2 電磁粒子模擬算法

網(wǎng)格生成器用于生成一個非均勻的、用于時域有限差分計算的網(wǎng)格，該網(wǎng)格由多個連續(xù)排列的Yee元胞組成。它能生成直角坐標(biāo)系和圓柱坐標(biāo)系下的三維網(wǎng)格以及直角坐標(biāo)系、圓柱坐標(biāo)系和極坐標(biāo)系下的二維網(wǎng)格。在進(jìn)行模擬之前，只需指定大致的網(wǎng)格步長，并標(biāo)注那些關(guān)鍵點必須位于全網(wǎng)格點位置的區(qū)域，網(wǎng)格生成器就會自動生成非均勻網(wǎng)格。如果某兩個區(qū)域所指定的網(wǎng)格步長不同，則會在兩個區(qū)域連接處的幾個Yee元胞上采用漸變的網(wǎng)格步長，以保證平滑的過渡。

時域電磁場求解器包括中心差分求解器和時偏求解器，用于在每個時間步計算和更新整個模擬區(qū)域內(nèi)的電磁場值，并處理各種電磁場邊界條件。在對電磁場邊界條件的處理上，也采用模塊化設(shè)計的思想，針對每一種邊界都提供功能實現(xiàn)部分和接口調(diào)用部分，并且盡量考慮到今后擴(kuò)展的需要，形成功能實現(xiàn)和接口調(diào)用的文件規(guī)范。另外，作為時域電磁場求解器的輔助部分，還提供了一些小的求解器:泊松求解器和靜磁場求解器用于提供電磁場初始條件，頻域求解器用于諧振腔本征模分析。

粒子求解器用于在每個時間步先計算帶電粒子的作用場，再推動粒子運動，然后計算粒子引起的電荷密度場和電流密度場并處理各種帶電粒子邊界條件。在對帶電粒子邊界條件的處理上，也采用模塊化設(shè)計的思想，形成文件規(guī)范。粒子計算部分涉及多種帶電粒子，粒子既可以通過初始放置產(chǎn)生，又可以通過各種發(fā)射過程產(chǎn)生，其中后者是軟件實現(xiàn)的一個重點。

3 模型及參數(shù)介紹

圖2所示為介質(zhì)壁加速器剖面圖。它由3個HGI組成，中間利用金屬法蘭分隔，并饋入梯形電壓(如圖3所示)，上升沿和下降沿分別為2 ns，半高寬10 ns，頂寬8 ns，輸入電壓300 kV;入口段采用束發(fā)射，且輸入電流 30 mA，束壓40 keV，發(fā)射半徑5 mm，粒子束橫向均勻分布，法蘭內(nèi)半徑8 mm，厚度為2 mm;電介質(zhì)材料(即HGI)內(nèi)半徑17.5 mm，外半徑25 mm，HGI厚度為1.5 cm，介質(zhì)材料均勻;質(zhì)子束注入后21 ns時開始饋入電壓;模擬總時間設(shè)為60 ns，徑向網(wǎng)格0.5 mm，角向60°，軸向1 mm，加速段前端4 cm，后端50 cm。

圖2 介質(zhì)壁加速器建模圖Fig.2 The modeling diagram of DWA

圖3 單個HGI電壓波形Fig.3 A single HGI voltage waveform

4 模擬結(jié)果及分析

圖4為加速段電場分布，由圖4(a)可知，HGI所在區(qū)域軸向電場大體是勻強(qiáng)電場。由圖4(b)可知，HGI所在區(qū)域徑向電場:對于質(zhì)子而言，加速段入口處為橫向聚焦效應(yīng)，出口處呈橫向發(fā)散。

圖5為59 ns時質(zhì)子平均能量沿Z方向的演化，由圖5可知，隨著時間的推移，受到3個HGI同時作用的質(zhì)子束團(tuán)能量能達(dá)到0.935 MeV，且該束團(tuán)逐漸追趕未被加速和小于3個HGI同時作用的質(zhì)子束團(tuán)。

圖6為HGI出口5 cm處面電流大小隨時間演化，由圖6可知，在加速器出口5 cm，31.6 ns時電流值最大，達(dá)到1.2 A。

圖4 加速器加速段電場分布Fig.4 Electric field distribution of accelerator section in DWA

圖5 59 ns時質(zhì)子平均能量沿Z方向的演化Fig.5 Average energy of proton evolutions along the Z direction at 59 ns

圖6 HGI出口5 cm處面電流隨時間演化Fig.6 The plane current with time evolution at 5 cm of HGI export

表1為其它參數(shù)不變且保持加速梯度20 MV/m的情況下，加速器輸出束流大小和能量隨HGI厚度的變化關(guān)系。由表1可知，HGI越厚，加速電壓峰值越大，得到的束流能量峰值越高，束流值越大。這需要綜合考慮加速器尺寸和所需的束流和束能大小，以選取合適的單個HGI厚度。

表2為其它參數(shù)不變，僅通過調(diào)節(jié)HGI厚度來改變加速梯度，分析輸出束流大小隨加速梯度的變化關(guān)系。由表2可知，峰值電流隨加速梯度的減小而增大，即在特定要求，可以通過損失部分加速梯度得到更大的輸出電流。

表3為其它參數(shù)不變，改變輸入束流大小，分析輸出束流與輸入束流大小的關(guān)系。由表3可知，在距離加速段后5 cm位置處測得電流遵循如下規(guī)律:輸出束流峰值隨輸入電流大小線性增長，且輸入束流每增加3 mA，對應(yīng)輸出峰值增加0.12 A，增長比率為1:40。

表4為其它參數(shù)不變，改變輸入束流半徑，分析輸出束流與輸入束流半徑的關(guān)系。由表4可知，在距離加速段后5 cm位置處測得電流遵循如下規(guī)律:輸出束流峰值隨輸入束流半徑的增加而減少。

表5為其它參數(shù)不變，改變3個HGI的分隔法蘭內(nèi)半徑，研究輸出束流與法蘭內(nèi)半徑的變化。在距離加速段后5 cm位置處測得電流遵循如下規(guī)律:輸出束流峰值隨法蘭內(nèi)半徑的增大而增加。

表6為其它參數(shù)不變，根據(jù)圖6顯示的HGI軸向勻強(qiáng)電場的特性，將3個HGI兩端觸發(fā)電壓延時設(shè)為2:1，分析束流隨脈沖時序的關(guān)系。由表6可知，在距離加速段后5 cm位置處測得電流遵循如下規(guī)律:隨著延時間距的增大，束流先增大后減小，最佳延時為21，21.5，21.7 ns，且電流值為不考慮延時的1.73倍。

表1 改變HGI厚度束流和束能變化Table 1 The changes of the beam current and energy by changing the thickness of HGI

表2 改變HGI的加速梯度束流變化Table 2 The changes of the beam current by changing the accelerating gradients of HGI

表3 改變輸入束流大小輸出束流變化Table 3 The changes of the output beam by changing the input beam size

表4 改變HGI內(nèi)半徑束流變化Table 4 The changes of the output beam by changing the inner radius of HGI

表5 改變輸入束流半徑輸出束流變化Table 5 The changes of the output beam by changing the input beam radius

表6 改變3個HGI脈沖時序束流變化Table 6 The changes of the beam current by changing three pulse sequences of HGI

5 結(jié)論

經(jīng)模擬分析得到了1 MeV介質(zhì)壁質(zhì)子加速器主要參數(shù)的優(yōu)化方案:(1)在保持加速梯度20 MV/m的情況下，HGI越厚，得到的束流能量峰值越高，束流值越大，需綜合考慮加速器尺寸和所需的束流和束能，以選取合適的HGI厚度;(2)峰值電流隨加速梯度的減小而增大;(3)輸出束流峰值隨輸入電流大小線性增長，且輸入束流每增加3mA，對應(yīng)輸出峰值增加0.12 A，增長比率為1:40;(4)輸出束流峰值隨輸入束流半徑的增加而減少，隨法蘭內(nèi)半徑的增大而增加;(5)隨著延時間距的增大，束流先增大后減小，最佳延時為21，21.5，21.7 ns，且電流值為不考慮延時的1.73倍。

［1］ MALCOLM B.Some future directions for repetitive pulsed power［J］.IEEE Trans.Plasma Science，2002，30: 262－266.

［2］ NUNALLY W C.Critical component requirement for compact pulse power system architectures［J］.IEEE Trans.Plasma Science，2005，33(4):1262－1267.

［3］ 施衛(wèi)，馬湘蓉，薛紅.半絕緣GaAs光電導(dǎo)開關(guān)的瞬態(tài)熱效應(yīng)［J］.物理學(xué)報，2010，59(8):5700－5705.

［4］ SAMPAYAN S，CAPORASO G，CHEN Y J，et al.Development of a compact radiography accelerator using dielectric accelerator technology［C］.Proc of the Particle Accelerator Conference，2005.716－718.

［5］ 楊超.基于負(fù)氫離子源的全三維PIC/MCC數(shù)值模擬算法研究［D］.成都:電子科技大學(xué)，2012.

［6］ 楊超，劉大剛，周俊，等.一種新型徑向三腔同軸虛陰極振蕩器全三維粒子模擬研究［J］.物理學(xué)報，2011，60(8):084102.

［7］ 王衛(wèi)，鄧建，夏連勝，等.基于大功率激光二極管的光導(dǎo)開關(guān)導(dǎo)通特性［J］.強(qiáng)激光與粒子束，2014，26 (4):045102.

［8］ 劉毅，夏連勝，王衛(wèi)，等.玻璃陶瓷Blumlein線的傳輸特性［J］.強(qiáng)激光與粒子束，2014，26(4):045018.

［9］ 楊超，龍繼東，王平，等.潘寧源放電的全三維電磁粒子模擬/蒙特卡羅碰撞數(shù)值算法研究［J］.物理學(xué)報，2013，62(20): 205207.

Numerical Simulation and Optimization of 1 MeV Dielectric Wall Proton Liner Accelerator

WANG Shu－qing，LIAO Peng，PENG Li，YANG Chao，YIN Mao－wei

(School of National Defense Science and Technology，Southwest University of Science and Technology，Mianyang 621010，Sichuan，China)

In this paper，1 MeV dielectric wall proton accelerator which is developed by China Academy of Engineering Physics was numerically simulated by the self－developed electromagnetic particle simulation software.A reliable optimization scheme was obtained after the system analyzing of its major parameters.The simulation results are the following:With the accelerating gradient of 20 MV/m，the thicker the HGI is，the higher beam energy peak and beam value become.Accelerator’s size，the needed beam and beam energy should be considered synthetically in order to choose a suitable thick of HGI.The peak current increases with the accelerating gradient decreasing，which means larger output current would be obtained by losing part of accelerating gradient in the specific requirements.The output beam peak grows linearly with the input current and when the input current increases by every 3 mA，the output beam peak increases by 0.12 A，in which growth rate is 1:40.The output beam peak decreases with the input beam radius increasing and increases with the flange inner radius increasing.With the rise of delay interval，the beam current increases first and then decreases.The optimal delay are 21，21.5，21.7 ns and the current is 1.73 times the value of the one without considering the delay.

Dielectric wall accelerator;Particle－in－cell simulation;Proton therapy

TL53

A

1671－8755(2015)03－0094－05

2015－03－27

國家自然科學(xué)基金重點基金(11035004);中物院科學(xué)技術(shù)發(fā)展基金(2013A0402018);國家自然科學(xué)基金青年基金(51407169);西南科技大學(xué)博士基金項目(13ZX7106);浙江省科技計劃項目(2013C33073)。

王樹青(1992—)，男，本科生。通信作者:楊超，男，博士，研究方向為模擬粒子理論及算法。E－mail:ychao1983@126.com