先进核反应堆的设计原理

先进核动力反应堆的设计原理

核反应堆是将核子链式反应所造成的质量亏损转化为粒子动能、辐射能等形式能量的装置，继而通过外部的冷却循环系统及汽轮发电机组最终转化为我们给苹果手机充电所需要的电能。

从1941年到现在，人们设计出了各种各样的反应堆，有高温气冷堆、轻水堆、重水堆、钠冷快堆等，它的设计思路主要是基于核临界理论，在化石能源供应日趋紧张的今天，发展一种高效、安全、干净的核能已经成为历史发展的必然选择，目前，我国也建造了中国实验快堆(CFER)、中国先进研究堆（CARR）,并在先进核动力堆型及燃料后处理的领域上走在了世界前列，所采取的主要计算方法有CITATION程序、MATLAB方法等。下面，我们就结合理论与实际来介绍一种典型的核反应堆的设计过程。 一、 圆柱型轻水堆

一座核电站的设计图可以说是相当复杂，足够上万个人忙活几年的，但其核心就在于反应堆的临界参数，我们先采用较简单的双群理论来阐述一座理论堆型的设计原理。

取堆芯为浓缩度为3%的二氧化铀，反射层为H2O，快群与慢群能量分界按常规取1ev,相关中子反应截面如下:

1,ra0,1,rs3.22,ra0.022,2,rs3.451,ca0.04,1,cs3.752,ca0.5416,2,cs0.372

121,cf(E1E2)0.18

示意图如图一：

设堆芯与反射层的厚度分别为R、T，分别用角标c和r表示各项参数，如堆芯的中子慢化面积为Mc，反射层的中子慢化面积则为Mr，快群与热群的中子通量密度表示为：

1(r)2(r)



E0ECEC0

(r,E)dE1

(r,E)dE(2)



根据分析各种核子的中子反应截面我们可以看出，在轻水堆内，快群中子主要是由热中子引起的裂变产生的，

它又通过慢化吸收和泄

漏而消失，而热群中子则来源于快群中子的慢化，并主要由于吸收和泄漏而消失，很显然，在反应堆稳定工况时，各项参数都是恒定的，由中子平衡可以建立反应堆稳态时芯部的快群及热群的中子扩散方程如下：

D1,C1,C(r)r,c(r)

2

2

1Keff

[(f)1,c1,c(r)(f)2,c2,c(r)]3

D2,c2,c(r)a2,c2,c(r)12,c1,c(r)(4)

相应的，我们可以根据中子平衡关系写出反射层内中子扩散方程：

D2,r2,C(r)a2,r2,r(r)12,r1,r(r)

2

D1,r1,r(r)r,r1,r(r)0

2

为简便起见，我们通常可以将上述方程写成如下形式：

1,r(r)k1,r1,r(r)0(5)2,r(r)k2,r(r)

2

2

1,r

2

2

12,rD2,r

1,r(r)0(6)

k1,rk2,r

2

2

1,rD1,ra2,rD2,r

1Lr

2

其中

由(3)式可以得到快群中子通量密度表达式为：

1,c(r)

112,c

[D2,c2,c(r)a2,c2,c(r)]

2

将其代入（3）中，便得到只包含热群通量密度微分方程：

2,c

(r)的四阶偏



4

2,c(r)(

1

c



1L

2c

)2,c(r)

2

k1

'

cLc

2,c(r)0

它可以用因式分解的方法求解，将其改写成：

()()2,c0

2

2

2

2

（7）

同理有： 2X(r)2X(r)0

()()1,c0（8）

2

2

2

2

式中



2

1212

[(

1

c



1L

2c

)



2

[(

11Lc

2

c

)

方程(7)、(8)具有如下形式的解

X(r)X(r)0(9)Y(r)Y(r)0(10)

2

2

2

2

因其解X及Y为两个独立的函数，因而的解可以表示成两个函数的线性叠加：

1,cAX(r)CY(r)2,cAY(r)CY(r)

'

'

在这里A、C、A'、C'为四个待定常数，但从（4）式可以得到：

D2,cAX(r)a,2,cAX(r)12,cAX(r)

令S1=A/A=

’

2''

D2,ca,2,c

12,c

2

,相应地S2=C/C=

’

a,2,cD2,c

12,c

2

所以我们可以将快群、热群中子的解写成

1,cAX(r)CY(r)(11)2,c

AY(r)s1

'



CY(r)s2

'

(12)

其中A/C为待定常数，它们的值受反射层尺寸性质的约束。现在我们的任务就是根据（9）、（10）来定出X（r）、Y(r)的解。 先解（9）式，我们把X(r)写成R(r)Z(z)，将X(r)按圆柱坐标系展开得到：

X(r,z)r

2

2



1X(r,z)r

r



X(r,z)z

2

2

R(r,z)0

2

将含有径向与纵向变量的方程分离得到：

1R(r)

2

2

[

R(r)r

2

2



2

1R(r)r

r

]Bz(13)

22

Z(z)z

Bz(14)

其中BZ为待定常数。事实上，通过经验法求解（14），我们容易得到：

B(

2

z



H

),Z(z)cosBZZ

2

H为圆柱反应堆的高。 再来求解（13）

，令x见的零阶贝塞尔方程式：

x

2

,将其代入（13）中去，便得到常

dR(x)dx

2

2

x

dR(x)dx

xR(x)0

2

其普遍解为：

R(x)

AJ0)EY0)

此处J、Y分为第一类及第二类贝塞尔系数。它们在坐标轴上

曲线如下：

根据边界条件，当r=R时，R(x)=0,且当r=0时，R(x)不得为无穷大，看出只有J曲线符合要求，因此我们得到

R(x)AJ0)

所以将R(r)及Z(z)的解带入X(r)中，就得到：

_

X(r)J0(r)cosBZZ,其中

_



B

2

z

(



H

)

2

。用类似的方式我们可以求出（10）式

的解为

_

Y(r)

I0(r)cosBZZ

, 

_

其中/的值由反应堆材料

的性质定出。然后分别把径向与纵向两个方向的中子通量相叠加，我们就得出了堆芯双群的中子通量方程：

1,c(r,z)[AJ0)CI0(r)]cosBZZ(15)2,c(r,z)[AJ0)/s1CI0(r)/s2]cosBZZ(16)

_

_

再来看反射层中的中子通量方程：

1,r(r)k1,r1,r(r)0(5)2,r(r)k2,r(r)

2

21,r

2

2

12,rD2,r

1,r(r)0(6)

（5）为二阶齐次线性方程，具有类似于（9）式的解： 

1,r

(r)FZ1(r) 其中F为待定系数

（6）为非齐次方程，其解可写为通解与特解的线性组合： 2,r(r,z)GZ2(r)s31,r(r)GZ2(r)s3FZ1(r)（17）

其中S3为反射层的耦合系数，将（17）代入（6）中可得：

s31,r(r)s3k2,r1,r(r)

12,r

2

22

12,rD2,r

1

2

1,r(r)0

整理后得到：s3



D2,rk2,rk1,r

故而在反射层中，快群、热群的中子通量密度方程分别为：

1,r(r,z)FK0(k1,rr)cosBzz(18)2,c(r,z)GK0(k2,rr)cosBzz(19)

这样，我们就得到了堆芯及反射层里中子通量密度的

1,c,2,c,1,r,2,r分别的解，它们被（15）、（16）、（18）、（19）四个方程

所规定。其中包含四个未知数：A、C、F、G,下面就要利用边界条件来确定这四个待定系数的值，根据方程：我们得到如下方程组：

AXCYFZ10

s1AXs2CYs3FZ1GZ20AX'CY'1FZ1'0

s1AX's2CY's32FZ1'G2Z2'0

以上方程组是关于A、C、F、G的齐次方程组，根据克莱姆法则，为使其有非零解，即使A、C、F、G不同时为零，则其系数的判别式必须等于零，即：

X

s1XX's1X'

Ys2YY's2Y'

Z1s3Z11Z1'2s3Z1'

0Z202Z2'

0



通过求解方程，我们可以得到系数A、C、F、G之间的关系，然而该方程组为齐次线性方程组，它的解具有一待定常数因子，这可以理解为反应堆可以在任意功率下达到临界稳态运行，我们可以根据实际反应堆运行功率水平来定出这个常数因子，从而确定4个系数的具体值，然后根据有效增值系数Keff的表达式：

keff

k

(1LcB)(1cB)

2

2

2

其中反应堆几何曲率B的表达式为：

BAC

2

2

2

是两个曲率的线性组合，因为A、C具有一常数因子，在这里

规定A+C=1。所以当我们解出方程，也就得到反应堆临界判据的特征值Keff,而要解出方程组，则需要利用数值解法工具：PDQ、CITATION，在这里，我们经计算得到当R=0.27m Z=0.59m T

0.5

m时，我们设计

的理论反应堆可以达到临界条件，有关数值解法的问题，我们在下一讲中详细介绍。

二、为了使读者对堆芯设计得到一个更直观的概念，下面用较简单的单群理论介绍一种立方体石墨反应堆的设计过程：

首先定参数，设我们计划用浓缩度为20%的铀金属（楼主你太浪费了有木有啊）混合等体积的分析纯石墨制造一个立方体的反应堆。

按照中子通量公式： 2B20 规定

Bm

2

k1Lc

2

为材料曲率,因X,Y,Z为独立变量，故设：

c(x,y,z)X(x)Y(y)Z(z)

2222

=-B2

2

单独考虑X轴：

X(x)AcoBsx(C)xX(x)Acos(Bxx)

X(x)Bxx

2

，故

X(x)

有通解为：

C=0,故有：

，又因为立方体具有对称性，故Bxsxin(

系数A待定，同理可推得

c

(x,y,z)Acos(Bxx)cos(Byy)cos(Bzz)

设在Z方向设有反射层，在反射层内有：2M20

r(x,y,z)X(x)Y(y)Z(z)MXBYYBZZ

2

2

2

2

2

有

r(x,y,z)Ce(B

D

13s

Y

BZM)x

e

By

e

Bz

,其中M2=Lr2+为中子扩散面积，得通解：

(x,y,z)Cer

(BYBZM)x

e

By

e

Bz

，

C为待定系数。

因堆芯与反射层交界处两边中子通量及泄漏量必相等，有边界条件： 还是先考虑X方向c(a)r(a)

Dcc(x)

'

xa

Drr(x)

'

xa

带入可得：

Acos(Ba)Ce

(BYBZM)a

C

Ma

Acos(Ba)e

(BYBZM)a

DcABsin(Ba)DrMCe

DcABsin(Ba)DrMAcos(Ba)

tan(Ba)

DrMBD5



13s

将M、B的值被堆芯材料的性质唯一确定，D，通过查铀235

及石墨的热中子散射截面就可以得到临界边长a.

(2) 由（1）知： 易得：maxV=abc,

A

c

(x,y,z)Acos(Bxx)cos(Byy)cos(Bzz)

考虑到Bx值是堆芯的曲率，故有：

By/2Bz/2

ydycoszdzVdVABx/2cosBxdxBy/2cosBBz/2Bzxy

Bx/2

A2(

Bx



)2(

By



)2(

Bz



)ABxByBz(

AabcA(2

2





)

3

所以KH



max

1V





3



V

dv



)abc

3

8

不均匀系数KH为3.876

先进核动力反应堆的设计原理

核反应堆是将核子链式反应所造成的质量亏损转化为粒子动能、辐射能等形式能量的装置，继而通过外部的冷却循环系统及汽轮发电机组最终转化为我们给苹果手机充电所需要的电能。

从1941年到现在，人们设计出了各种各样的反应堆，有高温气冷堆、轻水堆、重水堆、钠冷快堆等，它的设计思路主要是基于核临界理论，在化石能源供应日趋紧张的今天，发展一种高效、安全、干净的核能已经成为历史发展的必然选择，目前，我国也建造了中国实验快堆(CFER)、中国先进研究堆（CARR）,并在先进核动力堆型及燃料后处理的领域上走在了世界前列，所采取的主要计算方法有CITATION程序、MATLAB方法等。下面，我们就结合理论与实际来介绍一种典型的核反应堆的设计过程。 二、 圆柱型轻水堆

一座核电站的设计图可以说是相当复杂，足够上万个人忙活几年的，但其核心就在于反应堆的临界参数，我们先采用较简单的双群理论来阐述一座理论堆型的设计原理。

取堆芯为质量分数为3%的二氧化铀，其余成分为铀238，反射层为H2O，快群与慢群能量分界按常规取1ev,相关中子反应截面如下:

1,ra0,1,rs3.22,ra0.022,2,rs3.451,ca0.04,1,cs3.752,ca0.5416,2,cs0.372

121,cf(E1E2)0.18

示意图如图一：

设堆芯与反射层的厚度分别为R、T，分别用角标c和r表示各项参数，如堆芯的中子慢化面积为Mc，反射层的中子慢化面积则为Mr，快群与热群的中子通量密度表示为：

1(r)2(r)



E0ECEC0

(r,E)dE1

(r,E)dE(2)



根据分析各种核子的中子反应截面我们可以看出，在轻水堆内，快群中子主要是由热中子引起的裂变产生的，

它又通过慢化吸收和泄

漏而消失，而热群中子则来源于快群中子的慢化，并主要由于吸收和泄漏而消失，很显然，在反应堆稳定工况时，各项参数都是恒定的，由中子平衡可以建立反应堆稳态时芯部的快群及热群的中子扩散方程如下：

D1,C1,C(r)r,c(r)

2

2

1Keff

[(f)1,c1,c(r)(f)2,c2,c(r)]3

D2,c2,c(r)a2,c2,c(r)12,c1,c(r)(4)

相应的，我们可以根据中子平衡关系写出反射层内中子扩散方程：

D2,r2,C(r)a2,r2,r(r)12,r1,r(r)

2

D1,r1,r(r)r,r1,r(r)0

2

为简便起见，我们通常可以将上述方程写成如下形式：

1,r(r)k1,r1,r(r)0(5)2,r(r)k2,r(r)

2

2

1,r

2

2

12,rD2,r

1,r(r)0(6)

k1,rk2,r

2

2

1,rD1,ra2,rD2,r

1Lr

2

其中

由(3)式可以得到快群中子通量密度表达式为：

1,c(r)

112,c

[D2,c2,c(r)a2,c2,c(r)]

2

将其代入（3）中，便得到只包含热群通量密度微分方程：

2,c

(r)的四阶偏



4

2,c(r)(

1

c



1L

2c

)2,c(r)

2

k1

'

cLc

2,c(r)0

它可以用因式分解的方法求解，将其改写成：

()()2,c0

2

2

2

2

（7）

同理有： 2X(r)2X(r)0

()()1,c0（8）

2

2

2

2

式中



2

1212

[(

1

c



1L

2c

)



2

[(

11Lc

2

c

)

方程(7)、(8)具有如下形式的解

X(r)X(r)0(9)Y(r)Y(r)0(10)

2

2

2

2

因其解X及Y为两个独立的函数，因而的解可以表示成两个函数的线性叠加：

1,cAX(r)CY(r)2,cAY(r)CY(r)

'

'

在这里A、C、A'、C'为四个待定常数，但从（4）式可以得到：

D2,cAX(r)a,2,cAX(r)12,cAX(r)

令S1=A/A=

’

2''

D2,ca,2,c

12,c

2

,相应地S2=C/C=

’

a,2,cD2,c

12,c

2

所以我们可以将快群、热群中子的解写成

1,cAX(r)CY(r)(11)2,c

AY(r)s1

'



CY(r)s2

'

(12)

其中A/C为待定常数，它们的值受反射层尺寸性质的约束。现在我们的任务就是根据（9）、（10）来定出X（r）、Y(r)的解。 先解（9）式，我们把X(r)写成R(r)Z(z)，将X(r)按圆柱坐标系展开得到：

X(r,z)r

2

2



1X(r,z)r

r



X(r,z)z

2

2

R(r,z)0

2

将含有径向与纵向变量的方程分离得到：

1R(r)

2

2

[

R(r)r

2

2



2

1R(r)r

r

]Bz(13)

22

Z(z)z

Bz(14)

其中BZ为待定常数。事实上，通过经验法求解（14），我们容易得到：

B(

2

z



H

),Z(z)cosBZZ

2

H为圆柱反应堆的高。 再来求解（13）

，令x见的零阶贝塞尔方程式：

x

2

,将其代入（13）中去，便得到常

dR(x)dx

2

2

x

dR(x)dx

xR(x)0

2

其普遍解为：

R(x)

AJ0)EY0)

此处J、Y分为第一类及第二类贝塞尔系数。它们在坐标轴上

曲线如下：

根据边界条件，当r=R时，R(x)=0,且当r=0时，R(x)不得为无穷大，看出只有J曲线符合要求，因此我们得到

R(x)AJ0)

所以将R(r)及Z(z)的解带入X(r)中，就得到：

_

X(r)J0(r)cosBZZ,其中

_



B

2

z

(



H

)

2

。用类似的方式我们可以求出（10）式

的解为

_

Y(r)

I0(r)cosBZZ

, 

_

其中/的值由反应堆材料

的性质定出。然后分别把径向与纵向两个方向的中子通量相叠加，我们就得出了堆芯双群的中子通量方程：

1,c(r,z)[AJ0)CI0(r)]cosBZZ(15)2,c(r,z)[AJ0)/s1CI0(r)/s2]cosBZZ(16)

_

_

再来看反射层中的中子通量方程：

1,r(r)k1,r1,r(r)0(5)2,r(r)k2,r(r)

2

21,r

2

2

12,rD2,r

1,r(r)0(6)

（5）为二阶齐次线性方程，具有类似于（9）式的解： 

1,r

(r)FZ1(r) 其中F为待定系数

（6）为非齐次方程，其解可写为通解与特解的线性组合： 2,r(r,z)GZ2(r)s31,r(r)GZ2(r)s3FZ1(r)（17）

其中S3为反射层的耦合系数，将（17）代入（6）中可得：

s31,r(r)s3k2,r1,r(r)

12,r

2

22

12,rD2,r

1

2

1,r(r)0

整理后得到：s3



D2,rk2,rk1,r

故而在反射层中，快群、热群的中子通量密度方程分别为：

1,r(r,z)FK0(k1,rr)cosBzz(18)2,c(r,z)GK0(k2,rr)cosBzz(19)

这样，我们就得到了堆芯及反射层里中子通量密度的

1,c,2,c,1,r,2,r分别的解，它们被（15）、（16）、（18）、（19）四个方程

所规定。其中包含四个未知数：A、C、F、G,下面就要利用边界条件来确定这四个待定系数的值，根据方程：我们得到如下方程组：

AXCYFZ10

s1AXs2CYs3FZ1GZ20AX'CY'1FZ1'0

s1AX's2CY's32FZ1'G2Z2'0

以上方程组是关于A、C、F、G的齐次方程组，根据克莱姆法则，为使其有非零解，即使A、C、F、G不同时为零，则其系数的判别式必须等于零，即：

X

s1XX's1X'

Ys2YY's2Y'

Z1s3Z11Z1'2s3Z1'

0Z202Z2'

0



通过求解方程，我们可以得到系数A、C、F、G之间的关系，然而该方程组为齐次线性方程组，它的解具有一待定常数因子，这可以理解为反应堆可以在任意功率下达到临界稳态运行，我们可以根据实际反应堆运行功率水平来定出这个常数因子，从而确定4个系数的具体值，然后根据有效增值系数Keff的表达式：

keff

k

(1LcB)(1cB)

2

2

2

其中反应堆几何曲率B的表达式为：

BAC

2

2

2

是两个曲率的线性组合，因为A、C具有一常数因子，在这里

规定A+C=1。所以当我们解出方程，也就得到反应堆临界判据的特征值Keff,而要解出方程组，则需要利用数值解法工具：PDQ、CITATION，在这里，我们经计算得到当R=0.27m Z=0.59m T

0.5

m时，我们设计

的理论反应堆可以达到临界条件，有关数值解法的问题，我们在下一讲中详细介绍。

二、为了使读者对堆芯设计得到一个更直观的概念，下面用较简单的单群理论介绍一种立方体石墨反应堆的设计过程：

首先定参数，设我们计划用浓缩度为20%的铀金属（楼主你太浪费了有木有啊）混合等体积的分析纯石墨制造一个立方体的反应堆。

按照中子通量公式： 2B20 规定

Bm

2

k1Lc

2

为材料曲率,因X,Y,Z为独立变量，故设：

c(x,y,z)X(x)Y(y)Z(z)

2222

=-B2

2

单独考虑X轴：

X(x)AcoBsx(C)xX(x)Acos(Bxx)

X(x)Bxx

2

，故

X(x)

有通解为：

C=0,故有：

，又因为立方体具有对称性，故Bxsxin(

系数A待定，同理可推得

c

(x,y,z)Acos(Bxx)cos(Byy)cos(Bzz)

设在Z方向设有反射层，在反射层内有：2M20

r(x,y,z)X(x)Y(y)Z(z)MXBYYBZZ

2

2

2

2

2

有

r(x,y,z)Ce(B

D

13s

Y

BZM)x

e

By

e

Bz

,其中M2=Lr2+为中子扩散面积，得通解：

(x,y,z)Cer

(BYBZM)x

e

By

e

Bz

，

C为待定系数。

因堆芯与反射层交界处两边中子通量及泄漏量必相等，有边界条件： 还是先考虑X方向c(a)r(a)

Dcc(x)

'

xa

Drr(x)

'

xa

带入可得：

Acos(Ba)Ce

(BYBZM)a

C

Ma

Acos(Ba)e

(BYBZM)a

DcABsin(Ba)DrMCe

DcABsin(Ba)DrMAcos(Ba)

tan(Ba)

DrMBD5



13s

将M、B的值被堆芯材料的性质唯一确定，D，通过查铀235

及石墨的热中子散射截面就可以得到临界边长a.

(2) 由（1）知： 易得：maxV=abc,

A

c

(x,y,z)Acos(Bxx)cos(Byy)cos(Bzz)

考虑到Bx值是堆芯的曲率，故有：

By/2Bz/2

ydycoszdzVdVABx/2cosBxdxBy/2cosBBz/2Bzxy

Bx/2

A2(

Bx



)2(

By



)2(

Bz



)ABxByBz(

AabcA(2

2





)

3

所以KH



max

1V





3



V

dv



)abc

3

8

不均匀系数KH为3.876