Development of an MCNP-tally based burnup code and validation through PWR benchmark exercises

Annals of Nuclear Energy 36 (2009) 626–633

Contents lists available at ScienceDirect

Annals of Nuclear Energy journal homepage: www.elsevier.com/locate/anucene

Development of an MCNP-tally based burnup code and validation through PWR benchmark exercises B. El Bakkari a,*, T. El Bardouni a, O. Merroun a, Ch. El Younoussi a, Y. Boulaich a, E. Chakir b a b

ERSN-LMR, Department of physics, Faculty of Sciences P.O.Box 2121, Tetuan, Morocco EPTN-LPMR, Faculty of Sciences Kenitra, Morocco

a r t i c l e

i n f o

Article history: Received 7 September 2008 Received in revised form 21 December 2008 Accepted 29 December 2008 Available online 12 February 2009

a b s t r a c t The aim of this study is to evaluate the capabilities of a newly developed burnup code called BUCAL1. The code provides the full capabilities of the Monte Carlo code MCNP5, through the use of the MCNP tally information. BUCAL1 uses the fourth order Runge Kutta method with the predictor–corrector approach as the integration method to determine the fuel composition at a desired burnup step. Validation of BUCAL1 was done by code vs. code comparison. Results of two different kinds of codes are employed. The first one is CASMO-4, a deterministic multi-group two-dimensional transport code. The second kind is MCODE and MOCUP, a link MCNP–ORIGEN codes. These codes use different burnup algorithms to solve the depletion equations system. Eigenvalue and isotope concentrations were compared for two PWR uranium and thorium benchmark exercises at cold (300 K) and hot (900 K) conditions, respectively. The eigenvalue comparison between BUCAL1 and the aforementioned two kinds of codes shows a good prediction of the systems’ k-inf values during the entire burnup history, and the maximum difference is within 2%. The differences between the BUCAL1 isotope concentrations and the predictions of CASMO4, MCODE and MOCUP are generally better, and only for a few sets of isotopes these differences exceed 10%. Ó 2009 Elsevier Ltd. All rights reserved.

1. Introduction Burnup and depletion codes have been developed and used in the nuclear industry since the introduction of digital computing. These codes solve the diffusion equation in one to three dimensions using few neutron energy groups. Only few of the major fission products were included in the calculation. During last years, the enhancement in digital computing capabilities and the amelioration of neutron cross section evaluations have led to the development of more sophisticated numerical techniques, such as Monte Carlo method, to look for system eigenvalues. In burnup codes, the highest fidelity approach uses neutron absorption and fission reaction information generated via neutronics codes to determine the nuclide composition at a next time step. This kind of model allows the integration of all the neutron flux information into the calculation without post-processing and additional manipulation of neutron flux and cross-sections set (Parma, 2002). Neutron absorption and fission reaction data for individual nuclide are available as output from the Monte Carlo codes like MCNP through the use of tallies. The only requirement is that a

* Corresponding author. Tel.: +21251235092; fax: +21239994500. E-mail address: [email protected] (B. El Bakkari). 0306-4549/$ - see front matter Ó 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.anucene.2008.12.025

pointwise energy-dependent cross section set which is available for each nuclide of interest at required temperature. BUCAL1 was developed to incorporate the neutron absorption tally/reaction information generated directly by MCNP5 in the calculation of fissioned or neutron-transmuted isotopes for different regions. Using MCNP tally information directly in the computation allows for almost all data to be used in the analysis. Also, using MCNP and the continuous – energy cross-sections allows for cross section manipulations to be avoided. Accomplishing this goal would allow performing straightforward and accurate calculation without having to use the calculated group fluxes to perform transmutation analysis in a separate code. The algorithm chosen for BUCAL1 to solve the depletion equation system is based on the use of the fourth order Runge Kutta method with the predictor–corrector approach. This method consists on dividing the time of irradiation into multi-time steps (hi) and to perform burnup calculation at each time step (i) using the initial conditions. The implementation of the predictor–corrector approach into BUCAL1 numerical algorithm allows the correction of flux variations for large time steps. So, this approach permits the use of small number of time steps in burnup calculations. The BUCAL1 strategy consists of using the nuclide inventory, MCNP tally information, power density, and other data to

B. El Bakkari et al. / Annals of Nuclear Energy 36 (2009) 626–633

determine the new nuclide inventory for a given region of the core at a new time step. Then the new inventories are automatically placed back into MCNP input file and the case run for a new subsequent time step. The outline of this paper is as follows: in the first section, the mathematical model and related approximation techniques adopted in BUCAL1 are presented. The second section is devoted to describe the procedures to generate the neutron cross section libraries employed by BUCAL1. In the penultimate section, a physical and geometrical description for the benchmarks used is presented. The results of the burnup study obtained with BUCAL1 and other codes used are presented and analyzed in the final section.

2. Mathematical model 2.1. Overview In routine reactor design burnup calculations, the key objective is to determine the time-dependent fuel material compositions as well as the eigenvalues as a function of burnup. Two basic mechanisms of fuel depletion are under consideration: (i) various nuclear reactions such as nuclear fissions, neutron captures, etc. and (ii) the decay of radioactive isotopes. Once material compositions are known, eigenvalues can then be calculated efficiently using MCNP code for specified geometry. Mathematically, the material balance process can be described at any time by the following depletion equation:

X dNi X cji rf ;j Nj / þ rc;k!i Nk / ¼ dt j k X kl!i Nl  ðrf ;i Ni / þ ra;i N i / þ ki Ni Þ þ

ð1Þ

l

Where P

dNi dt

= time rate of change in concentration of isotope i,

unit volume of isotope i from fisj cji rf ;j N j / = production rate perP sion of all fissionable nuclides, k rc;k!i N k / = production rate per

unit volume of isotope i from neutron transmutation of all isotopes P including (n, c), (n, 2n), etc., l kl!i N l = production rate per unit volume of isotope i from decay of all isotopes including b, b+, a, c, etc., rf ;i Ni / = removal rate per unit volume of isotope i by fission, ra;i Ni / = removal rate per unit volume of isotope i by neutron absorption (excluding fission), ki N i = removal rate per unit of isotope i by decay. Eq. (1) can also be written in a compact vector form:

dNi ðtÞ ¼ Ai ðtÞNi dt

ð2Þ

where the matrix A is called the transition matrix.

627

where i

h ¼ t i  t i1 i

i1 i1 i1 K a;j ¼ h fj ðti1 ; N i1 1 ; N2 ; . . . ; Nj ; . . . ; Nn Þ i

K b;j ¼ h fj t i1 þ i

K c;j ¼ h fj ti1 þ i

i

!

i

!

h K a;1 i1 K a;2 K a;j K a;n þ ; N i1 þ ; N2 þ ; . . . ; N i1 ; . . . ; N i1 j n þ 2 1 2 2 2 2 K b;j h K b;1 i1 K b;2 K b;n þ ; Ni1 þ ; N2 þ ; . . . ; N i1 ; . . . ; N i1 j n þ 2 1 2 2 2 2 i

i1 i1 K d;j ¼ h fj ðt i1 þ h ; N i1 þ K c;j ; . . . ; N i1 1 þ K c;1 ; N 2 þ K c;2 ; . . . ; N j n þ K c;n Þ

ð4Þ

The depletion algorithm of BUCAL1 assumes that the flux spectrum of reactor is constant during the entire burnup step. This assumption, however, may not be completely correct if one uses too large burnup steps since the spectrum of a reactor changes during such steps. Hence, the reaction rates and the power distribution calculated at the beginning of the step might not adequately account for the changes during the entire burnup step. In this case, BUCAL1 uses a more accurate approach for the depletion calculation known as ‘‘the predictor–corrector depletion approach”. This approach involves the following multistep process:  A burnup calculation is completed in BUCAL1 to the final time step [ti ? tf] (Predictor step).  Fluxes and reaction rates are recalculated in a steady- state MCNP calculation at the final time step ‘‘tf”.  Then the recalculated fluxes and reaction rates are used to burn over the full time step [ti ? tf] (Corrector step).  The average atom densities from these two calculations are taken as the end-of-time step material compositions. This approximation is true only if the flux shape between the two time steps varies linearly; this approach is usually an acceptable approximation. Implementing this approach allows the user to burn a system using fewer burnup steps than if no approximation were made on the average flux behavior. 2.3. Total reaction rate calculation Various reaction rates and the one-group flux in each individual active cell are provided by MCNP flux tallies as:

(

/i ¼

R

Rijk ¼

/i ðEÞdE R j rk ðEÞ/i ðEÞdE

ð5Þ

where rjk ðEÞ = Microscopic cross section of reaction type k for isotope j, /i ðEÞ = The region averaged one-group flux in cell i, Rijk = Reaction rates of type k with nuclide j in cell i. These tallies in MCNP come from track length estimation of cell flux and reaction rates (tally type F4 in MCNP). In the current version of BUCAL1, the two groups of nuclides under consideration are:

2.2. Solution technique The majority of burnup codes used by the nuclear community use the matrix exponential method to solve the depletion equation in its compact vector form (2). In BUCAL1 a solution technique based on the fourth order Runge Kutta method is proposed. Using this method, the solution of Eq. (1) can be transformed into a system of linear Eqs. (3) and (4) those are easy and fast to solve.

8 i 1 > N1 ¼ Ni1 > 1 þ 6 ðK a;1 þ 2K b;1 þ 2K c;1 þ K d;1 Þ > > < i i1 N 2 ¼ N2 þ 16 ðK a;2 þ 2K b;2 þ 2K c;2 þ K d;2 Þ > > ::: > > : i 1 N n ¼ Ni1 n þ 6 ðK a;n þ 2K b;n þ 2K c;n þ K d;n Þ

ð3Þ

 actinides (ACT) that contain heavy metal nuclides with atomic number Z P 90 and their decay daughters;  fission products (FP) produced by fissions and their decay/capture daughters. Specifically, the calculated reaction rates of fission products and actinides are shown on Table 1. Only the neutron capture cross section is considered for fission products since the neutron absorption of fission products is primarily via (n, c) reaction. For actinides, four types of reactions are considered including capture, fission, (n, 2n) and (n, 3n) because of their important fission products generating and higher-mass actinides evolving. Note that it is not practical to calculate MCNP reaction rate for all nuclides in BUCAL1 (910 nuclides) due to the excessive CPU

628

B. El Bakkari et al. / Annals of Nuclear Energy 36 (2009) 626–633

distributed with the MCNP5 code, based on the ENDF/B-VI.8 evaluation.

Table 1 MCNP-tallied reactions. Reaction type

MCNP reaction identifier

Actinides

(n, y) (n, f) (n, 2n) (n, 3n)

102 6 16 17

Fission products

(n, y)

102

4. Benchmarks and computational tools

time needed and the unavailability of many MCNP cross-sections. Only a limited set of important ones (e.g., strong neutron absorption contributors), typically 45 actinides and 102 fission products are considered in MCNP calculations. For the less capturing nuclides, BUCAL1 uses approximate model based on Eq. (6), to take into account their capture reaction rates (Tan Dat, 1996):

"

r

c i/

¼ /M

c i;2200

r

# rffiffiffiffirffiffiffiffiffiffiffiffiffi p 293 þ rIi T 4

ð6Þ

Where, rci / = microscopic capture rate, /M = Maxwellien flux, rci;2200 = thermal neutron capture cross section at 293 K, T = temperature of neutrons causing fission, r = ratio of epithermal flux per unit of lethargy by the Maxwellien thermal flux, Ii = Resonance integral capture cross-section on isotope i. 2.4. Flux normalization parameters MCNP generates reaction rates for the BUCAL1 depletion calculation. These generated reaction rates are normalized to fissionsource-neutron, which need to be multiplied by a constant factor to take into account the entire reactor power. The typical way of calculating this constant factor (CF) as recommended in the MCNP manual (X-5 Monte Carlo Team, 2003) is:

CF ¼ P  t=Q:Keff

ð7Þ

where, P = Total power of the entire system (watts), t = Average fission neutrons per fission event, Q = Average recoverable energy per fission event (J/fission), Keff = Eigenvalue of the system. MCNP calculates the system-averaged t and Q values. The power level used for the constant factor is the total power of the system and it is a user input. The Q value used in the calculation is an estimation of the recoverable energy per fission event. BUCAL1 uses the prompt recoverable energy per fission (Q-prompt) multiplied by the constant normalization factor (1.111) recommended by the MCNPX2.6 developers (Hendrichs, 2007). This constant normalization factor is an estimate of the delayed fission energy and capture gamma energy contribution. Fig. 1 shows the resumed flow diagram of BUCAL1 burnup code. 3. Library generation The NJOY99 code (MacFarlane, 2000) with its latest update file up259 has been used to process the source evaluated nuclear data files into libraries suitable for use with the MCNP code. Fig. 2 shows the procedures used to generate the library through NJOY data processing. The probability tables for considering the selfshielding effects in the unresolved energy range can be used in the MCNP5. These tables were generated using the PURR module of NJOY and introduced into the MCNP libraries. For benchmark calculations, the ENDF/B-VII (Chadwick et al., 2006) evaluated neutron reaction data released from the Brookhaven National Laboratory, and JEFF-3.1 (The JEFF Team, 2006) elaborated by the NEA data bank were used. These libraries were prepared at the required benchmark temperatures (Tables 2a and 3a. For the S(a, b) treatment we have used the standard libraries

In this study, two pin-cell benchmark exercises described by Xu (2003) and Weaver (2000) were studied to test the validity BUCAL1 code. These benchmarks use UO2 (Xu, 2003) and ThO2–UO2 (Weaver, 2000) as fuel material, respectively. Fig. 2 shows the pin-cell model representing the unit cell of a Westinghouse PWR fuel bundle. Both burnup calculations described in this study are based on this model. The results obtained are benchmarked against CASMO-4, MCODE (first benchmark case) and MOCUP (second benchmark case). The first case under study is a 2-D single pin-cell model of a standard Westinghouse 17  17 PWR assembly. Four different regions are considered: the fuel pellet (entire fuel region is taken as one lumped cell), the gap, the cladding and the associated coolant. The fuel is UO2 with very high U-235 enrichment (9.75w/o) which allows achieving high burnup. Boundaries of the cell are set to be reflecting, both in radial and axial directions in MCNP calculation. The geometrical and operational parameters as well as the initial compositions are presented on Tables 2a and 2b. The second case under study involves the analysis of a PWR pincell excised from a standard 17  17 pin assembly typical of large Westinghouse PWRs. The usual all-UO2 fuel pellets were replaced by a ThO2–UO2 mixture at 94% of theoretical density consisting of 75w/o Th, 25 w/o U on a heavy metal basis, with the latter enriched to 19.5 w/o U-235, to give an overall enrichment of 4.869 w/o U-235 in total heavy metal. Boundaries of the cell are set to be reflecting, both in radial and axial directions. The benchmark parameters are presented on Tables 3a and 3b (see Fig. 3). CASMO-4 (Weaver, 2000) is a deterministic multi-group twodimensional transport code for standard LWR burnup calculations from Studsvik. It is based on the evaluated data files JEF-2.2 and ENDF/B-VI, which were developed at the OECD/NEA Data Bank and Brookhaven National Nuclear Data Center, respectively. MCODE (Xu, 2003) is an MCNP4C–ORIGEN2.1 linkage code from MIT (Massachusetts Institute of Technology) that uses the matrix exponential method with the predictor–corrector algorithm as integration method. In the transport calculations, the cross sections for most of the isotopes are taken from ENDF/B-V and ENDF/B-VI, while for the others, data from libraries evaluated at different laboratories are taken. In the burnup calculations, the cross sections not provided by MCNP are taken from the ORIGEN one-group cross section library PWRUE.LIB (3-cycle PWR library, thermal spectrum) and the decay data are taken from DECAY.LIB (which comprises a total of 1307 different nuclides, including 129 actinides and 879 fission products). MOCUP (Weaver, 2000) is the MCNP–ORIGEN2 Coupled Utility Program. It employs the MCNP4B generalized-geometry Monte Carlo transport code to produce the neutronics solution and the ORIGEN2 code to compute the time-dependent compositions of the individually selected MCNP cells. All data communication between the two codes is accomplished through the MCNP and ORIGEN2 input/output files. This allows a general material (target, fuel, control, etc.) to be depleted in a neutral particle field, with the accuracy of a transport neutronics solution. Since the MCNP version 4B library does not contain temperature-dependent neutron cross sections of most actinides, a number of libraries from the UTXS compilation were imported. Also for some fission products, the evaluated data files elaborated at Los Alamos National Laboratory were imported via INEEL (Idaho National Engineering & Environmental Laboratory, US). The main features of each code are summarized in Table 4.

B. El Bakkari et al. / Annals of Nuclear Energy 36 (2009) 626–633

629

Run the program and read BUCAL1 input

Pre-process initial MCNP input and run MCNP

Extract beginning of time step reaction rates and flux values Calculate the capture reaction rates for the less capturing fission products

Loop on all time steps

Loop on burnup cycles

Run BUCAL1 depletion calculations for all active cells (predictor) Update MCNP input based on BUCAL1 output material composition and run MCNP Extract end of time step reaction rates and flux values Calculate the capture reaction rates for the less capturing fission products Run BUCAL1 depletion calculation for all active cells (corrector) Average the predictor and corrector material compositions, update MCNP input and re-run MCNP

Cooling test ?

No

Yes Loop on cooling time steps

Run BUCAL1 cooling for all active cells

No

Finish Cooling time steps ? Yes END

Fig. 1. Flow diagram of BUCAL1.

Evaluated Data

ACE-Foramt Library

RECONR

ACER

BROADR

PENDF

HEATR

PURR

GASPR

THERMR

Fig. 2. Flow diagram of NJOY99 processing for MCNP5 library.

In the validation procedure of BUCAL1 code, 102 fission products and 45 actinides are considered in neutron transport calcula-

Table 2a UO2 pin-cell model parameters. Parameters

Values

Fuel pellet radius (cm) Cladding inner radius (cm) Cladding outer radius (cm) Pin pitch (cm) Fuel density (g/cm3) Fuel temperature (K) Cladding density (g/cm3) Cladding temperature (K) Coolant density (g/cm3) Coolant temperature (K) Power density (kW/liter core) Specific power (W/gU)

0.4096 0.4178 0.4750 1.26 10.3 300 6.550 300 0.997 300 104.5 34.6679

tions. These nuclides have been chosen to account for more than 99% of total mass and neutron absorptions. 5. Results and discussion MCNP–BUCAL1 calculations are performed on an Intel core 2 Duo Personal Computer with CPU 2.2 GHz and a cache size of

630

B. El Bakkari et al. / Annals of Nuclear Energy 36 (2009) 626–633

Table 2b Initial compositions (cold conditions at 300 K).

Table 3b Initial compositions (at Hot Full Power conditions).

Nuclide

Weight percent (w/o)

Atom density (1/cm3)

Fuel (9.75 w/o UO2)

U-234 U-235 U-238 O-16

0.0688 8.5946 79.4866 11.8500

1.82239E+19 2.26826E+21 2.07128E+22 4.59686E+22

Cladding

Zircaloy-4 H-1

100 11.19

4.344182E+22 6 66295E+22

Cladding

Zircaloy-4

100

4.31438E+22

Coolant

O-16

88.81

3.33339E+22

Coolant

H-1 O-16

11.19 88.81

4.71053E+22 2.35662E+22

Fuel (ThO2–UO2)

Clad

Nuclide

Weight percent (w/o)

Atom density (1/cm3)

Th-232 U-234 U-235 U-238 O-16

65.909 0.034 4 291 17.740 12.026

1.61215E+22 8.24518E+18 1.03615E+21 4.2295 7E+21 4.26835E+22

MCODE) that lead to an eigenvalue statistical error of about 130 pcm. The entire MCNP–BUCAL1 calculations take about 15 h. 5.1. PWR UO2 benchmark

Fuel

H2O Fig. 3. Pin-Cell Model.

Table 3a ThO2–UO2 pin-cell model parameters. Parameters

Values

Fuel temperature (K) Power density (KW/kgHM) Power density (KW/liter cell) Fuel density (g/cm3) Cladding temperature (K) Cladding density (g/cm3) Coolant pressure (bars) Coolant temperature (K) Coolant density (g/cm3) Fuel pellet radius (cm) Cladding inner radius (cm) Cladding outer radius (cm) Pin pitch (cm)

900 38.1347 107.284 9.424 621.1 6.505 155.13 583.1 0.705 0.41274 0.41896 0.47609 1.2626

2 MB, under 64-bits Linux system. Each MCNP run is done for 225.000 neutron histories (same neutron histories required as

Fig. 4 shows the comparison of eigenvalue history among BUCAL1 using two neutron cross section libraries ENDF/B-VII and JEFF-3.1, CASMO-4 and MCODE. It can be seen from this figure that before 40 MWD/kgIHM there is nearly a constant difference Dk (300 pcm) between BUCAL1, CASMO-4 and MCODE. Then it grows almost in a linear way. From 0 to 90 MWD/kgIHM, the BUCAL1 system eigenvalues are closer to MCODE than to CASMO-4. At 100 MWD/kgIHM burnup, the eigenvalue difference for BUCAL1 using ENDF/B-VII is about 0.1% from CASMO-4 and 1.6% from MCODE. When using JEFF-3.1 library, the eigenvalue obtained by BUCAL1 differs by 0.1% and 1.7% from CASMO-4 and MCODE, respectively. Several causes might contribute to the eigenvalue difference remarked between the three codes, such as neutronic/burnup coupling algorithms and statistical error propagation. The ENDF/B-VII library provides an eigenvalue slightly higher than that obtained when using JEFF-3.1 library. In addition to the eigenvalue comparison, an isotope composition comparison at 100 MWD/kgIHM is presented in Table 5. From this table it appears that most of calculated material compositions agree well with CASMO-4 and MCODE. Taking into consideration the high burnup the differences are acceptable. Table 5 also confirms that overall there are more fissile species present from MCODE which explains the high reactivity produced by this code. The Am and Sm isotopes concentrations differences between BUCAL1 and CASMO-4 are significant, which need further attention. 5.2. PWR ThO2–UO2 benchmark Fig. 5 shows the comparison of eigenvalue history among CASMO-4, MOCUP and BUCAL1 using ENDF/B-VII and JEFF-3.1 nuclear data libraries. It can be seen that BUCAL1 results obtained for the two neutron libraries are in good agreement with those obtained by CASMO-4 and MOCUP. Considering that the point of major con-

Table 4 Summary of benchmarking codes.

Cross section libraries Code developer Transport treatment Resonance treatment Number of energy groups Burnup algorithm Actinide representation Fission products

MCODE

MOCUP

CASMO-4

BUCAL1

ENDF/B-V ENDF/B-VI + other evaluated libraries MIT Monte Carlo Monte Carlo Continuous Standard predictor– corrector 39 100

ENDF/B-V ENDF/B-VI

ENDF/B-VI, JEF2.2

ENDF/B-VII JEFF-3.1

INEEL Monte Carlo Monte Carlo Continuous Beginning-of-timestep representation 17 41

Studsvik KRAM characteristics Collision probability 70 Standard predictor– corrector Th231thru Cf252 200

ERSN-LMR Monte Carlo Monte Carlo Continuous Standard predictor– corrector 45 102

631

B. El Bakkari et al. / Annals of Nuclear Energy 36 (2009) 626–633

1.6 MCODE [Ref.11]

1.5

CASMO-4 [Ref.11] Bucal1/JEFF-3.1

Eigenvalue, k-inf

1.4 Bucal1/ENDFB-VII

1.3

1.2

1.1

1.0

0.9 0

10

20

30

40

50

60

70

80

90

100

Bu(MWD/kgIHM) Fig. 4. Eigenvalue comparison between MCODE, CASMO-4, and BUCAL1.

Table 5 Fractional Differencea in nuclide concentration at 100 MWd/kgIHM from CASMO-4. Isotopes

CASMO-4

MCODE

BUCAL1 ENDFB-VII

BUCAL1 JEFF-3.1

Mo-95 Tc-99 Ru-101 Rh103 Ag109 Cs-133 Cs-135 Nd-143 Nd-145 Sm-147 Sm-149 Sm-150 Sm-151 Sm-152 Eu-153 U-234 U-235 U-238 Np-237 Pu-238 Pu-239 Pu-240 Pu-241 Pu-242 Am-241 Am-242m Am-243 Total actinides Total fissile Total fertile

1.2228E+20 1.1686E+20 1.1927E+20 4.6015E+19 6.9910E+18 1.1452E+20 6.9820E+19 7.4246E+19 7.1091E+19 9.5715E+18 1.2455E+17 2.6757E+19 7.6817E+17 9.3945E+18 1.1838E+19 6.7125E+18 2.5952E+20 1.9672E+22 3.4234E+19 1.9665E+19 1.4767E+20 6.3106E+19 4.2801E+19 2.6228E+19 2.3505E+18 3.3827E+16 6.2320E+18 2.0280E+22 4.9281E+20 1.9788E+22

0.002 0.045 0.003 0.033 0.146 0.082 0.004 0.003 0.001 0.141 0.058 0.081 0.103 0.162 0.114 0.012 0.021 0.002 0.089 0.083 0.056 0.088 0.051 0.031 0.100 1.030 0.232 0.001 0.030 0.002

0.019 0.006 0.004 0.043 0.090 0.068 0.134 0.030 0.024 0.418 0.230 0.006 0.257 0.175 0.113 0.075 0.086 0.001 0.025 0.006 0.005 0.055 0.033 0.020 0.592 0.480 0.282 0.002 0.045 0.001

0.019 0.019 0.006 0.057 0.016 0.064 0.119 0.046 0.005 0.399 0.253 0.011 0.238 0.156 0.103 0.069 0.082 0.001 0.016 0.041 0.005 0.072 0.036 0.047 0.427 0.496 0.210 0.001 0.042 0.0003

a

Fractional diff. = (N  Ncasmo-4)/Ncasmo-4, where N is nuclide concentration (at./b-cm).

cern is the burnup value where reactivity reaches 0.03, (which is representative of an n-batch core-average end of cycle value, with allowance of 3% reactivity loss for leakage), this eignevalue comparison shows almost no difference at that point between BUCAL1/JEFF-3.1, CASMO-4 and MOCUP. Whereas, BUCAL1/ENDFBVII overestimates the last three codes results by approximately 0.8%. Generally, this is encouraging because one must achieve better accuracy for thorium fuelled cores than for all-uranium fuelling to achieve equal accuracy in cycle length estimates (Weaver, 2000)

Up to 40 MWD/kgIHM, BUCAL1 results using the two libraries under study become much closer to MOCUP than to CASMO-4. At 60 MWD/kgIHM (corresponding to projected end-of-life core-average burnup), BUCAL1/ENDFB-VII produces an eigenvalue difference of about +2% from CASMO-4 and +0.5% from MOCUP. BUCAL1/JEFF3.1 eigenvalue difference at this burnup value is about +1.2% and 0.05% from CASMO-4 and MOCUP, respectively. The concentrations of the 17 actinides whose information is provided in the benchmark exercise at 60.749 MWD/kgIHM (which

632

B. El Bakkari et al. / Annals of Nuclear Energy 36 (2009) 626–633

1.30 CASMO-4 [Ref. 4]

1.25

MOCUP [Ref. 4] Bucal1/ENDFB-VII

1.20

Bucal1/JEFF-3.1

1.15

k-inf

1.10 1.05 1.00 0.95 0.90 0.85 0.80

B1, where ρ = 0.03

0

20

40

60

80

Bu(MWD/kgIHM) Fig. 5. Eigenvalue comparison as a function of burnup.

Table 6 Fractional difference in isotope concentration at 60.749 MWd/kg. Isotopes

CASMO-4

MOCUP

BUCAL1 ENDFB-VII

BUCAL1 JEFF-3.1

Th-232 Pa-231 Pa-233 U-232 U-233 U-234 U-235 U-236 U-238 Np-237 Np-238 Np-239 Pu-238 Pu-239 Pu-240 Pu-241 Pu-242 Total fissile Total Actinide Depleteda RatioofTh232to U238Depletionb

1.53769E+22 1.70440E+18 1.95229E+19 1.56006E+18 2.74202E+20 5.15172E+19 1.78104E+20 1.39420E+20 3.88419E+21 1.82660E+19 5.46096E+16 7.61806E+17 8.90932E+18 5.37090E+19 1.82233E+19 1.90707E+19 9.96772E+18 7.54683E+20 0.062601217 2.15589

0.003 0.048 0.035 0.034 0.040 0.176 0.021 0.054 0.004 0.058 0.037 0.043 0.026 0.071 0.032 0.024 0.036 0.024 +0.010 +0.107

0.003 0.004 0.238 0.057 0.061 0.172 0.034 0.045 0.004 0.044 0.338 0.312 0.059 0.050 0.031 0.052 0.045 0.002 +0.022 +0.109

0.003 0.107 0.236 0.063 0.067 0.188 0.031 0.035 0.005 0.014 0.317 0.316 0.017 0.056 0.030 0.055 0.014 0.002 +0.020 +0.120

a b

Total Actinide Depleted = (NActinide,t  NActinide,0)/NActinide,0, where NActinide,t is the total amount of actinides at time t; NActinide,0 is the total amount of actinides at time 0. Ratio = (Th-232 depleted)/(U-238 depleted).

is the upper limit of discharge burnup if a 3-batch core refuelling scheme is considered) are provided on Table 6. One observation is that, the inventory prediction obtained by BUCAL1 for thorium and uranium actinide chains agrees well with CASMO-4 and MOCUP. Only few isotope concentration differences exceed the 10%, especially Pa-233, U-234, Np-238, and Np-239. Total end-of-life heavy metal destruction is about 2% higher in BUCAL1 using the two libraries ENDF/B-VII and JEFF-3.1. Another point of interest is the large difference in U-234 concentration, which merits further attention, even though this nuclide has a small effect on eigenvalue variations. 6. Conclusion New burnup code utility called BUCAL1 was developed. The code uses the neutron absorption and fission tally information provided by the Monte Carlo neutronics code MCNP (version 5) di-

rectly into burnup/production calculative model. This allows for a more direct solution technique to be employed without the use of the spatial and energy-dependent neutron flux results. Compared to other linkage codes, BUCAL1 excels in its unique way of coupling with the MCNP code. The results of verification study show that the solution technique used in BUCAL1 is accurate and robust. From the code – vs. – code benchmarking of the high burnup UO2 lattice under cold conditions, it is concluded that BUCAL1 is suitable and ready to be used in burnup calculations. Even for high burnup cases, the material composition predictions are also acceptable compared to the comprehensive uranium benchmark reported by OECD (DeHart, 1996). Based on the results of the intercomparison done between BUCAL1, CASMO-4 and MOCUP, it appears that BUCAL1 can do also thorium related calculations with an acceptable agreement. From an engineering point of view, BUCAL1 is sufficiently accurate and gives consistent results.

B. El Bakkari et al. / Annals of Nuclear Energy 36 (2009) 626–633

References Chadwick, M.B., et al., 2006. ENDF/B-VII.0: Next Generation Evaluated Nuclear Data Library For Nuclear Science And Technology, Nuclear data Sheets, vol. 107, pp. 2931–3060. http://www.nndc.bnl.gov/exfor7/endf00.htm. DeHart, M.D., Brady, M.C., Parks, C.V., 1996. OECD/NEA Burnup Credit Calculational Criticality Benchmark Phase 1-B Results, ORNL-6901. Oak Ridge National Laboratory, June 1996. Hendrichs, Jh.S., Mckinney, Gr.W., 2007. MCNPX, Version 26E-, LA-UR-07-6632, Los Alamos National Laboratory (November 17, 2007). MacFarlane, R.E., 2000. NJOY-99. Nuclear Data Processing System. . Parma, Ed.J., 2002. BURNCAL: A Nuclear Reactor Burnup code Using MCNP Tallies, Sandia National Laboratories, SAND2002-3868 (November 2002).

633

Tan Dat, H., 1996. Mise au point de la bibliothèque physique et développement du système PEPABAC(PROGICIEL FAKIR) à l’aide des codes PEPIN et APOLLO1, thèse du doctorat CEA Saclay – Note CEA-N-2815 (December 1996). The JEFF Team, 2006. The JEFF-3.1 library, JEFF Report 21, NEA/OECD No. 6190, 2006, ISBN92-64-02314-3. . Weaver, K.D., Zhao, X., Pilat, E.E., Hejzlar, P., 2000. A PWR thorium pin cell burnup benchmark – INEEL/CON-2000-00227,In: Advances in reactor physics and mathematics and computation into the next millennium (PHYSOR 2000), May 7–11, 2000, IDAHO National Engineering and Enviromental Laboratory. X-5 Monte Carlo Team, 2003. MCNP- A General Monte Carlo N-Particle Transport Code, Version 5-Los Alamos National Laboratory. April 24, 2003. Xu, Zh., 2003. Design Strategies for optimizing High Burnup Fuel in Pressurized Water Reactors – Doctor Thesis, Massachusetts Institute of Technology. January 2003.