Overview

Serpent is a three-dimensional continuous-energy Monte Carlo reactor physics burnup calculation code, developed at VTT Technical Research Centre of Finland since 2004. The code is specialized in two-dimensional lattice physics calculations, but the universe-based geometry description allows the modeling of complicated three-dimensional geometries as well. The suggested applications of Serpent include:

The main features of the code, including the major limitations, are introduced below.

Geometry and neutron tracking

Serpent uses a universe-based geometry model (similar to MCNP and Keno-VI), which allows the description of practically any two- or three-dimensional fuel or reactor configuration. The geometry consists of material cells, defined by elementary and derived surface types.

The code also provides some additional geometry features specifically for fuel design. These features include simplified definition of cylindrical fuel pins and spherical fuel particles, square and hexagonal lattices for LWR geometries and circular cluster arrays for CANDU fuels. The random dispersion of microscopic fuel particles in high-temperature gas-cooled reactor fuels and pebble distributions in pebble-bed type HTGR cores can be modeled using geometry routines specifically designed for the task (see below).

Woodcock delta-tracking method

An analog Monte Carlo game and the k-eigenvalue criticality source method are used for simulating a self-sustaining chain reaction, and an external source mode is available for sub-critical and non-multiplying systems. Neutron transport is based on a combination of conventional surface-to-surface ray-tracing and the Woodcock delta-tracking method (Woodcock, 1965). The geometry routine has proven efficient and well suited for lattice calculations. The combination of two tracking methods overcomes the efficiency problems normally encountered with delta-tracking in the presence of localized heavy absorbers (Leppänen, 2010b).

The main drawback of delta-tracking is that the track-length estimate of neutron flux is not available and reaction rates have to be calculated using the potentially less-efficient collision estimator. This is usually not a problem in lattice calculations when reaction rates are scored in regions of high collision density. However, the efficiency of the collision estimator becomes poor in small or optically thin volumes located far or isolated from the active source. For this reason the Serpent code is not particularly well suited for shielding and detector calculations.

Particle and pebble-bed fuel models for HTGR calculations

High-temperature gas-cooled reactor geometries differ significantly from conventional light water reactors. The fissile material is encapsulated inside microscopic fuel particles (TRISO), randomly dispersed in fuel compacts or pebbles made of graphite. In pebble-bed type geometries the spherical fuel pebbles are piled inside the reactor core, which brings another level of random heterogeneity in the calculation.

The Serpent code has two geometry models for this type of fuels. The first (implicit) model samples the positions of fuel particles inside the dispersed media during the transport simulation. The routine works by preserving the particle dimensions and packing fractions (Leppänen, 2007b).

The second (explicit) model reads the positions of fuel particles or pebbles from a separate input file. The particle/pebble distribution is covered by a 3D mesh used by the search algorithm in the geometry routine. The explicit model works on several levels (particles inside a pebble and pebbles inside the core) and it has been tested in double-heterogeneous geometries consisting of more than 1,000,000 randomly positioned units. The routine makes no approximations and increases the calculation time only by 20-50% compared to a regular-lattice configuration. The routine also enables the calculation of pebble-wise power distributions, written in a separate output file.

To simplify the construction of HTGR geometries, the code has a separate command line routine that generates random particle distribution files for the explicit geometry model.

Interaction physics

Serpent reads continuous-energy cross sections from ACE format data libraries. The interaction physics is based on classical collision kinematics and ENDF reaction laws. Unresolved resonance cross sections can be sampled from probability tables. The most recent version has the option to apply Doppler-broadening rejection correction on the free-gas scattering kernel (see Becker, 2009 for the theory behind DBRC).

ACE format cross section libraries based on JEF-2.2, JEFF-3.1, JEFF-3.1.1, ENDF/B-VI.8 and ENDFB/B-VII evaluated data files are included in the Serpent installation package. Data is available for 432 nuclides at 6 temperatures between 300 and 1800K. Thermal bound-atom scattering data is included for light and heavy water and graphite. The data format is shared with the widely-used Los Alamos MCNP code and any continuous-energy MCNP data library can be used with Serpent as well.

Unionized energy grid format

The continuous-energy cross section data is reconstructed on a unionized energy grid, used for all reaction modes. The use of a single energy grid results in a major speed-up in calculation, as the number of time-consuming grid search iterations is reduced to minimum. If all the original grid points are preserved, there is no loss of data or accuracy in the process.

The drawback of the approach is that computer memory is wasted for storing redundant data points. The grid size may become prohibitively large in burnup calculation, often involving over 200 actinide and fission product nuclides. To overcome this problem the code uses two methods for reducing the memory demand (Leppänen, 2009b).

Built-in Doppler broadening routine

A built-in Doppler broadening preprocessor routine allows the conversion of ACE format cross sections into a higher temperature. This capability results in a more accurate description of the interaction physics in temperature-sensitive applications, as the data in the cross section libraries is available only in 300K intervals. The method has been validated with good results and the routine works efficiently without a major increase in the overall calculation time (Viitanen, 2009). The Doppler broadening capability was added in update 1.1.2.

Burnup calculation

Serpent is capable of running burnup calculation both as a completely stand-alone application and as a part of an externally-coupled sequence. The number of depletion zones is not restricted, although memory usage may become a limiting factor in large systems.

Fission products and actinide daughter nuclides are selected for the calculation without additional user effort and fuel rods and burnable absorber pins can be divided into annular sub-regions to account for the rim-effects. The irradiation history is defined in units of time or burnup. Reaction rates are normalized to total power, specific power density, flux or fission rate and the normalization can be changed by dividing the irradiation cycle into discrete intervals. Volumes and masses needed for the normalization are calculated automatically for simple geometries, such as fuel pin lattices, or the values can be defined by the user. Predictor-corrector calculation is optional and used by default for each burnup step.

Radioactive decay and fission yield data used in the calculation is read from standard ENDF format data libraries. The decay libraries may contain data for almost 4000 nuclides and meta-stable states, all of which is available for the calculation. The total number of different nuclides produced from fission, transmutation and decay reactions is generally lower, in the order of 1500. Energy-dependent fission yields are available for all main actinides (31 nuclides in ENDF/B-VII data).

Integral one-group transmutation cross sections are calculated either within the transport cycle, or after the cycle from a flux spectrum build on the unionized energy grid structure. The indirect method speeds up the calculation by a factor of 3-4, and due to the high energy resolution of the flux spectrum, the errors in the results are practically negligible. The main drawback is that all information on the statistical accuracy is lost. Similar methodology has been successfully used with coupled Monte Carlo burnup calculation codes (Haeck, 2007; Fridman, 2008).

Serpent has two fundamentally different options for solving the Bateman depletion equations. The first method is the Transmutation Trajectory Analysis (TTA) method (Cetnar, 2006), based on the analytical solution of linearized depletion chains. The second option is the Chebyshev Rational Approximation Method (CRAM), an advanced matrix exponential solution developed at VTT (Pusa, 2010). The two methods have shown to yield consistent results, both when used with Serpent (Leppänen, 2010a) and in separate methodological studies (Isotalo, 2010).

Parallelization

Serpent uses the Message Passing Interface (MPI) for parallel calculation. Parallelization is implemented by dividing the neutron histories to several tasks and combining the results after the transport cycle. This approach is simple and efficient, but it lacks error tolerance and dynamic load sharing. The overall calculation time is dependent on the slowest task, which is why the method is best applied in a symmetric parallel environment.

In addition to the neutronics simulation, parallelization in the burnup calculation mode divides also the preprocessing and depletion routines into several tasks. Each task handles a different irradiated material and the speed-up may become significant if the number of materials is large. The extensive memory demand due to the use of the unionized energy grid may become a limiting factor in parallel burnup calculation, if several tasks are sharing the same memory space.

Output

The Serpent code can be used for producing homogenized multi-group constants for deterministic reactor simulator calculations. The standard output includes:

All result estimates are accompanied by the associated relative statistical errors.

Homogenization can be carried out for multiple universes simultaneously, which allows the production of group constant data for several fuel assemblies within a single run.

The effective delayed neutron parameters are calculated using a forward Monte Carlo method developed at the NRG (Meulekamp, 2006). The few-energy group structure for group constant generation is arbitrary and defined by the user. The code also calculates fission source entropies for convergence studies. The total entropy is divided into x-, y- and z-components to monitor source convergence separately in each direction.

User-defined detectors (tallies) can be set up for calculating various integral reaction rates. The spatial integration domain can be defined by a combination of cells, universes, lattices and materials, or using a three-dimensional super-imposed mesh. The number and structure of detector energy bins is unrestricted. Various response functions are available for the calculation, including material total and isotopic reaction cross sections and ACE format dosimetry data.

Output for burnup calculation consists of isotopic compositions, transmutation cross sections, activities and decay heat data. The results are given both as material-wise and total values. Group constants and all the other output parameters are calculated and printed for each burnup step.

All numerical output is written in Matlab m-format files to simplify the simultaneous post-processing of several calculation cases. The code also has a geometry plotter feature and a reaction rate plotter, which is convenient for visualizing the neutronics in thermal systems (see the gallery for examples).

Beta features

There are various experimental features in the Serpent code that are still under development and not thoroughly tested or documented. Some of these features are briefly introduced in the following.

External source simulation

The development of an external source simulation mode began in version 1.1.11. The mode replaces the criticality source method and extends the simulation capabilities to non-multiplying systems. The source routines are still under development and the methodology requires extensive validation.

Equilibrium xenon calculation

Version 1.1.9 of the Serpent code includes a calculation routine that iterates the equilibrium concentration of fission product poison Xe-135 in one or several fuel regions. The calculation is based on the decay of the isotope and its precursors and the production and absorption rates normalized to total power or power density. The iteration can be used in both transport and burnup calculation modes. When used with burnup calculation, the concentration of Xe-135 is handled separate from the other nuclides.

Criticality iteration

In addition to the conventional k-eigenvalue criticality source calculation, Serpent has various options to iterate k_eff to a user-given value. Most commonly these methods are used to iterate criticality.

Coolant with soluble absorber can be defined as a mixture of two separate material compositions, most commoly water and boric acid. The absorber concentration can be treated as a free variable and used for k_eff iteration.

The second option is the alpha eigenvalue mode, which is equivalent to performing a time-dependent simulation (Cullen, 2003). This method is based on introducing a time absorption or multiplication reaction to balance the neutron generation and loss rates. The technique is commonly known and used in various deterministic and Monte Carlo transport calculation codes.

Experimental leakage models

Developing a Monte Carlo based leakage correction to account for the non-zero boundary currents in infinite lattice calculations is one of the long-term goals in Serpent development. The code has two built-in leakage models, which should be considered experimental, rather than standard Monte Carlo techniques since the results are generally not satisfactory when compared to deterministic transport codes. The first option is to use albedo boundary conditions on reflective or periodic infinite-lattice boundaries. The value of the albedo coefficient can be iterated to yield the desired k_eff. The second model is an attempt to perform B₁ fundamental mode calculation during the transport simulation. The routine was originally developed for the PSG code (see Sec. 9.5 in Leppänen, 2007c), but not transferred to Serpent until update 1.1.5. The method is very similar to the alpha eigenvalue simulation, based on leakage absorption or multiplication to balance the source and loss rates.

Critical spectrum calculation

As an intermediate solution to the leakage problem, Serpent version 1.1.14 has the option to homogenize cross sections by performing a separate critical spectrum calculation. The method follows a two-step procedure, in which homogenized micro-group cross sections from the Monte Carlo simulation are used to form B₁ equations, which are then solved outside the main transport loop. The procedure is very similar to the methods used in various deterministic codes, and preliminary comparisons to HELIOS (Fridman, 2011) have shown some very promising results. So far the methodology is used for group constant generation only, but extending the spectrum correction to burnup calculation is under development.

Results and validity

Serpent has been extensively validated in LWR lattice calculations (see validation reports for example results). Effective multiplication factors and homogenized few-group reaction cross sections are within the statistical accuracy from reference MCNP results, when the same ACE libraries are used in the calculations. Differences to other Monte Carlo codes (Keno-VI) are small, but statistically significant discrepancies can be observed in some cases. Differences to deterministic lattice codes are generally larger, mainly due to the fundamental differences between the calculation methods.

The burnup calculation capability is new and the routines have not been thoroughly tested. The most extensive LWR validation case carried out so far was presented at the M&C 2009 conference in May 2009 (Leppänen, 2009a). In that study, a comparison to deterministic CASMO-4E results shows consistent evolution of k-eff and the depletion of U-235 and gadolinium isotopes in a PWR assembly burnup calculation. Some small discrepancies are observed in the buildup rates of Pu-239 and the main fission product poisons. The differences are most likely to originate from several factors, such as the microscopic neutron capture cross section of U-238 and the special Xe-135 treatment used by the CASMO code.

Probability table sampling for unresolved resonance cross sections was added in update 1.1.4 and revised in update 1.1.5. Preliminary tests show that the method works as expected. What complicates the validation, however, is the fact that the treatment has only a minor impact in the results in most calculation cases. The sampling routine slows down the transport cycle to some extent, depending on the flux spectrum and the number of nuclides and materials with probability table data. The Serpent code has a few options for optimizing the calculation.

The current code version lacks any statistical tests performed on the calculated values. All Monte Carlo results are essentially random variables, assumed to be normally distributed. This assumption is based on the central limit theorem, and if the conditions of the theorem are not valid for a sequence of batch values, the confidence intervals of the normal distribution may not be applicable. This is not considered a major problem in lattice physics applications, and analyses performed on data produced in a series of test calculations has confirmed the validity of the assumptions.

Performance

Comparison to a similar calculation using Monteburns suggests that Serpent may run over 80 times faster than coupled Monte Carlo burnup calculation codes based on MCNP.  It should be noted, however, that the difference does not result so much from the efficiency of the code, but rather from the fact that MCNP tends to slow down significantly as the number of reaction rate tallies becomes large. What is more important is that the Serpent code can run full-scale assembly burnup calculations similar to deterministic transport codes, and the overall calculation time is counted in hours or days, rather than weeks or months.

NOTE: Some of the information in this website is outdated -- see the Serpent discussion forum (link below) for latest news

Serpent Discussion Forum

Serpent User's Manual
(March 5, 2012)

Base version:

Serpent 1.1.7

Current update:

1.1.17 (March 5, 2012)

Important updates: