劉紫靜等

摘 要：發(fā)展新一代反應堆物理數(shù)值計算理論和方法的目標是實現(xiàn)三維全堆芯Pin-by-Pin輸運計算，MOC、SN、SPN是其中三種最具有代表意義的全堆芯中子輸運計算方法。該研究介紹了這三種全堆芯中子輸運計算方法的概念及特點，通過調(diào)研國內(nèi)外研究現(xiàn)狀比較分析了三種計算方法的計算精度、計算效率及幾何適應性，總結了三種全堆芯輸運方法的優(yōu)缺點。最終歸納了兩種具備工程應用可行性的全堆芯中子輸運計算方案：（1）二維MOC耦合一維輸運/節(jié)塊法直接進行全堆芯輸運計算；（2）采用二維MOC進行柵格計算，預制參數(shù)截面庫，再運用三維SPN方法進行全堆芯輸運計算。

關鍵詞：全堆芯中子輸運 MOC SN SPN

中圖分類號：TL329+.2 文獻標識碼：A 文章編號：1672-3791（2015）07（c）-0001-04

近年來，我國大力發(fā)展第四代反應堆，多樣的堆芯環(huán)境給核設計帶來了更復雜的中子學特性和非均勻性，并對中子輸運方法的堆芯分析能力提出了更高的要求，主要體現(xiàn)在計算精度、計算效率、處理幾何復雜度等方面[1]。

新一代反應堆物理數(shù)值計算理論和方法的目標是實現(xiàn)三維全堆芯Pin-by-Pin 輸運計算，精確求解任意幾何結構的中子輸運問題，避免任何均勻化近似；具備強大的幾何處理能力，能適應復雜中子輸運問題，可以并行計算[2]。目前，國內(nèi)外開展了大量研究工作，并開發(fā)了多個全堆芯輸運計算程序。主要形成了以下三種最具代表意義的全堆芯Pin-by-Pin 輸運計算方法：①以非均勻柵元為基礎的三維特征線法（Method of Characteristics，MOC）；②以均勻柵元為基礎的三維離散縱標法（Discrete Ordinate Method，SN）。③以均勻柵元為基礎的三維簡化球諧函數(shù)法（Simplified Spherical Harmonics Method，SPN）。

1 特征線法（MOC）

特征線法（MOC）從中子輸運方程的微分形式出發(fā)，沿著中子運行軌跡（特征線）進行積分求解。從理論上說，其整個求解過程不受幾何形狀的限制，無需做均勻化處理，能準確處理強各向異性問題，是全堆芯Pin-by-Pin 輸運計算的理想方法[3]。目前，MOC已逐漸成為許多程序包中的標準求解器。表1列出了部分基于MOC開發(fā)的程序。

MOC方法采用特征線掃描方式求解中子輸運方程，通常采用平源近似，在求解過程中無需保存角中子通量和有關系數(shù)矩陣，比較節(jié)省內(nèi)存[4]。幾何處理大量采用模塊化特征線追蹤方法，重復結構共用一套特征線信息，從而節(jié)省大量存儲空間，提高計算效率。因空間、能群、角度離散，各條特征線的求解相對獨立，所以MOC具有天生的并行性。

但是MOC為取得高計算精度，必需采用細致的特征線和平源近似區(qū)劃分，導致求解復雜問題的速度緩慢，在三維計算中更為突出，需要研究加速并行。其次，MOC在處理三維模型的反射邊界條件時，存在入射線和反射線空間或角度錯位的問題[5]。再次，模塊化特征線追蹤方法僅適用于直線規(guī)則邊界問題，限制了求解問題的任意性。

直接采用MOC方法進行三維復雜中子輸運問題計算，速度緩慢且占用大量存儲空間?？紤]堆芯的非均勻性主要來源于堆芯徑向布置，軸向布置相對均勻，國際上提出了一種“退化”的三維MOC方法，即徑向二維模塊化MOC耦合軸向一維擴散或輸運方法。此方法大大提升了計算效率，同時可以保持較高的計算精度。需要提到的是，若軸向一維方向采用擴散方法，在計算控制棒移動等工況時會出現(xiàn)一定誤差?；诖朔椒ㄩ_發(fā)的程序主要有DeCART[6]、Tiger-3D[7]、MOCHA_2D1D[8]等。

2 離散縱標法（SN）

離散縱標法（SN）將方向自變量離散，每個方向?qū)欢ǖ牧Ⅲw角，所有立體角之和為4π，在特定方向上求解中子輸運方程，并假定立體角內(nèi)的角通量是常數(shù)，采用求積代替積分，得到中子標通量。SN在反應堆屏蔽計算的中子及光子深穿透問題、壓力殼輻照損傷分析中得到了廣泛應用[3]。表2列出了一些基于SN開發(fā)的程序。

SN方法每個離散方向的方程相對獨立，應用迭代法求解時，數(shù)值過程比較簡單，便于編寫通用程序。相比于擴散近似、球諧函數(shù)等近似方法，SN對邊界條件的離散精確自洽，便于處理真空邊界。

但是為得到足夠的精度，需要足夠細密的離散點，這導致計算耗用大量內(nèi)存和時間，尤其對于離散方向數(shù)較大的三維問題，挑戰(zhàn)更為嚴峻。其次，SN通常采用簡化方法來處理幾何和源模型，使其在處理復雜幾何及精確源模擬上存在嚴重缺陷[3]。再次，SN在進行多維問題計算時，由于角度離散本身固有的缺陷，容易出現(xiàn)“射線效應”并引起假散射，導致計算誤差。SN通常采用菱形差分計算模型，這往往導致計算出現(xiàn)與物理不符的負中子通量密度，需要采用置零修正，這也進一步影響了計算精度。

3 簡化球諧函數(shù)法（SPN）

球諧函數(shù)方法（PN）針對輸運方程中的角度自變量做近似，將中子通量密度用球諧函數(shù)展開成N階級數(shù)，將中子輸運方程化成相互耦合微分方程組，并確定級數(shù)中的系數(shù)，再采用迭代策略進行求解。但是當階數(shù)N增大時，尤其對于多維度復雜幾何，PN方法將變得非常復雜和困難[9]。簡化球諧函數(shù)法（SPN）在一維平板幾何下導出球諧函數(shù)方程組，將一維空間微分算子直接換成三維空間微分算子得到三維直角坐標系下的SPN方程組。此方法由PN方法適當化簡而來，減小了方程的復雜性，保持了較高的精度[10]。表3列出了一些SPN程序。

SPN方法僅對空間位置離散，對角度變量的處理是連續(xù)的，具有角度旋轉(zhuǎn)不變性，不存在射線效應；可以采用有限元方法，特別適用于不規(guī)則幾何的非結構網(wǎng)格計算[11]。相比高階PN和SN方法，SPN更容易在滿足精度的前提下節(jié)省計算時間和資源。

但是，SPN方法在柵元均勻化的基礎上實施，需要進行計算網(wǎng)格離散，通過SPH因子修正均勻化截面，預制參數(shù)截面庫，這些均給全堆芯pin-by-pin計算引入了誤差。其次，SPN方法作為低階輸運近似本身具有一定誤差。

4 MOC、SN、SPN對比分析

三種全堆芯輸運計算方法各有優(yōu)勢，但又都存在不同程度的缺點。國際上針對三種方法開展了一些對比研究，主要分析了計算精度、計算效率、幾何適應性等。

日本名古屋大學[12]通過BWR組件計算比較了MOC、SPN、擴散方法的計算精度和計算時間。結果表明：在采用細致網(wǎng)格進行均勻柵元的組件計算時，經(jīng)SPH因子修正后的SPN在計算精度方面與MOC計相當，比擴散方法高2倍；在計算時間方面，比MOC快很多，僅為擴散方法的1.5倍。SPN與基于非均勻柵元的MOC相比，Keff的誤差約為0.02%，反應率的誤差約為0.5%。另一方面，考慮到第四代反應堆中包含不同形式的組件，SPN方法必須具備處理交錯網(wǎng)格的能力，才有工程上的實用性。但目前尚未有人開展此方面的工作。

加拿大蒙特利爾大學以OECD/NEA發(fā)布的NEA3D-TAB-2007中子輸運基準題為例，分析比較了DRAGON程序中三維MOC方法和SN方法對空間參數(shù)劇烈變化問題的計算性能，結果表明：MOC和SN方法均具備較好的計算精度，但SN方法在非均勻性較強的介質(zhì)中具有明顯的射線效應；而MOC方法在平源近似情況下，角度離散和網(wǎng)格劃分數(shù)量的不足將導致計算誤差急劇增加，且在同等的計算條件下，MOC方法的計算時間遠超于SN方法。

5 MOC、SN、SPN聯(lián)合研究

考慮到三種方法各有優(yōu)劣，因此可以取長補短，采用兩種方法聯(lián)合進行全堆芯pin-by-pin計算。目前國際上開展的相關主要研究有：

日本核工程公司、名古屋大學和核燃料工業(yè)公司合作開發(fā)了PWR堆芯燃料管理系統(tǒng)程序包AEGIS/SCOPE2[12]，其中AEGIS為柵格輸運計算程序，采用MOC方法計算得到柵元均勻化截面并進行SPH因子修正；SCOPE2為堆芯程序，采用近似輸運的節(jié)塊SPN方法。AEGIS/SCOPE2采用MOC加SPN的“兩步法”進行全堆芯輸運計算。計算的Keff與實驗值相差152 pcm，棒功率相差1%以內(nèi)。

2005年，歐盟提出了NURESIM計劃[8]，旨在研究更先進、精確的堆芯計算方法，以滿足四代反應堆的設計計算要求。該計劃提出了兩條技術路線：先進蒙卡方法和確定論方法研究。其中先進確定論方法的研究路線為：采用二維MOC進行組件柵格計算，開發(fā)了APOLLO2、HELIOS-2等程序；采用三維節(jié)塊SPN進行堆芯計算，研制了DYN3D程序等；并以蒙卡程序Serpent為基準，驗證了HELIOS-2/DYN3D系統(tǒng)計算分析鈉冷快堆的準確性[10]。

從以上研究情況和國際發(fā)展趨勢來看，采用MOC加SPN的“兩步法”進行全堆芯輸運計算可以得取得足夠的計算精度，計算代價也在可接受范圍內(nèi)。因此，該方案可作為工程應用的選擇之一。

6 結語

該研究介紹了MOC、SN、SPN三種全堆芯中子輸運計算方法的概念及特點，通過國內(nèi)外研究現(xiàn)狀分析對比了三種方法的全堆芯輸運計算能力，得到以下結論：①計算精度：MOC不存在均勻化，精度最高。SN存在射線效應以及計算出負，計算精度略低。SPN需柵元均勻化，在細致網(wǎng)格條件下精度接近輸運計算，但在非均勻性強的區(qū)域誤差偏大。②計算效率：MOC計算時間最長，占用大量內(nèi)存，采用模塊化MOC可提高效率，節(jié)省存儲空間。SN同樣需耗用較大的計算內(nèi)存和時間。SPN計算時間最短，計算效率最高。③幾何適應性：MOC理論上適用于任意幾何，若采用模塊化MOC則限制了幾何適應能力。SN和有限元SPN方法同樣具備很強的幾何適應性。但由于計算耗時、耗內(nèi)存，MOC、SN還多限于二維問題計算。

三種方法取長補短，已形成兩套可行的全堆芯pin-by-pin輸運計算方案，可為今后全堆芯輸運計算方法的工程應用提供參考依據(jù)：①二維MOC耦合一維輸運/節(jié)塊法直接進行全堆芯計算，如韓國DeCART程序。②二維MOC進行柵格計算，預制參數(shù)截面庫，再采用三維SPN進行全堆芯計算，如日本AEGIS/SCOPE2、歐盟APOLLP2/DYN3D程序。

參考文獻

[1] 陳珍平，王電喜，何桃，等.基于 CAD 技術的特征線中子輸運計算程序開發(fā)[J].核科學與工程，2013，32（4）：354-359.

[2] 湯春桃.中子輸運方程特征線解法及嵌入式組件均勻化方法的研究[D].上海：上海交通大學，2009.

[3] 韓靜茹.三維蒙特卡羅—離散縱標雙向耦合屏蔽計算方法研究[D].北京：華北電力大學，2012.

[4] 湯春桃.線性源近似的中子輸運方程特征線解法[J].原子能科學技術，2011，45（12）：1414 -1420.

[5] 柴曉明，姚棟，王侃.基于特征線方法的三維中子輸運程序（Ⅰ）—邊界條件的插值處理[J].核動力工程，2010（2）：11-15.

[6] Cho J Y，Kim K S，Lee C C，et al.Axial SPN and radial MOC coupled whole core transport calculation[J].Journal of Nuclear Science and Technology，2007，44（9）：1156-1171.

[7] 吳文斌，李慶，王侃.基于模塊化射線追蹤的矩陣 MOC方法 （1）—理論研究[J].核動力工程，2014，35（3）：129-133.

[8] 梁亮，吳宏春，鄭友琦，等.基于模塊化特征線方法的二維/一維耦合輸運程序開發(fā)[J].核動力工程，2014（S2）：131-134.

[9] 馬文娟.基于輻射傳輸方程高階球諧近似模型的時域 DOT/FDOT 成像方法研究[D].天津：天津大學，2012.

[10] 曹良志，吳宏春.非結構網(wǎng)格中子輸運方程的球諧函數(shù)解法研究[J].核動力工程，2004，25（5）： 395-398.

[11] 黃細珍.輻射傳遞的球諧函數(shù)法模擬[D].哈爾濱：哈爾濱工業(yè)大學，2010.

[12] Tatsumi M，Ohoka Y，Endo T，et al.sVerification and Validation of AEGIS/SCOPE2：The State-of-the-Art In-Core Fuel Management System for PWRs[C].Xian：18th International Conference on Nuclear Engineering，2010：9-15.