Dynamics of Fluid Fuel Reactors in the Presence of Periodic Perturbations

Abstract

The appearance of perturbations characterized by a periodic time behavior in fluid fuel reactors is connected to the possible precipitation of fissile compounds which are moved within the primary circuit by the fuel motion. In this paper the time-dependent response of a critical fluid fuel system to periodic perturbations is analyzed, solving the full neutronic model and comparing the results with approximate methods, such as point kinetics. A fundamental eigenvalue of the problem is defined, characterizing the trend of divergence of the power. Parametric studies on the reactivity insertion, the fuel velocity and the recirculation time are performed, evidencing the sensitivity of the eigenvalue on typical design parameters. Non-linear calculations in the presence of a negative feedback term are then performed, in order to assess the possibility to control a fluid fuel system when periodic reactivity perturbations are involved.

1. Introduction

Circulating fluid fuel reactors have been studied with alternate fortune during the last forty years. The interest in the development of this reactor design lies in the promising perspectives concerning the management of the fuel cycle and in the possibility to perform the transmutation of long-life radioactive wastes. For these reasons, the molten salt reactor concept meets Generation IV requirements in terms of sustainability [1].

In the years 1965–1969, an extensive experimental activity on fluid fuel reactors was performed at Oak Ridge Laboratories [2, 3]. The Molten Salt Reactor Experiment (MSRE) proved the feasibility of the design and allowed to study the peculiar physical phenomena characterizing this kind of systems. Experience in the field of fuel performances, material irradiation and corrosion, operation and maintenance was also acquired.

The experimental results and conclusions drawn from the MSRE program have served as a starting point for the research activities developed in recent years in different countries in Europe. The MOST Project [4], financed by the European Community within the V Framework Program, has taken profit of the MSRE results performing an assessment of the state of the art of the technology and carrying out a benchmark activity on the available computational tools for fluid fuel systems, tested against experimental MSRE data [5]. A Coordinated Research Project supported by IAEA, “Studies of Advanced Reactor Technology Options for Effective Incineration of Radioactive Waste" [6], has developed a research task devoted to the safety assessment of a novel design of molten salt reactor, the MOSART concept [7], performing dynamic calculations [8].

The physics of fluid fuel systems is characterized by some peculiar phenomena, connected to the movement of the multiplying material through the core and in the primary circuit. The motion of the fuel affects the delayed neutron precursors, which are dragged within the system by the fluid flow; as a consequence, the emission of delayed neutrons takes place in different spatial positions with respect to the corresponding prompt emissions. This effect results in a spatial redistribution of delayed neutron precursors, depending on their lifetime, and a general reduction of their role in the fission process, since a part of them is swept outside the reactor core and decays in the primary circuit.

The reduction of importance of delayed emissions associated to the motion of the fluid fuel affects both the static and dynamic behavior of the reactor. The critical condition depends on the velocity regime and the effective fraction of delayed neutrons is reduced, with obvious effects on the time-dependent response of the system.

The physics of such systems also allows the appearance of peculiar perturbations of the reactor material composition which constitutes an interesting aspect of neutron dynamics to be studied. In this paper, the problem of the possible appearance of perturbations with a periodic behavior in time is approached: this aspect is strongly connected to some physical phenomena occurring in fluid fuel reactors.

The operation of a molten salt reactor and the consequent modification of the fuel composition due to burnup can cause the precipitation of fissile compounds. Since the isotopic vector of the fissile components is modified by fissions, the addition of fissile material to keep reactivity can lead to the solubility limit for some fissile nuclides, which as a result precipitate. The consequent lump of fissile material is driven through the system by the fuel motion and represents a localized perturbation of the multiplicativity of the medium with a periodic behavior in time. Analogously, the formation of helium gas bubbles within the reactor during operation constitutes an additional, localized perturbation of the material characteristics of the core which is swept throughout the whole system by the fuel velocity. The presence of gas bubbles was experimentally detected during the operation of MSRE [2] and represents an issue to be solved for the assessment of the reliability of the technology. The amount of reactivity associated to these phenomena can be rather relevant, amounting to hundreds of pcm [9].

The analysis of the time-dependent response of a molten salt reactor to this kind of perturbations can be treated by taking profit of previous studies on the peculiar phenomena associated to periodic reactivity insertions. In previous works [10–12] the problem of a periodic reactivity insertion has been approached by means of an analytical solution of the point kinetic model for both critical and subcritical systems. More recently, a preliminary analysis on the precipitation of fissile lumps and the formation of gas bubble in fluid fuel systems has been carried out [13].

In this paper, a comprehensive study of the physical phenomena associated to periodic perturbations in fluid fuel systems is performed. The time-dependent behavior of the system in response to a fissile lump moving within the primary circuit is simulated with point kinetics and solving the full space-time neutronic model, in order to evidence the role of spatial phenomena connected to the movement of a localized perturbation in the reactor. Results for the pure neutronic model are presented, together with calculations in the presence of thermal feedback, in order to determine the stability condition of the reactor.

2. Physical-Mathematical Model for the Neutronics of Fluid Fuel Reactors

The study of the neutronics of fluid fuel reactors requires the proper modification of the neutron and precursor balance equations, in order to take into account the presence of a velocity field within the reactor core. Prompt fission emissions are not affected, due to their different time scale, while a convective term is introduced in the equation for precursors. The general model for neutrons and precursors can be written in terms of neutron density 𝑛 and delayed emissivities ℰ_𝑖 as [14] 𝜕𝑛(𝐫,𝐸,𝛀,𝑡)𝑀𝜕𝑡=[𝐿(𝑡)+_𝑝(𝑡)]𝑛(𝐫,𝐸,𝛀,𝑡)+^𝑅_𝑖=1ℰ_𝑖1(𝐫,𝐸,𝑡)+𝑆(𝐫,𝐸,𝛀,𝑡),𝜆_𝑖𝜕ℰ_𝑖(𝐫,𝐸,𝑡)+1𝜕𝑡𝜆_𝑖∇⋅(𝐮(𝐫,𝑡)ℰ_𝑖=𝑀(𝐫,𝐸,𝑡))_𝑖(𝑡)𝑛(𝐫,𝐸,𝛀,𝑡)−ℰ_𝑖(𝐫,𝐸,𝑡),𝑖=1,2,…,𝑅. (1)

𝑅 families of delayed neutron precursors are considered and the delayed emissivity is related to the precursor density 𝐶_𝑖(𝐫,𝑡) by the expressionℰ_𝑖𝜒(𝐫,𝐸,𝑡)=^𝑖(𝐸)𝜆4𝜋^𝑖𝐶_𝑖(𝐫,𝑡). (2)

In model (1), the leakage 𝐿, prompt multiplication 𝑀_𝑝 and delayed multiplication 𝑀_𝑖 operators are introduced, while the velocity field, represented by the vector 𝐮(𝐫,𝑡), is taken as an input function. This system of equations requires initial and boundary conditions for both the neutron density and the delayed emissivities. Boundary conditions for delayed emissivities need to properly represent the correlation between the exiting flow of precursors at the outlet of the reactor core and the undecayed fraction of precursors re-entering the system after a certain time interval. The general mathematical formulation of this condition can be written asℰ_𝑖(𝐫,𝐸,𝑡)𝐮(𝐫,𝑡)⋅(−𝐧)=∫_𝒜outℰ_𝑖(𝐫^′,𝐸,𝑡−𝜏(𝐫^′→𝐫,𝑡))×𝑒^{−𝜆𝑖𝜏(𝐫′→𝐫,𝑡)}𝐮(𝐫^′,𝑡−𝜏(𝐫^′→𝐫,𝑡))⋅𝐧^′×𝔉(𝐫^′→𝐫,𝑡)𝑑𝒜^′,𝐫∈𝒜_in, (3)where the recirculation time in the external circuit 𝜏 is introduced. Condition (3) equates the incoming flow of precursor into the reactor core at time 𝑡, left-hand side of the equation, to the undecayed fraction of the outgoing flow of precursor evaluated at a preceding instant in time 𝑡−𝜏, related to the transit time in the recirculation loop. A geometrical redistribution function 𝔉 is also defined, expressing the probability for a precursor exiting at point 𝐫^′ to re-enter the core at point 𝐫 [15].

The solution of the full model (1) in the case of the precipitation of a fissile compound can be modeled as a localized perturbation of the multiplicativity of the medium moving at velocity 𝐮. This solution approach allows to correctly describe all the spatial and spectral effects associated to the presence of the moving perturbation, but it obviously requires a large computational effort.

An alternative procedure to simulate this kind of transients implies the adoption of a point-like kinetic model. Since the physical-mathematical model has been modified for molten salt reactors, also the corresponding approximate models, such as point kinetics, need to be extended and adapted. A consistent formulation of point kinetics can be obtained by making use of the standard factorization-projection technique, and a generalized formulation of the time-dependent amplitude equations is obtained [14]𝑑𝑃=𝑑𝑡𝜌−𝛽^{0𝑥002𝑑𝑐}ℓ_𝑝𝑃+^𝑅_𝑖=1𝜇_𝑖𝐺_𝑖+𝑆^{0𝑥002𝑑𝑐},𝜁_𝑖𝑑𝐺_𝑖=𝑑𝑡^{0𝑥002𝑑𝑐}𝛽_𝑖ℓ_𝑝𝑃−(𝜇_𝑖+𝜂_𝑖)𝐺_𝑖,𝑖=1,2,…,𝑅, (4)where 𝑃 and 𝐺_𝑖 are the amplitude functions for the neutron distribution and for the 𝑖th precursor family, respectively. The kinetic parameters in system (4) have a different definition from the standard formulation and depend on the shape of neutrons and precursors, since both unknowns are factorized and projected in the procedure. Additional terms, such as 𝜁_𝑖 and 𝜂_𝑖, appear as a consequence of the fuel motion, allowing this reformulated point model to consistently describe the dynamics of a fluid fuel reactor, although all spatial and spectral effects are neglected. In the point kinetic framework the simulation of the transit of a localized perturbation affects the value of the reactivity 𝜌 and the value of the effective delayed neutron fraction ^0𝑥0𝛽^02𝑑𝑐 when perturbations of multiplicativity are concerned.

The study of periodic perturbations in molten salt reactors is here performed solving the full model (1) in cylindrical (𝑟−𝑧) geometry and comparing it to point kinetic results. The perturbation is moved within the system along the axial direction from the core inlet to the core outlet in the time interval 𝜏_𝑐, while its residence time in the external circuit is 𝜏_𝑙. The time-dependence of the reactivity 𝜌 during 𝜏_𝑐 can be qualitatively assimilated to the function sin²(𝜋𝑡/𝜏_𝑐), according to elementary perturbation theory in diffusion approximation [16], while in the interval 𝜏_𝑙 the kinetic parameters of the reactor are left unperturbed. Also the effective delayed neutron fraction ^0𝑥0𝛽^02𝑑𝑐 experiences a perturbation with a periodic behavior in time.

The analysis is focused on the stability characteristics of the power oscillations appearing as a consequence of the periodic reactivity introduction and on the role of spatial effects, which are neglected when point-like calculations are performed. The model (1) is solved by numerical inversion of the full operator, while the point kinetic solution is produced adopting a semianalytical approach [12], trying to evidence the asymptotic trend in the power evolution. In a previous work [11] it has been proved that a periodic perturbation with zero mean value results in a divergence of the power in a critical system, due to the accumulation of precursors. In the case of the precipitation of a fissile lump a reactivity oscillation with a positive mean value is involved, and consequently the power is expected to be diverging.

At last, results in the presence of thermal feedback are presented. The point kinetic model is coupled to a zero-dimensional thermal model for the nuclear fuel, performing a thermal balance over the full core. The fuel temperature is updated on the basis of the power production, the feedback on neutronics is performed by means of an integral parameter 𝛼, such as𝜌_tot=𝜌_pert−𝛼(𝑇_fuel−𝑇_fuel,0), (5)and the total reactivity 𝜌_tot is inserted into (4). The implementation of a thermal module allows to perform a preliminary analysis of the effect of thermal feedback on the stability of the reactor in the presence of these perturbations, evidencing the influence of various physical parameters.

3. Semi-Analytical Solution of the Point Kinetic Model

The interest in the study of periodic perturbations relies in the possibility to evidence a peculiar behavior connected to the different time scales characterizing the evolution of the neutron and precursor populations. Moreover, in some specific and physically meaningful cases, like the problem of a square wave reactivity insertion, the point kinetic problem can be treated fully analytically evidencing the asymptotic trend of evolution of the point-like reactor, in both critical and subcritical configurations [11, 12].

In the case of a circulating lump or a gas bubble the associated reactivity has a time behavior which cannot allow a totally analytical approach. The localized modifications of the cross sections in the core during the time interval 𝜏_𝑐 can be approximated by a sin²(𝜋𝑡/𝜏_𝑐) function, while during 𝜏_𝑙 the inserted reactivity vanishes. In the calculation performed the movement of the perturbation has been simulated on the discrete spatial mesh in the cylindrical (𝑟−𝑧) reactor considered, and the supposition of a“sin²" behavior based on perturbation theory has been proved. Analogously, the effect on the value of ^0𝑥0𝛽^02𝑑𝑐 has been evaluated.

In the adopted semianalytical approach the full period of the oscillation 𝑇=𝜏_𝑐+𝜏_𝑙 is discretized in 𝑁 time steps Δ𝑡, retaining the kinetic parameters constant over each interval. The differential problem over the 𝑖th time interval in a critical source-free system can be written as𝑑𝐾𝑑𝑡|𝐗⟩=_𝑖|𝐗⟩,(𝑖−1)Δ𝑡≤𝑡≤𝑖Δ𝑡, (6)where 𝐾_𝑖 is the point kinetic matrix with kinetic parameters associated to the 𝑖th time step. In the following of this section only one family of delayed neutron precursors is considered to simplify notation, but the procedure is applicable to the general case of 𝑅 families. System (6) is discretized on the same time step Δ𝑡 using an implicit Euler algorithm for the time derivative. During the first time step, the system of algebraic equations which is obtained can be written as|𝐗(Δ𝑡)⟩=𝜗₁|𝐗(0)⟩, (7)where |𝐗(Δ𝑡)⟩ is the state of the system after the first time step, |𝐗(0)⟩ is the initial condition of the system, and the matrix 𝜗₁ is defined as 𝜗₁=𝐾₁Δ𝑡+𝐼̂_2×2. (8)

The solution of problem (7) can be obtained by analytical tools preliminarily evaluating the real eigenvalues (𝜖₁⁽¹⁾ and 𝜖₁⁽²⁾) and eigenvectors (|𝚪₁⁽¹⁾⟩ and |𝚪₁⁽²⁾⟩) of the matrix 𝜗₁ and of its adjoint 𝜗^†₁, ⟨𝚪₁⁽¹⁾| and ⟨𝚪₁⁽²⁾|, solution of the following homogeneous algebraic problems:𝜗₁|𝚪₁^(𝑖)⟩=𝜖₁^(𝑖)|𝚪₁^(𝑖)⟩,⟨𝚪₁^(𝑖)|𝜗^†₁=⟨𝚪₁^(𝑖)|𝜖₁^(𝑖). (9)

Due to equality (8), these eigenvectors turn out to be also the eigenvectors of the point kinetic matrix 𝐾₁; the established relation for the eigenvalues is(𝐾₁Δ𝑡+𝐼̂_2×2)|𝚪₁^(𝑖)⟩=𝜖₁^(𝑖)|𝚪₁^(𝑖)⟩,(𝜔₁^(𝑖)Δ𝑡+1)|𝚪₁^(𝑖)⟩=𝜖₁^(𝑖)|𝚪₁^(𝑖)𝜖⟩,_1(𝑖)=𝜔₁^(𝑖)Δ𝑡+1, (10)where 𝜔₁⁽¹⁾ and 𝜔₁⁽²⁾ are the eigenvalues of the point kinetic matrix 𝐾₁. The solution after the first time step can thus be expressed as |𝐗(Δ𝑡)⟩=⟨𝚪₁⁽¹⁾|𝐗(0)⟩𝜖₁⁽¹⁾|𝚪₁⁽¹⁾⟩+⟨𝚪₁⁽²⁾|𝐗(0)⟩𝜖₁⁽²⁾|𝚪₁⁽²⁾⟩. (11)

The solution of the problem over the whole period 𝑇=𝑁Δ𝑡 is obtained by applying formula (11) repetitively over all time steps; the heavy formulation of the final solution can be simplified by introducing some assumptions related to the well-known characteristics of the point kinetic eigenvalues 𝜔_𝑗^(𝑖).

The second eigenvalue of the point kinetic matrix 𝜔_𝑗⁽²⁾ can be approximated as (𝜌−𝛽^{0𝑥002𝑑𝑐})/ℓ_𝑝 and is largely negative; the corresponding eigenvalue of matrix 𝜗_𝑖 can be reduced below unity by satisfying the assumption2Δ𝑡<−𝜔_𝑗⁽²⁾,𝑗=1,2,…,𝑁. (12)

In the cases considered for this analysis 𝜔_𝑗⁽²⁾ is of the order of −10¹s⁻¹ and from condition (12) this circumstance implies Δ𝑡<10⁻¹s. Satisfying condition (12) allows to introduce the following simplification:|𝜖_𝑗⁽²⁾|<1,∀𝑗⇒^𝑁_𝑗=1𝜖_𝑗⁽²⁾=0. (13)

Other two approximations are introduced:⟨𝚪_𝑗⁽¹⁾|𝚪⁽²⁾_𝑗±1⟩=0,⟨𝚪_𝑗⁽²⁾|𝚪⁽¹⁾_𝑗±1⟩=0, (14)requiring a “pseudo-orthogonality" of the eigenvectors associated to subsequent time steps which is exactly verified for Δ𝑡 approaching zero. Consequently, the choice of a sufficiently small Δ𝑡 justifies all the assumptions (13) and (14) and allows to write the solution over the period 𝑇 as|𝐗(𝑁Δ𝑡)⟩=⟨𝚪₁⁽¹⁾|𝐗(0)⟩∏^𝑁_𝑗=1𝜖_𝑗⁽¹⁾|𝚪_𝑁⁽¹⁾⟩≡⟨𝚿|𝐗(0)⟩𝜉|𝚿⟩. (15)

In (15), 𝜉 is the fundamental eigenvalue, while |𝚿⟩ is the fundamental state of the system, defined as𝜉≡^𝑁_𝑗=1𝜖_𝑗⁽¹⁾,⟨𝚿|≡⟨𝚪₁⁽¹⁾|,|𝚿⟩≡|𝚪_𝑁⁽¹⁾⟩. (16)

Applying recursively (15), it is possible to describe the asymptotic evolution of the system for each transit of the fissile lump:|𝐗(𝑚𝑇)⟩=𝜉^𝑛⟨𝚿|𝐗(0)⟩|𝚿⟩, (17)where the exponent 𝑚 indicates the number of cycles.

Observing (17), it appears that the asymptotic evolution of a multiplying system excited through a periodic perturbation is governed by a fundamental eigenvalue that operates as an amplification factor for each cycle. The real eigenstate is obviously not the one calculated by means of point kinetics, but it can be evaluated in general by performing the inversion of the full space-angle-energy-time model.

When performing the inversion of the full model (1) an analytical approach cannot be followed, however the fundamental eigenvalue for the asymptotic evolution of the system can be defined as𝜉≡lim_𝑚→∞||𝑛(𝐫,𝐸,𝛀,𝑚Δ𝑇)||||𝑛(𝐫,𝐸,𝛀,(𝑚−1)Δ𝑇)||. (18)In real calculations, the evaluation of formula (18) is performed after few cycles owing to the large computational effort associated to the inversion of the full model, trying to evidence the appearance of a stable value of 𝜉.

The assumption of the existence of a “fundamental state” for the full model (1), appearing when the asymptotic regime is established, can be understood on the basis of well-known physical properties of linear systems in the presence of a periodic input. After an initial transient a dynamic oscillating equilibrium condition is reached for both the populations of neutrons and precursors over a full period; as a consequence, the asymptotic evolution can be described by an integral parameter as in formula (18).

4. Numerical Results

Calculations are performed using the numerical code DYNAMOSS [15], solving the multigroup diffusion equations for a fluid fuel system in cylindrical (𝑟−𝑧) geometry and imposing an axial velocity field without dependence on the radial variable. Two different reactor configurations are considered: the first one assumes geometrical dimensions and material data typical of MSRE with two energy groups and one family of precursors [3]; the second one reproduces the main characteristics of the MOSART design [7] with four energy groups and six precursor families. In both systems, the reactor core has homogenized cross sections and the MOSART-like geometry includes a radial reflector. In Table 1 the main data for the two systems are summarized.

Table 1 

Geometrical data for the systems analyzed.

The first set of results concerns the MSRE system in the presence of a fissile lump moving through the system and a parametric analysis is performed on physically meaningful quantities influencing the power evolution. In Figure 1, the flux perturbation associated to the movement of the fissile solid component is shown at different instants during the circulation period in the core 𝜏_𝑐. The maximum reactivity associated to this perturbation is about 250 pcm . The appearance of spatial and spectral effects is clearly visible: the increased multiplicativity in this localized region causes a larger absorption of thermal neutrons with a consequent additional production of fast neutrons by fission.

(a) 𝑡=0 seconds (unperturbed system)

(b) 𝑡=2 seconds

(c) 𝑡=4 seconds

(d) 𝑡=8 seconds

Figure 1 

Behavior of the neutronic flux during the transit of the fissile lump. Two energy groups are considered (higher-level flux corresponds to the first group) and the solution is plotted at different instants during the circulation period.

The power evolution following the transit of the fissile lump is then simulated using the point-like approach described in the previous section and it is compared to the results obtained by inversion of the full space-time operator, in order to compare the capability to predict the trend of the power divergence by simplified models, with respect to more accurate ones. In this case the reactivity associated to the localized perturbation amounts to 25 pcm at its maximum. In Figure 2(a) the power evolution following the transit of the solid fissile compound along the axis of the system is plotted during three periods (𝑇=20 seconds). In Figures 2(b) and 2(c), the time dependence of reactivity and effective delayed neutron fraction is reported, showing that the effect of the transit of the perturbation on the value of ^0𝑥0𝛽^02𝑑𝑐 is negligible.

(a) Power evolution

(b) Reactivity

(c) Effective delayed neutron fraction

Figure 2 

Behavior of total power (a), reactivity (b), and effective delayed neutron fraction (c) for a transient in a critical system involving the transit along the reactor axis of a fissile lump. In (a), the point kinetic solution (thin line) and the inversion of the full space-time model (bold line) are reported. Corresponding asymptotic evolutions (dash-dotted lines) are also presented.

The first half period is characterized by an insertion of reactivity into the system, with a consequent increase of the power, which is reproduced by both the dynamic models adopted. During the second part of the period the non-equilibrium condition between neutrons and precursors causes a decrease of the power which is predicted with large differences between the two options. The importance of spatial effects can also be clearly evidenced by observing the modification of the neutron flux along the axial coordinate during the transit of the lump. The spatial modifications introduced by the perturbation on the fission cross section are neglected by the point model, thus giving large discrepancies in the prediction of the power evolution.

Figure 2(a) also reports the asymptotic evolution of the power for the full and point model, according to formulae (17) and (18). The agreement of the curves for point kinetics is very good, while in the full case the asymptotic establishment of a “fundamental mode” is not reached during the first oscillation periods, thus giving larger discrepancies. In Table 2, the relation between maximum reactivity insertion and value of the fundamental eigenvalue 𝜉, evaluated through model (17), is reported. The power starts rising very quickly, since for each cycle it is amplified by a factor 𝜉. The last example regards a super-prompt-critical configuration, since the maximum reactivity insertion is larger than the value of ^0𝑥0𝛽^02𝑑𝑐.

Table 2 

Behavior of the fundamental eigenvalue 𝜉 as a function of the reactivity insertion.

A square-wave problem producing the same results in terms of fundamental eigenvalue and eigenstate is now considered. Such a problem allows a fully analytical treatment. In Figure 3 the previous graphs are compared to the response for the corresponding square-wave problem. The half-period of the oscillation with a positive insertion of reactivity for both 𝜌 and ^0𝑥0𝛽^02𝑑𝑐 corresponds to the time of transit of the lump in the core, while the amplitude is slightly larger (few percents) than the mean value of 𝜌(𝑡) and ^0𝑥0𝛽^02𝑑𝑐(𝑡) on the interval in which a positive reactivity is inserted. The behavior of the state of the system is quite different during the oscillation period, but the reactor is exactly in the same conditions at the end of each oscillation.

(a)

(b)

(c)

(d)

Figure 3 

Behavior of the total power (a), neutron delayed precursor concentration (b), reactivity (c), and effective fraction of delayed neutrons (d) during the transit a fissile of lump along the axis of a critical system. The different curves in the graphs are associated to the different models adopted: real reactivity injection (bold line), square wave approximation (thin line).

Some additional parametric studies on important quantities are also performed, solving the time-dependent problem with the point kinetic model in order to easily evidence the main integral parameters coming into play. At first, the influence of the recirculation time in the external circuit 𝜏_𝑙 is analyzed. The geometry of the reactor core is left unchanged in terms of composition and velocity field, while the value of 𝜏_𝑙 is modified. The critical condition is determined and then the transit of the fissile lump is simulated by perturbing locally the multiplicativity, keeping the relative variation with respect to the critical value constant in all calculations.

In Figure 4, the behaviors of the main parameters of the problem as functions of 𝜏_𝑙 are reported. The influence on the reactivity insertion is very small, as could be expected since the relative perturbation on cross sections is kept constant, while the increase of the recirculation time causes a corresponding increase of the fundamental eigenvalue. This results in a faster divergence of the power, due to the fact that, with an increased 𝜏_𝑙, the role of delayed neutrons is reduced and the reactivity insertion in dollars turns out to be larger. This effect is confirmed by observing Figures 4(c) and 4(d), reporting the mean value of ^0𝑥0𝛽^02𝑑𝑐 during the transit time in the core 𝜏_𝑐 (𝛽₊) and the value experienced during the recirculation time 𝜏_𝑙 (𝛽₋). It must be noted also that the effect on these two values of the effective delayed neutron fraction is comparable, meaning that the change of 𝜏_𝑙 does not affect the shape in time of the oscillation of ^0𝑥0𝛽^02𝑑𝑐. As a further comment, the appearance of a saturation value for all the integral parameters can be observed, as the increased recirculation time approaches the asymptotic condition in which no delayed neutron precursors re-enter the core [15].

(a)

(b)

(c)

(d)

Figure 4 

Parametric study of the effect of the recirculation time 𝜏𝑙 on the integral parameters regulating the divergence of the power. (a): mean reactivity inserted in the interval 𝜏𝑐; (b): value of the fundamental eigenvalue 𝜉; (c): mean value of 0𝑥0𝛽02𝑑𝑐 during 𝜏𝑐; (d): value of 0𝑥0𝛽02𝑑𝑐 during 𝜏𝑙.

An analogous set of calculations is performed changing parametrically the fuel velocity. The geometry of the system is left unchanged, thus an increase of the fuel velocity implies a reduction of both transit times 𝜏_𝑐 and 𝜏_𝑙. In Figure 5 the results obtained are presented. The value of the reactivity is again basically unperturbed, while some peculiarities can be observed concerning the values of 𝜉 and ^0𝑥0𝛽^02𝑑𝑐. As the fuel velocity increases the effective delayed neutron fraction is reduced, as it is already observed in the parametric study on 𝜏_𝑙 and as could be expected. On the other side, this reduction of ^0𝑥0𝛽^02𝑑𝑐 does not imply an increase of the fundamental eigenvalue, which on its turn decreases. This effect can be explained by the fact that, with the increased velocity, the time spent by the lump on each cycle is reduced as well as its effect in increasing the power level and accumulating precursors. The appearing of a saturation value for all the parameters as the velocity increases and approach to the “infinite velocity" condition can also be observed [15].

(a)

(b)

(c)

(d)

Figure 5 

Parametric study of the effect of the fuel velocity 𝑢 on the integral parameters regulating the divergence of the power. Top left: mean reactivity inserted in the interval 𝜏𝑐; top right: value of the fundamental eigenvalue 𝜉; bottom left: mean value of 0𝑥0𝛽02𝑑𝑐 during 𝜏𝑐; bottom right: value of 0𝑥0𝛽02𝑑𝑐 during 𝜏𝑙.

At last, results in the presence of a simplified zero-dimensional thermal model are presented. The aim of this analysis is to enlighten the possibility to stabilize the reactor behavior and prevent the power divergence. The feedback model adopted introduces an additional reactivity term in the point kinetic equations (4) connected to the modification of the fuel temperature, as in formula (5). The reactor design considered is of MOSART-type: the starting reactor power is 2400 MW, the core inlet temperature 𝑇_in is 600^°C and the average fuel temperature 𝑇_fuel,0 is 660^°C. Parametric studies on the value of the feedback coefficient 𝛼 and on the recirculation time 𝜏_𝑙 are performed.

The possibility to control the reactor with a negative feedback effect is connected to the saturation of the precursor concentration in the system. In Figure 6 the evolution of the precursor amplitude for the first family, 𝐺₁, is plotted adopting different values of the feedback coefficient 𝛼. Reference values for the velocity and the recirculation time are adopted, see Table 1. The stabilization of the precursor concentration is clearly visible, implying the stabilization of the system power which keeps oscillating with a constant mean value. The increase of the role of thermal feedback obviously reduces the amplitude of the excursion from the initial value to the asymptotic one. In Figure 7, the response of three different precursor families is plotted during the first oscillation period, showing different behaviors connected to their decay constants. The role of thermal feedback is also evidenced.

Figure 6 

Precursor concentration evolution in the presence of a moving fissile lump evaluated with point kinetics coupled with a thermal feedback model. The stabilizing role of feedback is clearly visible.

Figure 7 

Precursor concentration for the first (∘), third (□), and fifth (∇) families during the first period of oscillation. Two different values of feedback coefficient are considered: 𝛼=1.53 pcm/°C (solid lines), 𝛼=3.06 pcm/°C (dash-dotted lines).

The effect of thermal feedback on the system power and on the amplitude for all six precursor families is shown in Figure 8, where the maximum values of these quantities over each oscillation are plotted as functions of the number of periods 𝑚. The comparison of the fastest decaying family 𝐺₆ to the one with the longest mean life 𝐺₁ shows that the number of periods needed in order to reach a stable asymptotic level is higher for the families with smaller decay constants. Consequently, the number of oscillations required to establish the asymptotic equilibrium condition is regulated mainly by the first family.

Figure 8 

Saturation trend for precursors and system power when including feedback effects in the evaluation of a periodic reactivity insertion. Feedback coefficient 𝛼=1.53 pcm/°C. The maximum values of power and delayed amplitudes over each oscillation are plotted as functions of the number of periods 𝑚.

A parametric study on the recirculation time 𝜏_𝑙 is then performed, as shown in Figure 9, where the maximum and minimum value of the stable oscillation for power and precursors is plotted as function of the recirculation time. An increased value of 𝜏_𝑙 results in larger oscillations of both the precursor concentration and the power, since an increased recirculation time reduces the value of the effective delayed neutron fraction. This phenomenon affects mainly the maximum value reached during the oscillation; the importance of feedback effects in order to reduce the amplitude of the power oscillations is clearly visible.

(a)

(b)

Figure 9 

Characteristics of the stable power and precursor oscillation in the presence of feedback as a function of the recirculation time 𝜏𝑙. The solid line represents the maximum value during the oscillation, the dash-dotted line the minimum. Different values of the feedback coefficient are considered: 𝛼=1.53 pcm/°C (∘), 𝛼=3.06 pcm/°C (□), 𝛼=4.59 pcm/°C (∇).

5. Conclusions

The problem of the precipitation of solid fissile salts in circulating fuel reactors has been thoroughly analyzed, enlightening its effects on the dynamic behavior of the system. The power oscillations induced by a localized perturbation of multiplicativity moving through the reactor have been simulated by means of the solution of the full space-time model and using a simplified point-like approach, evidencing the appearance of an asymptotically diverging fundamental state. The differences in the evaluation of the time-dependent response of the system are discussed and show the importance of full space dynamic tools for a realistic prediction of this type of transients in molten salt reactors. The results of a parametric analysis on important physical parameters such as the fuel velocity and recirculation time have been discussed.

Calculations in the presence of a thermal feedback model have been carried out, enlightening the stabilizing effect of negative feedback and the establishment, after an initial transient, of a stable oscillation around a constant power level. These numerical evaluations constitute a preliminary assessment on the safety issues connected to periodic reactivity perturbations in a reactor operating at full power.

Further developments on this research subject include the implementation of more detailed thermal models in order to take into account spatial effects, as well as improvements in the neutronic description. The introduction of more complex fluid-dynamic models and the development of a fully coupled neutronic-thermal-fluid-dynamic module is envisaged, in order to be able to perform a comprehensive analysis of the reactor design.

Acknowledgment

The first author gratefully acknowledges the financial support provided by Associazione per lo Sviluppo Scientifico e Tecnologico del Piemonte (ASP) during the development of the present work.