本發(fā)明涉及核反應(yīng)堆堆芯設(shè)計(jì)和核安全技術(shù)領(lǐng)域，具體涉及一種精確計(jì)算核反應(yīng)堆內(nèi)時(shí)空中子分布的方法。

背景技術(shù)：

在核電廠的設(shè)計(jì)中，核反應(yīng)堆的設(shè)計(jì)一直處于最核心的地位，而在核反應(yīng)堆的設(shè)計(jì)中，精確模擬核反應(yīng)堆在瞬態(tài)過程中時(shí)在不同時(shí)間點(diǎn)以及不同位置處的中子分布一直是一個(gè)巨大的挑戰(zhàn)。目前數(shù)值模擬是核反應(yīng)堆設(shè)計(jì)的最主要的手段，因此為了精確模擬核反應(yīng)堆在瞬態(tài)過程中時(shí)在不同時(shí)間點(diǎn)以及不同位置處的中子分布，數(shù)值模擬計(jì)算的方法便顯得尤為重要。

近年來隨著計(jì)算機(jī)計(jì)算能力的大幅度提高，在數(shù)值模擬中采用高保真時(shí)空中子動(dòng)力學(xué)模型已經(jīng)成為國(guó)際的主流。高保真時(shí)空中子動(dòng)力學(xué)能夠保證理論模型和計(jì)算模型不失真，計(jì)算所得解為真解，計(jì)算分辨率和計(jì)算精度高，計(jì)算結(jié)果精確，能夠達(dá)到高度保真的要求。要實(shí)現(xiàn)高保真時(shí)空中子動(dòng)力學(xué)的計(jì)算，就必須采用一步法對(duì)核反應(yīng)堆直接進(jìn)行幾何建模，采用多能群結(jié)構(gòu)，通過精確求解中子輸運(yùn)方程，在計(jì)算過程中不引入幾何近似以及傳統(tǒng)兩步法均勻化帶來的誤差。通過高保真時(shí)空中子動(dòng)力學(xué)的模擬可大大提高核反應(yīng)堆設(shè)計(jì)的精確性，大幅度減少實(shí)驗(yàn)量，提高反應(yīng)堆的安全性和經(jīng)濟(jì)性，對(duì)核反應(yīng)堆改進(jìn)設(shè)計(jì)和新型反應(yīng)堆設(shè)計(jì)均具有重大的理論價(jià)值和應(yīng)用前景。

但巨額的計(jì)算時(shí)間成為了限制高保真時(shí)空中子動(dòng)力學(xué)計(jì)算發(fā)展的最主要瓶頸。國(guó)際上目前主要流行通過預(yù)估校正準(zhǔn)靜態(tài)方法處理時(shí)空中子動(dòng)力學(xué)中與時(shí)間相關(guān)的因素。在預(yù)估校正準(zhǔn)靜態(tài)方法中，通過在較大時(shí)間步長(zhǎng)內(nèi)以高計(jì)算代價(jià)直接求解時(shí)空中子輸運(yùn)方程獲得預(yù)估的中子分布，此為預(yù)估，再在此求解輸運(yùn)方程的步長(zhǎng)內(nèi)以較小步長(zhǎng)和較小計(jì)算代價(jià)求解點(diǎn)堆方程獲得校正的中子分布的幅值，最后在輸運(yùn)計(jì)算的步長(zhǎng)末通過校正的幅值去修正預(yù)估的中子分布的幅值，此為校正。但是預(yù)估校正準(zhǔn)靜態(tài)方法在處理極短時(shí)間內(nèi)核反應(yīng)堆內(nèi)中子迅速變化的情況時(shí)會(huì)出現(xiàn)較大的誤差，導(dǎo)致最終計(jì)算結(jié)果存在較大誤差，因此極大地限制了預(yù)估校正準(zhǔn)靜態(tài)方法在高保真時(shí)空中子動(dòng)力學(xué)中的應(yīng)用。

因此，充分利用預(yù)估校正準(zhǔn)靜態(tài)極大節(jié)省計(jì)算時(shí)間的優(yōu)勢(shì)，并通過一定的手段克服預(yù)估校正準(zhǔn)靜態(tài)在處理極短時(shí)間內(nèi)核反應(yīng)堆內(nèi)中子迅速變化時(shí)失準(zhǔn)的缺陷，才能使預(yù)估校正準(zhǔn)靜態(tài)方法在高保真時(shí)空中子動(dòng)力學(xué)的模擬中發(fā)揮真正的作用。

技術(shù)實(shí)現(xiàn)要素：

為了克服上述現(xiàn)有技術(shù)的問題，解決在核反應(yīng)堆功率急劇變化的情況下，預(yù)估校正準(zhǔn)靜態(tài)方法計(jì)算失準(zhǔn)的問題，本發(fā)明的目的是提出一種精確計(jì)算核反應(yīng)堆內(nèi)時(shí)空中子分布的方法，該方法在現(xiàn)有預(yù)估校正準(zhǔn)靜態(tài)方法的基礎(chǔ)上，通過改進(jìn)現(xiàn)有的預(yù)估校正準(zhǔn)靜態(tài)方法，在不顯著增加計(jì)算時(shí)間的前提下，能夠精確捕捉核反應(yīng)堆在極短時(shí)間內(nèi)中子的迅速變化，徹底克服傳統(tǒng)預(yù)估校正準(zhǔn)靜態(tài)方法的缺陷，使其在高保真時(shí)空中子動(dòng)力學(xué)計(jì)算中真正發(fā)揮作用。

為了達(dá)到上述目的，本發(fā)明采用如下技術(shù)方案：

一種精確計(jì)算核反應(yīng)堆時(shí)空中子分布的方法，包括如下步驟：

步驟1：進(jìn)行核反應(yīng)堆穩(wěn)態(tài)計(jì)算，采用一步法直接進(jìn)行輸運(yùn)計(jì)算，得到核反應(yīng)堆處于穩(wěn)態(tài)時(shí)每個(gè)能群、每個(gè)平源區(qū)的中子通量密度以及每個(gè)平源區(qū)內(nèi)的每組緩發(fā)中子先驅(qū)核初始的濃度，具體包括如下步驟：

1)從截面文件中讀取各個(gè)材料的原始多群宏觀截面信息與動(dòng)力學(xué)參數(shù)信息；

2)從輸入卡片中讀取核反應(yīng)堆的幾何信息與計(jì)算條件信息；

3)根據(jù)輸入卡片中讀取的幾何信息進(jìn)行直接的幾何建模：首先根據(jù)輸入卡片中的幾何描述得到核反應(yīng)堆的幾何布置；其次根據(jù)核反應(yīng)堆的幾何布置建立特征線方法計(jì)算所需的邊界條件以及內(nèi)部特征線的長(zhǎng)度信息，為中子輸運(yùn)計(jì)算提供模塊化特征線信息；

4)根據(jù)1)、2)、3)中得到的信息采用特征線方法進(jìn)行中子輸運(yùn)計(jì)算，得到各個(gè)平源區(qū)的中子通量密度分布，具體的計(jì)算公式如下：

式中：

ω——角度方向

——梯度算子

g——當(dāng)前能群編號(hào)

g'——非當(dāng)前能群編號(hào)

g——能群總數(shù)

r——空間位置

σt,g(r)——r處第g群的宏觀總截面

σs,g'→g(r)——r處g'能群到g能群的散射截面

χg(r)——r處第g能群的裂變譜

νσf,g——第g群的中子產(chǎn)生截面

——r,ω處第g能群中子角通量密度

φg(r)——r處第g能群的中子標(biāo)通量密度

φg'(r)——r處第g'能群的中子標(biāo)通量密度

sf(r)——r處裂變?cè)?/p>

keff——輸運(yùn)計(jì)算得到的有效增殖因子

由此得到各個(gè)平源區(qū)的中子通量密度；

5)根據(jù)4)中計(jì)算所得到的各個(gè)平源區(qū)的中子通量密度以及1)中讀取的動(dòng)力學(xué)參數(shù)信息，得到穩(wěn)態(tài)下的各組緩發(fā)中子先驅(qū)核密度，具體的公式如下：

式中：

r——空間位置

k——緩發(fā)中子先驅(qū)核編號(hào)

g——能群編號(hào)

g——能群總數(shù)

keff——輸運(yùn)計(jì)算得到的有效增殖因子

ck(r)——r處臨界狀態(tài)第k組緩發(fā)中子先驅(qū)核密度

βk(r)——r處第k組緩發(fā)中子份額

νσf,g——第g群的中子產(chǎn)生截面

λk(r)——r處第k組緩發(fā)中子先驅(qū)核的衰變常數(shù)

φg(r)——r處臨界狀態(tài)第g群的中子通量密度

sf(r)——r處裂變?cè)?/p>

6)根據(jù)1)、2)、3)得到的信息，對(duì)粗網(wǎng)有限差分cmfd方程進(jìn)行中子共軛計(jì)算，得到各個(gè)粗網(wǎng)的共軛中子通量密度，具體的計(jì)算公式如下：

式中：

u——坐標(biāo)軸方向標(biāo)號(hào)

x,y,z——直角坐標(biāo)系坐標(biāo)軸方向

hu——粗網(wǎng)在u方向的高度

——粗網(wǎng)在u方向右邊界的凈中子流

——粗網(wǎng)在u方向左邊界的凈中子流

g——當(dāng)前能群標(biāo)號(hào)

g'——非當(dāng)前能群標(biāo)號(hào)

g——總能群數(shù)

σrg——第g能群的移出截面

σs,g→g'——第g群到第g'群的散射截面

vσf,g——第g能群的中子產(chǎn)生截面

χg'——第g'能群的裂變譜

keff——輸運(yùn)計(jì)算得到的有效增殖因子

φg——粗網(wǎng)第g能群的共軛中子通量密度

φg'——粗網(wǎng)第g'能群的共軛中子通量密度

7)將每個(gè)平源區(qū)的中子通量密度利用因子分解分解成為幅值函數(shù)與形狀函數(shù)的乘積，此時(shí)幅值函數(shù)為1，利用4)求解得到的臨界狀態(tài)下的中子通量密度求得初始時(shí)刻的中子形狀函數(shù)，具體的計(jì)算公式如下：

式中：

g——能群標(biāo)號(hào)

r——空間位置

φg(r,0)——穩(wěn)態(tài)時(shí)r處第g能群的中子通量密度

ψg(r,0)——穩(wěn)態(tài)時(shí)r處第g能群的中子通量密度形狀函數(shù)

n(0)——穩(wěn)態(tài)時(shí)中子通量密度的幅值函數(shù)

步驟2：讀取輸入卡片擾動(dòng)信息，根據(jù)得到的擾動(dòng)信息執(zhí)行截面擾動(dòng)，打破核反應(yīng)堆的穩(wěn)態(tài)狀態(tài)，從而開始核反應(yīng)堆時(shí)空中子動(dòng)力學(xué)計(jì)算；

步驟3：執(zhí)行時(shí)空中子動(dòng)力學(xué)計(jì)算的第一個(gè)時(shí)間步，以0.25ms為輸運(yùn)計(jì)算時(shí)間步長(zhǎng)，進(jìn)行中子輸運(yùn)形式的固定源計(jì)算，得到0.25ms時(shí)刻的中子通量密度以及中子通量密度形狀函數(shù)，并由中子通量密度以及中子通量密度形狀函數(shù)得到中子通量密度幅值，具體包括如下步驟：

1)以0.25ms為輸運(yùn)計(jì)算時(shí)間步長(zhǎng)進(jìn)行中子輸運(yùn)形式的固定源計(jì)算，得到0.25ms時(shí)刻的中子通量密度，具體公式如下：

式中：

ω——角度方向

——梯度算子

g——當(dāng)前能群標(biāo)號(hào)

g'——非當(dāng)前能群標(biāo)號(hào)

g——能群總數(shù)

r——空間位置

n——第n個(gè)輸運(yùn)計(jì)算步

k——緩發(fā)中子標(biāo)號(hào)

——n時(shí)刻r處第g群的宏觀總截面

——n時(shí)刻r處g'能群到g能群的散射截面

——n時(shí)刻r處第g能群的裂變譜

——n時(shí)刻第g群的中子產(chǎn)生截面

——n時(shí)刻r,ω處第g能群中子角通量密度

——n時(shí)刻r處第g能群的中子標(biāo)通量密度

——n時(shí)刻r處第g'能群的中子標(biāo)通量密度

——n-1時(shí)刻r處第g能群的中子標(biāo)通量密度

keff——輸運(yùn)計(jì)算得到的有效增殖因子

——n時(shí)刻r處第g能群的固定源

——n時(shí)刻r處的裂變?cè)?/p>

ag(r)——r處第g能群的固定源系數(shù)

bg(r)——r處第g能群的固定源系數(shù)

cg(r)——r處第g能群的固定源系數(shù)

vg——第g能群的中子速度

δtn——第n個(gè)輸運(yùn)步的步長(zhǎng)

——n時(shí)刻第g能群的等效緩發(fā)份額

χdk,g——第k組緩發(fā)中子在第g能群的緩發(fā)裂變譜

βk——第k組緩發(fā)中子份額

——n-1時(shí)刻r處第g能群的等效緩發(fā)源項(xiàng)

2)計(jì)算歸一化常數(shù)，歸一化常數(shù)是求解得到0.25ms時(shí)刻處中子通量密度形狀函數(shù)的必須要素，具體公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間變量

φ^*——粗網(wǎng)的共軛中子通量密度

v——中子速度

——t時(shí)刻的中子通量密度形狀函數(shù)

——穩(wěn)態(tài)下的中子通量密度形狀函數(shù)

φi,g(0)——穩(wěn)態(tài)下i平源區(qū)的第g能群的中子通量密度

c——?dú)w一化常數(shù)

3)利用1)和2)中計(jì)算得到的中子通量密度與歸一化常數(shù)，得到0.25ms時(shí)刻的中子通量密度形狀函數(shù)，具體計(jì)算公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間變量

——t時(shí)刻i平源區(qū)的第g能群中子通量密度形狀函數(shù)

φi,g(t)——t時(shí)刻i平源區(qū)的第g能群的中子通量密度

n(t)——t時(shí)刻的幅值函數(shù)

φ(t)——t時(shí)刻的中子通量密度

c——?dú)w一化常數(shù)

v——中子速度

φ^*——粗網(wǎng)的共軛中子通量密度

——t時(shí)刻中子通量密度的形狀函數(shù)

<·>——對(duì)全能量、全空間積分

4)由1)計(jì)算得到的中子通量密度與3)中計(jì)算得到的中子通量密度形狀函數(shù)得到0.25ms時(shí)刻的中子通量密度幅值，具體公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間變量

n(t)——t時(shí)刻的中子通量密度幅值函數(shù)

φi,g(t)——t時(shí)刻i平源區(qū)的第g能群的中子通量密度

——t時(shí)刻i平源區(qū)的第g能群的中子通量密度形狀函數(shù)

步驟4：在tn-1,transport至tn,transport時(shí)間間隔內(nèi)以遠(yuǎn)大于步驟3中0.25ms的時(shí)間步長(zhǎng)進(jìn)行求解輸運(yùn)形式的固定源方程，n＝1,n，其中n為整個(gè)計(jì)算過程劃分的輸運(yùn)計(jì)算步數(shù)，求得tn,transport時(shí)刻的中子通量密度，此步驟為預(yù)估校正準(zhǔn)靜態(tài)中的預(yù)估計(jì)算步，具體的計(jì)算公式如下：

式中：

ω——角度方向

——梯度算子

g——當(dāng)前能群標(biāo)號(hào)

g'——非當(dāng)前能群標(biāo)號(hào)

g——能群總數(shù)

r——空間位置

n——第n個(gè)輸運(yùn)計(jì)算步

k——緩發(fā)中子標(biāo)號(hào)

——n時(shí)刻r處第g群的宏觀總截面

——n時(shí)刻r處g'能群到g能群的散射截面

——n時(shí)刻r處第g能群的裂變譜

——n時(shí)刻第g群的中子產(chǎn)生截面

——n時(shí)刻r,ω處第g能群中子角通量密度

——n時(shí)刻r處第g能群的中子標(biāo)通量密度

——n時(shí)刻r處第g'能群的中子標(biāo)通量密度

——n-1時(shí)刻r處第g能群的中子標(biāo)通量密度

keff——輸運(yùn)計(jì)算得到的有效增殖因子

——n時(shí)刻r處第g能群的固定源

——n時(shí)刻r處的裂變?cè)?/p>

ag(r)——r處第g能群的固定源系數(shù)

bg(r)——r處第g能群的固定源系數(shù)

cg(r)——r處第g能群的固定源系數(shù)

vg——第g能群的中子速度

δtn——第n個(gè)輸運(yùn)步的步長(zhǎng)

——n時(shí)刻第g能群的等效緩發(fā)份額

χdk,g——第k組緩發(fā)中子在第g能群的緩發(fā)裂變譜

βk——第k組緩發(fā)中子份額

——n-1時(shí)刻r處第g能群的等效緩發(fā)源項(xiàng)

步驟5：將步驟4中得到的中子通量密度通過因子分解分解為中子通量密度幅值函數(shù)與中子通量密度形狀函數(shù)的乘積，在tn-1,transport與tn,transport兩個(gè)時(shí)間點(diǎn)計(jì)算點(diǎn)堆參數(shù)，在tn-1,transport與tn,transport之間對(duì)點(diǎn)堆參數(shù)進(jìn)行線性插值，在tn-1,transport與tn,transport之間求解點(diǎn)堆方程組，得到tn,transport時(shí)刻的點(diǎn)堆幅值，利用tn,transport時(shí)刻的點(diǎn)堆幅值對(duì)tn,transport時(shí)刻的中子通量密度進(jìn)行幅值的校正，此為預(yù)估校正準(zhǔn)靜態(tài)中的校正步，具體包含如下步驟：

1)將步驟4中得到的中子通量密度通過因子分解分解為中子通量密度幅值函數(shù)與中子通量密度形狀函數(shù)的乘積，進(jìn)而利用步驟1中求得的粗網(wǎng)的共軛中子通量密度，求解出tn,transport的中子通量密度形狀函數(shù)，具體公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間變量

——t時(shí)刻i平源區(qū)的第g能群的預(yù)估中子通量密度

n^p(t)——t時(shí)刻預(yù)估的中子通量密度幅值函數(shù)

n(t)——t時(shí)刻的幅值函數(shù)

——t時(shí)刻i平源區(qū)的第g能群的中子通量密度形狀函數(shù)

c——?dú)w一化常數(shù)

v——中子速度

φ^*——粗網(wǎng)的共軛中子通量密度

——t時(shí)刻的中子通量密度形狀函數(shù)

——穩(wěn)態(tài)下的中子通量密度形狀函數(shù)

φi,g(0)——穩(wěn)態(tài)下i平源區(qū)的第g能群的中子通量密度

<·>——對(duì)全能量、全空間積分

2)分別計(jì)算tn-1,transport與tn,transport時(shí)刻的精確點(diǎn)堆參數(shù)，具體公式如下:

f(t)＝<φ^*(r)χg(r)sf(r,t)>

式中：

t——時(shí)間

g——當(dāng)前能群標(biāo)號(hào)

g'——非當(dāng)前能群標(biāo)號(hào)

g——總能群數(shù)

r——空間位置

k——緩發(fā)中子標(biāo)號(hào)

u——坐標(biāo)軸方向標(biāo)號(hào)

x,y,z——直角坐標(biāo)系坐標(biāo)軸方向

hu——粗網(wǎng)在u方向的高度

——粗網(wǎng)第g群在u方向右邊界的凈中子流

——粗網(wǎng)第g群在u方向左邊界的凈中子流

σrg——第g能群的移出截面

σs,g→g'——第g群到第g'群的散射截面

vσf,g——第g能群的中子產(chǎn)生截面

χg(r)——r處第g能群的裂變譜

keff——輸運(yùn)計(jì)算得到的有效增殖因子

φ^*(r)——r處粗網(wǎng)的共軛中子通量密度

——t時(shí)刻r處的中子通量密度形狀函數(shù)

sf(r,t)——t時(shí)刻r處的裂變?cè)?/p>

keff——輸運(yùn)計(jì)算得到的有效增殖因子

f(t)——t時(shí)刻精確點(diǎn)堆參數(shù)分母

ρ(t)——t時(shí)刻的反應(yīng)性

βk(r)——r處第k組緩發(fā)中子份額

χdk,g(r)——r處第k組緩發(fā)中子在第g能群的緩發(fā)裂變譜

——t時(shí)刻第k組緩發(fā)中子份額

——t時(shí)刻緩發(fā)中子份額總和

v(r)——r處中子速度

λ(t)——t時(shí)刻等效中子代時(shí)間

λk(t)——t時(shí)刻第k組緩發(fā)中子衰變常數(shù)

ck(r,t)——t時(shí)刻r處臨界狀態(tài)第k組緩發(fā)中子先驅(qū)核密度

λk(r)——r處第k組緩發(fā)中子先驅(qū)核的衰變常數(shù)

ck(t)——t時(shí)刻第k組緩發(fā)中子先驅(qū)核濃度

<·>——在全能量、全空間進(jìn)行積分

3)在tn-1,transport與tn,transport之間以δtpk為時(shí)間步長(zhǎng)等距劃分點(diǎn)堆計(jì)算所需的間隔，在點(diǎn)堆所計(jì)算的時(shí)間間隔tn-1,pk至tn,pk內(nèi)用tn-1,transport與tn,transport時(shí)間點(diǎn)上的點(diǎn)堆參數(shù)進(jìn)行插值，其中進(jìn)而求解點(diǎn)堆方程組，得到tn,pk時(shí)刻的幅值n^c，具體公式如下：

式中：

t——時(shí)間變量

i——緩發(fā)中子標(biāo)號(hào)

n(t)——t時(shí)刻幅值

ci(t)——t時(shí)刻第i組緩發(fā)中子先驅(qū)核濃度

ρ(t)——t時(shí)刻反應(yīng)性

βi(t)——t時(shí)刻第i組緩發(fā)中子份額

β(t)——t時(shí)刻緩發(fā)中子份額總和

λ(t)——t時(shí)刻等效中子代時(shí)間

λi——第i組緩發(fā)中子衰變常數(shù)

4)當(dāng)點(diǎn)堆計(jì)算的時(shí)間點(diǎn)達(dá)到tn,transport時(shí)，校正tn,transport時(shí)刻的中子通量密度，具體公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間

——t時(shí)刻i平源區(qū)第g能群校正后的中子通量密度

n^c(t)——t時(shí)刻校正的幅值函數(shù)

n^p(t)——t時(shí)刻預(yù)估的幅值函數(shù)

——t時(shí)刻i平源區(qū)第g能群的中子通量密度形狀函數(shù)

步驟6：重復(fù)執(zhí)行步驟4與步驟5，直到時(shí)空中子動(dòng)力學(xué)計(jì)算結(jié)束。

與現(xiàn)有技術(shù)相比，本發(fā)明具有如下突出優(yōu)點(diǎn)：

1)徹底消除了傳統(tǒng)預(yù)估校正準(zhǔn)靜態(tài)方法在處理核反應(yīng)堆功率急劇變化時(shí)計(jì)算結(jié)果偏差較大的缺陷，在不顯著增加計(jì)算時(shí)間的前提下大幅度提高計(jì)算精度，對(duì)預(yù)估校正準(zhǔn)靜態(tài)方法在高保真時(shí)空中子動(dòng)力學(xué)計(jì)算中得到真正應(yīng)用具有極重要作用；

2)能夠?qū)崿F(xiàn)變步長(zhǎng)的計(jì)算策略，能夠?qū)崿F(xiàn)在核反應(yīng)堆功率急劇變化的時(shí)間點(diǎn)采取小時(shí)間步長(zhǎng)進(jìn)行計(jì)算，在核反應(yīng)堆功率變化平緩的時(shí)間點(diǎn)采取大時(shí)間步長(zhǎng)進(jìn)行計(jì)算，能夠靈活充分應(yīng)用計(jì)算機(jī)資源，提高計(jì)算的經(jīng)濟(jì)性。

附圖說明

圖1為堆芯穩(wěn)態(tài)計(jì)算流程。

圖2為輸運(yùn)修正點(diǎn)堆幅值流程圖。

圖3為點(diǎn)堆修正輸運(yùn)幅值流程圖。

圖4為c5g7‐td基準(zhǔn)題功率變化對(duì)比圖。

圖5為c5g7‐td基準(zhǔn)題誤差對(duì)比圖。

具體實(shí)施方式

下面結(jié)合附圖和具體實(shí)施方式對(duì)本發(fā)明進(jìn)行詳細(xì)的說明：

步驟1：如圖1所示，進(jìn)行核反應(yīng)堆穩(wěn)態(tài)計(jì)算，采用一步法直接進(jìn)行輸運(yùn)計(jì)算，得到核反應(yīng)堆處于穩(wěn)態(tài)時(shí)每個(gè)能群、每個(gè)平源區(qū)的中子通量密度以及每個(gè)平源區(qū)內(nèi)的每組緩發(fā)中子先驅(qū)核初始的濃度，具體包括如下步驟：

1)從截面文件中讀取各個(gè)材料的原始多群宏觀截面信息與動(dòng)力學(xué)參數(shù)信息；

2)從輸入卡片中讀取核反應(yīng)堆的幾何信息與計(jì)算條件信息；

3)根據(jù)輸入卡片中讀取的幾何信息進(jìn)行直接的幾何建模：首先根據(jù)輸入卡片中的幾何描述得到核反應(yīng)堆的幾何布置；其次根據(jù)核反應(yīng)堆的幾何布置建立特征線方法計(jì)算所需的邊界條件以及內(nèi)部特征線的長(zhǎng)度信息，為中子輸運(yùn)計(jì)算提供模塊化特征線信息；

4)根據(jù)1)、2)、3)中得到的信息采用特征線方法進(jìn)行中子輸運(yùn)計(jì)算，得到各個(gè)平源區(qū)的中子通量密度分布，具體的計(jì)算公式如下：

式中：

ω——角度方向

——梯度算子

g——當(dāng)前能群編號(hào)

g'——非當(dāng)前能群編號(hào)

g——能群總數(shù)

r——空間位置

σt,g(r)——r處第g群的宏觀總截面

σs,g'→g(r)——r處g'能群到g能群的散射截面

χg(r)——r處第g能群的裂變譜

νσf,g——第g群的中子產(chǎn)生截面

——r,ω處第g能群中子角通量密度

φg(r)——r處第g能群的中子標(biāo)通量密度

φg'(r)——r處第g'能群的中子標(biāo)通量密度

sf(r)——r處裂變?cè)?/p>

keff——輸運(yùn)計(jì)算得到的有效增殖因子

由此得到各個(gè)平源區(qū)的中子通量密度；

5)根據(jù)4)中計(jì)算所得到的各個(gè)平源區(qū)的中子通量密度以及1)中讀取的動(dòng)力學(xué)參數(shù)信息，得到穩(wěn)態(tài)下的各組緩發(fā)中子先驅(qū)核密度，具體的公式如下：

式中：

r——空間位置

k——緩發(fā)中子先驅(qū)核編號(hào)

g——能群編號(hào)

g——能群總數(shù)

keff——輸運(yùn)計(jì)算得到的有效增殖因子

ck(r)——r處臨界狀態(tài)第k組緩發(fā)中子先驅(qū)核密度

βk(r)——r處第k組緩發(fā)中子份額

νσf,g——第g群的中子產(chǎn)生截面

λk(r)——r處第k組緩發(fā)中子先驅(qū)核的衰變常數(shù)

φg(r)——r處臨界狀態(tài)第g群的中子通量密度

sf(r)——r處裂變?cè)?/p>

6)根據(jù)1)、2)、3)得到的信息，對(duì)粗網(wǎng)有限差分cmfd方程進(jìn)行中子共軛計(jì)算，得到各個(gè)粗網(wǎng)的共軛中子通量密度，具體的計(jì)算公式如下：

式中：

u——坐標(biāo)軸方向標(biāo)號(hào)

x,y,z——直角坐標(biāo)系坐標(biāo)軸方向

hu——粗網(wǎng)在u方向的高度

——粗網(wǎng)在u方向右邊界的凈中子流

——粗網(wǎng)在u方向左邊界的凈中子流

g——當(dāng)前能群標(biāo)號(hào)

g'——非當(dāng)前能群標(biāo)號(hào)

g——總能群數(shù)

σrg——第g能群的移出截面

σs,g→g'——第g群到第g'群的散射截面

vσf,g——第g能群的中子產(chǎn)生截面

χg'——第g'能群的裂變譜

keff——輸運(yùn)計(jì)算得到的有效增殖因子

φg——粗網(wǎng)第g能群的共軛中子通量密度

φg'——粗網(wǎng)第g'能群的共軛中子通量密度

7)將每個(gè)平源區(qū)的中子通量密度利用因子分解分解成為幅值函數(shù)與形狀函數(shù)的乘積，此時(shí)幅值函數(shù)為1，利用4)求解得到的臨界狀態(tài)下的中子通量密度求得初始時(shí)刻的中子形狀函數(shù)，具體的計(jì)算公式如下：

式中：

g——能群標(biāo)號(hào)

r——空間位置

φg(r,0)——穩(wěn)態(tài)時(shí)r處第g能群的中子通量密度

ψg(r,0)——穩(wěn)態(tài)時(shí)r處第g能群的中子通量密度形狀函數(shù)

n(0)——穩(wěn)態(tài)時(shí)中子通量密度的幅值函數(shù)

步驟2：讀取輸入卡片擾動(dòng)信息，根據(jù)得到的擾動(dòng)信息執(zhí)行截面擾動(dòng)，打破核反應(yīng)堆的穩(wěn)態(tài)狀態(tài)，從而開始核反應(yīng)堆時(shí)空中子動(dòng)力學(xué)計(jì)算；

步驟3：如圖2所示，執(zhí)行時(shí)空中子動(dòng)力學(xué)計(jì)算的第一個(gè)時(shí)間步，以0.25ms為輸運(yùn)計(jì)算時(shí)間步長(zhǎng)，進(jìn)行中子輸運(yùn)形式的固定源計(jì)算，得到0.25ms時(shí)刻的中子通量密度以及中子通量密度形狀函數(shù)，并由中子通量密度以及中子通量密度形狀函數(shù)得到中子通量密度幅值，具體包括如下步驟：

1)以0.25ms為輸運(yùn)計(jì)算時(shí)間步長(zhǎng)進(jìn)行中子輸運(yùn)形式的固定源計(jì)算，得到0.25ms時(shí)刻的中子通量密度，具體公式如下：

式中：

ω——角度方向

——梯度算子

g——當(dāng)前能群標(biāo)號(hào)

g'——非當(dāng)前能群標(biāo)號(hào)

g——能群總數(shù)

r——空間位置

n——第n個(gè)輸運(yùn)計(jì)算步

k——緩發(fā)中子標(biāo)號(hào)

——n時(shí)刻r處第g群的宏觀總截面

——n時(shí)刻r處g'能群到g能群的散射截面

——n時(shí)刻r處第g能群的裂變譜

——n時(shí)刻第g群的中子產(chǎn)生截面

——n時(shí)刻r,ω處第g能群中子角通量密度

——n時(shí)刻r處第g能群的中子標(biāo)通量密度

——n時(shí)刻r處第g'能群的中子標(biāo)通量密度

——n-1時(shí)刻r處第g能群的中子標(biāo)通量密度

keff——輸運(yùn)計(jì)算得到的有效增殖因子

——n時(shí)刻r處第g能群的固定源

——n時(shí)刻r處的裂變?cè)?/p>

ag(r)——r處第g能群的固定源系數(shù)

bg(r)——r處第g能群的固定源系數(shù)

cg(r)——r處第g能群的固定源系數(shù)

vg——第g能群的中子速度

δtn——第n個(gè)輸運(yùn)步的步長(zhǎng)

——n時(shí)刻第g能群的等效緩發(fā)份額

χdk,g——第k組緩發(fā)中子在第g能群的緩發(fā)裂變譜

βk——第k組緩發(fā)中子份額

——n-1時(shí)刻r處第g能群的等效緩發(fā)源項(xiàng)

2)計(jì)算歸一化常數(shù)，歸一化常數(shù)是求解得到0.25ms時(shí)刻處中子通量密度形狀函數(shù)的必須要素，具體公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間變量

φ^*——粗網(wǎng)的共軛中子通量密度

v——中子速度

——t時(shí)刻的中子通量密度形狀函數(shù)

——穩(wěn)態(tài)下的中子通量密度形狀函數(shù)

φi,g(0)——穩(wěn)態(tài)下i平源區(qū)的第g能群的中子通量密度

c——?dú)w一化常數(shù)

3)利用1)和2)中計(jì)算得到的中子通量密度與歸一化常數(shù)，得到0.25ms時(shí)刻的中子通量密度形狀函數(shù)，具體計(jì)算公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間變量

——t時(shí)刻i平源區(qū)的第g能群中子通量密度形狀函數(shù)

φi,g(t)——t時(shí)刻i平源區(qū)的第g能群的中子通量密度

n(t)——t時(shí)刻的幅值函數(shù)

φ(t)——t時(shí)刻的中子通量密度

c——?dú)w一化常數(shù)

v——中子速度

φ^*——粗網(wǎng)的共軛中子通量密度

——t時(shí)刻中子通量密度的形狀函數(shù)

<·>——對(duì)全能量、全空間積分

4)由1)計(jì)算得到的中子通量密度與3)中計(jì)算得到的中子通量密度形狀函數(shù)得到0.25ms時(shí)刻的中子通量密度幅值，具體公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間變量

n(t)——t時(shí)刻的中子通量密度幅值函數(shù)

φi,g(t)——t時(shí)刻i平源區(qū)的第g能群的中子通量密度

——t時(shí)刻i平源區(qū)的第g能群的中子通量密度形狀函數(shù)

步驟4：如圖3所示，在tn-1,transport至tn,transport時(shí)間間隔內(nèi)以遠(yuǎn)大于步驟3中0.25ms的時(shí)間步長(zhǎng)進(jìn)行求解輸運(yùn)形式的固定源方程，n＝1,n，其中n為整個(gè)計(jì)算過程劃分的輸運(yùn)計(jì)算步數(shù)，求得tn,transport時(shí)刻的中子通量密度，此步驟為預(yù)估校正準(zhǔn)靜態(tài)中的預(yù)估計(jì)算步，具體的計(jì)算公式如下：

式中：

ω——角度方向

——梯度算子

g——當(dāng)前能群標(biāo)號(hào)

g'——非當(dāng)前能群標(biāo)號(hào)

g——能群總數(shù)

r——空間位置

n——第n個(gè)輸運(yùn)計(jì)算步

k——緩發(fā)中子標(biāo)號(hào)

——n時(shí)刻r處第g群的宏觀總截面

——n時(shí)刻r處g'能群到g能群的散射截面

——n時(shí)刻r處第g能群的裂變譜

——n時(shí)刻第g群的中子產(chǎn)生截面

——n時(shí)刻r,ω處第g能群中子角通量密度

——n時(shí)刻r處第g能群的中子標(biāo)通量密度

——n時(shí)刻r處第g'能群的中子標(biāo)通量密度

——n-1時(shí)刻r處第g能群的中子標(biāo)通量密度

keff——輸運(yùn)計(jì)算得到的有效增殖因子

——n時(shí)刻r處第g能群的固定源

——n時(shí)刻r處的裂變?cè)?/p>

ag(r)——r處第g能群的固定源系數(shù)

bg(r)——r處第g能群的固定源系數(shù)

cg(r)——r處第g能群的固定源系數(shù)

vg——第g能群的中子速度

δtn——第n個(gè)輸運(yùn)步的步長(zhǎng)

——n時(shí)刻第g能群的等效緩發(fā)份額

χdk,g——第k組緩發(fā)中子在第g能群的緩發(fā)裂變譜

βk——第k組緩發(fā)中子份額

——n-1時(shí)刻r處第g能群的等效緩發(fā)源項(xiàng)

步驟5：將步驟4中得到的中子通量密度通過因子分解分解為中子通量密度幅值函數(shù)與中子通量密度形狀函數(shù)的乘積，在tn-1,transport與tn,transport兩個(gè)時(shí)間點(diǎn)計(jì)算點(diǎn)堆參數(shù)，在tn-1,transport與tn,transport之間對(duì)點(diǎn)堆參數(shù)進(jìn)行線性插值，在tn-1,transport與tn,transport之間求解點(diǎn)堆方程組，得到tn,transport時(shí)刻的點(diǎn)堆幅值，利用tn,transport時(shí)刻的點(diǎn)堆幅值對(duì)tn,transport時(shí)刻的中子通量密度進(jìn)行幅值的校正，此為預(yù)估校正準(zhǔn)靜態(tài)中的校正步，具體包含如下步驟：

1)將步驟4中得到的中子通量密度通過因子分解分解為中子通量密度幅值函數(shù)與中子通量密度形狀函數(shù)的乘積，進(jìn)而利用步驟1中求得的粗網(wǎng)的共軛中子通量密度，求解出tn,transport的中子通量密度形狀函數(shù)，具體公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間變量

——t時(shí)刻i平源區(qū)的第g能群的預(yù)估中子通量密度

n^p(t)——t時(shí)刻預(yù)估的中子通量密度幅值函數(shù)

n(t)——t時(shí)刻的幅值函數(shù)

——t時(shí)刻i平源區(qū)的第g能群的中子通量密度形狀函數(shù)

c——?dú)w一化常數(shù)

v——中子速度

φ^*——粗網(wǎng)的共軛中子通量密度

——t時(shí)刻的中子通量密度形狀函數(shù)

——穩(wěn)態(tài)下的中子通量密度形狀函數(shù)

φi,g(0)——穩(wěn)態(tài)下i平源區(qū)的第g能群的中子通量密度

<·>——對(duì)全能量、全空間積分

2)分別計(jì)算tn-1,transport與tn,transport時(shí)刻的精確點(diǎn)堆參數(shù)，具體公式如下:

f(t)＝<φ^*(r)χg(r)sf(r,t)>

式中：

t——時(shí)間

g——當(dāng)前能群標(biāo)號(hào)

g'——非當(dāng)前能群標(biāo)號(hào)

g——總能群數(shù)

r——空間位置

k——緩發(fā)中子標(biāo)號(hào)

u——坐標(biāo)軸方向標(biāo)號(hào)

x,y,z——直角坐標(biāo)系坐標(biāo)軸方向

hu——粗網(wǎng)在u方向的高度

——粗網(wǎng)第g群在u方向右邊界的凈中子流

——粗網(wǎng)第g群在u方向左邊界的凈中子流

σrg——第g能群的移出截面

σs,g→g'——第g群到第g'群的散射截面

vσf,g——第g能群的中子產(chǎn)生截面

χg(r)——r處第g能群的裂變譜

keff——輸運(yùn)計(jì)算得到的有效增殖因子

φ^*(r)——r處粗網(wǎng)的共軛中子通量密度

——t時(shí)刻r處的中子通量密度形狀函數(shù)

sf(r,t)——t時(shí)刻r處的裂變?cè)?/p>

keff——輸運(yùn)計(jì)算得到的有效增殖因子

f(t)——t時(shí)刻精確點(diǎn)堆參數(shù)分母

ρ(t)——t時(shí)刻的反應(yīng)性

βk(r)——r處第k組緩發(fā)中子份額

χdk,g(r)——r處第k組緩發(fā)中子在第g能群的緩發(fā)裂變譜

——t時(shí)刻第k組緩發(fā)中子份額

——t時(shí)刻緩發(fā)中子份額總和

v(r)——r處中子速度

λ(t)——t時(shí)刻等效中子代時(shí)間

λk(t)——t時(shí)刻第k組緩發(fā)中子衰變常數(shù)

ck(r,t)——t時(shí)刻r處臨界狀態(tài)第k組緩發(fā)中子先驅(qū)核密度

λk(r)——r處第k組緩發(fā)中子先驅(qū)核的衰變常數(shù)

ck(t)——t時(shí)刻第k組緩發(fā)中子先驅(qū)核濃度

<·>——在全能量、全空間進(jìn)行積分

3)在tn-1,transport與tn,transport之間以δtpk為時(shí)間步長(zhǎng)等距劃分點(diǎn)堆計(jì)算所需的間隔，在點(diǎn)堆所計(jì)算的時(shí)間間隔tn-1,pk至tn,pk內(nèi)用tn-1,transport與tn,transport時(shí)間點(diǎn)上的點(diǎn)堆參數(shù)進(jìn)行插值，其中進(jìn)而求解點(diǎn)堆方程組，得到tn,pk時(shí)刻的幅值n^c，具體公式如下：

式中：

t——時(shí)間變量

i——緩發(fā)中子標(biāo)號(hào)

n(t)——t時(shí)刻幅值

ci(t)——t時(shí)刻第i組緩發(fā)中子先驅(qū)核濃度

ρ(t)——t時(shí)刻反應(yīng)性

βi(t)——t時(shí)刻第i組緩發(fā)中子份額

β(t)——t時(shí)刻緩發(fā)中子份額總和

λ(t)——t時(shí)刻等效中子代時(shí)間

λi——第i組緩發(fā)中子衰變常數(shù)

4)當(dāng)點(diǎn)堆計(jì)算的時(shí)間點(diǎn)達(dá)到tn,transport時(shí)，校正tn,transport時(shí)刻的中子通量密度，具體公式如下：

式中：

i——平源區(qū)標(biāo)號(hào)

g——能群標(biāo)號(hào)

t——時(shí)間

——t時(shí)刻i平源區(qū)第g能群校正后的中子通量密度

n^c(t)——t時(shí)刻校正的幅值函數(shù)

n^p(t)——t時(shí)刻預(yù)估的幅值函數(shù)

——t時(shí)刻i平源區(qū)第g能群的中子通量密度形狀函數(shù)

步驟6：重復(fù)執(zhí)行步驟4與步驟5，直到時(shí)空中子動(dòng)力學(xué)計(jì)算結(jié)束。

下面以時(shí)空中子動(dòng)力學(xué)基準(zhǔn)題c5g7‐td系列的td0‐5的計(jì)算結(jié)果為例說明本發(fā)明的效果：c5g7‐td基準(zhǔn)題系列的td0‐5是一個(gè)二維七群非均勻的時(shí)空中子動(dòng)力學(xué)基準(zhǔn)題，基準(zhǔn)解為以0.25ms為步長(zhǎng)的全隱式向后差分計(jì)算所得，以2ms為時(shí)間步長(zhǎng)進(jìn)行傳統(tǒng)預(yù)估校正準(zhǔn)靜態(tài)計(jì)算與本發(fā)明修正的預(yù)估校正準(zhǔn)靜態(tài)計(jì)算，總功率變化對(duì)比如圖4所示，傳統(tǒng)預(yù)估校正準(zhǔn)靜態(tài)與本發(fā)明的修正預(yù)估校正準(zhǔn)靜態(tài)和基準(zhǔn)解的誤差如圖5所示。從計(jì)算結(jié)果可以看出，本發(fā)明在基本不增加計(jì)算時(shí)間的前提下極大地改善了傳統(tǒng)預(yù)估校正準(zhǔn)靜態(tài)的精度，此計(jì)算結(jié)果證明本發(fā)明具有相當(dāng)高的實(shí)用價(jià)值以及創(chuàng)新性。