基于高品質(zhì)因數(shù)波導(dǎo)型縮比例船形腔的多物理場耦合仿真研究

劉景源,裴士倫,張?zhí)炀?殷治國,汪 洋

(中國原子能科學(xué)研究院 核技術(shù)綜合研究所,北京 102413)

高能強流質(zhì)子加速器在核物理與粒子物理等前沿研究領(lǐng)域、大眾健康和先進能源等國民經(jīng)濟領(lǐng)域,均有著十分廣泛而重要的應(yīng)用[1-2]。近年來,中國原子能科學(xué)研究院提出1套可產(chǎn)生2 GeV、6 MW連續(xù)波質(zhì)子束的強流回旋加速器組合方案,其中包含100 MeV直邊分離扇回旋加速器、800 MeV螺旋分離扇回旋加速器和2 GeV 連續(xù)波(CW)固定場交變梯度加速器(FFA),3臺加速器的高頻腔均工作在44.4 MHz,最終可實現(xiàn)加速器組合全流程連續(xù)束、等時性高效率加速和高平均質(zhì)子束流功率輸出[3]。在分圈式等時性加速器中,當(dāng)束流加速到引出階段,提高螺旋軌道圈的間距,加大留給引出裝置的安裝空間,可有效降低束流被引出裝置阻擋,提高引出效率,實現(xiàn)單圈引出。在回旋式等時性加速器中,在運動方向上,束流的縱向聚焦力較弱,縱向空間電荷效應(yīng)引起的流強閾值與束流圈能量增益的3次方呈正比[4-5],較橫向空間電荷效應(yīng)引起的流強閾值低很多。根據(jù)PSI實驗室高功率質(zhì)子加速器的長年運行調(diào)試經(jīng)驗,受空間電荷影響,束流流強與圈能量增益的3次方呈正比,該實驗室590 MeV回旋加速器先后使用的矩形腔體或歐米伽腔體已較難滿足更高的電磁儲能需求[6-8],經(jīng)過對矩形、歐米伽形和跑道形等多種形狀的波導(dǎo)型高頻腔進行系統(tǒng)研究,提出全流線形的船形高頻腔設(shè)計方案[9]。高頻腔徑向變軌范圍約2.8 m,加速電壓可達到2 MV,諧振頻率44.4 MHz、Q值高達90 000。

為掌握船形高頻腔的加工工藝,采用ANSYS HFSS和ANSYS workbench設(shè)計工作在177.6 MHz的1∶4縮比例船形腔,開展腔體頻率調(diào)諧的仿真研究工作和腔體的樣機研制,同時利用現(xiàn)有的高頻功率源開展相關(guān)高頻功率實驗研究。本文將采用穩(wěn)態(tài)熱分析和Fluent兩種熱仿真軟件對該腔體樣機進行多物理場耦合仿真計算,研究船形腔溫度分布、結(jié)構(gòu)變形和頻率變化。

1 設(shè)計指標(biāo)要求

該船形腔經(jīng)過適當(dāng)調(diào)整,可在不對腔體的高頻性能造成任何影響的情況下,改造為一臺引出電子束能量5～6 MeV、引出靶上平均流強16～20 mA、引出靶上功率90～100 kW的高能電子加速器,圖1所示為加速器基本結(jié)構(gòu)與原理示意圖[10]。

圖1 高能電子加速器基本結(jié)構(gòu)與原理示意圖

綜合考慮2 GeV CW FFA高頻腔的性能指標(biāo)和1∶4縮比例船形腔樣機在電子回旋加速器工程應(yīng)用中的要求,縮比例腔需滿足以下要求:1) 縮比例腔體諧振頻率調(diào)諧范圍要求能達到±170 kHz左右;2) 高功率運行時,損耗在腔體上的RF熱功率約為100 kW,冷卻水若能將這部分熱量全部帶走,冷卻水管出口相對35 ℃入口水溫的溫升要求<10 ℃;3) 在計入大氣壓力、電動缸拉力或壓力、RF熱功率、重力及各種位移約束等的情況下,機械結(jié)構(gòu)變形導(dǎo)致的等效應(yīng)力應(yīng)盡可能低,以方便選用合適形態(tài)的無氧銅板材及不銹鋼材。圖2示出縮比例腔機械設(shè)計總裝結(jié)構(gòu),該腔體主要由主腔體、支撐組件及調(diào)諧電動缸組件等組成,頻率調(diào)諧采用電動缸上下拉壓船形腔體復(fù)合板外壁的方法實現(xiàn)。

圖2 縮比例腔機械設(shè)計總裝結(jié)構(gòu)

2 多物理場耦合仿真理論

常用的波導(dǎo)型高頻腔有矩形波導(dǎo)型和圓柱形波導(dǎo)型,兩者在電磁場分布和特性參數(shù)計算上類似[11]。2 GeV FFA中所使用的高頻腔為矩形波導(dǎo)型,圖3為高頻腔的示意圖,此種類型高頻腔中的模式有兩種,分別為TEmnp模和TMmnp模,m、n及p分別對應(yīng)x、y及z方向上出現(xiàn)電磁場極大值的個數(shù),a、b及d分別為腔體在x、y及z方向的長度。TEmnp?；騎Mmnp模的諧振頻率可表示為:

圖3 矩形波導(dǎo)型高頻腔示意圖

(1)

若束流沿z軸穿過高頻腔并獲得加速,腔體需在其運動方向上提供加速電場?？紤]到高頻腔中用基模加速,則矩形波導(dǎo)型高頻腔的工作模只能選TM110模,即a>b>d或b>a>d,圖4示出TM110模的電磁場分布形式。為使束流能無阻礙穿過高頻腔,在腔上沿x軸(即圓形加速器的半徑方向)開設(shè)了長條形束流通道。腔內(nèi)電場Ez沿x軸和y軸均呈半正弦分布,在靠近腔體沿x軸兩端的電場Ez太低,不能用于加速,由此束流通道在沿x軸方向上的長度要小于a,既要滿足設(shè)計需求又要小于a。由式(1)可知,矩形波導(dǎo)型高頻腔TM110模的頻率與d無關(guān),因此該類型的高頻腔可設(shè)計成長窄型(即d較小),這樣非常有利于增加加速器在粒子束運動方向上安裝空間,實際高頻腔與圖4在形狀上會有所差異,需借助3維計算機軟件來對其進行計算和優(yōu)化[9]。

圖4 矩形波導(dǎo)型高頻腔中TM110模場分布

矩形波導(dǎo)型高頻腔在高功率運行時,在腔體內(nèi)表面會產(chǎn)生大量熱損耗,這些熱量主要通過熱傳導(dǎo)和熱對流被帶走。根據(jù)運行經(jīng)驗,冷卻水出口溫度相對入口的溫升應(yīng)≤10 ℃,船形腔冷水流速一般<3 m/s。流速過高容易引起機械震動,流速過低無法及時帶走腔體產(chǎn)生的熱量,經(jīng)綜合考慮后縮比例腔采用2 m/s[12-13]。整個腔體本體外部布置TU1水管,要求冷卻水能帶走的熱耗散功率約為100 kW,計算分析按照平均35 ℃水的對流換熱系數(shù)來計算。選擇10 mm×10 mm、壁厚1 mm的TU1水冷管進行計算。

傳熱仿真計算中,首先需要確定水冷管的流動狀態(tài),流動狀態(tài)根據(jù)式(2)、(3)計算得到Re。

(2)

(3)

式中:S為非圓斷面的過流截面面積;L為過流截面上流體與水冷管壁接觸的周長;v為介質(zhì)流速;ρ為介質(zhì)密度;η為動力黏度系數(shù)。物性參數(shù)列于表1,可計算出Re=19 875,水冷管處于湍流狀態(tài)。由于104

表1 35 ℃水的物性參數(shù)

(4)

式中:Prw為普朗特數(shù);κw為冷卻介質(zhì)的導(dǎo)熱率。式(4)的使用條件是:當(dāng)流體被加熱時n=0.4,而流體被冷卻時n=0.3;管長與直徑比≥60;壁面溫度和冷卻水主流溫度之差在20～30 ℃[14],對于本文計算腔體,屬于流體被加熱過程,n取0.4。計算得對流換熱系數(shù)約為8 100 W/(m2·K)。

通常矩形波導(dǎo)型高頻腔采用鋁或銅等金屬板加工而成,腔體由薄壁金屬板焊接成型[12]。腔體周圍需設(shè)計固定結(jié)構(gòu)和支撐結(jié)構(gòu),防止矩形波導(dǎo)型高頻腔在受到重力、大氣壓力、調(diào)諧過程中產(chǎn)生的拉和壓力產(chǎn)生變形。2 GeV FFA中所使用的高頻腔采用Y/2態(tài)TU1母材板電子束焊加工。Y/2態(tài)母材的抗拉強度、屈服強度和伸長率分別約為220 MPa、103 MPa和65%[15-16];Y/2態(tài)電子束焊拉伸試樣的斷裂均位于焊縫處,抗拉強度分別約為235 MPa和225 MPa。其一次和二次許用應(yīng)力限值分別為68.7 MPa和206 MPa。

3 多物理場耦合仿真

圖5示出了在ANSYS Workbench中開展“RF-熱-結(jié)構(gòu)-RF”全自洽閉環(huán)多物理場有限元分析的具體流程,使用Steady-State Thermal模塊開展固熱分析。導(dǎo)入到Geometry中的Solidworks機械模型包含有兩種類型的實體部件:金屬和非金屬實體部件。在HFSS中,與RF場相關(guān)的非金屬實體和金屬實體與RF場計算相關(guān),只計算這些區(qū)域以使整個模擬計算可順利進行并效率更高。為了使計算分析更接近實際情況,RF計算完成后需調(diào)節(jié)腔體內(nèi)部儲能到4 J,由此可使高頻功率損耗達到設(shè)計要求的100 kW。在固熱分析中,RF損耗將其設(shè)置為熱負載。在這種情況下,將對流系數(shù)和壁面冷卻水溫度作為邊界條件,計算得到腔體溫度分布。在靜態(tài)結(jié)構(gòu)分析中,將腔體溫度分布結(jié)果作為負載導(dǎo)入靜態(tài)結(jié)構(gòu)分析模塊,設(shè)置重力、大氣壓力、固定支撐等約束,計算得到整個腔體機械模型的變形及應(yīng)力分布。為了評估結(jié)構(gòu)變形對RF場的影響,選擇對應(yīng)的實體設(shè)置并設(shè)置結(jié)構(gòu)變形反饋功能。將計算得到的變形反饋給HFSS,然后再評估機械變形對腔體RF性能的影響。多次迭代上述過程,可獲得穩(wěn)定的計算結(jié)果。但在大多數(shù)情況下,一次或兩次迭代即可。

圖5 基于Steady-State Thermal模塊的船形腔多物理場仿真流程圖

固熱分析的方法在設(shè)置邊界條件時,對水冷管各處做了溫度分布均勻、對流換熱系數(shù)相等的假設(shè)。使用Fluent開展流固熱分析時,通過設(shè)置運動參數(shù)、物性參數(shù)和計算方程,仿真過程中逐步計算得到水冷管溫度,這樣的結(jié)果更符合實際情況。圖6示出了基于Fluent的多物理場仿真分析流程,使用流體仿真模塊Fluent模塊代替Steady-State Thermal模塊,該方法較固熱分析略微復(fù)雜,對模型的要求較高,在不影響整體分析結(jié)果的基礎(chǔ)上對導(dǎo)入的模型進行一定的調(diào)整,修補不適用于Fluent的模型缺陷。綜合考慮計算機性能和計算時間等條件,在進行網(wǎng)格劃分和流體傳熱仿真時,不考慮船形腔的支撐組件,將HFSS計算得到的腔體內(nèi)壁功率設(shè)置為Fluent的熱負載。設(shè)置冷卻水進口流速和溫度,根據(jù)式(3)得到Re,求解模型采用k-ε模型。

圖6 基于Fluent的船形腔多物理場仿真流程圖

4 計算結(jié)果及討論

表2列出HFSS仿真計算的腔體功率損耗分布,在諧振頻率177.848 MHz下腔體無載Q值約為44 157.7。

表2 功率損耗分布

4.1 溫度分布特性

圖7分別示出穩(wěn)態(tài)熱分析和Fluent仿真的腔體溫度分布,穩(wěn)態(tài)熱分析顯示腔體最高溫度為77 ℃,Fluent仿真計算得到腔體最高溫度為82.6 ℃。兩種方法均顯示腔體最高溫度分布出現(xiàn)在鼻錐兩端根部區(qū)域。對比可發(fā)現(xiàn),穩(wěn)態(tài)熱分析結(jié)果中腔體溫度沿y-z平面對稱分布,而Fluent仿真計算結(jié)果顯示,腔體出口附近溫度比進口附近高約4.7 ℃。腔體熱量可被冷卻水帶走,且冷卻水進出口溫升小于10 ℃,滿足設(shè)計指標(biāo)要求。

a——穩(wěn)態(tài)熱仿真;b——Fluent仿真

在Fluent計算結(jié)果中,腔體總流量為Q=6.713 L/s,進出口平均溫升為ΔT=3.6 ℃,計算得到腔體總功率為100.89 kW。在Fluent的計算結(jié)果中提取水冷管壁面溫度Tw、熱流密度q及管內(nèi)主流溫度Tf,水冷管壁面溫度根據(jù)牛頓冷卻公式(式(5))計算得到腔體管道對流換熱系數(shù)為7 532～8 573 W/(m2·K),與式(2)～(4)計算的理論值基本一致。

q=h(Tw-Tf)

(5)

4.2 腔體變形及頻率調(diào)節(jié)

為了獲得±170 kHz的頻率調(diào)諧范圍,分別在船形腔上部和下部各設(shè)置3個電動缸。利用電動缸對船形腔體銅壁施加拉力或壓力后,會使其內(nèi)部真空部分的體積發(fā)生變化。因為拉力或壓力均施加在船形腔內(nèi)RF磁場較強、電場較弱的區(qū)域附近,所以拉腔體會使其頻率降低、壓腔體會使其頻率升高。圖8、9分別示出了利用調(diào)諧裝置將船形腔頻率調(diào)變-164 kHz和+176 kHz時得到的機械變形及應(yīng)力分布,最高應(yīng)力出現(xiàn)在電動缸復(fù)合板位置。當(dāng)腔體頻率調(diào)變-164 kHz時,主腔體最大變形為2.64 mm,最高應(yīng)力為187.22 MPa(二次應(yīng)力),此時船形腔上部和下部的電動缸對應(yīng)的總拉力均為30 000 N。當(dāng)船形腔頻率調(diào)變+176 kHz時,主腔體最大變形為2.96 mm,最高應(yīng)力為187.16 MPa(二次應(yīng)力),此時船形腔上部和下部電動缸對應(yīng)的總壓力均為14 000 N。最高應(yīng)力低于TU1二次應(yīng)力3[P]=R態(tài)204 MPa。

a——機械變形分布;b——應(yīng)力分布

a——機械變形分布;b——應(yīng)力分布

在將船形腔諧振頻率調(diào)高或調(diào)低幾乎相同的量時,總壓力僅約為總拉力的1/2,這是因為船形腔正常工作時內(nèi)部為真空狀態(tài),主腔體銅壁外表面受到1個大氣壓的壓力作用。無論是利用電動缸拉還是壓腔體,所造成的變形都會在拉力或壓力施加點附近造成較高的應(yīng)力,實際加工腔體時的焊縫應(yīng)盡最大可能遠離這些區(qū)域。

圖10、11分別示出將電動缸設(shè)置為拉力30 000 N和壓力14 000 N,基于Fluent的多物理場有限元仿真得到的溫度分布腔體機械變形及應(yīng)力分布,最大變形和最高應(yīng)力都分布在上下電動缸復(fù)合板區(qū)域。電動缸總拉力為30 000 N時,船形腔主腔體最大形變?yōu)?.61 mm,最高應(yīng)力為187.31 MPa。電動缸總壓力為14 000 N時,船形腔主腔體最大形變?yōu)?.50 mm,最高應(yīng)力為187.18 MPa,腔體調(diào)諧范圍為329 kHz。

a——機械變形分布;b——應(yīng)力分布

對比圖8和圖11可得:相同參數(shù)下,兩種仿真方式得到的腔體頻率變化范圍相差11 kHz;基于Fluent的多物理場仿真計算得到的變形小于基于穩(wěn)態(tài)熱分析的結(jié)果;前者計算得到的腔體受到的最大應(yīng)力大于后者的結(jié)果。表3、4分別給出最大變形和最高應(yīng)力。造成這些差異的原因除了計算誤差和溫度分布條件不一致外,還有靜態(tài)結(jié)構(gòu)分析約束條件的差異。受計算機性能限制,基于Fluent的仿真流程需要減少對支撐組件的靜力學(xué)計算,只能以腔體固定板端面作為固定面。而約束條件的改變,造成原本可以被腔體、固定板和支撐組件共同分攤的應(yīng)力,只能由腔體來承擔(dān)。

表3 腔體最大變形

表4 腔體最高應(yīng)力

a——機械變形分布;b——應(yīng)力分布

5 縮比例船形腔冷測

圖12示出腔體在15 ℃的環(huán)境溫度下,通過分子泵將腔體真空度提高到8.1×10-8mbar,電機驅(qū)動位置從-2.5 mm變化到2.5 mm測得的腔體頻率變化曲線。腔體頻率變化范圍為177.561～177.935 MHz,腔體加工滿足設(shè)計需求。在大氣狀態(tài)下,利用2個耦合度極小的小耦合環(huán),通過S21掃頻曲線對無載品質(zhì)因數(shù)進行測量,測得Qs≈42 314,約為模擬計算值44 157.7的95.8%。

圖12 腔體冷測頻率曲線

6 結(jié)論

本文利用ANSYS workbench開展了基于穩(wěn)態(tài)和Fluent熱分析的船形腔多物理耦合仿真計算。結(jié)果表明,腔體最高溫度分布在鼻錐兩端根部區(qū)域,主腔體表面溫度分布均勻,水冷管布置合理,可帶走腔體耗散能量。兩種熱分析結(jié)果對比顯示,Fluent仿真結(jié)果更接近實際情況。靜態(tài)結(jié)構(gòu)分析表明應(yīng)力主要集中在上下電動缸復(fù)合板區(qū)域,最高應(yīng)力低于TU1二次應(yīng)力的許用要求?；诜€(wěn)態(tài)熱分析的多物理耦合仿真和實際冷測均顯示腔體頻率變化范圍滿足設(shè)計需求。受計算機算力影響,基于Fluent的多物理耦合仿真對模型優(yōu)化,無法完全模擬腔體的受力約束情況,導(dǎo)致與穩(wěn)態(tài)熱分析相比,腔體頻率變化范圍的結(jié)果相差11 kHz。