Merge pull request #483 from samuelshaner/velocity-tally

Inverse Velocity Score
This commit is contained in:
Paul Romano 2015-10-26 04:27:33 -05:00
commit 9ab91e7df7
10 changed files with 124 additions and 21 deletions

View file

@ -1392,13 +1392,15 @@ The ``<tally>`` element accepts the following sub-elements:
options are "flux", "total", "scatter", "absorption", "fission",
"nu-fission", "delayed-nu-fission", "kappa-fission", "nu-scatter",
"scatter-N", "scatter-PN", "scatter-YN", "nu-scatter-N", "nu-scatter-PN",
"nu-scatter-YN", "flux-YN", "total-YN", "current", and "events". These
correspond to the following physical quantities:
"nu-scatter-YN", "flux-YN", "total-YN", "current", "inverse-velocity" and
"events". These correspond to the following physical quantities:
:flux:
Total flux in particle-cm per source particle. Note: The ``analog``
estimator is actually identical to the ``collision`` estimator for the
flux score.
Total flux in particle-cm per source particle.
.. note::
The ``analog`` estimator is actually identical to the ``collision``
estimator for the flux score.
:total:
Total reaction rate in reactions per source particle.
@ -1484,6 +1486,14 @@ The ``<tally>`` element accepts the following sub-elements:
specified. Furthermore, it may not be used in conjunction with any
other score.
:inverse-velocity:
The flux-weighted inverse velocity where the velocity is in units of
meters per second.
.. note::
The ``analog`` estimator is actually identical to the ``collision``
estimator for the inverse-velocity score.
:events:
Number of scoring events. Units are events per source particle.

View file

@ -60,20 +60,22 @@ module constants
! Values here are from the Committee on Data for Science and Technology
! (CODATA) 2010 recommendation (doi:10.1103/RevModPhys.84.1527).
real(8), parameter :: &
PI = 3.1415926535898_8, & ! pi
MASS_NEUTRON = 1.008664916_8, & ! mass of a neutron in amu
MASS_PROTON = 1.007276466812_8, & ! mass of a proton in amu
AMU = 1.660538921e-27_8, & ! 1 amu in kg
N_AVOGADRO = 0.602214129_8, & ! Avogadro's number in 10^24/mol
K_BOLTZMANN = 8.6173324e-11_8, & ! Boltzmann constant in MeV/K
INFINITY = huge(0.0_8), & ! positive infinity
ZERO = 0.0_8, &
HALF = 0.5_8, &
ONE = 1.0_8, &
TWO = 2.0_8, &
THREE = 3.0_8, &
FOUR = 4.0_8
real(8), parameter :: &
PI = 3.1415926535898_8, & ! pi
MASS_NEUTRON = 1.008664916_8, & ! mass of a neutron in amu
MASS_NEUTRON_MEV = 939.565379_8, & ! mass of a neutron in MeV/c^2
MASS_PROTON = 1.007276466812_8, & ! mass of a proton in amu
AMU = 1.660538921e-27_8, & ! 1 amu in kg
C_LIGHT = 2.99792458e8_8, & ! speed of light in m/s
N_AVOGADRO = 0.602214129_8, & ! Avogadro's number in 10^24/mol
K_BOLTZMANN = 8.6173324e-11_8, & ! Boltzmann constant in MeV/K
INFINITY = huge(0.0_8), & ! positive infinity
ZERO = 0.0_8, &
HALF = 0.5_8, &
ONE = 1.0_8, &
TWO = 2.0_8, &
THREE = 3.0_8, &
FOUR = 4.0_8
! ============================================================================
! GEOMETRY-RELATED CONSTANTS
@ -257,7 +259,7 @@ module constants
EVENT_ABSORB = 2
! Tally score type
integer, parameter :: N_SCORE_TYPES = 21
integer, parameter :: N_SCORE_TYPES = 22
integer, parameter :: &
SCORE_FLUX = -1, & ! flux
SCORE_TOTAL = -2, & ! total reaction rate
@ -279,7 +281,8 @@ module constants
SCORE_SCATTER_YN = -18, & ! angular flux-weighted scattering moment (0:N)
SCORE_NU_SCATTER_YN = -19, & ! angular flux-weighted nu-scattering moment (0:N)
SCORE_EVENTS = -20, & ! number of events
SCORE_DELAYED_NU_FISSION = -21 ! delayed neutron production rate
SCORE_DELAYED_NU_FISSION = -21, & ! delayed neutron production rate
SCORE_INVERSE_VELOCITY = -22 ! flux-weighted inverse velocity
! Maximum scattering order supported
integer, parameter :: MAX_ANG_ORDER = 10

View file

@ -3128,6 +3128,8 @@ contains
end if
case ('kappa-fission')
t % score_bins(j) = SCORE_KAPPA_FISSION
case ('inverse-velocity')
t % score_bins(j) = SCORE_INVERSE_VELOCITY
case ('current')
t % score_bins(j) = SCORE_CURRENT
t % type = TALLY_SURFACE_CURRENT

View file

@ -985,6 +985,7 @@ contains
score_names(abs(SCORE_NU_SCATTER_PN)) = "Scattering Prod. Rate Moment"
score_names(abs(SCORE_NU_SCATTER_YN)) = "Scattering Prod. Rate Moment"
score_names(abs(SCORE_DELAYED_NU_FISSION)) = "Delayed-Nu-Fission Rate"
score_names(abs(SCORE_INVERSE_VELOCITY)) = "Flux-Weighted Inverse Velocity"
! Create filename for tally output
filename = trim(path_output) // "tallies.out"

View file

@ -353,6 +353,8 @@ contains
str_array(j) = "nu-scatter-yn"
case (SCORE_EVENTS)
str_array(j) = "events"
case (SCORE_INVERSE_VELOCITY)
str_array(j) = "inverse-velocity"
case default
str_array(j) = reaction_name(tally%score_bins(j))
end select

View file

@ -654,6 +654,8 @@ contains
str_array(j) = "nu-scatter-yn"
case (SCORE_EVENTS)
str_array(j) = "events"
case (SCORE_INVERSE_VELOCITY)
str_array(j) = "inverse-velocity"
case default
str_array(j) = reaction_name(t%score_bins(j))
end select

View file

@ -127,6 +127,27 @@ contains
end if
case (SCORE_INVERSE_VELOCITY)
if (t % estimator == ESTIMATOR_ANALOG) then
! All events score to an inverse velocity bin. We actually use a
! collision estimator in place of an analog one since there is no way
! to count 'events' exactly for the inverse velocity
if (survival_biasing) then
! We need to account for the fact that some weight was already
! absorbed
score = p % last_wgt + p % absorb_wgt
else
score = p % last_wgt
end if
score = score / material_xs % total &
/ (sqrt(TWO * p % E / (MASS_NEUTRON_MEV)) * C_LIGHT)
else
! For inverse velocity, we need no cross section
score = flux / (sqrt(TWO * p % E / (MASS_NEUTRON_MEV)) * C_LIGHT)
end if
case (SCORE_SCATTER, SCORE_SCATTER_N)
if (t % estimator == ESTIMATOR_ANALOG) then
! Skip any event where the particle didn't scatter

View file

@ -0,0 +1 @@
ba1010f940c50314d61aae9f729b7bb476a6b35c5556fbc689ba3fde70ff50b0ba5dad0db3b38349a1562d67817047090ac450a64d895c8f3393163fe7914763

View file

@ -0,0 +1,29 @@
k-combined:
9.903196E-01 4.279617E-02
tally 1:
1.049628E-03
2.261930E-07
4.056389E-04
3.411247E-08
2.243766E-03
1.069671E-06
6.354432E-04
8.370608E-08
tally 2:
1.048031E-03
2.263789E-07
4.276093E-04
4.136358E-08
2.667795E-03
1.528407E-06
6.199140E-04
7.840402E-08
tally 3:
1.048031E-03
2.263789E-07
4.276093E-04
4.136358E-08
2.667795E-03
1.528407E-06
6.199140E-04
7.840402E-08

View file

@ -0,0 +1,32 @@
#!/usr/bin/env python
import os
import sys
sys.path.insert(0, os.pardir)
from testing_harness import TestHarness, PyAPITestHarness
import openmc
class ScoreInverseVelocityTestHarness(PyAPITestHarness):
def _build_inputs(self):
filt = openmc.Filter(type='cell', bins=(21, 22, 23, 27))
tallies = [openmc.Tally(tally_id=i) for i in range(1, 4)]
[t.add_filter(filt) for t in tallies]
[t.add_score('inverse-velocity') for t in tallies]
tallies[0].estimator = 'tracklength'
tallies[1].estimator = 'analog'
tallies[2].estimator = 'collision'
self._input_set.tallies = openmc.TalliesFile()
[self._input_set.tallies.add_tally(t) for t in tallies]
super(ScoreInverseVelocityTestHarness, self)._build_inputs()
def _cleanup(self):
super(ScoreInverseVelocityTestHarness, self)._cleanup()
f = os.path.join(os.getcwd(), 'tallies.xml')
if os.path.exists(f): os.remove(f)
if __name__ == '__main__':
harness = ScoreInverseVelocityTestHarness('statepoint.10.*', True)
harness.main()