Why don't nuclear reactors go kaboom? A reactor kinetics primer part - 2

Advanced test reactor Foto: Matt Howard, Source: Wikimedia Licens: Creative Commons Attribution-Share Alike 2.0 generic

Its time for some more fun with reactor kinetics, in the last post we ended by looking at the point kinetics equation with one group of delayed neutrons. In this post as I promised we will talk about reactivity feedbacks. To brush up your memory, reactivity is defined as:

\Large \rho = \frac{k_{\infty}-1}{k_{\infty}}

It basically tells you how much you are deviating from a perfectly critical system where each neutron gets consumed and produces exactly one new neutron. If reactivity is positive the system will have a increasing neutron population and vice versa. From an earlier post we can recall that k_{\infty} can be defined through the four factor formula as:

\Large k_{\infty} = \eta_T f p \epsilon

Where the different factors are:

  • \Large\eta_{T} is the number of neutrons emitted by fission per thermal neutron absorbed in the fuel. It can be expressed as \large \eta_{T} = \nu \frac{\Sigma_f}{\Sigma_a} Where \nu is the number of neutrons produced per fission, \Sigma_f is the fission cross section in the fuel and \Sigma_a is the total absorption cross section in the fuel (absorption includes fission as well as other reactions that doesn't cause fission).
  • \Large f which is the thermal utilization factor. It gives the probability that a neutron gets absorbed in the fuel instead of in non fuel materials in the material mixture. It is defined as \large f = \frac{\Sigma_{aF}}{\Sigma_{aF}+\Sigma_{aM}} Where \Sigma_{aF} is the absorption cross section of fuel and \Sigma_{aM} is absorption cross section of all other materials.
  • \Large p is the probability that a neutron will get from fast energies to thermal energies without being absorbed in resonances (more on those later) while slowing down.
  • \Large\epsilon is the fast fission factor and it gives the ratio between the number of neutrons produced by all fission (including those caused by fast neutrons) and neutrons produced by fission in thermal energies.

Lets write out some of the factors explicitly in the four factor formula:

\Large k_{\infty} = \eta_T f p \epsilon= \nu \frac{\Sigma_f}{\Sigma_a}\frac{\Sigma_{aF}}{\Sigma_{aF}+\Sigma_{aM}}p \epsilon

The first two terms is just the probability of a neutron getting absorbed in the fuel multiplied with the number of neutrons produced per neutron absorbed in the fuel. What we are interested in now is how this changes based on what happens to the reactor.

Everything that happens in a nuclear reactor comes down to neutron cross sections. The neutron cross section of a nuclide is basically a measure of how probable it is that a neutron flying by will interact with the nuclide. Think of the neutron as a dart and the cross section as the area of the board you are trying to hit. The thing with cross sections are that they are not just a representation of the size of the nucleus. In our ordinary world we are used to the fact that the area of an object doesn't change, a dart board certainly doesn't massively shrink or expand due to outside temperature or how much its raining, but neutron interactions are quantum mechanical phenomena so our gut feelings doesn't apply very well.

In the quantum mechanical world energy is everything. Simplified one can say that in a nuclide the protons and neutrons (collectively called nucleons) form concentric shells.  When the nuclide is in its ground state all the nucleons are in the lowest shells they can occupy (all of them can not be in the same shell due to the Pauli exclusion principle), if you where to kick on this nuclide violently one or a few nucleons might get bumped into a higher shell, this is called an excited state. Each excited state has a specific energy level and a nuclide can not have any energy level not corresponding to excited states. Now imagine that a neutron with energy E is hurdling towards a nuclide, if the nuclide catches the neutron a new nuclide (with one more neutron obviously) will be created with an energy equal to the neutrons incoming energy (plus some binding energy etc). If this energy happens to match the energy of an excited state of the new nuclide it is much more probable that the neutron will be captured. This is called an absorption resonance. The further away from matching an excited state the less the probability of absorption.

There are, as always, more complications to the picture but that gives a basic insight into why cross sections changes based on the energy of the neutrons and we don't really need to go into any fine detail on nuclear reactions to understand reactor behavior.

The cross section vs neutron energy can look something like this.

We can see all the resonance peaks in the energies from about 1 eV (electron-volt) up to around 5000 eV. Again simplified one can say that above 5000 eV there are no excited states anymore and below 1 eV it all goes to the ground state. The ground state is the most stable state (with a few exceptions in some nuclides like Americium-242) and it makes sense that the cross section would be largest for reactions going into the ground state because nature always strives for stability. Now we understand how cross sections, and in turn reactivity, can change. Just change the energy of the neutrons! Lets turn to what determines neutron energies and how it impacts reactivity.

When a neutron is born by fission it has a high energy, around 2 MeV. As we can see from the above plot cross sections are very small at those energies, about a thousand times smaller than at energy around 0.1 eV. We can remember that reaction rates(and in other words power production) is proportional to the product of the flux, the cross section and material density. If everything else stays the same then a smaller cross section means a higher demand for fissile materials. This cost something and most of the time we want to build reactors where we slow down the neutrons. Slowing down is done by a process called moderation and that is just a fancy way of saying that the neutrons bounce against atoms (called moderator) for a while until they have slowed down. When the neutrons have an energy about equal to the thermal energy (the viggling and vibrations due to heat, about 0.025 eV) of the moderator they can not slow down anymore and one say the neutrons have become thermal. Moderation is best done by light atoms because a neutron loses more energy when colliding with a light atom compared to a heavy one. Think of billiard balls, if you shoot one billiard ball on another ball it is possible that there will be a complete energy transfer and the first ball stops while the second goes on. If you throw a billiard ball on a wall however it will bounce back with pretty much the same energy. The  lightest atom and thus most efficient moderator is hydrogen, water has plenty of hydrogen (H2O remember) so water is an excellent moderator. The only reason hydrogen it is not the best moderator is because hydrogen has a tendency to absorb neutrons, that is why one instead sometimes want to use heavy water (less absorption but worse moderation) or graphite (much worse moderation but less absorption) as moderators.

Moderator density feedback

Lets consider primarily water moderated reactors (they constitute the overwhelming majority of the reactor fleet). The density of water can change quite a lot depending on its temperature. Heat it up and it expands, heat it up even more and it starts to boil and form steam and steam has much lower density than water. How does this impact reactivity? If we look at the four factor formula again:

\Large k_{\infty} = \eta_T f p \epsilon= \nu \frac{\Sigma_f}{\Sigma_a}\frac{\Sigma_{aF}}{\Sigma_{aF}+\Sigma_{aM}}p \epsilon

The first factor is only related to fuel and its not influenced by moderator density. The second factor is related to the the density of moderator \Sigma_{aM} = \rho_m\sigma_{aM}, where \rho_m is moderator density. If moderator density goes down clearly f grows, so the thermal utilization is increased by lowering of density and this contributes positively to reactivity.

The resonance espace probability, p, is influenced by how well the neutrons are moderated. Resonances are very complicated, especially because they overlap and interfere with each other, so I won't go into any mathematical treatment of resonances. If we consider a normal reactor with fuel rods submerged in water an effect is taking placed called resonance self shielding. The water has no resonances that can snag up neutrons, the resonance absorption takes place solely in the fuel. We notice from the figure above that the resonances are situated between the fast region where neutrons are born and the thermal region. If we can take the neutrons out of the fuel the moment they are born, let them thermalize in the moderator and then return them to the fuel we in effect avoid the resonances completely, this is the self shielding and it increases reactivity. If we decrease water density it becomes more likely that the neutrons will bounce back into the fuel before they have become thermal and thus reactivity decreases from a lowering of the density.

The two different effects on reactivity compete and depending on how the reactor is designed it decides if the end result is negative or positive reactivity from a density change. If a reactor has more water than what is strictly necessary for moderation then the effect of a density decrease will be positive, otherwise negative. All water moderated reactors are designed so that the moderator density feedback is negative. The Chernobyl design on the other hand had positive feedback from a lowering of density. The moderation was primarily done by graphite and water was used as coolant, when water started boiling less neutrons got absorbed in the water but the graphite was still moderating just as efficiently.

The fourth parameter, the fast fission factor \epsilon, increases slightly if moderator density goes down, but the effect is small and can be neglected in this treatment.

Moderator temperature feedback.

As mentioned above there is a limit to how slow the neutrons can get when bouncing against the moderator atoms and that is the speed at which the moderator atoms moves by themselves. This opens up for another feedback, if the moderator heats up then the thermal neutrons will get more energetic and this brings about moderator temperature feedback. If water is the moderator the density and temperature changes can not be separated because a temperature change causes a density change. In a graphite moderated reactor there is hardly any moderator density changes (moderator is big blocks of graphite and graphite doesn't expand all that much) and there moderator temperature is most important. When the moderator heats up and the thermal neutrons gets more energetic they can get pushed up into the lower resonance region and the resonance escape parameter is effected, depending on the fuel composition this can have a negative or positive influence on reactivity. If most of the fuel has low lying fission resonances then reactivity will go up and vice versa of the low lying resonances are capture resonances.

Doppler feedback

A third important feedback is what happens in the fuel itself, this is again an effect on the resonance escape probability. The fuel atoms also move due to its temperature. The neutron cross section doesn't only care about the neutrons energy, its actually determined by the collision energy. If the fuel atoms are moving towards (or away) from the neutrons then that changes the collision energy, this means that if we heat up or cool down the fuel that effects the cross section and we have a feedback from that as well. This feedback is called doppler feedback because it is analogues to the doppler effect in acoustics. The effect of the doppler feedback is that it makes the resonances broader and flatter, this allows the resonance to interact with more neutrons and the resonance escape probability is changed. If the resonance happens to be a capture resonance it will lower reactivity, if it is a fission resonance it will increase reactivity and total doppler feedback depends on what kind of resonance dominates in the particular nuclide of interest.

Doppler feedback is the most rapid feedback, heating of the fuel is practically instantaneous when power increases while moderator feedbacks are delayed by the time it takes heat to transfer from the fuel into the moderator.

How do we squeeze all this into the point kinetics equation

There are more or less complicated ways to calculate the different reactivity contributions, we won't bother with that and later just settle for tabulated values. We simply defined a reactivity coefficient as:

\Large \alpha_x = \frac{d\rho}{dx} \approx \frac{\Delta\rho}{\Delta x}

Where x can be moderator temperature, moderator density, moderator temperature or whatever other factor you can think of. We then squeese this into the point kinetics equation like this:

\Large \frac{d\phi}{dt} = \frac{(\rho_0+\alpha_x*\Delta x(t) )-\beta}{\Lambda}k_{\infty}\phi + v\lambda C

Now we are in a bit of a dilemma, now we have introduced a lot of parameters material parameters into the point kinetics equation that we need to calculate somehow. We need moderator density, moderator temperature and fuel temperature. That means we need some kind of model to calculate the fluid flow through the reactor, we need some kind of geometry model to model heat transfer from fuel to water and so on. This will be the topic for the next blog post. But we can already now look at a simplified model, lets look at some event that is so quick that heat transfer to the moderator isn't all that important. Then we can assume that all energy produced simply goes into heating the fuel and we can ignore all feedbacks except the doppler feedback. We will continue to stick to one group of delayed neutrons.

To make things a bit more elegant we can replace neutron flux with neutron density \large \phi(t) = n(t)v and remembering that \large \Lambda = \frac{1}{\Sigma_{a}v}, then we get these equations:

\Large \frac{dn}{dt} = \frac{(\rho_0+\alpha_x*\Delta x(t) )-\beta}{\Lambda}k_{\infty}n(t) + \lambda C

\Large \frac{dC}{dt} = \frac{\beta k_{\infty}n(t)}{\Lambda} - \lambda C

We don't have any nice analytical expression for the time dependence of the reactivity and we have to dump analytical solutions and look to numerical methods. Let's do the simplest method, forward euler method. Forward euler is simply a taylor series with all higher order terms tossed into the trash can. With forward euler we have that:

\Large f_{u+1} = f_u + h\frac{df_u}{dt}

That means that the value of the function at time step u+1 is equal to the value of the function at the preceding time step plus the value of the first derivative of the function at time n multiplied by the time step size h. We can then write out equations as this:

\Large n_{u+1} = n_u + h[\frac{(\rho_0+\alpha_T*\Delta T_u )-\beta}{\Lambda}k_{\infty}n_u + \lambda C_u]

\Large C_{u+1} = C_u + h[\frac{\beta k_{\infty}n_u}{\Lambda} - \lambda C_u]

Where \alpha_T and T reefers to the fuel doppler coefficient and the fuel temperature respectively. Now we need some initial conditions and some way to calculate fuel temperature. We can say that initial neutron density,n_0 is 1 (100% normal operating power or whatnot). At time 0 the time derivative of C should be zero (i.e constant neutron precursor density) so that gives us this:

\Large \frac{dC}{dt}(t=0) = 0 = \frac{\beta k_{\infty}n_o}{\Lambda} - \lambda C_0\Leftrightarrow C_0 =\frac{\beta k_{\infty}n_o}{\Lambda\lambda}

We can also say that at time 0 the reactivity is zero, \rho_0 = 0, (system is critical).

That leaves us with the fuel temperature. At normal operating conditions fuel temperatures are somewhere in the neighborhood of 1000 degrees Celsius or 1300 Kelvin. What we are really interested in is the total energy in the fuel, the enthalpy. For that we need to know the heat capacity of the fuel, lets take uranium dioxide as fuel material. It has a heat capacity of about 300 \large \frac{J}{KgK}. In other words it takes 300 joule to heat one kg of uranium dioxide by 1 kelvin. Lets consider one PWR fuel assembly, it contains about 500 kg of uranium dioxide. If we assume 1300 Kelvin as operating temperature then the initial enthalpy of the assembly is:

 300*1300*500 = 195*10^6J

There are 157 in most PWR's and typical power of a reactor would  be 1000 MW electric which means about 3000 MW thermal power. Each assembly then produces about  3000/157 \approx 20 MW. The hottest assemblies in a reactor can produce perhaps 50% more power than the average assembly so we end up with about 30 MW as the initial power of the assembly. This means that the cooling is taking away about 30 MW of power for the assembly to be in steady state. We can as an approximation assume the heat transfer from the assembly will continue at this level during the transient.

During a time step of size h the enthalpy change will then simply be the power during the time step, multiplied by the time step size minus the energy removed:

\Large \Delta E_n = P_nh - 30*10^6h

Temperature change is enthalpy change divided by heat capacity:

\Large \Delta T_n = \frac{\Delta E_n}{300}

Now we have everything we need to plug it all into a script and calculate the behavior of the power given different reactivity insertions and different values of the doppler coefficient. The program I wrote for this can be downloaded from this link. Doppler.exe is a windows executable and the rest of the files are the source code, license terms etc. Read the README file if you are using linux and want to run it.

Lets start by some fairly typical values for all the parameters. The dopper reactivity coefficient is usually around a few parts in hundred thousand per Kelvin. In nuclear engineering 1*10^{-5} is usually abbrevated as pcm, so lets say doppler coefficient is 5 pcm/K. We can let inserted reactivity be 220 pcm as in the example from the last blog post. Let prompt neutron lifetime be 1*10^{-3}, delayed neutron lifetime 8.66 seconds and effective delayed neutron fraction (beta effective) 0.0065 (all of these values gets set if you click on the "Typical PWR" button). Remember from last blog post that these values (except for the doppler) meant that the power doubled every 14 seconds. What happens now? This happens:

Within about two seconds the doppler feedback returns the power down to normal again and the power peak goes up to about 125%. No doubling of power, no need for operators to do anything. 220 pcm of reactivity is still within the region where delayed neutron fraction slows down everything, if we put in a reactivity that is larger than the delayed neutron fraction we get prompt criticality. Lets see if the doppler can handle that, plug in 1000 pcm of reactivity instead of 220. Then we get this:

Now we get a more dramatic power spike, max power is 400%! But again within about 2 seconds the power has dropped down again to normal. Lets look at the fuel temperature as well.

Fuel temperature increases from 1300 K to about 1550 K, melting temperature for uranium dioxide is almost 3000 K so we are far away from fuel damage. The doppler feedback in other words manages to handle a pretty hefty reactivity spike before the operator even has a chance to spray his coffee on the screen. Now you might ask what would even cause a reactivity spike in the first place? Not much actually, but one possibility is for instance if a control rod drops out of a BWR or is ejected out of the top of a PWR. The reactivity worth of a PWR control rod is in the neighborhood of  700 pcm.

Lets look at positive doppler coefficient! If we plug in +1 pcm/K and push in 300 pcm of reactivity we get this:

Scheisse!! Fuel now fails (melts, bursts, whatever happens quickest) within 2 seconds and power goes to 8300% even though initial reactivity was well below prompt criticality, what the f*uck just happened? It's simple, for every kelvin that the fuel temperature increases reactivity goes up by one pcm. When the fuel has heated up by 350 kelvin (and 350 pcm of reactivity has been added on top of the 300 pcm of initial reactivity) we are suddenly prompt critical, right back to the nasty situation of a doubling time in the tenths of a second. The delayed neutrons buys some time for some kind of automated emergency procedure but they get overpowered. This is the reason it is illegal to build reactors with positive reactivity coefficient, it just goes to hell like Chernobyl did.

Enough fun with doppler coefficients, next installment of this series will deal with heat removal from the fuel and fluid properties, this will allows us to look at moderator feedbacks.