1 Introduction
The diffusion approximation of the neutron transport equation is widely used in nuclear reactor analysis allowing wholecore 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, assemblybyassembly coarsemesh 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 pinbypin 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 spaceangular 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 threedimensional neutron transport equation is transformed to a set of onedimensional 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 wholecore 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 onedimensional 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 wholecore 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 pinbypin 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 lightwater 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 steadystate conditions or dynamic behavior, some spectral problems are considered (Stacey, 2007; Bell and Glasstone, 1970). The steadystate condition is usually described by solution of a spectral problem (eigenvalue problem); the fundamental eigenvalue (the largest eigenvalue) is called keffective 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 timeeigenvalue 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 steadystate multigroup neutron transport problem. To solve spectral problems with nonsymmetrical matrices we use welldesigned algorithms and relevant free software including the library SLEPc (Scalable Library for Eigenvalue Problem Computations, http://slepc.upv.es/). We use a KrylovSchur algorithm, a variation of Arnoldi method, described in (Stewart, 2001).
The paper is organized as follows. The steadystate 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 twodimensional test problems for IAEA2D benchmark problem and HWR reactor using the twogroup 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 twodimensional or threedimensional 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 crosssection, — transport crosssection, — removal crosssection, — scattering crosssection, — spectra of neutrons, — generation crosssection, — 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 mtype, and
The Marshaktype 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 followswhere
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 (keffective). 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 timeindependent 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 nonstationary 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.
4 Numerical examples
To study the properties of the eigenvalues and eigenfunctions of dufferent types, several benchmarks are studied. The first benchmark is the IAEA2D two dimensional hexagonal problem of VVERtype 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 nonsymmetric effects on the eigenvalues and eigenfunctions (the changed IAEA2D test with reflector). Finally, we studied the complex eigenvalues and eigenfunctions for the heavy water hexagonal reactor test HWR (Chao and Shatilla, 1995).
The twogroup 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 KrylovSchur 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).
4.1 IAEA2D without reflector
The geometrical model of the IAEA2D 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.
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 
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 IAEA2D 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 
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 
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 
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 
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.
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 
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 
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.
4.2 IAEA2D 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 
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 
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 
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 
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.
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 
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 
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.
4.3 Azimutally nonsymmetric test IAEA2D with reflector
To investigate azimutally nonsymmetric geometry effects on the eigenfunction behaviour, we replaced two unrodded assemblies in the northeast part of the core by rodded ones (material 3, see Fig. 19).
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 
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 
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 
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 
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.
According to (12)(14), the prompt neutron generation time for the diffusion fundamental eigenvalue and for the fundamental eigenvalue. Thus, for the nonsymmetric 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 heavywater 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.1105805e2  8.16457e3  2.26216e3 
2  0.89752185  2.2306487e2  2.30623e2  
2  1  1.38255219  1.1174585e2  8.22378e3  2.22750e3 
2  0.89749043  2.2387609e2  2.26849e2  
3  1  1.37441741  1.0620368e2  8.08816e3  2.14281e3 
2  0.88836771  1.6946527e2  2.04887e2  
4  1  1.31197955  1.2687953e2  1.23115e2  0.0 
2  0.87991376  5.2900925e2  0.0  
6  1  1.38138909  1.056312e2  7.76568e3  2.39469e3 
2  0.90367052  2.190298e2  2.66211e2  
7  1  1.30599110  1.1731321e2  1.10975e2  0.0 
2  0.83725587  4.3330365e3  0.0  
8  1  1.29192957  1.1915316e2  1.15582e2  0.0 
2  0.81934103  3.0056488e4  0.0  
9  1  1.06509884  2.8346221e2  2.61980e2  0.0 
2  0.32282849  3.3348874e2  0.0 
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 
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.1645e05  0.983793 + 1.2072e05 
3  0.983594 1.1645e05  0.983793 1.2072e05 
4  0.964240 + 2.1564e05  0.964523 + 2.2337e05 
5  0.964240 2.1564e05  0.964523 2.2337e05 
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.5570e05  0.896317 + 3.6750e05 
10  0.895682 3.5570e05  0.896317 + 3.6750e05 
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 
Diffusion  SP  

1  42.263 + 0.0  41.380 + 0.0 
2  84.867 0.06130  83.821 0.06358 
3  84.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 
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 
Diffusion  SP  

1  0.04437 + 0.0  0.04395 + 0.0 
2  0.05755 1.15549e05  0.05735 1.22333e05 
3  0.05755 + 1.15549e05  0.05735 + 1.22333e05 
4  0.06807 6.35264e06  0.06798 6.66947e06 
5  0.06807 + 6.35264e06  0.06798 + 6.66947e06 
6  0.07219 + 0.0  0.07213 + 0.0 
7  0.07415 + 0.0  0.07412 + 0.0 
8  0.07453 + 0.0  0.07452 + 0.0 
9  0.07577 1.52484e06  0.07574 1.60360e06 
10  0.07577 + 1.52484e06  0.07574 + 1.60360e06 