Modelling dynamic processes in a nuclear reactor by state change modal method

Modelling of dynamic processes in nuclear reactors is carried out, mainly, using the multigroup neutron diffusion approximation. The basic model includes a multidimensional set of coupled parabolic equations and ordinary differential equations. Dynamic processes are modelled by a successive change of the reactor states. It is considered that the transition from one state to another occurs promptly. In the modal method the approximate solution is represented as eigenfunction expansion. The numerical-analytical method is based on the use of dominant time-eigenvalues of a group diffusion model taking into account delayed neutrons.


Introduction
To describe neutronic processes occurring in a nuclear reactor [6], a complex integrodifferential equation is involved, in which the neutron flux depends on time, energy, spatial and angular variables [11]. As a rule, the group diffusion approximation are used [7] in practical calculations. Currently diffusion models are derived and applied using sophisticated homogenization methodologies (see [9]).
In the theory and practice of neutronics calculations, fast methods to obtain approximate solutions are given great attention. A class of methods for modeling the nonstationary group neutron diffusion is associated with the multiplicative representation of the space-time factorization methods and the quasistatic method [5]. An approximate solution is represented as a product of the time-dependent amplitude the spatial shape function, which is often associated with the fundamental eigenfunction of neutron diffusion equations.
When using the quasistatic method, the problem is significantly simplified, thus, it is doubtful to obtain good accuracy for an approximate solution, in particular, in dynamic regimes with complex flux redistributions. The general approach of the modal method has been developed [10]: the solution is represented as a sum of dominant eigenvalues with time-dependent coefficients. Nonstationary processes can be described using the approximate solution expansion in timeeigenvalue of α-eigenvalue problem [12]. In a simpler model without delayed neutrons, the modal methods were used in [8]. The principal point is connected with the fact that when using this approach, we deal with an unbound system of the coefficient equations. It should also be noted that the eigenvalues are complex for both λand α-eigenvalue problem. To set the initial state, this leads to the need to solve the appropriate adjoint spectral problems.

Problem statement
The neutron flux is modelled in multigroup diffusion approximation. The neutron dynamics is considered in the bounded convex two-dimensional or three-dimensional area Ω (x = {x 1 , ..., x d } ∈ Ω, d = 2, 3) with boundary ∂Ω. The neutron diffusion is described by: .., c M } and matrices: where δ gg is the Kronecker symbol. We shall use the set of vectors φ, whose components satisfy the albedo boundary conditions: where n -outer normal to the boundary ∂Ω.
Here φ g (x, t) -neutron flux of g group at point x and time t, G -number of energy groups, M -number of types of delayed neutrons and other terms are standard.
The Cauchy problem is formulated for equations (1) when

State change modal method
The nuclear reactor is always non-stationary. We will use the following simplified description of the dynamic processes in a nuclear reactor. In selected time interval, the non-stationary neutron flux is determined by the nuclear reactor state. The reactor state is characterized by the constant coefficients of multigroup diffusion equations (1). Dynamic processes in a nuclear reactor can be considered as a change of states. At a certain time t = t s , s = 1, 2, ... an instantaneous change of state occurs. The state s is defined by the parameters in equations (1): Simulation of the reactor dynamic behavior is considered as a consistent solution of subtasks for the reactor states. The initial condition for the state s (at t = t s−1 ) is the final state of the reactor for the state s − 1.
An approximate description of the non-stationary process at a separate stage is based on modal approximation. An approximate solution is considered in the form of decomposition in eigenfunctions of time and α-eigenvalue problem. Finite-element approximation in space is used.
Let's u = {φ, c}. Rewrite the system of equations (1) as where I is the identity matrix. Initial conditions: After approximating over the space by the finite element method from (4), (5) we turn to the Cauchy problem for a linear system of ordinary differential equations with constant coefficients: where h is the discretization parameter. The main feature of the problems we are considering is that the matrices A h and B h are real and asymmetric. The modal approximation corresponds to the representation of the approximate solution (u h ≈ u N ) of problem (6), (7) in the following form where N is the number of dominant eigenvalues, w n (x) are eigenfunctions. Let us define eigenfunctions and eigenvalues as the solution of the α-eigenvalue problem: In general, eigenfunctions and eigenvalues (9) are complex.Then in the representation (8) we obtain A special attention should be paid to define the coefficients a n (t s−1 ) = b n , n = 1, 2, ..., N . In the case of real eigenvalues, we have This representation is not very suitable for practical use with modal approximation, when we work only with dominant eigenfunctions. The initial condition includes two components u s ). Dynamic behaviour of these components is due to different time-scale processes: slow processes for delayed neutrons and fast processes for neutron flux when the reactor state changes. We can model the slow phase of the reactor dynamics with modal approximation: only the function c s h (x) is approximated. The standard approach for the decomposition of the function u s h (x) over the system of non-orthogonal functions v n (x), n = 1, 2, ..., N h consists in using the biorthogonal system of functions [3]. Consider the spectral problem adjoint to (9) The eigenfunctions of problems (9) and (10) are orthogonal in the sense of the equality where (·, ·) means corresponding scalar product. In view of this, one can obtain In the approximate solution of problem (6), (7) only the first N coefficients b n in (11) are used (see (8)): where v n (x) = (φ n (x), c n (x)). In this case, the spectral problems (9), (10) are solved for N dominant eigenvalues. The solution of the adjoint spectral problem is involved only for calculating the initial condition coefficients. The state change modal method is based on the following calculating scheme. Off-line calculation. Calculation of the coefficients of the mathematical model of the multigroup diffusion approximation for the separate reactor states, performed in advance. It can be included calculated dominant eigenvalues and eigenfunctions of the α-eigenvalue problem (9). These data can be supplemented by dominant eigenvalues and eigenvalues of the conjugate eigenvalue problem (10). On-line calculation. Real-time modeling is performed using the modal solution of the problem (6), (7). The coefficients in the representation (12) are calculated from the initial condition using (11). The solution for other time intervals is determined according to (8).

The test: the dynamics of the VVER-1000 reactor during the transition from the supercritical mode to the subcritical mode
A test problem for a VVER-1000 reactor [4] is considered in the two-dimensional approximation using a two-group diffusion approximation with delayed neutrons.

Supercritical state: α-eigenvalue problem
Below the results of a numerical solution of the α-eigenvalue problem (9) are presented. The dominant eigenvalues α n = λ (α) n , n = 1, 2, ..., N are searched for at N h . Similar calculations of the eigenvalues for the VVER-1000 test problem without delayed neutrons can be found in [1].
In our example, the main eigenvalue is negative and therefore the major harmonic will increase, and all others will fade. This demonstrates the regular mode of the reactor operation. The value α = λ (α) 1 itself determines the neutron flux amplitude and is directly related to the reactor period in the regular regime.

Adjoint spectral problem
Analogous data were obtained for approximate solution of the adjoint spectral problem (10). The eigenvalues of the spectral problems (9) and (10) coincide. Their difference from each other is an indirect measure of the accuracy of the numerical solution. Data on the dominant eigenvalues, which are given in Table 2 (k = 96, p = 3), show that the eigenvalues of the main and adjoint spectral problems are close to each other with good accuracy.
The spectral problems under consideration are characterized by small imaginary parts of the eigenvalues. Therefore, we can expect that the eigenfunctions of problem (9) are close to orthogonal. For convenience of comparison, the eigenfunctions are normalized in L 2 (Ω). The maximum non-orthogonality (for (φ 1 )) does not exceed 10 %. The biorthogonality condition of the eigenfunctions of the fundamental functions (see (9)) and the adjoint spectral problems (see (10)) is valid with approximately the same accuracy. This observed error can be related to an approximate calculation of eigenvalues and eigenfunctions.
Within the modal method, we can not rely on high accuracy when considering a relatively small number of dominant eigenvalues. Therefore, in the example under consideration, we can assume that the eigenvalues are real, and the corresponding eigenfunctions are orthogonal. Instead of (11) coefficients are used to approximate the initial condition.

Subcritical state
In the supercritical mode, due to the sufficiently large magnitude of the main eigenvalue, the regular regime of the reactor is rapidly developing, where Here v 0 1 (x) is the first mode of the supercritical state. We consider the problem with the transition from this supercritical state at t 0 = 0 to the subcritical state.
The subcritical stage is characterized by a 15% increase in the coefficient Σ 2 for material 4 in the VVER-1000 test diffusion constants. The initial state is characterized by specifying the initial conditions at t 0 = 0 as The calculation results of the dominant eigenvalues for the subcritical state are presented in Tables 3. In this case even the first eigenvalues not significantly differ from each other. For an approximate solution, we use the following formulation: where the coefficients b n , n = 1, 2, ..., N are calculated according to the given initial condition (13). These coefficients for N = 50 are shown in Fig.1. As we see, an approximate solution can be described by first mode only. We distinguish two phases of the dynamic process: fast and slow. In the fast phase, the initial condition (14) is rearranged to the initial condition, which corresponds to (15): from the function u(x, 0) to the function u N (x, 0). The slow phase is associated with the solution evolution according to (15). Within the state change modal technology, the fast phase is not modeled at all.
The beginning and the end of the fast phase are illustrated through the calculational data shown in Fig. 2. The results were obtained with N = 10. We can note very small changes in the topology of the initial and reconstructed initial conditions. Let us pay attention to the substantial restructuring of the solution, which is illustrated by large changes in the neutron flux amplitudes for the first and second groups.

Comparison with the nonstationary problem solution
An approximate solution, obtained using modal approximation, can be compared with the dynamic problem solution (1). Fully implicit scheme on a uniform grid in time with a sufficiently small step τ = 0.0025 is used (see details in [2]). We consider the dynamics of the neutron power of the nuclear reactor P and the delayed neutrons source C at the initial stage during the transition from the critical state to the subcritical one. Here There is a rapid change in neutron power over a short period of time, while the delayed neutrons source changes rather slow. The dynamics of the slow phase is illustrated in Figs. 3. Here the solution of full equations (dynamic in Figs. 3), and the modal approximation solution (modal) are presented. We see that the integral characteristics of the reactor dynamics at slow stage are calculated with good accuracy.

Conclusions
The problem of reactor dynamic simulation is considered on the basis of multigroup neutron diffusion equations taken into account delayed neutrons. The modal approximation is used: Numerical simulation of reactor non-stationary processes is carried out on the basis of a successive change in the states of the reactor. These states are characterized by a set of constant parameters to describe the multigroup neutron flux behavior. The state change modal method was developed. The phase, which described fast transition to the approximate solution, is selected as a set of dominant modes. At a slow phase of the reactor dynamics, the solution is based on the evolution of dominant modes.
The VVER-1000 two-dimensional test problem is used to perform the transition from the supercritical state to the subcritical state. Calculations of dominant modes for a reactor supercritical state are performed. The major mode solution, which determines the reactor regular regime, is used as the initial condition for transition to the subcritical state. Comparison of the calculation results obtained by using the modal approximation and direct dynamics calculation, shows the acceptable accuracy in calculation of neutron power and delayed neutrons source.