Dipartimento di Energetica, Politecnico di Torino, Corso Duca degli Abruzzi, Torino 24-10129, Italy
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]
families of
delayed neutron precursors are considered and the delayed emissivity is related
to the precursor density by the
expression
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 aswhere 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]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 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 , 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 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 asand the total reactivity 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 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“" behavior
based on perturbation theory has been proved. Analogously, the effect on the
value of 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 aswhere 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
aswhere is the state of
the system after the first time step, is the initial
condition of the system, and the matrix is defined as
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:
Due to equality (8), these eigenvectors
turn out to be also the eigenvectors of the point kinetic matrix ; the established relation for the eigenvalues
iswhere and are the
eigenvalues of the point kinetic matrix . The solution after the first time step can thus be
expressed as
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 and is largely
negative; the corresponding eigenvalue of matrix can be reduced
below unity by satisfying the assumption
In the cases considered for this analysis is of the order
of and from condition (12) this circumstance
implies . Satisfying
condition (12) allows to introduce the following
simplification:
Other two approximations are
introduced: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
In (15), is the
fundamental eigenvalue, while is the
fundamental state of the system, defined as
Applying recursively (15), it is possible to
describe the asymptotic evolution of the system for each transit of the fissile
lump: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 asIn 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 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.
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 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 ( 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 is negligible.
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 .
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 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 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.
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 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 . 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].
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 during ; (d): value of 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 . 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 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].
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 during ; bottom right: value of 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 MW, the core
inlet temperature is 600°C and the
average fuel temperature 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: pcm/°C (solid lines), 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 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.
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: pcm/°C (), pcm/°C (), 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.