From 5c4dcff7cddfe8ce60ef07defc2ff2f9d150d465 Mon Sep 17 00:00:00 2001 From: Colin Josey Date: Tue, 5 Apr 2016 21:09:31 -0400 Subject: [PATCH] Update the WMP documentation This commit responds to comments on PR #618. It also clarifies some points about the multipole algorithm that are not well documented. --- docs/source/methods/cross_sections.rst | 61 +++++++++++++++++++------- 1 file changed, 44 insertions(+), 17 deletions(-) diff --git a/docs/source/methods/cross_sections.rst b/docs/source/methods/cross_sections.rst index 24bf9dab0..50a97d575 100644 --- a/docs/source/methods/cross_sections.rst +++ b/docs/source/methods/cross_sections.rst @@ -74,41 +74,67 @@ allows on-the-fly Doppler broadening to arbitrary temperature. The multipole method was introduced by [Hwang]_ and the faster windowed multipole method by [Josey]_. In the multipole format, cross section resonances -are represented by poles, `p_j`, and residues, `r_j`, in the complex plane. The -0 K cross sections in the resolved resonance region can be computed by summing -up a contribution from each pole: +are represented by poles, :math:`p_j`, and residues, :math:`r_j`, in the complex +plane. The 0K cross sections in the resolved resonance region can be computed +by summing up a contribution from each pole: .. math:: \sigma(E, T=0\text{K}) = \frac{1}{E} \sum_j \text{Re} \left[ \frac{i r_j}{\sqrt{E} - p_j} \right] Assuming free-gas thermal motion, cross sections in the multipole form can be -analytically Doppler broadened to give (after some approximation) the form: +analytically Doppler broadened to give the form: .. math:: \sigma(E, T) = \frac{1}{2 E \sqrt{\xi}} \sum_j \text{Re} \left[i r_j - \sqrt{\pi} W(z) \right] + \sqrt{\pi} W_i(z) - \frac{r_j}{\sqrt{\pi}} C \left(\frac{p_j}{\sqrt{\xi}}, + \frac{u}{2 \sqrt{\xi}}\right)\right] +.. math:: + W_i(z) = \frac{i}{\pi} \int_{-\infty}^\infty dt \frac{e^{-t^2}}{z - t} +.. math:: + C \left(\frac{p_j}{\sqrt{\xi}},\frac{u}{2 \sqrt{\xi}}\right) = + 2p_j \int_0^\infty du' \frac{e^{-(u + u')^2/4\xi}}{p_j^2 - u'^2} .. math:: z = \frac{\sqrt{E} - p_j}{2 \sqrt{\xi}} .. math:: \xi = \frac{k_B T}{4 A} +.. math:: + u = \sqrt{E} -where `T` is the temperature of the resonant scatterer, `k_B` is the Boltzmann -constant, `A` is the mass of the target nucleus, and `W` is the Faddeeva -function. +where :math:`T` is the temperature of the resonant scatterer, :math:`k_B` is the +Boltzmann constant, :math:`A` is the mass of the target nucleus. For +:math:`E \gg k_b T/A`, the :math:`C` integral is approximately zero, simplifying +the cross section to: .. math:: - W(z) = e^{-z^2} \text{erfc}(-iz) + \sigma(E, T) = \frac{1}{2 E \sqrt{\xi}} \sum_j \text{Re} \left[i r_j + \sqrt{\pi} W_i(z)\right] -It is prohibitively expensive to evaluate a Faddeeva function for all of the -poles every time a cross section lookup is performed in Monte Carlo. To -mitigate that computational cost, the WMP method only evaluates poles within a -certain energy "window" around the incident neutron energy and accounts for the -effect of resonances outside that window with a polynomial fit. +The :math:`W_i` integral simplifies down to an analytic form. We define the +Faddeeva function, :math:`W` as: -Note that the implementation of WMP in OpenMC assumes that inelastic scattering -does not occur in the resolved resonance region. +.. math:: + W(z) = e^{-z^2} \text{Erfc}(-iz) +Through this, the integral transforms as follows: + +.. math:: + \text{Im} (z) > 0 : W_i(z) = W(z) +.. math:: + \text{Im} (z) < 0 : W_i(z) = -W(z^*)^* + +There are freely available algorithms_ to evaluate the Faddeeva function. For +many nuclides, the Faddeeva function needs to be evaluated thousands of times to +calculate a cross section. To mitigate that computational cost, the WMP method +only evaluates poles within a certain energy "window" around the incident +neutron energy and accounts for the effect of resonances outside that window +with a polynomial fit. This polynomial fit is then broadened exactly. This +exact broadening can make up for the removal of the :math:`C` integral, as +typically at low energies, only curve fits are used. + +Note that the implementation of WMP in OpenMC currently assumes that inelastic +scattering does not occur in the resolved resonance region. This is usually, +but not always the case. Future library versions may eliminate this issue. .. only:: html @@ -117,7 +143,7 @@ does not occur in the resolved resonance region. .. [Brown] Forrest B. Brown, "New Hash-based Energy Lookup Algorithm for Monte Carlo codes," LA-UR-14-24530, Los Alamos National Laboratory (2014). -.. [Hwang] R. N. Hwag, "A Rigorous Pole Representation of Multilevel Cross +.. [Hwang] R. N. Hwang, "A Rigorous Pole Representation of Multilevel Cross Sections and Its Practical Application," *Nucl. Sci. Eng.*, **96**, 192-209 (1987). @@ -130,3 +156,4 @@ does not occur in the resolved resonance region. .. _NJOY: http://t2.lanl.gov/codes.shtml .. _ENDF/B data: http://www.nndc.bnl.gov/endf .. _Leppanen: http://dx.doi.org/10.1016/j.anucene.2009.03.019 +.. _algorithms: http://ab-initio.mit.edu/wiki/index.php/Faddeeva_Package