From c10859ae024dde5e5ad7efd086522cb49d82f764 Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Thu, 3 Nov 2016 17:21:12 -0400 Subject: [PATCH] Move diff tallies from MeV to eV --- src/cross_section.F90 | 28 ++++------ src/tally.F90 | 68 ++++++++++++------------ tests/test_diff_tally/inputs_true.dat | 2 +- tests/test_diff_tally/test_diff_tally.py | 2 +- 4 files changed, 45 insertions(+), 55 deletions(-) diff --git a/src/cross_section.F90 b/src/cross_section.F90 index 5aac3b837..39552a716 100644 --- a/src/cross_section.F90 +++ b/src/cross_section.F90 @@ -712,15 +712,14 @@ contains ! temperature. !=============================================================================== - subroutine multipole_deriv_eval(multipole, Emev, sqrtkT_, sigT, sigA, sigF) + subroutine multipole_deriv_eval(multipole, E, sqrtkT, sigT, sigA, sigF) type(MultipoleArray), intent(in) :: multipole ! The windowed multipole ! object to process. - real(8), intent(in) :: Emev ! The energy at which to + real(8), intent(in) :: E ! The energy at which to ! evaluate the cross section - ! in MeV - real(8), intent(in) :: sqrtkT_ ! The temperature in the form - ! sqrt(kT (in MeV)), at which - ! to evaluate the XS. + real(8), intent(in) :: sqrtkT ! The temperature in the form + ! sqrt(kT), at which to + ! evaluate the XS. real(8), intent(out) :: sigT ! Total cross section real(8), intent(out) :: sigA ! Absorption cross section real(8), intent(out) :: sigF ! Fission cross section @@ -730,8 +729,6 @@ contains real(8) :: sqrtE ! sqrt(E), eV real(8) :: invE ! 1/E, eV real(8) :: dopp ! sqrt(atomic weight ratio / kT) - real(8) :: E ! energy, eV - real(8) :: sqrtkT ! sqrt(kT (in eV)) integer :: i_pole ! index of pole integer :: i_window ! index of window integer :: startw ! window start pointer (for poles) @@ -741,15 +738,11 @@ contains ! ========================================================================== ! Bookkeeping - ! Convert to eV. - E = Emev * 1.0e6_8 - sqrtkT = sqrtkT_ * 1.0e3_8 - ! Define some frequently used variables. sqrtE = sqrt(E) invE = ONE / E dopp = multipole % sqrtAWR / sqrtkT - T = sqrtkT_**2 / K_BOLTZMANN + T = sqrtkT**2 / K_BOLTZMANN if (sqrtkT == ZERO) call fatal_error("Windowed multipole temperature & &derivatives are not implemented for 0 Kelvin cross sections.") @@ -795,12 +788,9 @@ contains sigF = sigF + real(multipole % data(RM_RF, i_pole) * w_val) end if end do - sigT = -HALF*multipole % sqrtAWR / sqrt(K_BOLTZMANN*1.0e6_8) * T**(-1.5)& - * sigT - sigA = -HALF*multipole % sqrtAWR / sqrt(K_BOLTZMANN*1.0e6_8) * T**(-1.5)& - * sigA - sigF = -HALF*multipole % sqrtAWR / sqrt(K_BOLTZMANN*1.0e6_8) * T**(-1.5)& - * sigF + sigT = -HALF*multipole % sqrtAWR / sqrt(K_BOLTZMANN) * T**(-1.5) * sigT + sigA = -HALF*multipole % sqrtAWR / sqrt(K_BOLTZMANN) * T**(-1.5) * sigA + sigF = -HALF*multipole % sqrtAWR / sqrt(K_BOLTZMANN) * T**(-1.5) * sigF end if end subroutine multipole_deriv_eval diff --git a/src/tally.F90 b/src/tally.F90 index e99d6072a..f3cd015f7 100644 --- a/src/tally.F90 +++ b/src/tally.F90 @@ -2984,8 +2984,8 @@ contains dsigT = ZERO associate (nuc => nuclides(p % event_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3011,8 +3011,8 @@ contains dsigA = ZERO associate (nuc => nuclides(p % event_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3037,8 +3037,8 @@ contains dsigA = ZERO associate (nuc => nuclides(p % event_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3062,8 +3062,8 @@ contains dsigF = ZERO associate (nuc => nuclides(p % event_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3087,8 +3087,8 @@ contains dsigF = ZERO associate (nuc => nuclides(p % event_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3123,8 +3123,8 @@ contains do l = 1, mat % n_nuclides associate (nuc => nuclides(mat % nuclide(l))) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8 .and. & + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E .and. & micro_xs(mat % nuclide(l)) % total > ZERO) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) @@ -3140,8 +3140,8 @@ contains dsigT = ZERO associate (nuc => nuclides(i_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3161,8 +3161,8 @@ contains do l = 1, mat % n_nuclides associate (nuc => nuclides(mat % nuclide(l))) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8 .and. & + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E .and. & (micro_xs(mat % nuclide(l)) % total & - micro_xs(mat % nuclide(l)) % absorption) > ZERO) then call multipole_deriv_eval(nuc % multipole, p % last_E, & @@ -3182,8 +3182,8 @@ contains dsigA = ZERO associate (nuc => nuclides(i_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3204,8 +3204,8 @@ contains do l = 1, mat % n_nuclides associate (nuc => nuclides(mat % nuclide(l))) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8 .and. & + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E .and. & micro_xs(mat % nuclide(l)) % absorption > ZERO) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) @@ -3221,8 +3221,8 @@ contains dsigA = ZERO associate (nuc => nuclides(i_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3242,8 +3242,8 @@ contains do l = 1, mat % n_nuclides associate (nuc => nuclides(mat % nuclide(l))) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8 .and. & + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E .and. & micro_xs(mat % nuclide(l)) % fission > ZERO) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) @@ -3259,8 +3259,8 @@ contains dsigF = ZERO associate (nuc => nuclides(i_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3280,8 +3280,8 @@ contains do l = 1, mat % n_nuclides associate (nuc => nuclides(mat % nuclide(l))) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8 .and. & + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E .and. & micro_xs(mat % nuclide(l)) % nu_fission > ZERO) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) @@ -3299,8 +3299,8 @@ contains dsigF = ZERO associate (nuc => nuclides(i_nuclide)) if (nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % last_E, & p % sqrtkT, dsigT, dsigA, dsigF) end if @@ -3371,8 +3371,8 @@ contains do l=1, mat % n_nuclides associate (nuc => nuclides(mat % nuclide(l))) if (nuc % mp_present .and. & - p % E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % E <= nuc % multipole % end_E/1.0e6_8) then + p % E >= nuc % multipole % start_E .and. & + p % E <= nuc % multipole % end_E) then call multipole_deriv_eval(nuc % multipole, p % E, & p % sqrtkT, dsigT, dsigA, dsigF) deriv % flux_deriv = deriv % flux_deriv & @@ -3444,8 +3444,8 @@ contains associate (nuc => nuclides(mat % nuclide(l))) if (mat % nuclide(l) == p % event_nuclide .and. & nuc % mp_present .and. & - p % last_E >= nuc % multipole % start_E/1.0e6_8 .and. & - p % last_E <= nuc % multipole % end_E/1.0e6_8) then + p % last_E >= nuc % multipole % start_E .and. & + p % last_E <= nuc % multipole % end_E) then ! phi = Sigma_MT ! (1 / phi) * (d_phi / d_T) = (d_Sigma_MT / d_T) / Sigma_MT ! (1 / phi) * (d_phi / d_T) = (d_sigma_MT / d_T) / sigma_MT diff --git a/tests/test_diff_tally/inputs_true.dat b/tests/test_diff_tally/inputs_true.dat index a82c8bca8..f82b5d8f0 100644 --- a/tests/test_diff_tally/inputs_true.dat +++ b/tests/test_diff_tally/inputs_true.dat @@ -1 +1 @@ -58a47bc1cce52a90cb7d0ed5e8143fe0ac2a3208fb7b4c050647ab2c2f3486272250138b7e548b2bee912085b4a70d818833760e1f90e9ad6aa713caa12473bf \ No newline at end of file +be155010540ca076356596aef8511edd58c91772a5e956e5af77437a39ce3d238b91b57df30c2d599fc6e01c887bf73897ee8ed91425c2d475d97202689c2b88 \ No newline at end of file diff --git a/tests/test_diff_tally/test_diff_tally.py b/tests/test_diff_tally/test_diff_tally.py index ce5b4a68e..5e34fd886 100644 --- a/tests/test_diff_tally/test_diff_tally.py +++ b/tests/test_diff_tally/test_diff_tally.py @@ -30,7 +30,7 @@ class DiffTallyTestHarness(PyAPITestHarness): self._input_set.tallies = openmc.Tallies() filt_mats = openmc.MaterialFilter((1, 3)) - filt_eout = openmc.EnergyoutFilter((0.0, 1.0, 20.0)) + filt_eout = openmc.EnergyoutFilter((0.0, 1.0e6, 20.0e6)) # We want density derivatives for both water and fuel to get coverage # for both fissile and non-fissile materials.