Overview

Serpent is a multi-purpose three-dimensional continuous-energy Monte Carlo particle transport code, developed at VTT Technical Research Centre of Finland, Ltd. The development started in 2004, and the code has been publicly distributed by the OECD/NEA Data Bank and RSICC since 2009. Serpent started out as a simplified reactor physics code, but the capabilities of the current development version, Serpent 2, extend well beyond reactor modeling. The applications can be roughly divided into three categories:

1) Traditional reactor physics applications, including spatial homogenization, criticality calculations, fuel cycle studies, research reactor modeling, validation of deterministic transport codes, etc.
2) Multi-physics simulations, i.e. coupled calculations with thermal hydraulics, CFD and fuel performance codes
3) Neutron and photon transport simulations for radiation dose rate calculations, shielding, fusion research and medical physics

The main features and capabilities of the code are introduced below. Support for users is provided at the Serpent Discussion Forum, and the Serpent Wiki acts as an on-line user manual.

The recommended publication for referencing Serpent is:

Leppänen, J., et al. (2015) "The Serpent Monte Carlo code: Status, development and applications in 2013." Ann. Nucl. Energy, 82 (2015) 142-150.

although it should be noted that this review article does not cover some of the latest features, including photon transport mode and advanced geometry types.

NOTE: Please report any scientific publications where you have extensively used Serpent to the developer team, so we can include them on the list of publication. This is important especially for theses, which are otherwise difficult to keep track of.

Geometry and particle tracking

Similar to other Monte Carlo codes the basic geometry description in Serpent relies on a universe-based constructive solid geometry (CSG) model, which allows the description of practically any two- or three-dimensional fuel or reactor configuration. The CSG geometry consists of homogeneous material cells, defined by elementary and derived surface types that are combined using Boolean operators (intersections, unions and complements). Serpent supports conventional square and hexagonal lattices, and provides special geometry types for CANDU and randomly-dispersed particle fuels. In addition to CSG-type universes Serpent has the option to import CAD and unstructured mesh based geometries.

Woodcock delta-tracking method

Particle transport in Serpent is based on the combination of conventional surface-tracking and the Woodcock delta-tracking method (Woodcock, 1965). The tracking routine has proven efficient and well suited for geometries where the neutron or photon mean-free-path is long compared to the dimensions. This is typically the case in reactor physics calculations involving fuel assemblies and especially HTGR micro-particle fuels. Complex geometries with highly-refined spatial detail are also encountered in various fusion applications. The traditional delta-tracking method is subject to certain efficiency problems related to localized heavy absorbers, which in Serpent are avoided by switching to surface-tracking when necessary (Leppänen, 2010b).

The main drawback of delta-tracking is that the track-length estimate of particle flux is not available, and reaction rates have to be calculated using the potentially less-efficient collision estimator. This is usually not a problem at all in reactor calculations when reaction rates are scored in regions of high collision density. When the collision rate is low, the efficiency of the estimator can be improved by introducing additional virtual collisions over the particle flight path. Serpent also provides a special detector type based on the track-length estimator for calculating reaction rates in small or optically thin volumes, in which the efficiency of the collision estimator is often poor.

Explicit particle and pebble-bed fuel model for HTGR calculations

High-temperature gas-cooled reactor geometries differ significantly from conventional light water reactors. The fissile material is encapsulated inside microscopic tristructural-isotropic (TRISO) fuel particles, 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 explicit HTGR geometry model in Serpent reads the coordinates of fuel particles or pebbles from a separate input file, and generates the geometry as it is defined, without any approximations. The model works on several levels (particles inside a pebble and pebbles inside the core) and it has been tested in realistic double-heterogeneous reactor configurations consisting of over 60 million randomly positioned units (Suikkanen, 2010; Rintala, 2015) The computational overhead from handling the unstructured configuration is moderate compared to a similar regular-lattice model. The routine also enables the calculation of pebble-wise power distributions over the reactor core without defining additional tallies.

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

CAD and unstructured mesh-based geometry types

Serpent has the capability to import CAD and unstructured mesh based geometries as part of the universe structure. CAD models are read in the stereolitography (STL) format, in which the surfaces of geometry bodies are represented by a mesh of flat triangles. The STL file format is widely used for 3D printing, and therefore supported by most CAD tools. The geometry type was introduced in Serpent 2 in 2014 (Leppänen, 2014f), and it has been used for modeling complicated structures in fusion (Leppänen, 2015a; 2016a; Sirén, 2016) and fission (Talamo, 2015; 2016a; 2016b) reactor geometries.

The unstructured mesh based geometry type (Leppänen, 2014d) was developed as a by-product of the multi-physics interface used for coupling Serpent 2 to OpenFOAM CFD calculations (see below). The geometry consists of a three-dimensional tetra-, hexa- or polyhedral mesh read in the standard OpenFOAM mesh file format.

Interaction physics

Serpent reads continuous-energy cross sections from ACE format data libraries. The interaction physics is based on classical collision kinematics, ENDF reaction laws and probability table sampling in the unresolved resonance region. Improved treatment for the free-gas scattering kernel near resonances is also available, based on the DBRC Doppler-broadening rejection correction method (Becker, 2009).

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 installation package of Serpent 1. Interaction 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. Since the data format is shared with MCNP, any continuous-energy ACE format data library generated for MCNP can be used with Serpent as well. The data format determines the "laws of physics" for neutron interactions, and the results from Serpent calculations can be expected to agree with MCNP to within statistics.

Unionized energy grid format

Continuous-energy cross sections read from the library files are reconstructed on a unionized energy grid, used for all reaction modes (Leppänen, 2009b). The use of a single energy grid results in a major speed-up in calculation, as the number of CPU time consuming grid search iterations is reduced to minimum. Macroscopic cross sections for each material are pre-generated before the transport simulation. Instead of calculating the cross sections by summing over the constituent nuclides during tracking, the values are read from pre-generated tables, which is another effective means to improve the performance.

The drawback of the unionized energy grid approach is that computer memory is wasted for storing redundant data points. The grid size may become prohibitively large in burnup calculations, often involving more than 250 actinide and fission product nuclides. To overcome this issue, Serpent 2 provides different optimization modes for small and large burnup calculation problems, in which the the unionized energy grid approach is used selectively (Leppänen, 2012a). The lowest optimization modes allow running large burnup calculation problems with tens or hundreds of thousands of depletion zones, while the higher modes provide considerable speed-up in assembly-level calculations.

Doppler-broadening of cross sections

A built-in Doppler-broadening preprocessor routine (Viitanen, 2009) allows adjusting the temperatures of ACE format cross sections. 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 significant computational overhead.

Another option for adjusting the material temperatures is the target motion sampling (TMS) on-the-fly temperature treatment routine (Viitanen, 2012a; 2012b; 2013a; 2014a; 2015a; 2015b; 2015c). Instead of averaging the cross sections over the Maxwellian distribution (actual Doppler-broadening), the TMS method accounts for thermal motion explicitly, by making a coordinate transformation in the target-at-rest frame before handling the collision physics. The method was developed especially for the purpose of multi-physics coupling (see below), where it can be used for modeling a wide variation of material temperatures or even continous temperature distributions.

The adjustment of neutron thermal scattering cross section and energy-angle distributions in S(α,β) format is based on interpolation between tabular data. This can be performed in a pre-processing stage or on-the-fly during the transport simulation (Viitanen, 2016). The most significant limitation in the temperature treatment routines is that the adjustments cannot be currently applied to unresolved resonance probability table sampling.

Photon transport mode

Photon physics routines were implemented in Serpent 2 in 2015 (Kaltiaisenaho, 2016). The physics model currently covers the basic interactions (Rayleigh and Compton scattering, photoelectric effect and electron-positron pair production) for photon energies ranging from 1 keV to 100 MeV. Secondary photons are produced by atomic relaxation and bremsstrahlung, handled using the thick-target bremsstrahlung (TTB) approximation. The physics model is comparable to the methods used in other Monte Carlo transport codes (e.g., MCNP6, PENELOPE, Geant4, EGS5, EGSnrc, FLUKA). In addition to the standard ACE format cross section libraries Serpent reads photon interaction data from supplementary data files, which is why the physics model is not fully compatible with that used in MCNP.

The source distribution for photon transport simulations can be obtained from a radioactive decay source. In this source mode Serpent combines the compositions of activated materials into photon emission spectra read from ENDF format radioactive decay data files. Source generation can be combined with a burnup or activation calculation performed using built-in automated calculation routines (see below). The methodology has been applied, for example, to spent fuel transport cask calculations and the the evaluation of shut-down dose rates following a plasma shot in the ITER fusion reactor (Leppänen, 2016a).

The original incentive for developing a photon transport mode was to account for gamma heating in coupled multi-physics simulations. Modeling of accurate heat deposition in coolant and structural materials requires accounting for the direct contribution of neutrons and fission and capture gammas, which in turn requires a coupled neutron-photon transport mode. Such calculation mode is currently under development. The implementation of photon physics routines has also allowed broadening the scope of Serpent applications from traditional reactor physics calculations to radiation transport and shielding.

Burnup calculation

The burnup calculation capability in Serpent was established early on, and is entirely based on built-in calculation routines, without coupling to any external solvers. The number of depletion zones is not restricted, although memory usage may require reducing the optimization when the number of burnable materials is large.

Fission and activation products and actinide daughter nuclides are selected for the calculation without additional user effort, and burnable materials can be sub-divided into depletion zones automatically. The irradiation history is defined in units of time or burnup. Reaction rates are normalized to total power, specific power density, flux, fission or source rate, and the normalization can be changed by dividing the irradiation cycle into a number of separate depletion intervals. A restart features allows performing fuel shuffling or applying any modifications in the input by dividing the calculation into several parts. Volumes and masses needed for the normalization are calculated automatically for simple geometries, such as 2D fuel pin lattices. The values can also be obtained from a Monte Carlo based volume calculation routine or entered manually.

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. The concentrations of all included nuclides with decay data are tracked in the burnup calculation, and the number of nuclides with cross sections typically ranges from 200 to 300. Energy-dependent fission yields are available for all main actinides (31 nuclides in ENDF/B-VII data). Isomeric branching ratios for neutron reactions are not included in the ACE format data libraries. Serpent uses fixed ratios for the most important nuclides (e.g Am-241 and Pm-147) by default, and has the option to read energy-dependent data from ENDF format files.

Flux-volume-averaged one-group transmutation cross sections are calculated either during the transport simulation, or by collapsing the continuous-energy reaction cross sections after the calculation has been completed using a flux spectrum collected on the unionized energy grid. The spectrum collapse 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. Similar methodology has been used with other coupled Monte Carlo burnup calculation codes (Haeck, 2007; Fridman, 2008a; 2008b).

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 for Serpent at VTT (Pusa, 2010; 2011; 2012;. 2013a;. 2013b;. 2013c;. 2013d;. 2014). The two methods have shown to yield consistent results, both when used with Serpent (Leppänen, 2009a) and in separate methodological studies (Isotalo, 2011a).

Burnup algorithms include the conventional explicit Euler and predictor-corrector methods, but Serpent 2 also offers various higher-order methods and sub-step solutions for burnup calculation (Isotalo, 2011b; 2011c; 2013b; 2015a; 2015b). The stability of 3D burnup calculations can be improved by implicit algorithms (Dufek, 2014).

Fission product poisons Xe-135 and Sm-149 can be handled separately from the other nuclides, and iterated to their equilibrium concentration during the transport simulation. The equilibrium calculation is independent of the depletion routine, and the iteration can also be performed in transport mode without burnup calculation.

Coupled multi-physics simulations

Two-way coupling to thermal hydraulics, CFD and fuel performance codes has been a major topic in Serpent development for the past several years. The multi-physics coupling scheme in Serpent 2 is designed to operate on two levels:

1) Internal coupling to built-in solvers for fuel behavior and thermal hydraulics
2) External coupling via a universal multi-physics interface

The built-in solvers are integrated to the transport simulation at source code level, and designed to provide solutions to coupled problems at a relatively low computational cost. The solvers include FINIX (Ikonen, 2013a; 2013b; 2015; 2016; Valtavirta, 2014b) -- a thermo-mechanical fuel behavior module for the modeling of temperature feedback inside fuel pins in steady-state and transient conditions, and COSY -- a three-dimensional system/component scale thermal hydraulics solver based on a porous medium three-field flow model. The FINIX solver is available by request, but the development of COSY is still under way.

The multi-physics interface is designed for coupling Serpent to external solvers. The interface has several structured mesh formats for coupling to thermal hydraulics codes, and additional interface types are available for fuel performance (Valtavirta, 2013b) and CFD (Leppänen, 2014e) code coupling. The CFD interface is based on an unstructured tetra-, hexa- or polyhedral mesh read in the standard OpenFOAM mesh file format. The purpose of the multi-physics interface is essentially to separate the state variables from the geometry description, which allows handing all data flow between the coupled codes without modifications in the main input files. The methodology relies heavily on the capability to model continuously-varying density distributions (Leppänen, 2013b) and the on-the-fly temperature treatment routines (see description above).

Work on coupled multi-physics applications continues. Coupled calculations have been carried out in steady-state, transient and burnup modes. The transient capability in Serpent allows the modeling of both prompt super-critical reactivity excursions and slow transients below prompt criticality. A delayed neutron model (Valtavirta, 2016) allows the tracking of precursor concentrations over long time periods.

In addition to fission reactor applications, Serpent has also been coupled to plasma scenario simulations to provide a realistic source distribution for fusion neutronics calculations (Sirén, 2016). This work is an essential part of expanding the use of Serpent to fusion research.

Variance reduction

When Serpent started out as a reactor physics code, obtaining sufficient statistics for the results was just a matter of running a sufficient number of neutron histories. The uniform fission site method was later implemented to improve the statistical accuracy in full-core calculations, in which the outermost fuel pins in assemblies located at the core-reflector boundary typically receive a low number of scores. Similar methodology is also used in the MC21 code (Kelly, 2012).

The implementation of more efficient general-purpose variance reduction techniques was started fairly recently, along with the effort to broaden the scope of applications to radiation transport and shielding calculations. The methodology relies on a conventional super-imposed weight-window mesh. The importances used for obtaining the weight-window boundaries can be produced by state-of-the-art calculation tools, such as ADVANTG or MAVRIC, or using a built-in light-weight solver based on the response-matrix method. Serpent can read standard MCNP WWINP format files, although the methodology is still under development and subject to several limitations.

Parallelization

Serpent can be run in parallel in computer clusters and multi-core workstations. Parallelization at core level is handled by thread-based OpenMP, which has the advantage that all CPU cores within the computational node are accessing the same memory space. Calculations can be divided into several nodes by distributed-memory MPI parallelization.

In addition to the particle transport simulation, parallelization in the burnup calculation mode divides also the preprocessing and depletion routines between several CPU's.

Results and output

Spatial homogenization was the main intended application for Serpent when the project was started in 2004. The group constant generation capability in Serpent 2 currently covers (Leppänen, 2016b):

Homogenized few-group reaction cross sections
Scattering and scattering production matrices
Transport cross sections and diffusion coefficients calculated using several methods
Discontinuity factors
Form factors for pin-power reconstruction
Albedos and partial albedos
Poison cross sections for Xe-135 and Sm-149 and their precursors
Effective delayed neutron fractions
Transmutation cross sections for micro-depletion

Homogenization can be performed in infinite spectrum, or using a leakage correction based on the deterministic solution of the B1 equations. Serpent also has a built-in homogeneous diffusion flux solver for calculating discontinuity factors in regions where the net current over the boundaries is not reduced to zero by reflective boundary conditions. This is the case, for example, in reflectors and assembly colorset configurations. The calculation of homogenized group constants is fully automated, and Serpent provides an automated burnup sequence capable of performing branch calculations for state-point variation.

User-defined detectors (tallies) can be set up for calculating various integral reaction rates for neutrons and photons. 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 results can be divided into an arbitrary number of energy and time bins. The standard tally types are based on the collision estimate of particle flux, and special detectors are available for calculating surface currents and reaction rates inside simple surface types using the track-length estimator. Photon tallies also include a pulse-height detector.

Various response functions are available for the calculation of integral reaction rates, including material-wise macroscopic and isotopic microscopic cross sections, ACE format dosimetry data and user-defined functions. Built-in mass-energy attenuation coefficients are available for calculating photon dose rates. Energy deposition and radiation dose can also be evaluated using analog detector types, which in some case provide more physical results compared to the use of flux-to-dose conversion factors.

Serpent calculates adjoint-weighted point kinetics parameters and effective delayed neutron fractions using the iterated fission probability (IFP) method (Leppänen, 2014b), relying on an implementation similar to that developed for the MCNP code (Kiedrowski, 2011).

Output for burnup calculation consists of isotopic compositions, activities, spontaneous fission rates, decay heat and radiotoxicity 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 post-processing of the results. The code also has a geometry plotter feature and a reaction rate plotter, which is convenient for visualizing the neutronics and tally results (see the gallery for examples).

Validation

Each Serpent update is checked by comparison to MCNP by running a standard set of assembly calculation problems. Effective multiplication factors and homogenized few-group reaction cross sections are within the statistical accuracy from the reference results, when the same ACE libraries are used in the calculations. Validation against MCNP has also been carried out with equally good results for calculations involving individual nuclides, by comparing the flux spectra produced by the two codes. 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.

Validation of burnup calculation routines is considerably more difficult, due to the lack of a perfect reference code. In addition to discrepancies in the transport simulation, there are additional factors related to decay, fission yield and isomeric branching ratio data, formulation of transmutation and decay chains, depletion algorithms, and so on. Comparison of effective multiplication factors and other integral parameters shows generally good agreement between different calculation codes, but significant discrepancies can be found in the concentrations of individual nuclides. Burnup calculations are extremely sensitive to differences in the fundamental physics data, and changing one decay library to another can alone result in orders of magnitude discrepancies for some radionuclides. The lack of high-quality experimental reference data is a common problem in the validation of all burnup calculation codes.

Systematic validation for criticality safety analyses using experimental configurations and data from the International Handbook of Evaluated Criticality Safety Benchmark Experiments is currently under way. Similar validation project is planned for radiation shielding calculations.

Serpent-generated group constants have been used as the input data for various nodal diffusion codes, such as PARCS, DYN3D, and several in-house codes developed at VTT. The best way to validate the Monte Carlo based code sequence is to compare the results to reference Serpent full-core calculations. This approach has been applied, for example, for the Serpent-ARES code sequence, using the MIT BEAVRS Benchmark as the test case. The study involved initial core HZP and HFP calculations and fuel cycle simulations, with comparisons to 3D Monte Carlo calculations and experimental results (Leppänen, 2014b;. 2016b). The calculations also demonstrated the practical feasibility of using the continuous-energy Monte Carlo method for producing the full set of group constants for full-scale PWR fuel cycle simulations.

 
Updates on website
February 15, 2017
List of publications updated.
February 6, 2017
New code version, website content brought up to date.
November 8, 2016
New publications.
October 27, 2016
New code version, website for the 6th International Serpent UGM set up


Serpent Discussion Forum




Serpent Wiki



Serpent 1 User's Manual
(March 6, 2013)

NOTE: This manual covers only the basic features in Serpent 1. The capabilities of Serpent 2 are described at the on-line Wiki and the discussion forum.


Serpent 1 base version:

1.1.7 (NEA / RSICC)

Current update:

1.1.19 (April 2, 2013)

NOTE: Serpent 1 is no longer actively maintained, and all users are strongly encouraged to switch to Serpent 2 instead.


Serpent 2 base version:

2.1.0 (available by request)

Current update:

2.1.28 (Feb. 6, 2017)


Some recent and upcoming events
October 10-14, 2016

Serpent workshop at the 26th Symposium of AER on VVER Reactor Physics and Reactor Safety in Helsinki, Finland.
September 26-29, 2016

6th International Serpent User Group Meeting in Milan, Italy.
May 1, 2016

Serpent multi-physics workshop at the PHYSOR 2016 conference in Sun Valley, Idaho, May 1-5, 2016. Presentations available for download.
December 7, 2015

Presentation on the past, present and future challenges of developing the Serpent code in an ANS Seminar at MIT
October 21, 2015

Serpent workshop at the 7th International Conference on Modeling and Simulation in Nuclear Science and Engineering in Ottawa, ON, Canada
October 13-16, 2015

5th International Serpent User Group meeting in Knoxville, TN, USA
June 11-12, 2015

Serpent and fusion neutronics workshop at the University of Cambridge, UK
May 8, 2015

Tuomas Viitanen defended his Doctoral Thesis: Development of a stochastic temperature treatment technique for Monte Carlo neutron tracking at Aalto University.
February 26-27, 2015

Serpent and multi-physics workshop at LPSC - Grenoble, France
September 28 - October 3, 2014

Serpent workshop at the PHYSOR-2014 conference in Kyoto, Japan
September 17-19, 2014

4th International Serpent User Group meeting in Cambridge, UK
November 6-8, 2013

The Third International Serpent User Group Meeting in Berkeley, California, USA, organized by the University of California, Berkeley. (see meeting website)
October 27-31, 2013

Serpent contribution in a Monte Carlo codes invited session and several Serpent-related presentations at the Joint International Conference on Supercomputing in Nuclear Applications + Monte Carlo 2013 (SNA+MC 2013), Paris, France.
May 24, 2013

Maria Pusa defended her Doctoral Thesis on Numerical methods for nuclear fuel burnup calculations at Aalto University.
September 19-21, 2012

The Second International Serpent User Group Meeting in Madrid, Spain, organized by the Universidad Politécnica de Madrid (see meeting website)
April 30 - May 2, 2012

Serpent workshop at UC Berkeley, USA
April 15-20, 2012

Several Serpent-related papers at the PHYSOR-2012 conference in Knoxville, TN, USA
January 31, 2012

Beta-testing phase of Serpent 2 started
September 15-16, 2011

2011 Serpent International User Group Meeting, Dresden, Germany (also see the topic at the discussion forum and the meeting website)
February 7-8, 2011

Serpent presentations at the Workshop on Recent Developments and Advanced Applications in the Monte Carlo Method, UNIST, Ulsan, Korea
October 4, 2010

Presentation on Serpent development in an ANS Seminar at MIT
March, 2010

Serpent 1.1.7 available at RSICC  (Code Number C00757)
January 15, 2010

Serpent cross section libraries released as a separate NEA package (Package-ID NEA-1854)
January 6, 2010

NEA Base version upgraded to 1.1.7 (Package-ID NEA-1840)
October 13, 2009

M.Sc. Thesis: "Implementing a Doppler-Preprocessor of Cross Sections in Reactor Physics Code Serpent" completed at Helsinki University of Technology (Viitanen, 2009)
May 26, 2009

Serpent 1.1.0 available at the OECD  / NEA Data Bank (Package-ID NEA-1840)
April 8, 2009

Serpent 1.1.0 submitted to the OECD / NEA Data Bank for public distribution