Numerical modeling of neutron transport in SP3 approximation by finite element method

03/27/2019 ∙ by Alexander V. Avvakumov, et al. ∙ rambler IBRAE 0

The SP3 approximation of the neutron transport equation allows improving the accuracy for both static and transient simulations for reactor core analysis compared with the neutron diffusion theory. Besides, the SP3 calculation costs are much less than higher order transport methods (SN or PN). Another advantage of the SP3 approximation is a similar structure of equations that is used in the diffusion method. Therefore, there is no difficulty to implement the SP3 solution option to the multi-group neutron diffusion codes. In this work, the application of the SP3 methodology based on solution of the λ- and α-spectral problems has been tested for the IAEA-2D and HWR reactor benchmark tests. The FEM is chosen to achieve the 3D geometrical generality, using GMSH as a generic mesh generator. The results calculated with the diffusion and SP3 methods are compared with the reference transport calculation results. It was found for the HWR reactor test that some eigenvalues are complex when calculating using both diffusion and SP3 options.

READ FULL TEXT VIEW PDF
POST COMMENT

Comments

There are no comments yet.

Authors

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

1 Introduction

The diffusion approximation of the neutron transport equation is widely used in nuclear reactor analysis allowing whole-core calculations with reasonable accuracy. The main feature of neutron diffusion equation is following: it is assumed that the neutron current is proportional to the neutron flux gradient (the Fick’s law). There are also three assumptions: neutron absorption much less likely than scattering, linear spatial variation of the neutron distribution and isotropic scattering (Stacey, 2007). To provide the validity of diffusion theory, the modern diffusion codes use, as a rule, assembly-by-assembly coarse-mesh calculation scheme with effective homogenized cross sections, prepared by more accurate transport approximations. To improve the diffusion code restrictions related with limitations on mesh spacing, different approaches are used including nodal and finite element methods (Avvakumov et al., 2017; Lawrence, 1986).

For many situations of interest (for instance, the pin-by-pin calculation taking account strongly absorbing control rods), the applicability of neutron diffusion theory are limited. Therefore, a more rigorous approximation for the neutron transport is required.

The solution of the neutron transport equation is very complicated problem because of seven independent variables: five for space-angular description, one for energy and one for time. To simplify the transport problem, different approaches are used such as the spherical harmonics () approximation (Azmy and Sartori, 2010). The approximation of the neutron transport equation is derived by expansion of the angular dependence of the neutron flux in the N spherical harmonics. During the last time, the simplest version of the method, namely the simplified approximation became widespread (McClarren, 2010). The major feature of the method is following: the three-dimensional neutron transport equation is transformed to a set of one-dimensional equations. The number of the trial functions is equal to 2(N+1) compared with the method which uses (N+1) trial functions. This leads to significant reduce in the computation time for typical whole-core calculations.

The approximation was first derived by Gelbard (Gelbard, 1960, 1961, 1962) in the early 1960s. He replaced the spatial derivatives with Laplacian and divergence operators in a one-dimensional planar geometry. The resulting equations are elliptic, for example, the

equations consist of two equations of diffusion type with two unknown fluxes: the scalar flux and the second angular flux moment. More rigorous theoretical foundation of the

methodology has been derived by Brantley and Larsen (Brantley and Larsen, 2000) on the basis of variational methods.

The method, as expected, can provide accuracy improvement compared with the common used diffusion method. Besides, implementation of the equations into the diffusion code is not difficult because of the similar structure of the and diffusion equations. For this reason the method was adopted in different whole-core calculation codes, such as DYN3D (Beckert and Grundmann, 2008), PARCS (Downar et al., 2010) and others. According to (Tada et al., 2008), application of the theory to the pin-by-pin calculation for BWR geometry resulted in remarkable improvement in the calculation accuracy compared with the diffusion method. Another report (Brewster, 2018) shows the comparison of the diffusion and methods to calculate the control rod reactivity in a light-water reactor. As compared with the Monte Carlo reference calculation, the method gives twice as accurate result compared with the diffusion method. Besides, as it turned out, the computation time using the method is only 1.5 times longer than that using the diffusion method (Tada et al., 2008).

Thus, the method can be considered as an improved approximation of the neutron transport equation compared with the diffusion method. In this regard, it will be very useful to compare the spectral parameters, calculated by both the diffusion and methods. To characterize the reactor steady-state conditions or dynamic behavior, some spectral problems are considered (Stacey, 2007; Bell and Glasstone, 1970). The steady-state condition is usually described by solution of a spectral problem (-eigenvalue problem); the fundamental eigenvalue (the largest eigenvalue) is called k-effective of the reactor core (Stacey, 2007; Bell and Glasstone, 1970). The reactor dynamic behavior can naturally be described on the basis of the approximate solution expansion in time-eigenvalue of -eigenvalue problem (Ginestar et al., 2002; Verdu et al., 2010; Verdu and Ginestar, 2014). At large times, one can talk about the asymptotic behavior of a neutron flux, whose amplitude is exp(t). Because the operator matrix is nonsymmetric, some of its eigenvalues can be complex (Aragones et al., 2007)

. Previously the complex eigenvalues and eigenfunctions were found in the spectral problems for some numerical tests

(Avvakumov et al., 2017).

In this paper we consider the approximation for the steady-state multi-group neutron transport problem. To solve spectral problems with nonsymmetrical matrices we use well-designed algorithms and relevant free software including the library SLEPc (Scalable Library for Eigenvalue Problem Computations, http://slepc.upv.es/). We use a Krylov-Schur algorithm, a variation of Arnoldi method, described in (Stewart, 2001).

The paper is organized as follows. The steady-state and dynamic models of a nuclear reactor based on the multigroup equations are given in Section 2. In Section 3 we discuss various spectral problems. Some numerical examples of calculation of spectral characteristics of two-dimensional test problems for IAEA-2D benchmark problem and HWR reactor using the two-group system of diffusion and equations is discussed in Section 4. The results of the work are summarized in Section 5.

2 Problem statement

Let’s consider the symmetric form of the equation for the neutron flux (Ryu and Joo, 2013). The neutron dynamics is considered in the limited convex two-dimensional or three-dimensional area () with boundary . The neutron transport is described by the system of equations

(1)

where

Here — number of energy groups, — scalar flux, — pseudo 0th moment of angular flux, — second moment of angular flux, — total cross-section, — transport cross-section, — removal cross-section, — scattering cross-section, — spectra of neutrons, — generation cross-section, — density of sources of delayed neutrons, — decay constant of sources of delayed neutrons, — number of types of delayed neutrons.

The density of sources of delayed neutrons is described by the equations

(2)

where is the fraction of delayed neutrons of m-type, and

The Marshak-type conditions are set at the boundary of the area :

(3)

System of equations (1) and (2) is supplemented with boundary conditions (3) and corresponding initial conditions:

(4)

Let’s write the boundary problem (1)–(4

) in operator form. The vectors

, , and matrices are defined as follows

where

is the Kronecker symbol. We shall use the set of vectors , whose components satisfy the boundary conditions (3). Using the set definitions, the system of equations (1) and (2) can be written as following

(5)

Without taking into account delayed neutrons (all neutrons are considered as prompt), we have

(6)

The Cauchy problem is formulated for equations (6) when

(7)

where , . Для уравнений (5) задается также начальное условие

(8)

where .

3 Spectral problems

To characterize the reactor dynamic processes described by Cauchy problem (5)-(7), let’s consider some spectral problems (Bell and Glasstone, 1970; Stacey, 2007).

The spectral problem, which is known as the -spectral problem, is usually considered. For the system of equations (6), (7), we have

(9)

where

The minimal eigenvalue is used for characterisation of neutron field, thus

is the effective multiplication factor (k-effective). The value is related to the critical state of the reactor, and the corresponding eigenfunction is the stationary solution of the Eq (5), (6). At , one can speak about supercriticality, at — about subcriticality.

The spectral problem (9) cannot directly be connected with the dynamic processes in a nuclear reactor. The eigenvalues of the multiplication factor of the reactor and the corresponding eigenfunctions do not depend on the time delay for the emission of delayed neutrons. The reason is that the problem (9) on eigenvalues is the problem of finding time-independent solutions of the neutron transport equation, and the term describing the contribution of fission to the neutron balance is equal to the total number of fission neutrons, both instantaneous and delayed divided by . At the best, we can get only the limiting case — the stationary critical state. The more acceptable spectral characteristics for the non-stationary equation (5) are related the spectral problem

(10)

The -eigenvalue problem without delayed neutrons (6) can be written as follows

(11)

where

The fundamental eigenvalue

is called (Bell and Glasstone, 1970) the -eigenvalue or the period eigenvalue, because it is inversely related to the reactor period. The problem of the period eigenvalues essentially takes into account the contribution of delayed neutrons. In particular, the long lifetime of the predecessors of delayed neutrons makes a large contribution to the slowly decreasing eigenfunctions of the reactor period, and this does not occur when only instantaneous neutrons are taken into account.

The asymptotic behaviour of Cauchy problem solution (5)-(7) at large times can be connected with the eigenvalue . In this regular mode, the reactor behaviour is described by the function . If , then the reactor is critical; if , then we get the neutron flux decreasing (subcritical state), and if , then we get the neutron flux increasing (supercritical state).

There is a simple approximate relationship between the dominant eigenvalues of the -spectral (-spectral) and -spectral problems without delayed neutrons (Verdu et al., 2010):

(12)

where — the prompt neutron generation time and — the prompt neutron lifetime. The prompt neutron generation time is determined as the neutron flux functional using adjoint fundamental -eigenfunction rather then the -eigenfunction (Verdu et al., 2010).

For the spectral problems with delayed neutrons, the relationship is as folows

(13)

which corresponds to inhour’s equation. It is obvious that if we consider small perturbations of reactivity (), then the dominant -eigenvalue does not practically depend on and is related only with delayed neutrons parameters and . Under such assumptions one can expect that the shapes of the fundamental eigenfunctions of different spectral problems are similar and, therefore, the relationships (12) and (13) are close to exact solution. In this case, we can derive the following relationship between and .

(14)

where the second term in the rihght part of the equation (14) can be neglected.

It is supposed that the relationships (12)-(14) are valid for the rest eigenvalues of the -spectral and -spectral problems (Verdu et al., 2010).

4 Numerical examples

To study the properties of the eigenvalues and eigenfunctions of dufferent types, several benchmarks are studied. The first benchmark is the IAEA-2D two dimensional hexagonal problem of VVER-type without reflector (Chao and Shatilla, 1995). The second benchmark is the similar test but with radial reflector (Chao and Shatilla, 1995). The aim of the third test is to investigate azimutally non-symmetric effects on the eigenvalues and eigenfunctions (the changed IAEA-2D test with reflector). Finally, we studied the complex eigenvalues and eigenfunctions for the heavy water hexagonal reactor test HWR (Chao and Shatilla, 1995).

The two-group model () is used in all tests. The method of finite elements (Brenner and Scott, 2008; Quarteroni and Valli, 2008) on triangular calculation grids is used for the approximate solution of the spectral problem. The standard Lagrangian finite elements are used. The software has been developed using the engineering and scientific calculation library FEniCS (Logg et al., 2012). SLEPc has been used for numerical solution of the spectral problems. We used a Krylov-Schur algorithm with an accuracy of . The following parameters were varied in the calculations:

  • — the number of triangles per one assembly (Fig. 1);

  • — the order of finite element.

In this work we compare the calculations with the previous diffusion model calculations (Avvakumov et al., 2014, 2017).

Figure 1: Discretization of assembly into 6, 24 and 96 finite elements.

4.1 IAEA-2D without reflector

The geometrical model of the IAEA-2D reactor core (Chao and Shatilla, 1995) consists of a set of hexagonal assemblies and is presented in Fig. 2, where the assemblies of various types are marked with various digits. The total size of assembly equals 20 cm.


Figure 2: Geometrcial model of the IAEA-2D reactor core without reflector.

Diffusion neutronics constants in the common units are given in Table 1. The following delayed neutrons parameters are used: one group of delayed neutrons with effective fraction and decay constant s. Neutron velocity cm/s and cm/s.

Material 1 2 3 4
1.5 1.5 1.5 1.5
0.4 0.4 0.4 0.4
0.01 0.01 0.01 0.0
0.08 0.085 0.13 0.01
0.02 0.02 0.02 0.04
0.1922222 0.1922222 0.1922222 0.1822222
0.7533333 0.7483333 0.7033333 0.8233333
0.00 0.00 0.00 0.00
0.135 0.135 0.135 0.00
Table 1: Diffusion neutronics constants for IAEA-2D.

4.1.1 Solution of Lambda Modes spectral problem

As a reference solution for the diffusion model, we used the previous results (Avvakumov et al., 2014); for the model — the solution obtained using very fine mesh (). The maximum difference in assembly power between two models is about 2 percent for the rodded assemblies (material 3, see Fig. 2).

The results of the solution of the effective multiplication factor for test IAEA-2D without a reflector are shown in Table 2. Hereinafter, for -spectral problems, the following notation is used: — effective multiplication factor for the diffusion model; — effective multiplication factor for the model; — absolute deviation from the reference value in pcm ();

— the standard deviation of the relative power in percent. These data demonstrate the convergence of the computed eigenvalues with refinement of the calculation mesh and increase in polynomial degree.

The results of the first 10 eigenvalues for are shown in Table 3. The power and error distributions for the diffusion and models are presented in Figs 3 and 4 for . Hereinafter, for each assembly the following data are given: the reference solution (the diffusion or model), the solution for and the relative error from the reference solution.

1 0.97335 473 3.80 0.97445 490 4.02
6 2 0.97760 48 0.45 0.97881 54 0.52
3 0.97801 7 0.07 0.97925 10 0.09
1 0.97654 154 1.28 0.97772 163 1.38
24 2 0.97799 9 0.08 0.97923 12 0.11
3 0.97807 1 0.01 0.97934 1 0.02
1 0.97765 43 0.36 0.97888 47 0.40
96 2 0.97807 1 0.02 0.97933 2 0.02
3 0.97808 0 0.01 0.97935
Ref. 0.97808 0.97935
Table 2: The effective multiplication factor.
Diffusion SP
1 0.97808 + 0.0 0.979351 + 0.0
2 0.96318 + 0.0 0.964604 + 0.0
3 0.96318 + 0.0 0.964604 + 0.0
4 0.93844 + 0.0 0.940253 + 0.0
5 0.93844 + 0.0 0.940253 + 0.0
6 0.91966 + 0.0 0.921844 + 0.0
7 0.90220 + 0.0 0.904467 + 0.0
8 0.87141 + 0.0 0.874997 + 0.0
9 0.84957 + 0.0 0.853155 + 0.0
10 0.84957 + 0.0 0.853154 + 0.0
Table 3: The eigenvalues for .
Figure 3: Power and error distributions using the diffusion model.
Figure 4: Power and error distributions using the model.

4.1.2 Solution of -spectral problem without delayed neutrons

As a reference solution for both the diffusion and transport models, we used the solutions obtained using very fine mesh (). Hereinafter, for -spectral problems, the following notation is used: -eigenvalue by diffusion model; -eigenvalue by model; — absolute deviation from the reference value.

The calculation results for the -spectral problem without delayed neutrons using the different meshes and the finite element orders are shown in Table 4. These data demonstrate the convergence of approximate computed eigenvalues with refinement of the calculation mesh and increase in polynomial degree.

1 556.3 100.8 532.7 104.1
6 2 465.6 10.1 440.0 11.4
3 457.0 1.5 430.7 2.1
1 488.1 32.6 463.0 34.4
24 2 457.4 1.9 431.0 2.4
3 455.7 0.2 428.9 0.3
1 464.6 9.1 438.4 9.8
96 2 455.8 0.3 428.9 0.3
3 455.5 428.6
Ref. 455.5 428.6
Table 4: The - eigenvalues.

The spectral problem results for the first 10 eigenvalues are shown in Table 5. The eigenvalues are well separated. In this example, the fundamental eigenvalue is less compared the rest and therefore the main harmonic will attenuate more slowly. A regular mode of the reactor is thereby defined. The value determines the amplitude of neutron flux and is connected directly with reactor period in the regular mode.

Diffusion SP
1 455.540 + 0.0 428.561 + 0.0
2 760.532 + 0.0 730.398 + 0.0
3 760.543 + 0.0 730.408 + 0.0
4 1267.192 + 0.0 1228.835 + 0.0
5 1267.192 + 0.0 1228.836 + 0.0
6 1647.145 + 0.0 1601.437 + 0.0
7 2083.289 + 0.0 2031.778 + 0.0
8 2696.887 + 0.0 2616.862 + 0.0
9 3188.356 + 0.0 3092.715 + 0.0
10 3188.363 + 0.0 3092.722 + 0.0
Table 5: The eigenvalues for .

The eigenfunctions for fundamental eigenvalue () of the -spectral problem without delayed neutron are shown in Fig. 5. Due to the fact that the reactor state is close to critical (), the fundamental eigenfunctions of the -spectral problem are close to the fundamental eigenfunctions of the -spectral problem. The eigenfunctions are shown in Fig. 6, Fig. 7.

Figure 5: Eigenfunctions , using model.
Figure 6: Eigenfunctions , using model.
Figure 7: Eigenfunctions , using model.

4.1.3 Solution of -spectral problem with delayed neutrons

As a reference solution for both the diffusion and transport models, we used the solutions obtained using very fine mesh (). The -spectral problem results with delayed neutrons using the different meshes and finite element orders are shown in Table 6. Compared with the previous case without delayed neutrons, these data demonstrate the similar convergence of the computed eigenvalues.

1 0.06465 0.00264 0.06410 0.00295
6 2 0.06232 0.00031 0.06153 0.00035
3 0.06206 0.00005 0.06122 0.00007
1 0.06296 0.00095 0.06224 0.00109
24 2 0.06207 0.00005 0.06123 0.00008
3 0.06202 0.00001 0.06116 0.00001
1 0.06228 0.00027 0.06147 0.00032
96 2 0.06202 0.00001 0.06116 0.00001
3 0.06201 0.06115
Ref. 0.06201 0.06115
Table 6: The -eigenvalues.
Diffusion SP
1 0.06201 + 0.0 0.06115 + 0.0
2 0.06837 + 0.0 0.06796 + 0.0
3 0.06837 + 0.0 0.06796 + 0.0
4 0.07279 + 0.0 0.07258 + 0.0
5 0.07279 + 0.0 0.07258 + 0.0
6 0.07446 + 0.0 0.07430 + 0.0
7 0.07547 + 0.0 0.07536 + 0.0
8 0.07662 + 0.0 0.07652 + 0.0
9 0.07717 + 0.0 0.07709 + 0.0
10 0.07721 + 0.0 0.07711 + 0.0
Table 7: The eigenvalues for .

The first 10 spectral problem eigenvalues are shown in Table 7. Due to the contribution of delayed neutrons, the fundamental eigenvalue is much smaller compared with the case without delayed neutrons.

The eigenfunctions for fundamental eigenvalue () of the -spectral problem with delayed neutron are shown in Fig. 8. The eigenfunctions are shown in Fig. 9, Fig. 10. The eigenfunctions of the problems without and with delayed neutrons are close to each other in topology.

Figure 8: Eigenfunctions , .
Figure 9: Eigenfunctions , .
Figure 10: Eigenfunctions , .

According to (12)-(14), the prompt neutron generation time for both the diffusion and fundamental eigenvalues. Assuming the same value of for both -spectral problems, the uncertainty in the fundamental eigenvalues is less than 0.1 percent.

4.2 IAEA-2D with reflector

This test differs from the previous one only additional external row of reflector assemblies (material 4, see Table 1).

4.2.1 Solution of Lambda Modes spectral problem

As a reference solution for the diffusion model, we used the previous results (Avvakumov et al., 2015); for the model — the solution obtained using the MCNP4C code (Bahabadi et al., 2016). As well as in the previous benchmark calculations, the maximum difference in assembly power between two models is about 2 percent for the rodded assemblies (material 3, see Fig. 2).

The comparison of the calculated effective multiplication factors is shown in Table 8. The results of the first 10 eigenvalues for are presented in Table 9. One can see that there are several eigenvalues of multiplicity two.

The power distributions and calculation errors for using the diffusion model are shown in Fig 11 and for using the model are shown Fig 12.

1 1.01041 490 13.29 1.01159 536 14.14
6 2 1.00623 72 1.88 1.00711 88 2.19
3 1.00558 7 0.22 1.00636 13 0.35
1 1.00699 148 4.54 1.00792 169 4.96
24 2 1.00561 10 0.30 1.00640 17 0.42
3 1.00551 0 0.02 1.00626 3 0.17
1 1.00591 36 1.28 1.00671 48 1.42
96 2 1.00552 1 0.04 1.00626 3 0.18
3 1.00551 0 0.01 1.00625 2 0.18
Ref. 1.00551 1.00623
Table 8: The effective multiplication factor.
Diffusion SP
1 1.005510 + 0.0 1.006245 + 0.0
2 0.996490 + 0.0 0.997254 + 0.0
3 0.996490 + 0.0 0.997254 + 0.0
4 0.976791 + 0.0 0.977759 + 0.0
5 0.976791 + 0.0 0.977759 + 0.0
6 0.958684 + 0.0 0.959895 + 0.0
7 0.928980 + 0.0 0.930969 + 0.0
8 0.924186 + 0.0 0.925931 + 0.0
9 0.904788 + 0.0 0.907349 + 0.0
10 0.904788 + 0.0 0.907349 + 0.0
Table 9: The eigenvalues for .
Figure 11: Power and error distributions using the diffusion model.
Figure 12: Power and error distributions using the model.

4.2.2 Solution of -spectral problem without delayed neutrons

As a reference solution, we use the fine mesh solutions obtained using the diffusion or transport model (). The -spectral problem results are shown in Table 10.

1 184.95 84.14 205.92 91.32
6 2 113.58 12.77 130.02 15.42
3 101.98 1.17 116.72 2.12
1 126.66 25.85 143.85 29.25
24 2 102.58 1.77 117.31 2.71
3 100.88 0.07 114.83 0.23
1 107.82 7.01 122.84 8.24
96 2 100.97 0.16 114.94 0.34
3 100.81 114.60
Ref. 100.81 114.60
Table 10: The -eigenvalues.

The results of the first 10 eigenvalues for are presented in Table 9. As before, the eigenvalues are well separated. In this example, the fundamental eigenvalue is negative and therefore the main harmonic will increase, while all others will attenuate.

Diffusion SP
1 100.81 + 0.0 114.60 + 0.0
2 62.93 + 0.0 49.42 + 0.0
3 62.93 + 0.0 49.42 + 0.0
4 405.31 + 0.0 390.15 + 0.0
5 405.31 + 0.0 390.15 + 0.0
6 710.64 + 0.0 693.47 + 0.0
7 1141.43 + 0.0 1118.67 + 0.0
8 1469.68 + 0.0 1438.31 + 0.0
9 1494.37 + 0.0 1468.54 + 0.0
10 1494.37 + 0.0 1468.54 + 0.0
Table 11: The eigenvalues for .

The eigenfunctions for fundamental eigenvalue () of the -spectral problem without delayed neutrons are shown in Fig. 13. Due to the fact that a state of the reactor is close to critical (), the fundamental eigenfunctions of the -spectral problem are close to the fundamental eigenfunctions of the -spectral problem. The eigenfunctions are shown in Fig. 14, Fig. 15.

Figure 13: Eigenfuncions , .
Figure 14: Eigenfuncions , .
Figure 15: Eigenfuncions , .

4.2.3 Solution of -spectral problem with delayed neutrons

As a reference solution, we use the fine mesh solutions obtained using the diffusion or transport model (). The -spectral problem results are shown in Table 12.

1 68.2268 67.8084 88.9461 87.6086
6 2 1.2810 0.8626 11.1554 9.8179
3 0.4506 0.0322 1.8063 0.4688
1 9.0267 8.6083 25.1658 23.8283
24 2 0.4686 0.0502 1.9832 0.6457
3 0.4202 0.0018 1.3787 0.0412
1 0.7018 0.2834 4.9794 3.6419
96 2 0.4225 0.0041 1.3994 0.0619
3 0.4184 1.3375
Ref. 0.4184 1.3375
Table 12: The -eigenvalues.
Diffusion SP
1 0.4184 + 0.0 1.3375 + 0.0
2 0.0281 + 0.0 0.0238 + 0.0
3 0.0281 + 0.0 0.0238 + 0.0
4 0.0628 + 0.0 0.0622 + 0.0
5 0.0628 + 0.0 0.0622 + 0.0
6 0.0695 + 0.0 0.0692 + 0.0
7 0.0737 + 0.0 0.0736 + 0.0
8 0.0741 + 0.0 0.0740 + 0.0
9 0.0754 + 0.0 0.0752 + 0.0
10 0.0763 + 0.0 0.0762 + 0.0
Table 13: The eigenvalues for .

The calculation results for the first 10 eigenvalues are shown in Table 13. Due to the contribution of delayed neutrons, the fundamental eigenvalue is much smaller than in the case without delayed neutrons. Again the fundamental eigenvalue is negative and therefore the main harmonic will increase, while all others will attenuate.

The eigenfunctions for fundamental eigenvalue () of the -spectral problem with delayed neutrons are presented in Fig. 16. The eigenfunctions are shown in Fig. 17, Fig. 18. The eigenfunctions of the problems without and with delayed neutrons are close to each other in topology.

Figure 16: Eigenfunctions , .
Figure 17: Eigenfunctions , .
Figure 18: Eigenfunctions , .

According to (12)-(14), the prompt neutron generation time for the diffusion fundamental eigenvalue and for the fundamental eigenvalue. Assuming the same value of for both -spectral problems, the uncertainty in the fundamental eigenvalues is less than 0.1 percent.

4.3 Azimutally non-symmetric test IAEA-2D with reflector

To investigate azimutally non-symmetric geometry effects on the eigenfunction behaviour, we replaced two unrodded assemblies in the north-east part of the core by rodded ones (material 3, see Fig. 19).

Figure 19: Geometrcial model of the non-symmetric test IAEA-2D.

4.3.1 Solution of Lambda Modes spectral problem

As a reference solution for both the diffusion and transport models, we used the solutions obtained using very fine mesh (). The effective multiplication factors are shown in Table 14. The results of the first 10 eigenvalues for are presented in Table 15. As can be see from Table 15 compared with Table 9 for unperturbed case, all eigenvalues became well separated (pairs of the eigenvalues are vanished).

1 1.00809 509 1.00931 556
6 2 1.00374 74 1.00465 90
3 1.00306 6 1.00387 12
1 1.00454 154 1.00550 175
24 2 1.00310 10 1.00391 16
3 1.00300 0 1.00376 1
1 1.00341 41 1.00424 49
96 2 1.00300 0 1.00377 2
3 1.00300 1.00375
Ref. 1.00300 1.00375
Table 14: The effective multiplication factor.
diffusion SP
1 1.002996 + 0.0 1.003751 + 0.0
2 0.994571 + 0.0 0.995365 + 0.0
3 0.986297 + 0.0 0.987243 + 0.0
4 0.970315 + 0.0 0.971407 + 0.0
5 0.968980 + 0.0 0.970207 + 0.0
6 0.945551 + 0.0 0.947166 + 0.0
7 0.928439 + 0.0 0.930441 + 0.0
8 0.923863 + 0.0 0.925611 + 0.0
9 0.903265 + 0.0 0.905868 + 0.0
10 0.901593 + 0.0 0.904253 + 0.0
Table 15: The eigenvalues for .

4.3.2 Solution of -spectral problem without delayed neutrons

As a reference solution for both the diffusion and models, we used the solutions obtained using very fine mesh (). The -spectral problem results are shown in Table 16. The results of the first 10 eigenvalues for are presented in Table 17.

1 143.12 88.55 164.63 96.08
6 2 67.99 13.42 84.75 16.20
3 55.82 1.25 70.78 2.23
1 81.93 27.36 99.47 30.92
24 2 56.45 1.88 71.41 2.86
3 54.65 0.08 68.80 0.25
1 62.00 7.43 77.28 8.73
96 2 54.74 0.17 68.91 0.36
3 54.57 68.55
Ref. 54.57 68.55
Table 16: The -eigenvalues.
Diffusion SP
1 -54.57 + 0.0 -68.55 + 0.0
2 97.07 + 0.0 83.18 + 0.0
3 242.22 + 0.0 226.42 + 0.0
4 513.07 + 0.0 496.61 + 0.0
5 530.98 + 0.0 512.74 + 0.0
6 898.88 + 0.0 878.37 + 0.0
7 1148.46 + 0.0 1125.66 + 0.0
8 1481.13 + 0.0 1449.58 + 0.0
9 1512.16 + 0.0 1486.05 + 0.0
10 1527.83 + 0.0 1501.83 + 0.0
Table 17: The eigenvalues for .

Let’s consider changes in the eigenfunctions due to the rodded assembly insertion. The eigenfunctions for fundamental eigenvalue () of the -spectral problem without delayed neutrons are shown in Fig. 20. The eigenfunctions are shown in Fig. 21, Fig. 22. As can be seen from Fig. 20 to Fig. 22, the overal structure of the eigenfunctions is preserved taking into account the neutron flux perturbations.

Figure 20: Eigenfunctions , .
Figure 21: Eigenfunctions , .
Figure 22: Eigenfunctions , .

According to (12)-(14), the prompt neutron generation time for the diffusion fundamental eigenvalue and for the fundamental eigenvalue. Thus, for the non-symmetric test we obtained similar neutronic properties compared with the symmetric test. The corresponding fundamental eigenvalues of the -spectral problem with delayed neutrons are calculated using (14): and .

4.4 HWR test problem

This benchmark is a model of large heavy-water reactor HWR (Chao and Shatilla, 1995). The geometry of the HWR test is presented in Fig. 23. Fuel assemblis (1, 2, 3 and 6 in Fig. 23) located in the central part of the core, are surrouded by the target zone and reflector layer (7 and 9 in Fig. 23). There are two types of rodded assemblies (4 and 8). The assembly size is equal to 17.78 cm.

Diffusion constants are given in Table 18. The following delayed neutrons parameters are used: one group of delayed neutrons with effective fraction and decay constant s. Neutron velocity cm/s and cm/s.

Material Group , cm , cm , cm , cm
1 1 1.38250058 1.1105805e-2 8.16457e-3 2.26216e-3
2 0.89752185 2.2306487e-2 2.30623e-2
2 1 1.38255219 1.1174585e-2 8.22378e-3 2.22750e-3
2 0.89749043 2.2387609e-2 2.26849e-2
3 1 1.37441741 1.0620368e-2 8.08816e-3 2.14281e-3
2 0.88836771 1.6946527e-2 2.04887e-2
4 1 1.31197955 1.2687953e-2 1.23115e-2 0.0
2 0.87991376 5.2900925e-2 0.0
6 1 1.38138909 1.056312e-2 7.76568e-3 2.39469e-3
2 0.90367052 2.190298e-2 2.66211e-2
7 1 1.30599110 1.1731321e-2 1.10975e-2 0.0
2 0.83725587 4.3330365e-3 0.0
8 1 1.29192957 1.1915316e-2 1.15582e-2 0.0
2 0.81934103 3.0056488e-4 0.0
9 1 1.06509884 2.8346221e-2 2.61980e-2 0.0
2 0.32282849 3.3348874e-2 0.0
Table 18: Diffusion constants for HWR test.
Figure 23: Geometrcial model of the HWR test.

4.4.1 Solution of Lambda Modes spectral problem

1 0.991985 2.0 1.16 0.992178 5.0 0.80
6 2 0.991989 2.4 0.31 0.992166 3.8 0.24
3 0.991964 0.1 0.08 0.992132 0.4 0.07
1 0.991983 1.8 0.05 0.992165 3.7 0.08
24 2 0.991965 0.0 0.01 0.992133 0.5 0.01
3 0.991963 0.2 0.01 0.992128 0.0 0.00
1 0.991969 0.4 0.08 0.992140 1.2 0.01
96 2 0.991963 0.2 0.02 0.992129 0.1 0.00
3 0.991963 0.2 0.01 0.992128
Ref. 0.991965 0.992128
Table 19: The effective multiplication factor.

As a reference solution for the diffusion model, we used the results obtained by (Chao and Shatilla, 1995); for the model — the solution obtained using very fine mesh ().

The effective multiplication factors for the HWR test are shown in Table 19.

diffusion SP
1 0.991963 + 0.0 0.992128 + 0.0
2 0.983594 + 1.1645e-05 0.983793 + 1.2072e-05
3 0.983594 1.1645e-05 0.983793 1.2072e-05
4 0.964240 + 2.1564e-05 0.964523 + 2.2337e-05
5 0.964240 2.1564e-05 0.964523 2.2337e-05
6 0.943290 + 0.0 0.943733 + 0.0
7 0.923872 + 0.0 0.924257 + 0.0
8 0.918657 + 0.0 0.918798 + 0.0
9 0.895682 + 3.5570e-05 0.896317 + 3.6750e-05
10 0.895682 3.5570e-05 0.896317 + 3.6750e-05
Table 20: The eigenvalues for .

The results of the first 10 eigenvalues for are presented in Table 20. The eigenvalues of the -spectral problem are the complex values with small imaginary parts, and the eigenvalues are the real values. One can see pairs of the complex conjugate values.

4.4.2 Solution of -spectral problem without delayed neutrons

As a reference solution for the diffusion and models we use the solutions obtained using very fine mesh ().

The -spectral problem results at different computational parameters are shown in Table 21. The results of the first 10 eigenvalues for are presented in Table 22. The eigenvalues are well separated. The eigenvalues , , of the -spectral problem, like for the -spectral problem, are the complex values with small imaginary parts, and the eigenvalues , are the real values.

The eigenfunctions for fundamental eigenvalue () of the -spectral problem are shown in Fig. 24. The real part of the eigenfunctions is shown in Fig. 25. Fig. 26 shows the imaginary part of these eigenfunctions. The eigenfunctions of the -spectral and -spectral problems are close to each other in topology.

1 42.281 0.018 41.246 0.134
6 2 42.135 0.128 41.190 0.190
3 42.259 0.004 41.362 0.018
1 42.196 0.067 41.228 0.152
24 2 42.253 0.010 41.354 0.026
3 42.263 0.000 41.379 0.001
1 42.241 0.022 41.330 0.050
96 2 42.262 0.001 41.377 0.003
3 42.263 41.380
Ref. 42.263 41.380
Table 21: The -eigenvalues.
Diffusion SP
1 042.263 + 0.0 41.380 + 0.0
2 084.867 0.06130 83.821 0.06358
3 084.867 + 0.06130 83.821 + 0.06358
4 182.914 0.11367 181.471 0.11805
5 182.914 + 0.11367 181.471 + 0.11805
6 293.017 + 0.0 290.940 + 0.0
7 371.528 + 0.0 369.374 + 0.0
8 515.465 0.16397 512.337 0.17197
9 515.465 + 0.16397 512.337 + 0.17197
10 518.670 + 0.0 517.975 + 0.0
Table 22: The eigenvalues for .
Figure 24: The eigenfunctions (left) and (right).
Figure 25: Real part of eigenfunctions (left) and (right).
Figure 26: Imaginary part of eigenfunctions (left) and (right).

4.4.3 Solution of -spectral problem with delayed neutrons

As a reference solution for the diffusion and models we use the solutions obtained using very fine mesh ().

The -spectral problem results are shown in Table 23. Due to the contribution of delayed neutrons, the fundamental eigenvalue is much smaller than in the case without delayed neutrons. The results of the first 10 eigenvalues for is shown in Table 24. The eigenvalues , , of the -spectral problem, like as before, are the complex values with small imaginary parts, and the eigenvalues , are the real values.

1 0.04431 0.00006 0.04383 0.00012
6 2 0.04430 0.00007 0.04386 0.00009
3 0.04437 0.00000 0.04394 0.00001
1 0.04432 0.00005 0.04386 0.00009
24 2 0.04436 0.00001 0.04394 0.00001
3 0.04437 0.00000 0.04395 0.00000
1 0.04435 0.00002 0.04392 0.00003
96 2 0.04437 0.00000 0.04395 0.00000
3 0.04437 0.04395
Ref. 0.04437 0.04395
Table 23: The -eigenvalues.
Diffusion SP
1 0.04437 + 0.0 0.04395 + 0.0
2 0.05755 1.15549e-05 0.05735 1.22333e-05
3 0.05755 + 1.15549e-05 0.05735 + 1.22333e-05
4 0.06807 6.35264e-06 0.06798 6.66947e-06
5 0.06807 + 6.35264e-06 0.06798 + 6.66947e-06
6 0.07219 + 0.0 0.07213 + 0.0
7 0.07415 + 0.0 0.07412 + 0.0
8 0.07453 + 0.0