cosRMC復(fù)雜曲面開(kāi)發(fā)及在聚變中子學(xué)的應(yīng)用

劉仕倡，申靖文，嚴(yán)伊蔓，杜 華，孫業(yè)帥，秦 瑤，鄭 俞， 盧 棚，李 夏，張 顯，劉松林，陳義學(xué)

(1.華北電力大學(xué) 核科學(xué)與工程學(xué)院，北京 102206； 2.國(guó)家電投集團(tuán)科學(xué)技術(shù)研究院有限公司 核電軟件技術(shù)中心，北京 102209； 3.中國(guó)科學(xué)院 等離子體物理研究所，安徽 合肥 230031)

中國(guó)聚變工程實(shí)驗(yàn)堆(CFETR)[1]是我國(guó)于2011年開(kāi)始設(shè)計(jì)研究的超導(dǎo)托卡馬克反應(yīng)堆，在DEMO級(jí)別聚變電站與ITER間起到了過(guò)渡銜接的作用。托卡馬克型的聚變堆堆芯主要包括：第一壁、包層、偏濾器、真空室、輻射屏蔽、生物屏蔽等。包層是聚變堆堆芯的一個(gè)重要組成部分，在強(qiáng)中子輻照下，包層主要有氚增殖、能量轉(zhuǎn)換及輻射屏蔽等作用，因此，包層中子學(xué)的設(shè)計(jì)與聚變堆中的大量問(wèn)題均相關(guān)。此外，托卡馬克裝置具有非常復(fù)雜的幾何結(jié)構(gòu)，包含各種結(jié)合體并且結(jié)構(gòu)十分精細(xì)，而蒙特卡羅方法幾何描述精細(xì)的特點(diǎn)無(wú)疑非常適宜于解決這樣的問(wèn)題[2]。但如果用于模擬計(jì)算的蒙特卡羅程序不具備高次曲面的描述能力，僅使用一階或二階曲面并不能滿足幾何建模精細(xì)化的需求。那么對(duì)于這類復(fù)雜幾何模型的建立就僅能采取幾何近似的方式，如可將托卡馬克裝置的環(huán)形真空室的圓環(huán)面近似為很多段直線，不僅耗時(shí)耗力，還會(huì)導(dǎo)致精度損失，不利于聚變堆的準(zhǔn)確模擬。國(guó)際上具有處理聚變堆復(fù)雜幾何能力的蒙特卡羅程序有美國(guó)Los Alamos實(shí)驗(yàn)室的MCNP[3]、日本原子能機(jī)構(gòu)的MVP[4]、法國(guó)原子能委員會(huì)的TRIPOLI[5]等。國(guó)內(nèi)方面，中國(guó)科學(xué)院核能安全技術(shù)研究所的SuperMC[6]在磁約束核聚變裝置的蒙特卡羅模擬方面進(jìn)行了大量研究，北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所的JMCT[7]也進(jìn)行了一些研究。其中，MCNP是應(yīng)用最為廣泛的經(jīng)典蒙特卡羅程序，其復(fù)雜曲面處理方法被廣泛應(yīng)用于國(guó)內(nèi)外蒙特卡羅程序，然而對(duì)該方法的數(shù)學(xué)模型和計(jì)算原理進(jìn)行系統(tǒng)推導(dǎo)和闡述的研究報(bào)道較少。

本文基于自主化蒙特卡羅程序cosRMC，研究及闡述蒙特卡羅復(fù)雜曲面建模的數(shù)學(xué)模型和計(jì)算方法，并對(duì)原有的幾何模塊進(jìn)行拓展，使其具備對(duì)復(fù)雜幾何模型的描述能力，最后將其應(yīng)用于CFETR三維模型的聚變中子學(xué)分析中，對(duì)聚變堆包層設(shè)計(jì)時(shí)的3個(gè)關(guān)鍵的中子學(xué)參數(shù)進(jìn)行計(jì)算分析。

1 cosRMC復(fù)雜幾何建模功能

1.1 cosRMC及其幾何建模功能介紹

cosRMC是由清華大學(xué)工程物理系反應(yīng)堆工程計(jì)算分析實(shí)驗(yàn)室(REAL團(tuán)隊(duì))與國(guó)家電投集團(tuán)科學(xué)技術(shù)研究院有限公司核電軟件技術(shù)中心聯(lián)合研發(fā)的以蒙特卡羅粒子輸運(yùn)方法為理論基礎(chǔ)的自主化蒙特卡羅程序[8]。主要用于核工程領(lǐng)域中固定源計(jì)算、臨界計(jì)算、燃耗計(jì)算、瞬態(tài)計(jì)算和均勻化并群計(jì)算等。

cosRMC采用文本建模方式，基于構(gòu)造幾何實(shí)體(CSG)對(duì)系統(tǒng)幾何結(jié)構(gòu)進(jìn)行描述。即通過(guò)對(duì)交界面不同指向的區(qū)域進(jìn)行布爾運(yùn)算，如交、并和非等運(yùn)算，構(gòu)成對(duì)體的描述。目前，cosRMC程序可進(jìn)行描述的面包括一般平面、任意球面、軸心位于坐標(biāo)軸上或平行于坐標(biāo)軸的圓柱面或圓錐面等標(biāo)準(zhǔn)形式的曲面，所有曲面的曲面方程均使用解析方程表示方式。在本文工作之前，cosRMC無(wú)法實(shí)現(xiàn)對(duì)一些復(fù)雜曲面的描述，如任意位置的圓柱面、圓錐面、橢圓面和圓環(huán)面等。

1.2 復(fù)雜幾何建模功能實(shí)現(xiàn)

粒子在蒙特卡羅程序幾何中的輸運(yùn)需確定粒子沿飛行方向與曲面的交點(diǎn)，以圓環(huán)面為例，通過(guò)解析方法簡(jiǎn)化為求解一個(gè)四次多項(xiàng)式方程解析解的問(wèn)題。

1) 中子軌跡與復(fù)雜曲面的聯(lián)立求解

(1)

進(jìn)一步在粒子輸運(yùn)問(wèn)題中，在射線簇{x=ξ+aX,y=η+BX,z=ζ+γX;-∞

{[(1-β2)+ρβ2]X2+[2(aξ+γζ)+2ρβη-

2ρβy0]+(ξ2+ζ2+ρη2-2ρηy0+B0)}2=

A0[(1-β2)X2+2(aξ+γζ)X+ξ2+ζ2]

(2)

式(2)可轉(zhuǎn)化為：

(GX2+MX+U)2=A0(FX2+LX+T)

(3)

其中，F(xiàn)=1-β2，G=F+ρβ2，L=2(aξ+γζ)，M=L+2ρβ(η-y0)，T=ξ2+ζ2，U=T+ρη(η-2y0)+B0。

所以，射線簇{x=ξ+aX,y=η+BX,z=ζ+γX;X>0}與圓環(huán)面求交點(diǎn)(x,y,z)的問(wèn)題可轉(zhuǎn)化成一個(gè)求解式(4)所示四次方程的正根X的問(wèn)題，僅需判斷射線簇是與圓環(huán)面的內(nèi)表面還是外表面存在交點(diǎn)即可。

X4+BX3+CX2+DX+E=0

(4)

2) 復(fù)雜曲面多項(xiàng)式方程的求解

由于舍入誤差的存在，使得利用計(jì)算機(jī)程序求解形如式(4)的四次多項(xiàng)式方程均非常困難。Maisonnier和Everett提出了一種求解“射線簇與三維或四維曲面與交點(diǎn)”的單精度四次算法[9]，該算法可分為兩步，第1步是將粒子穿過(guò)圓環(huán)面時(shí)四次多項(xiàng)式方程的求解轉(zhuǎn)化為三次方程的求解。

首先討論三階曲面的解析方程，可利用泰勒展開(kāi)求解：

f(x)=d+cx+bx2+x3=

(5)

對(duì)式(5)求導(dǎo)代換后可得到新的三次解析方程：

g(x)=y3+py+q

(6)

Δ3=-4(27)W

(7)

同樣對(duì)四階曲面的解析方程也進(jìn)行泰勒展開(kāi)：

F(X)=E+DX+CX2+BX3+X4=

(8)

通過(guò)推導(dǎo)可得出，四階曲面方程的判別式為Δ4=Δ3。根據(jù)判別式條件就變成了一個(gè)類似于三階曲面方程求解的過(guò)程。

第2步是根據(jù)多項(xiàng)式系數(shù)判定來(lái)確定是否存在任意正實(shí)根，主要包括3種情況：當(dāng)Δ3>0(W<0)時(shí)，有3個(gè)不同的實(shí)根；當(dāng)Δ3=0(W=0)時(shí)，有3個(gè)含有重復(fù)根的實(shí)根；當(dāng)Δ3<0(W>0)時(shí)，有單個(gè)實(shí)根和1對(duì)共軛復(fù)根。

根據(jù)以上理論方法，本文在cosRMC原有的幾何模塊中新添加了橢圓面、雙曲面、拋物面、圓環(huán)面等復(fù)雜曲面，依舊采用助記名的方式存儲(chǔ)，如表1所列。

1.3 功能測(cè)試

聚變堆概念設(shè)計(jì)堆PPCS(power plant conceptual study)是2000年于歐洲啟動(dòng)的最為系統(tǒng)的聚變堆研究計(jì)劃，目的是幫助評(píng)估聚變能形式和確立歐洲聚變計(jì)劃的可持續(xù)性和優(yōu)先權(quán)，屬于世界上目前較具規(guī)模和代表性的磁約束聚變動(dòng)力反應(yīng)堆概念設(shè)計(jì)系列之一，如圖1所示。作為典型的磁約束聚變堆概念設(shè)計(jì)系列之一，PPCS具有非常復(fù)雜的幾何機(jī)構(gòu)，在使用蒙特卡羅方法對(duì)其進(jìn)行聚變中子學(xué)分析時(shí)，需構(gòu)建諸多包含復(fù)雜曲面的幾何體，以及包層模塊重復(fù)結(jié)構(gòu)的描述，如圖1c中標(biāo)注所示。

表1 新增cosRMC復(fù)雜曲面卡Table 1 New cosRMC complex surface card

a——PPCS的cosRMC模型；b——PPCS的MCNP模型；c——PPCS模型中選取的柵元 圖1 PPCS的cosRMC及MCNP的模型示意圖Fig.1 Model diagram of PPCS with cosRMC and MCNP

利用cosRMC程序?qū)PCS進(jìn)行建模(圖1)后進(jìn)行模擬計(jì)算，中子源采用14 MeV的各向同性點(diǎn)源。MCNP和cosRMC模擬的總粒子數(shù)均為1×108。選取了18個(gè)有高次曲面的柵元，圖1c中標(biāo)注的數(shù)字即為選取柵元的柵元號(hào)。統(tǒng)計(jì)其體通量，與MCNP的計(jì)算結(jié)果進(jìn)行比較，令相對(duì)統(tǒng)計(jì)誤差為：

(9)

σMCNP和σcosRMC為MCNP和cosRMC的計(jì)數(shù)器相對(duì)統(tǒng)計(jì)誤差。柵元通量的對(duì)比如圖2所示?？煽闯觯琧osRMC程序?qū)λ袔缀误w上包含高次曲面的體通量的計(jì)算結(jié)果，與MCNP的計(jì)算結(jié)果吻合良好，相對(duì)偏差幾乎均小于10%。且cosRMC單核計(jì)算時(shí)間為250.980 6 min，MCNP單核計(jì)算時(shí)間為258.74 min，計(jì)算效率相當(dāng)。說(shuō)明了cosRMC程序中實(shí)現(xiàn)的復(fù)雜幾何建模功能的正確性。

圖2 PPCS模型cosRMC與MCNP模擬計(jì)算結(jié)果比較Fig.2 Comparison of simulation results of PPCS between cosRMC and MCNP

2 CFETR建模與聚變中子學(xué)分析

2.1 CFETR建模

CFETR的三維模型體積大、幾何結(jié)構(gòu)復(fù)雜，托卡馬克裝置中包含環(huán)形真空室、環(huán)向線圈、極向線圈等結(jié)構(gòu)部件使得在使用蒙特卡羅方法對(duì)其進(jìn)行聚變中子學(xué)分析時(shí)，需構(gòu)建諸多包含四階曲面的幾何體。使用拓展了復(fù)雜幾何建模功能的cosRMC程序?qū)FETR三維模型進(jìn)行精細(xì)化建模，如圖3所示，與MCNP程序構(gòu)建的模型圖進(jìn)行對(duì)比，如圖4所示。從模型圖可看出cosRMC程序描述的CFETR三維模型的正確性。

圖3 CFETR的cosRMC模型徑向與縱向剖面圖Fig.3 Radial and longitudinal profiles of cosRMC model of CFETR

圖4 CFETR的MCNP模型縱向剖面圖Fig.4 Longitudinal profile of MCNP model of CFETR

進(jìn)一步對(duì)CFETR的包層中子學(xué)中的3個(gè)重要參數(shù)進(jìn)行計(jì)算分析。CFETR的包層模塊如圖5所示，沿環(huán)向上和極向上進(jìn)行分割，分為內(nèi)包層模塊和外包層模塊。中子源采用氘氚聚變中子源，分布于5個(gè)柵元中，中子源強(qiáng)分別為1%、6%、13%、30%和50%。源中子能量為14 MeV的高斯聚變譜，方向各向同性。CFETR的氚增殖比(TBR)、中子壁載荷(NWL)和核熱沉積計(jì)算過(guò)程中，聚變功率采用第1階段的運(yùn)行功率200 MW。MCNP和cosRMC模擬的總粒子數(shù)均為6.9×106。MCNP和cosRMC均采用FENDL-2.1核截面數(shù)據(jù)庫(kù)進(jìn)行計(jì)算。

圖5 CFETR包層結(jié)構(gòu)編號(hào)及位置示意圖Fig.5 Diagram of CFETR blanket structure number and location

2.2 TBR

氘氚聚變反應(yīng)堆能否投入商業(yè)運(yùn)行的首要條件是能實(shí)現(xiàn)氚的自持。聚變中子學(xué)分析計(jì)算時(shí)，增殖材料中6Li、7Li的產(chǎn)氚率為：

T6=NLi?σ6Li(n,α)(r,E)φ(r,E)dEdr

(10)

T7=NLi?σ7Li(n,n′α)(r,E)φ(r,E)dEdr

(11)

其中，σ6Li(n,α)和σ7Li(n,n′α)分別為6Li的(n,α)反應(yīng)截面和7Li的(n,n′α)反應(yīng)截面。那么，堆芯等離子體平均1次氘氚聚變反應(yīng)在包層內(nèi)增殖的TBR為：

(12)

其中，Sn為聚變中子源強(qiáng)。利用cosRMC對(duì)其增殖包層各模塊的氚增殖比進(jìn)行計(jì)算，與MCNP的計(jì)算結(jié)果進(jìn)行比較(表2)?？煽闯?，兩程序?qū)Ω髟鲋嘲鼘与霸鲋潮鹊挠?jì)算結(jié)果吻合良好，相對(duì)偏差均小于1%。并且，通過(guò)模擬計(jì)算可發(fā)現(xiàn)，增殖包層總的TBR為1.2，能滿足氘氚聚變反應(yīng)堆運(yùn)行時(shí)的氚自持需求。

2.3 NWL

NWL即加載到第一壁表面(W-PFC)單位面積上的平均聚變中子功率。之所以要計(jì)算NWL，是因在托卡馬克裝置中第一壁是面向等離子體的部件，在包層中起著重要作用，它的安全直接影響包層模塊的穩(wěn)定性，其計(jì)算方法如式(13)所示，其中A為表面積。

(13)

表3列出利用cosRMC計(jì)算的CFETR第一壁表面的NWL分布，各包層模塊的結(jié)果與MCNP符合良好，相對(duì)偏差均在1%以內(nèi)。進(jìn)一步計(jì)算可得到第一壁表面的平均NWL為0.305 1 MW/m2，NWL峰值出現(xiàn)在位于赤道平面的BLK3包層。

表2 cosRMC和MCNP計(jì)算的CFETR增殖包層各模塊TBR的比較Table 2 Comparison of TBR in each module of CFETR blanket

表3 CFETR第一壁表面的NWL分布計(jì)算結(jié)果比較Table 3 Comparison of NWL in each module of CFETR first wall surface

2.4 核熱沉積

托卡馬克裝置中用于約束等離子體的磁鐵結(jié)構(gòu)通常是由帶有不銹鋼包層的銅制成的，以提供抵抗變形的強(qiáng)度。銅線圈具有極好的電導(dǎo)率，允許通過(guò)的電流強(qiáng)度高達(dá)數(shù)十MA。但同時(shí)會(huì)在銅線圈上產(chǎn)生較大的熱沉積。限制等離子體放電的因素之一即是當(dāng)高強(qiáng)度的電流通過(guò)磁鐵后，磁鐵溫度的升高，應(yīng)使核熱沉積最小化。因此，研究磁鐵線圈系統(tǒng)上產(chǎn)生的核熱沉積是十分必要的，其計(jì)算方法如式(14)所示。

(14)

其中：H(E)為熱數(shù)；σt(E)為總截面；ρ為核數(shù)密度；m為柵元質(zhì)量。

表4列出利用cosRMC程序計(jì)算的CFETR的BLK1包層內(nèi)各部件的核熱沉積，兩程序模擬計(jì)算結(jié)果吻合良好，相對(duì)偏差均小于5%。同一包層上不同部件的核熱沉積不同是因?yàn)楦鞑考幬恢貌煌?，中子通量密度不同，受到的中子輻照?qiáng)度也不同。

3 結(jié)論

本文研究及推導(dǎo)了蒙特卡羅復(fù)雜曲面建模的數(shù)學(xué)模型和計(jì)算方法，并在自主化蒙特卡羅程序cosRMC原有的幾何模塊上拓展了復(fù)雜曲面建模功能。通過(guò)構(gòu)建聚變堆概念設(shè)計(jì)PPCS模型并進(jìn)行計(jì)算，柵元通量與MCNP的相對(duì)偏差小于10%，說(shuō)明了該功能的正確性。

表4 CFETR的BLK1包層各部件核熱沉積計(jì)算結(jié)果比較Table 4 Comparison of nuclear thermal deposition in different components of BLK1 blanket

進(jìn)一步利用cosRMC對(duì)CFETR進(jìn)行三維精細(xì)化建模，并對(duì)包層模塊設(shè)計(jì)中的3個(gè)重要中子學(xué)參數(shù)進(jìn)行計(jì)算分析。結(jié)果表明，cosRMC與MCNP計(jì)算的各增殖包層TBR的相對(duì)偏差均小于1%。且增殖包層總的TBR為1.2，能滿足氘氚聚變堆運(yùn)行時(shí)的氚自持需求。cosRMC與MCNP計(jì)算的NWL相對(duì)偏差均在1%以內(nèi)，NWL最大值出現(xiàn)在位于赤道平面的BLK3包層。cosRMC與MCNP計(jì)算的BLK1包層內(nèi)各部件的核熱沉積的相對(duì)偏差均小于5%。研究結(jié)果證明了cosRMC應(yīng)用于聚變堆包層中子學(xué)分析的正確性和有效性。CFETR中子學(xué)參數(shù)的計(jì)算分析，也為其設(shè)計(jì)和優(yōu)化提供了參考。

感謝清華大學(xué)工程物理系核能科學(xué)與工程管理研究所反應(yīng)堆工程計(jì)算分析實(shí)驗(yàn)室(REAL團(tuán)隊(duì))在程序使用過(guò)程中提供的技術(shù)支持和幫助。