Prediction of Neutronic and Kinetic Parameters of Ghana Research Reactor 1 ( GHARR-1 ) after 19 Years of Operation Using Monte CarloN Particle ( MCNP ) Code

The Ghana Research Reactor-1 (GHARR-1) core was modified with an addition of a 9.0 mm layer of beryllium to the top shim tray to compensate for reactivity loss due to fuel depletion after 19 years of operation. Neutronic and kinetic parameters have been predicted using Monte Carlo N-Particle Code version 5 (MCNP5) to determine whether they were within acceptable operating margins. Excess reactivity, control rod worth, moderator reactivity coefficient, delayed neutron fraction and neutron generation time have been predicted as 3.86, 6.98, −0.1218 mk/ ̊C, 8.17507 × 10 ∆k/k, and 8.147 × 10 s respectively. These parameters compared favorably with those provided in the initial Safety Analysis Report.


Introduction
The Ghana Research Reactor-1 (GHARR-1) has undergone core conversion from the use of 90.2% high enriched uranium (HEU) fuel to the use of 12.5% low enriched uranium (LEU) fuel.The conversion of the reactor core was done under a project supported by the International Atomic Energy Agency (IAEA) with the aim of reducing proliferation risks associated with HEU fuel [1].Con-version of the core was of benefit to Ghana Atomic Energy Commission (GAEC), as the reactor core had been in operation for 19 years.
During the 19 years of operation, the excess reactivity of the HEU core dropped from 4.0 mk to about 2.3 mk due to fuel depletion and accumulation of fission products, some of which are neutron poisons.To restore the excess reactivity to the required 4.0 mk, the HEU core was modified with an addition of 9.0 mm of beryllium to the top shim tray.This modification was supported by an assessment of various safety parameters [2].The assessment of safety parameters and associated operational margins was necessary to provide an indication of whether the reactor would operate safely in all modes.Besides assessment of the safety parameters, it was also necessary to document the operating conditions of the HEU core after addition of beryllium to the top shim tray as well as to update the final safety analysis report (SAR) related to conversion of core as required by the regulatory authority.Data generated under this work will also be used as reference for operating conditions of the LEU cores since it was envisaged that the operating conditions of the LEU would be similar to those of the HEU [3].
The safety parameters that were updated in the final SAR for the HEU core are; thermal hydraulics, neutronics and kinetics respectively.Thermal hydraulic parameters are not discussed under this work but have been discussed in [2].This work assesses neutronic and kinetic parameters of GHARR-1's HEU core after 19 years of operation and subject to addition of 9.0 mm of beryllium to the top ship tray.

Description GHARR-1
The reactor was designed by the China Institute of Atomic Energy (CIAE) and was commissioned in 1986.It was designed to operate with a nominal power of 30 kW and a maximum neutron flux of 1 × 10 12 n/cm 2 •s.The core was designed to use HEU as fuel, beryllium as reflector, light water as a moderator as well as a coolant and natural convection for cooling.The reactor was used for educational purposes, human resource development in nuclear science and technology, neutron activation analysis, and production of short-lived radioisotopes for use in hospitals.The core had 354 lattice positions arranged in 10 concentric circles at a pitch distance of 10.95 mm.It had 344 uranium-aluminum (U-Al) alloy fuel elements, 6 dummy aluminium elements, 4 tie rods, 2 grid plates and a central guide tube for the cadmium control rod.It was surrounded by an annular beryllium reflector and rested on a block of beryllium reflector plate.The top part of the core had a shim tray which was used for adding layers of beryllium for increased neutron reflection and reactivity insertion into the core.The top shim tray design allowed for reactivity compensation due to fuel depletion and accumulation of fission products as the reactor operated.Figure 1 shows a schematic diagram of the vertical cross-section of GHARR-1's HEU core.
The reactor was under moderated with a high negative moderator temperature reactivity coefficient (approximately −0.l mk/˚C for the temperature range World Journal of Nuclear Science and Technology

Theory
Neutron interactions with reactor materials cause different phenomena such as nuclear fission, induced radioactivity and associated decay heating.These interactions play a central role in creating the power distributions, within nuclear materials, that drive the heat transfer process in a reactor.They also help in the determination of reactor safety, core properties, reactivity control, reactor kinetics, xenon stability, fuel depletion, and are also important for isotope production [6].Reactor properties can deviate from nominal conditions due to planned or unplanned changes in reactivity and abnormal or accident conditions.The deviation from nominal conditions has an effect on the neutron flux in the reactor core.
An understanding of the time-dependent behavior of neutrons in a nuclear reactor is important in the analysis of nuclear reactor safety.The transient neutron flux behaviour resulting from a departure from criticality may arise during operational changes such as control rods movement, environmental changes such as those of boron concentration, or accidental disturbances in steady-state operating condition.This phenomenon, among others, may be classified as either reactor kinetics or reactor dynamics.
Reactor kinetics also known as reactor kinetics without feedback, is the study of the neutron flux's time-dependence for postulated changes in the macroscopic cross sections.Reactor dynamics or reactor kinetics with feedback and with spatial effects, is the study of the time-dependence of a neutron flux on the macroscopic cross sections as a function of the neutron flux level [7].The kinetic parameters of importance in safety evaluation of a nuclear reactor are the effective delayed neutron fraction, prompt neutron lifetime and neutron generation time [8].Neutronics of a nuclear reactor is the study of the distribution and multiplication of neutrons in the reactor and its effect on reactor power distribution.Diffusion codes are important in the study of neutronics and also in the determination of the neutron flux from the solution of the neutron diffusion equation.The neutron flux is important in computing the power distribution [9].It also looks at how neutrons multiply in a reactor.The neutronic parameters that are important for the dynamic reactor behaviour are the excess reactivity, shutdown margin, control rod worth, power peaking factors and the various reactivity coefficients [10].
The multiplication of neutronsis is important if a reactor has to sustain a chain reaction.In order to have a sustainable chain reaction, neutrons emitted during fission should be able to cause further fissions in other fissile or fissionable nuclei.The fission chain reaction is quantitatively expressed in terms of the effective neutron multiplication factor (K eff ).The neutron multiplication factor, expressed by Equation ( 1), is a ratio of the number of neutrons in the current generation and the number of neutrons in the preceding generation.This factor is dependent on various parameters that also vary with temperature and are linked to reactor safety.

Number of neutrons in one generation Number of neutrons in the preceding generation
When K eff is equal to unit, the reactor is said to be critical, when above or below unit, the reactor is said to be super critical or subcritical respectively.The measure of a reactor's departure from criticality is called reactivity (ρ) and is defined as a fractional change in neutron population per generation.Equation (2) shows the mathematical expression of reactivity.
The ability to predict correctly the effective multiplication factor (K eff ) and computation of reactivity (ρ) is fundamental to solutions of criticality problems.
The K eff and ρ are important in the computation of such parameters as; the delayed neutron fraction, control rod worth, shutdown margin, neutron life time, and neutron generation time (Λ).
Correct computation of the prompt neutron generation time for a nuclear reactor under design or whose configuration has been modified is important for prediction and optimization of various reactor parameters.The neutron generation time has a bearing on various interactions and environmental conditions within a reactor.The prompt neutron generation time (Λ) is the average time from birth of a neutron to a point it reproduces itself through fission.The prompt neutron generation time can also be defined as the inverse of the neutron production rate as expressed by Equation (3).This is applicable when considering a one-group diffusion-theory approximation.where V in the average neutron velocity, ν is the average number of neutrons produced per fission and ∑ f in the fission macroscopic cross section.
There is, however, a difference between neutron generation time and neutron life time.Neutron life time is the time period between the birth of a neutron and its absorption.The neutron life time can be used to determine the neutron generation time using Equation (3).In light water reactors, the prompt neutron generation time increases with fuel burn up.This may be caused by a reduction in fissile materials as well as absorption in the moderator which makes the prompt neutron live much longer [11].
One of the commonly used method for determining the neutron generation time is the 1/V absorber method.This method applies the K eff and ρ to compute the neutron generation time.In this method, all reactor materials in the core model are doped with small amounts of a 1/V neutron absorber material.The 1/V neutron absorber material inserts a negative reactivity expressed by Equation (4).
where Κ pt is the perturbed multiplication factor due to the 1/V neutron absorber material and Κ u is the steady state unperturbed multiplication factor.According to Bretscher [12] [13] the negative reactivity inserted is equal to the product of the neutron lifetime (l p ), the atomic density (n) of the 1/V neutron absorber material (atoms/b-cm), the speed of thermal neutrons (V = 2.2E5 cm/s), and the absorption cross section (σ = 3837b) of the 1/V neutron absorber material.This relation is expressed by Equation (5).
The neutron generation time is given as the limit of l p as the n approaches zero.Equation ( 6) gives the expression for the neutron generation time.The individual probabilistic events that comprise a process are simulated sequentially.The probability distributions governing these events are statistically sampled to describe the total phenomenon.It consists of following each of the many particles from a source throughout its life to its death in some terminal category (absorption, escape, etc.).Probability distributions are randomly sampled using transport data to determine the outcome at each step of its life.In this work, MCNP was used because of its ability to model complex three-dimensional geometries and its extensive validation [14].

Modification to GHARR-1's MCNP Model
The Reactor Burn up System (REBUS) code has been used to generate an inventory of isotopes produced due to fission products and fuel burn up after 19 years of GHARR-1's operation.The REBUS code was simulated under GHARR-1's operation scheme of 2.5 hours per day, 4 days a week and 52 weeks a year.This was done to produce an inventory of isotopes that was reflective of GHARR-1's operating period.This inventory of isotopes was used to update the material card's fuel composition in GHARR-1's MCNP5 model.The material card is used to define the composition of different materials in the reactor model.A 9.0 mm layer of beryllium was added to the top shim tray of the reactor core as can be seen in Figure 1.The addition of beryllium added a reactivity worth of 1.3 mk which was sufficient to compensate for reactivity loss in the core due to fuel burn up and accumulation of fission products after 19 years of operation.

MCNP Simulation Criticality Problem
Criticality problems in MCNP are simulated using the KCODE card which requires an input of the nominal number of source histories (N) per cycle, an initial guess of K eff , the number of source cycles (I c ) to be skipped before K eff can accumulate, the number of active cycles (I A ), and the total number of cycles (I t ) [15].
The control rod worth, delayed neutron fraction, neutron generation time, shutdown margin, and power peaking factor cases have been simulated with the following specifications in the model; N was set at 300,000, K eff was given an initial guess of 1.004, I t was set as 250 cycles with I c set at 20 cycles and the remaining I A set 230 cycles being active.These values were adequate to provide sufficient statistical accuracy in the results.The value for the number of cycles to be skipped was chosen to provide enough cycles for K eff to statistically converge around a certain solution, allowing for the most accurate results to be obtained.
The initial guess for K eff of 1.004 was based on the K eff for GHARR-1's clean core value, thus this value was sufficient to allow for good statistical accuracy in the prediction of K eff .A total of 69 million particles were used for each simulation.

Excess Reactivity
Simulation for K eff was performed with the control rod fully withdrawn and the KCODE card parameters set as described in Section 3.1.The K eff obtained was used to compute the excess reactivity using Equation (7).

Shutdown Margin
The core has been simulated with the control rod fully inserted into the reactor core.The K effin obtained has been used to compute ρ sdm using Equation (8).

Control Rod Worth
Two cases have been simulated to obtain the eigenvalues needed to compute the control rod worth.For the first case, the control rod was fully withdrawn, giving a multiplication factor K effout , and in the second case, the control rod was fully inserted into the reactor core giving the multiplication factor K effin .Equation ( 9) was used to compute the control rod worth (ρ crw ).

Reactivity Coefficients
Reactivity coefficients have been simulated with the following KCODE card input; N was set at 200,000, K eff was given an initial guess of 1.004, I t was set at 220 cycles with the first 20 cycles I c being inactive and the 200 cycles I A being active.
A total of 40 million particle histories have been used.

Moderator Temperature Coefficient
To obtain multiplication factors required to determine the moderator temperature coefficient, changes were made to the TMP card and all Cell cards that contain water in the MCNP model for GHARR-1.The TMP data card allows for input of equilibrium energy (MeV) of atoms in the fuel at a particular temperature.The TMP card and all Cell cards containing water in the model were respectively updated with the equilibrium temperature (MeV) and density corresponding to the temperature (˚C) that was being simulated for.The multiplication factors were predicted for the temperature range 20˚C to 100˚C.Table 1 shows the temperatures that have been used for predicting K eff and their corresponding equilibrium temperatures, densities and percent voiding.Corresponding reactivities for the temperature ranges have been computed using Equation (8).Each moderator reactivity computed has been plotted against the corresponding temperature (˚C) to obtain the relationship between temperature (˚C) and reactivity.The reactivity response to change in temperature has been determined using the 2 nd order differential equation.Equation (10) where α x is the reactivity coefficient for parameter (x).

Moderator Void Coefficient
The MCNP5 code does not have a card to input percent voiding in water due to changes in temperature.An indirect method has been used where densities corresponding to percent void as presented in Table 1 for the temperature range 20˚C to 100˚C were used.This indirect method is an accurate way for water voiding in MCNP Daniel [16].The temperature range 20˚C to 100˚C was chosen to cover both normal and higher temperatures.During normal operations of the reactor, the coolant temperature does not change significantly, however higher temperatures were chosen to study the behaviour of the moderator void coefficient at higher temperature.All Cell cards containing water and the TMP card were respectively updated with densities and corresponding temperatures (˚C) for the range 20˚C to 100˚C.Respective multiplication factors have been predicted for each temperature.Corresponding reactivities for the temperature range 20˚C to 100˚C were computed using Equation (8).Table 1 shows the densities and their relative percent voids and temperatures that were used.The relationship between temperature and reactivity was determined using the 4 th order polynomial equation.Equation (10) was used to compute the moderator void coefficient.

Fuel Temperature Coefficient
The TMP data card was updated with temperatures for the range 30˚C to 600˚C.This range was adequate to cover the operating range for the U-Al alloy fuel whose melting temperature is 640˚C.It was also sufficient to give an indication of the variation of the fuel temperature coefficient under accident conditions.The multiplication factors were predicted for the temperature range 30˚C to 600˚C.Corresponding reactivities for the temperature range 30˚C to 600˚C were computed using Equation (8).The relationship between temperature and reactivity was determined using the 4 th order polynomial equation.Equation (10) World Journal of Nuclear Science and Technology was used to compute the fuel temperature coefficient.

Effective Delayed Neutron Fraction
The effective delayed neutron fraction (β eff ) has been determined using the Prompt method which requires prediction of the multiplication factor for two cases [17].The first case involves prediction of the multiplication factor for prompt neutrons (Κ p ) only without the delayed neutrons.The second case involves prediction of the multiplication factor (Κ eff ) which combines both prompt and delayed neutrons.
In this method, the TOTNU data card with entry NO has been used to obtain the multiplication factor for prompt neutrons only.The TOTNU data card with entry NO prevents the influence of delayed neutrons on multiplicity per fission.Thus, the multiplicity per fission is due to the average prompt neutron and subsequently Κ p is predicted for all fissionable nuclides with available prompt values.
The TOTNU card without any entry has been used to predict the average number of neutrons from fission.The TOTNU card without any entry considers both prompt and delayed neutrons and the resultant effective multiplication factor expressed by Equation ( 1) is given as Κ eff .
For both cases, the KCODE card parameters were set as described in Section 3.1.The values for Κ p and the Κ eff have been used to compute the delayed neutron fraction using Equation (11).

Neutron Generation Time
The neutron generation time has been determined using the 1/V absorber method described in Section 2. Boron-10 has been used as a doping material because of its 1/V neutron absorption property.Eleven cases have been used to predict multiplication factors required for the determination of the neutron generation time (Λ).The first case was used to predict the unperturbed multiplication factor (Κ u ) an equivalent of the multiplication factor Κ eff described by Equation ( 1) for a reactor without any doping materials added.The remaining 10 cases have been used to predict the perturbed multiplication factors (Κ pt ) for10 different concentrations of boron-10 in the range 6.0 × 10 −8 -15 × 10 −8 atoms/b-cm.Concentrations of boron lower than 6 × 10 −8 atoms/b-cm have not been used since they do not provide significant change in reactivity as it has been reported by Hanson A. L. and Diamond D. J. [18].Table 2 shows the concentrations of boron-10 that were used to dope reactor core materials.For each of the 11 cases, the KCODE card has been used to predict multiplication factors with parameters set as described in Section 3. The corresponding negative reactivity inserted and the prompt neutron lives were computed using Equations ( 4) and ( 5).The various neutron life time (l p ) obtained using Equation ( 5) have been Table 2. Boron-10 concentrations for reactor core material perturbation.plotted against their corresponding boron concentrations.The neutron generation time (Λ) has been determined as the intercept on the Y-axis based on the graphical application of Equation ( 6).

Results and Discussion
Table 3 shows a comparison of neutronic and kinetic parameters from this work and those provided for in the initial GHARR-1 Safety Analysis Report.Neutronic parameters for GHARR-1 after 19 years of operation were observed to be in good agreement with those provided in the initial SAR.The excess reactivity was predicted to be within the expected margin of 3.5 mk to 4.0 mk and it was also determined to be less than the limiting value of 1/2β eff as required if the reactor is to operate normally.If the excess reactivity exceeds 1/2β eff , the reactor would tend to prompt criticality with a power excursion within seconds rendering it impossible to control.The control rod worth was predicted to be slightly higher compared to the value reported in the SAR even though the two values were fairly in agreement.The higher value from this work could be attributed to the increased margin between reactivities ρ out and ρ in described by Equation (9).It is expected that the reactivity ρ in inserted when the control rod is fully inserted into the reactor core would be lower for the core after 19 years of operation compared to the value provided in the initial SAR due to increased negative reactivity resulting from fuel depletion and accumulation of neutron poisons.This causes a wider margin between the two reactivities making the control rod worth slightly higher.The control rod worth is an important design parameter for the safety of a reactor.The shutdown margin (ρ sdm ) was predicted to be slightly higher compared to the value provided in the initial SAR.This could be attributed to a reduced value of the multiplication factor Κ effin for the control rod fully inserted into the reactor core.From Equation ( 8), the reduced Κ effin implies an increased value of ρ sdm .The reduction in the value of multiplication factor Κ effin compared to that in the SAR is influenced by increased neutron absorption due to accumulated neutron poisons in the fuel after 19 years.Despite the addition of a layer of beryllium to the top shim tray, the neutron poisons accumulated in the fuel still add a small amount of negative reactivity into the reactor core.This negative reactivity causes the ρ sdm to be slightly higher compared to that reported in the SAR.
Figure 2 compares the axial power density distribution of the reactor core after 19 years and that of the core provided in the initial SAR.The axial power density distribution has its maximum value at the center of the core.This could be attributed to the presence of water as a moderator in the control rod guide tube when the control rod is fully withdrawn, and the presence of the annular beryllium reflector around the center that causes a dense population of neutrons due to low neutron leakage in the center compared to the top and bottom parts of the reactor core.As can be observed, the power density starts decreasing axially from the center of the core.However, towards the top shim tray and the bottom beryllium block, it starts increasing again.The increase in power density towards the top and bottom beryllium blocks is influenced by increased neutron population as a result of neutrons being reflected by beryllium.
Both power distributions for the core after 19 years and that of the core provided in the initial SAR tilt lower towards the top of the reactor core.This implies the power at the top of the core in both cases is lower compared to the bottom of the core.This could be attributed to high neutron leakage at the top of the core compared to the bottom parts of both cores.The top part of the initial core is bare and that of the core after 19 years is reflected by 9 mm of beryllium compared to 50 mm of beryllium for the bottom parts of both cores.
Both the core after 19 years and that provided in the initial SAR show a similar power density distribution.However, the distribution for the core after 19 years is lower compared to that of the core provided in the initial SAR for the bottom part of the power density distribution.The lowered power density distribution observed near the bottom for the core after 19 years compared to the World Journal of Nuclear Science and Technology The maximum power peaking factor for the core after 19 years was predicted to be 1.3522 and that of the core provided in the initial SAR was predicted to be 1.3525.Both values are within the expected range or 1.3 -1.5 for light water moderated research reactors.
Figure 3 compares the integral control rod worth of the core after 19 years with that of the core provided in the initial SAR.The plot shows that the control rod worth for the core after 19 years is lower compared to that of the core provided in the initial SAR for the length 0.0 mm to 215.0 mm from the bottom of the core.Beyond 215.0 mm, the control rod worth for the core after 19 years is higher compared to that of the core provided in the initial SAR.The lower neutron population due to fuel depletion and presence of neutron poisons in the bottom part of the core after 19 years causes the control rod to absorb fewer neutrons.This is observed through the lower values of the control rod worth for the core after 19 years compared to that of the core provided in the initial SAR.
The integral control rod is an important design feature which shows how efficient a control rod is in absorbing neutrons.
Figure 4 shows a plot of the neutron lifetime (l p ) against boron-10 concentration.The neutron generation time (Λ) was determined as the intercept on the World Journal of Nuclear Science and Technology  neutron-lifetime axis.The mathematical computation is elaborated in Equations ( 3) to (6).The value for the neutron generation time for the core after 19 years of operation was determined to be slightly higher compared to that of the core provided in the initial SAR.This could be attribute, in part, to the reduced average number of neutrons produced per fission due to fuel depletion and accumulation of neutron poisons in the fuel despite the addition of 9.0 mm of beryllium to the top shim tray.In light water reactors, the prompt neutron generation time increases with fuel burnup.This may be caused by a reduction in fissile materials as well as absorption in the moderator which makes the prompt neutron live much longer.Mathematically, this could be explained using Equation (3) from which the variable is the number of neutrons produced per fission.Reduction of this value implies an increase in the prompt neutron generation time.

Moderator Void Coefficient
The moderator void coefficient for the temperature range 20˚C to 100˚C was predicted to be negative.Figure 5 shows the relationship between water voiding and reactivity.The correlation between water voiding and reactivity was determined to be negatively strong with a value of −0.9839.This entails, an increase in water voiding due to an increase in temperature, would cause a decrease in reactivity.This is a safety design feature of the reactor which renders it safe in case of an accidental insertion of abnormal reactivity that could cause boiling and voiding.In case boiling occurred, the reactivity would reduce thus preventing power excursion.Equation (12) shows the moderator void coefficient as a It can be observed from Table 4 and Figure 5 that the temperature feedback becomes more negative for greater temperatures of water.Thus the reactor will not have any detrimental effects in case of any abnormal operations resulting in temperature increase.

Conclusion
Neutronics and kinetic parameters of GHARR-1's HEU reactor core after 19 years subject to addition of 9.0 mm of beryllium to the top shim tray have been observed to be within the acceptable operating margins.The parameters were in good agreement with those of the core provided in the initial SAR.The accumulation of fission products and neutron poisons in the fuel meat had some effects on various parameters despite addition of beryllium.This caused some parameters under this work to either be higher or lower compared to those provided in the initial SAR.However, neutronic and kinetic parameters were observed to be within the allowed margins and as such the reactor was safe to operate with the addition of 9.0 mm of beryllium which compensated for the reactivity loss.The results obtained under this work are fit to be used as a basis for comparing operating conditions of the LEU cores and also for future works to be carried out on the LEU core model.

B
. M. Mweetwa et al.DOI: 10.4236/wjnst.2018.84014164 World Journal of Nuclear Science and Technology Carlo N-Particle Transport (MCNP) is a general-purpose, continuousenergy, generalized-geometry, time-dependent, coupled neutron/photon/ electron Monte Carlo Transport code.The MCNP code is useful for complex problems that cannot be modeled by computer codes that use deterministic methods.Various cases can be simulated in several transport modes: neutron only, photon only, electron only, combined neutron/photon transport where the photons are produced by neutron interactions, neutron/photon/electron, photon/electron, or electron/photon.The neutron energy regime is from 10 −11 MeV to 20 MeV, and the photon and electron energy regimes are from 1 keV to 1000 MeV.The capa-World Journal of Nuclear Science and Technology bility to calculate K eff eigenvalues for fissile systems is also a standard feature.

Figure 2 .
Figure 2. Comparison of axial power distribution of the ırradiated and clean cores.

Figure 3 .
Figure 3.Comparison of the Integral control rod for the irradiated and clean cores.

Figure 4 .
Figure 4. Plot of neutron life time against boron concentration.

Table 1 .
Relationship between temperature, equilibrium energy density and percent void.
has been used to compute the moderator temperature coefficients.B. M. Mweetwa et al.DOI: 10.4236/wjnst.2018.84014167 World Journal of Nuclear Science and Technology

Table 3 .
Comparison of neutronic and kinetic parameters from this work and the SAR.