From 72d84bc286b946723e278fdcc6474d2e93c26a77 Mon Sep 17 00:00:00 2001
From: Alex Lindsay <alexlindsay239@gmail.com>
Date: Mon, 7 May 2018 14:42:09 -0600
Subject: [PATCH] Use consistent normalization for Zernikes.

---
 openmc/filter_expansion.py | 13 +++++++------
 src/math.F90               | 12 +++---------
 2 files changed, 10 insertions(+), 15 deletions(-)

diff --git a/openmc/filter_expansion.py b/openmc/filter_expansion.py
index 872eaac9f5..b092410828 100644
--- a/openmc/filter_expansion.py
+++ b/openmc/filter_expansion.py
@@ -318,15 +318,15 @@ class ZernikeFilter(ExpansionFilter):
 
     This filter allows scores to be multiplied by Zernike polynomials of the
     particle's position normalized to a given unit circle, up to a
-    user-specified order. The Zernike polynomials follow the definition by `Noll
-    <https://doi.org/10.1364/JOSA.66.000207>`_ and are defined as
+    user-specified order. The standard Zernike polynomials follow the definition by
+    Born and Wolf, *Principles of Optics* and are defined as
 
     .. math::
-        Z_n^m(\rho, \theta) = \sqrt{2n + 2} R_n^m(\rho) \cos (m\theta), \quad m > 0
+        Z_n^m(\rho, \theta) = R_n^m(\rho) \cos (m\theta), \quad m > 0
 
-        Z_n^{m}(\rho, \theta) = \sqrt{2n + 2} R_n^{m}(\rho) \sin (m\theta), \quad m < 0
+        Z_n^{m}(\rho, \theta) = R_n^{m}(\rho) \sin (m\theta), \quad m < 0
 
-        Z_n^{m}(\rho, \theta) = \sqrt{n + 1} R_n^{m}(\rho), \quad m = 0
+        Z_n^{m}(\rho, \theta) = R_n^{m}(\rho), \quad m = 0
 
     where the radial polynomials are
 
@@ -335,7 +335,8 @@ class ZernikeFilter(ExpansionFilter):
         \frac{n+m}{2} - k)! (\frac{n-m}{2} - k)!} \rho^{n-2k}.
 
     With this definition, the integral of :math:`(Z_n^m)^2` over the unit disk
-    is exactly :math:`\pi` for each polynomial.
+    is :math:`\frac{\epsilon_m\pi}{2n+2}` for each polynomial where :math:`\epsilon_m` is
+    2 if :math:`m` equals 0 and 1 otherwise.
 
     Specifying a filter with order N tallies moments for all :math:`n` from 0 to
     N and each value of :math:`m`. The ordering of the Zernike polynomial
diff --git a/src/math.F90 b/src/math.F90
index 20bacf88cf..c633131645 100644
--- a/src/math.F90
+++ b/src/math.F90
@@ -600,12 +600,6 @@ contains
     real(8) :: k1, k2, k3, k4          ! Variables for R_m_n calculation
     integer :: i,p,q                   ! Loop counters
 
-    real(8), parameter :: SQRT_N_1(0:10) = [&
-         sqrt(1.0_8), sqrt(2.0_8), sqrt(3.0_8), sqrt(4.0_8), &
-         sqrt(5.0_8), sqrt(6.0_8), sqrt(7.0_8), sqrt(8.0_8), &
-         sqrt(9.0_8), sqrt(10.0_8), sqrt(11.0_8)]
-    real(8), parameter :: SQRT_2N_2(0:10) = SQRT_N_1*sqrt(2.0_8)
-
     ! n == radial degree
     ! m == azimuthal frequency
 
@@ -669,11 +663,11 @@ contains
     do p = 0, n
       do q = -p, p, 2
         if (q < 0) then
-          zn(i) = zn_mat(p+1, abs(q)+1) * sin_phi_vec(abs(q)) * SQRT_2N_2(p)
+          zn(i) = zn_mat(p+1, abs(q)+1) * sin_phi_vec(abs(q))
         else if (q == 0) then
-          zn(i) = zn_mat(p+1, q+1) * SQRT_N_1(p)
+          zn(i) = zn_mat(p+1, q+1)
         else
-          zn(i) = zn_mat(p+1, q+1) * cos_phi_vec(abs(q)+1) * SQRT_2N_2(p)
+          zn(i) = zn_mat(p+1, q+1) * cos_phi_vec(abs(q)+1)
         end if
         i = i + 1
       end do