Difference between revisions of "Delta- and surface-tracking"

From Serpent Wiki
Jump to: navigation, search
(Advantages and limitations of delta-tracking)
m
(73 intermediate revisions by one other user not shown)
Line 6: Line 6:
 
* [[input syntax manual#set forcedt|set forcedt]] - enforces the use of delta-tracking in a given list of materials
 
* [[input syntax manual#set forcedt|set forcedt]] - enforces the use of delta-tracking in a given list of materials
 
* [[input syntax manual#set blockdt|set blockdt]] - enforces the use of surface-tracking in a given list of materials
 
* [[input syntax manual#set blockdt|set blockdt]] - enforces the use of surface-tracking in a given list of materials
* [[input syntax manual#set minxs|set minxs]] - definse the mean-free-path of collisions used to score the collision flux estimator
+
* [[input syntax manual#set cfe|set cfe]] - definse the mean-free-path of collisions used to score the collision flux estimator
  
 
The output parameters are:
 
The output parameters are:
Line 12: Line 12:
  
  
== Transport algorithm in Monte Carlo simulation ==
+
== Transport algorithms in Monte Carlo simulation ==
  
 
The Monte Carlo simulation consists of a large number particle histories, in which the random
 
The Monte Carlo simulation consists of a large number particle histories, in which the random
Line 21: Line 21:
 
The transport simulation follows a random walk from one interaction to the next. The procedure can be described as follows:
 
The transport simulation follows a random walk from one interaction to the next. The procedure can be described as follows:
 
#Sample path length (distance to next collision)
 
#Sample path length (distance to next collision)
#Transport neutron to the collision point
+
#Transport particle to the collision point
 
#Sample interaction
 
#Sample interaction
 
If the sampled interaction is scattering, the procedure restarts from beginning by sampling the
 
If the sampled interaction is scattering, the procedure restarts from beginning by sampling the
Line 29: Line 29:
 
cross section (denoted here as <math>\Sigma</math>).  If it is assumed that the particle travels through an infinite
 
cross section (denoted here as <math>\Sigma</math>).  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:
 
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:
 +
  
 
<math>
 
<math>
 
l = -\frac{1}{\Sigma}\log\xi
 
l = -\frac{1}{\Sigma}\log\xi
 
</math>
 
</math>
 +
  
 
where <math>\xi</math> is a uniformly distributed random variable on the unit interval.
 
where <math>\xi</math> is a uniformly distributed random variable on the unit interval.
Line 39: Line 41:
 
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.
 
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- and delta-tracking ==
+
=== Surface-tracking ===
  
If it is known that the particle makes it to the next material boundary and interacts somewhere
+
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.
beyond the other side, the point of crossing can be taken as the starting point of a new path. This
+
is the general idea in the surface-tracking algorithm, in which the particle track is stopped at each
+
boundary crossing, and a new path sampled using the cross section of the next material.
+
  
The algorithm requires calculating the distance to the nearest boundary surface. The only way
+
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.  
 
to accomplish this is to loop over all candidate surfaces and pick the shortest value.  
  
Line 53: Line 52:
 
complex geometries:
 
complex geometries:
 
*Determining the distance to the nearest boundary can become computationally expensive if the cells are comprised of a large number of surfaces
 
*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 the dimensions
+
*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
 
An alternative to surface-tracking is the Woodcock delta-tracking
Line 61: Line 62:
 
of the transported particle.
 
of the transported particle.
  
Since virtual collisions do not change the  
+
Since virtual collisions do not change the statistics of the
 
random walk in any way, the material total cross section <math>\Sigma</math>
 
random walk in any way, the material total cross section <math>\Sigma</math>
 
can be adjusted with an arbitrary virtual collision cross section <math>\Sigma_0</math>:
 
can be adjusted with an arbitrary virtual collision cross section <math>\Sigma_0</math>:
 +
  
 
<math>
 
<math>
\Sigma'(\mathbf{r}, E) = \Sigma(\mathbg{r},E) + \Sigma_0(\mathbf{r},E)
+
\Sigma'(\mathbf{r}, E) = \Sigma(\mathbf{r},E) + \Sigma_0(\mathbf{r},E)
 
</math>
 
</math>
 +
  
 
without changing the outcome of the simulation.
 
without changing the outcome of the simulation.
It is then possible to adjust the cross sections of all  
+
It is therefore possible to adjust the cross sections of all  
material regions 1, 2, ... in the
+
material regions 1, 2, 3, ... in the
 
system such that:
 
system such that:
 +
  
 
<math>
 
<math>
 
\Sigma_1'(E) = \Sigma_2'(E) = \Sigma_3'(E) \dots = \Sigma_\mathrm{m}(E)
 
\Sigma_1'(E) = \Sigma_2'(E) = \Sigma_3'(E) \dots = \Sigma_\mathrm{m}(E)
\</math>
+
</math>
 +
 
  
 
where <math>\Sigma_\mathrm{m}</math> is called the ''majorant'' cross section.  
 
where <math>\Sigma_\mathrm{m}</math> is called the ''majorant'' cross section.  
Line 83: Line 88:
 
to define the virtual collision cross sections at all if the majorant is simply
 
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:
 
taken as the maximum of all material totals at each energy point:
 +
  
 
<math>
 
<math>
\Sigma_\mathr{m}(E) = \max{\big[\Sigma(\mathbf{r},E)\big]}
+
\Sigma_\mathrm{m}(E) = \max{\big[\Sigma(\mathbf{r},E)\big]}
 
</math>
 
</math>
 +
  
 
Unlike the physical total cross section, which depends on the material  
 
Unlike the physical total cross section, which depends on the material  
 
located at the particle position, the majorant cross section is completely
 
located at the particle position, the majorant cross section is completely
 
independent of the spatial coordinates.
 
independent of the spatial coordinates.
 
 
The point of having a macroscopic cross section that is uniform
 
The point of having a macroscopic cross section that is uniform
throughout the geometry is that when used for sampling path lengths:
+
throughout the geometry is that when it is used for sampling path lengths:
  
 
<math>
 
<math>
l = -\log(\xi)/\Sigma_\mathrm{m}
+
l = -\frac{1}{\Sigma}_\mathrm{m}\log\xi
 
</math>
 
</math>
 +
  
 
the values are statistically valid regardless of the number of material
 
the values are statistically valid regardless of the number of material
 
boundaries crossed.  
 
boundaries crossed.  
 
 
At the end point of the sampled path the tracking routine
 
At the end point of the sampled path the tracking routine
 
performs rejection sampling. The probability to accept the collision is given
 
performs rejection sampling. The probability to accept the collision is given
 
by ratio of the physical total cross section to the majorant:
 
by ratio of the physical total cross section to the majorant:
 +
  
 
<math>
 
<math>
Line 110: Line 117:
 
</math>
 
</math>
  
If the
+
 
collision is rejected, a new path length is sampled using \Sigma_\mathrm{m}
+
If the collision is rejected, a new path length is sampled using <math>\Sigma_\mathrm{m}</math>
and the neutron is moved to the next collision site candidate.
+
and the particle is moved to the next tentative collision site.
  
 
Since the majorant cross section is always larger than or equal to the
 
Since the majorant cross section is always larger than or equal to the
total cross section, the path lengths sampled in delta-tracking
+
total cross section, the path lengths sampled using delta-tracking
 
are, on the average, shorter than those sampled with surface-tracking.
 
are, on the average, shorter than those sampled with surface-tracking.
 
The average physical distance between two collisions is preserved, as some
 
The average physical distance between two collisions is preserved, as some
Line 125: Line 132:
 
calculate the surface distances or stop the particle at the material boundaries. This becomes
 
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
 
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 fuel geometries]] and advanced
+
compared to dimensions. The difference between surface- and delta-tracking is emphasized in HTGR [[particle fuels]] and advanced
 
[[CAD-based geometry type|CAD]] and an [[Unstructured mesh-based geometry type|mesh]] based geometry types.  
 
[[CAD-based geometry type|CAD]] and an [[Unstructured mesh-based geometry type|mesh]] based geometry types.  
  
 
Since the majorant cross section used for sampling the path lengths does not depend on the spatial coordinates
 
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 modelling
+
and the material total is needed only at discrete locations, variations of delta-tracking can be used for transporting particles through
inhomogeneous material compositions. This capability is utilized in the [[TMS on-the-fly temperature treatment routine]] and [[multi-physics interface]] in Serpent 2.
+
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
 
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
 
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 in the geometry.
+
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 [[Result estimators#Track-lenght estimator (TLE)|track-length estimate (TLE)]] of
+
 
flux, and reaction rate estimates need to be calculated using the
+
Another drawback is that delta-tracking rules out the use of the [[Result estimators#Implicit estimators|track-length estimate (TLE)]] of particle
potentially less efficient [[Result estimators#collision estimator (CFE)|collision estimator (CFE)]].
+
flux, and integral reaction rate estimates need to be calculated using the
 +
potentially less efficient [[Result estimators#Implicit estimators|collision estimator (CFE)]].
 +
 
 +
 
 +
<div><ul>
 +
<li style="display: inline-block;"> [[File:Dtxs1.png]] </li>
 +
<li style="display: inline-block;"> [[File:Dtxs2.png]] </li>
 +
''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 <sup>155</sup>Gd and <sup>157</sup> 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.''
 +
</ul></div>
  
 
== Hybrid method used in Serpent ==
 
== Hybrid method used in Serpent ==
Line 144: Line 159:
  
 
<references/>
 
<references/>
 +
 +
[[Category:Theory]]

Revision as of 14:22, 29 November 2017

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,[1] and a mathematical verification was derived by Coleman in 1968.[2] Delta-tracking is well described in a text book by Lux and Koblinger,[3] and a description of the methodology used in Serpent is found in an article in Annals of Nuclear Energy from 2010.[4]

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).


  • Dtxs1.png
  • Dtxs2.png
  • 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 155Gd and 157 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.

Hybrid method used in Serpent

References

  1. ^ Woodcock, E. R., Murphy, T., Hemmings, P. J., and Longworth, T. C. "Techniques used in the GEM code for Monte Carlo neutronics calculations in reactors and other systems of complex geometry." ANL-7050, Argonne National Laboratory, 1965.
  2. ^ Coleman, W. A. "Mathematical verification of a certain Monte Carlo sampling technique and applications of the technique to radiation transport problems." Nucl. Sci. Eng., 31 (1968) 76–81.
  3. ^ Lux, I. and Koblinger, L. "Monte Carlo Particle Transport Methods: Neutron and Photon Calculations." CRC Press, Inc. (1991).
  4. ^ Leppänen, J. "Performance of Woodcock delta-tracking in lattice physics applications using the Serpent Monte Carlo reactor physics burnup calculation code." Ann. Nucl. Energy 37 (2010) 715–722.