Delta- and surface-tracking

This a brief description on the delta-tracking based transport routine used in Serpent. The method is also used by other Monte Carlo codes, most notably in the HOLE geometry package in MONK and MCBEND. The original delta-tracking algorithm was introduced by Woodcock in 1965, and a mathematical verification was derived by Coleman in 1968. Delta-tracking is well described in a text book by Lux and Koblinger, and a description of the methodology used in Serpent is found in an article in Annals of Nuclear Energy from 2010.

The input parameters related to delta-tracking are:

• set dt - sets the probability threshold used for selecting between surface- and delta-tracking
• set forcedt - enforces the use of delta-tracking in a given list of materials
• set blockdt - enforces the use of surface-tracking in a given list of materials
• set cfe - definse the mean-free-path of collisions used to score the collision flux estimator

The output parameters are:

• TODO

Transport algorithms in Monte Carlo simulation

The Monte Carlo simulation consists of a large number particle histories, in which the random walk of an individual particle is followed, or tracked, through the geometry from its birth to eventual absorption or escape. Tracking is most typically carried out in a constructive solid geometry (CSG), composed of homogeneous material cells, which are defined by combinations of elementary and derived surface types. Serpent 2 also has two advanced geometry types, based on STL format CAD models and an unstructured polyhedral mesh. 

The transport simulation follows a random walk from one interaction to the next. The procedure can be described as follows:

1. Sample path length (distance to next collision)
2. Transport particle to the collision point
3. Sample interaction

If the sampled interaction is scattering, the procedure restarts from beginning by sampling the distance to the next collision. The direction and energy are changed in each scattering event.

By definition, the interaction probability per traveled path length is given by the macroscopic total cross section (denoted here as $\Sigma$). If it is assumed that the particle travels through an infinite homogeneous medium characterized by constant total cross section, it can be shown that the free path length follows an exponential distribution. This distribution can be sampled using the inverse method, and the distance to the next collision site is given by:

$l = -\frac{1}{\Sigma}\log\xi$

where $\xi$ is a uniformly distributed random variable on the unit interval.

The prerequisite of using this simple formula for sampling the distance to the next collision site is that the material is infinite and homogeneous. This is the case when the particle remains within a single material region, but not if its track crosses the boundary between two materials.

Surface-tracking

Since the interaction probability per traveled path length does not depend on the events that occured in the particle's history, any point in the path can be considered the stating point of a new path. If the sampled path crosses the boundary between two materials, the track can be stopped at the crossing point and a new path length sampled using the cross section of the next material. This is the general idea in the surface-tracking algorithm.

In order to stop the particle at the boundary between two materials, the algorithm needs to calculate the distance to the nearest boundary surface. The only way to accomplish this is to loop over all candidate surfaces and pick the shortest value.

Surface tracking is considered the standard tracking algorithm and it is used by virtually every Monte Carlo transport calculation code. The method has a few drawbacks related to its efficiency in complex geometries:

• Determining the distance to the nearest boundary can become computationally expensive if the cells are comprised of a large number of surfaces
• The fact that the particle track has to be stopped at each boundary crossing becomes a computational bottleneck when the mean-free-path is long compared to geometry dimensions

Delta-tracking

An alternative to surface-tracking is the Woodcock delta-tracking algorithm, which is based on the rejection sampling of particle path lengths. The procedure relies on the concept of a virtual collision, which is a fictive interaction that preserves the energy and direction of the transported particle.

Since virtual collisions do not change the statistics of the random walk in any way, the material total cross section $\Sigma$ can be adjusted with an arbitrary virtual collision cross section $\Sigma_0$:

$\Sigma'(\mathbf{r}, E) = \Sigma(\mathbf{r},E) + \Sigma_0(\mathbf{r},E)$

without changing the outcome of the simulation. It is therefore possible to adjust the cross sections of all material regions 1, 2, 3, ... in the system such that:

$\Sigma_1'(E) = \Sigma_2'(E) = \Sigma_3'(E) \dots = \Sigma_\mathrm{m}(E)$

where $\Sigma_\mathrm{m}$ is called the majorant cross section.

In practice, it is not necessary to define the virtual collision cross sections at all if the majorant is simply taken as the maximum of all material totals at each energy point:

$\Sigma_\mathrm{m}(E) = \max{\big[\Sigma(\mathbf{r},E)\big]}$

Unlike the physical total cross section, which depends on the material located at the particle position, the majorant cross section is completely independent of the spatial coordinates. The point of having a macroscopic cross section that is uniform throughout the geometry is that when it is used for sampling path lengths:

$l = -\frac{1}{\Sigma}_\mathrm{m}\log\xi$

the values are statistically valid regardless of the number of material boundaries crossed. At the end point of the sampled path the tracking routine performs rejection sampling. The probability to accept the collision is given by ratio of the physical total cross section to the majorant:

$P = \frac{\Sigma(\mathbf{r},E)}{\Sigma_\mathrm{m}(E)}$

If the collision is rejected, a new path length is sampled using $\Sigma_\mathrm{m}$ and the particle is moved to the next tentative collision site.

Since the majorant cross section is always larger than or equal to the total cross section, the path lengths sampled using delta-tracking are, on the average, shorter than those sampled with surface-tracking. The average physical distance between two collisions is preserved, as some paths are extended over multiple virtual collisions.

Advantages and limitations of delta-tracking

The advantage of delta-tracking over the surface-tracking algorithm is that there is no need to calculate the surface distances or stop the particle at the material boundaries. This becomes significant for computational performance in geometries where the particle mean-free-path is long compared to dimensions. The difference between surface- and delta-tracking is emphasized in HTGR particle fuels and advanced CAD and an mesh based geometry types.

Since the majorant cross section used for sampling the path lengths does not depend on the spatial coordinates and the material total is needed only at discrete locations, variations of delta-tracking can be used for transporting particles through inhomogeneous materials. This capability is utilized in the TMS on-the-fly temperature treatment routine and the multi-physics interface in Serpent 2.

Delta-tracking also has its drawbacks. Since the majorant cross section reflects the largest interaction probability within the system, the efficiency of the rejection sampling loop may become poor in the presence of localized heavy absorbers (control rods, burnable absorber pins, etc.) that dominate the majorant cross section, but occupy a relatively small volume of the geometry (see example below).

Another drawback is that delta-tracking rules out the use of the track-length estimate (TLE) of particle flux, and integral reaction rate estimates need to be calculated using the potentially less efficient collision estimator (CFE).

• 
• 

Example of the effect of localized heavy absorber on majorant cross section and rejection sampling efficiency. Left: Majorant and macroscopic cross sections in a system with localized heavy absorber (Gd-fuel pins in BWR assembly). The majorant is dominated by the high capture cross sections of ¹⁵⁵Gd and ¹⁵⁷ Gd, even though the burnable absorber pins occupy a relatively small volume of the geometry. Right: Rejection probability in coolant and moderator where neutrons spend most of their lifetime. The efficiency of the rejection sampling scheme becomes poor especially at low energy.