Added support for change-in-angle (mu) filters for scattering scores

This commit is contained in:
Adam Nelson 2015-08-25 19:42:19 -04:00
parent 3d1af11b6c
commit ec788085b3
9 changed files with 112 additions and 14 deletions

View file

@ -35,7 +35,8 @@ FILTER_TYPES = {1: 'universe',
6: 'mesh',
7: 'energy',
8: 'energyout',
9: 'distribcell'}
9: 'distribcell',
10:'mu'}
SCORE_TYPES = {-1: 'flux',
-2: 'total',

View file

@ -377,7 +377,7 @@ class StatePoint(object):
raise ValueError(msg)
# Read the bin values
if FILTER_TYPES[filter_type] in ['energy', 'energyout']:
if FILTER_TYPES[filter_type] in ['energy', 'energyout', 'mu']:
bins = self._get_double(
n_bins+1, path='{0}{1}/bins'.format(subbase, j))
@ -417,7 +417,6 @@ class StatePoint(object):
path='{0}{1}/n_score_bins'.format(base, tally_key))[0]
tally.num_score_bins = n_score_bins
scores = [SCORE_TYPES[j] for j in self._get_int(
n_score_bins, path='{0}{1}/score_bins'.format(base, tally_key))]
n_user_scores = self._get_int(

View file

@ -373,7 +373,7 @@ class Tally(object):
if score in self.scores:
return
else:
self._scores.append(score)
self._scores.append(score.strip())
@num_score_bins.setter
def num_score_bins(self, num_score_bins):
@ -896,7 +896,7 @@ class Tally(object):
bins = list(itertools.product(*xyz))
# Create list of 2-tuples for energy boundary bins
elif filter.type in ['energy', 'energyout']:
elif filter.type in ['energy', 'energyout', 'mu']:
bins = []
for k in range(filter.num_bins):
bins.append((filter.bins[k], filter.bins[k+1]))
@ -1251,6 +1251,22 @@ class Tally(object):
filter_bins = np.tile(filter_bins, tile_factor)
df[filter.type + ' [MeV]'] = filter_bins
elif 'mu' is filter.type:
bins = filter.bins
num_bins = filter.num_bins
# Create strings for
template = '{0:.2f} - {1:.2f}'
filter_bins = []
for i in range(num_bins):
filter_bins.append(template.format(bins[i], bins[i+1]))
# Tile the mu bins into a DataFrame column
filter_bins = np.repeat(filter_bins, filter.stride)
tile_factor = data_size / len(filter_bins)
filter_bins = np.tile(filter_bins, tile_factor)
df[filter.type] = filter_bins
# universe, material, surface, cell, and cellborn filters
else:
filter_bins = np.repeat(filter.bins, filter.stride)

View file

@ -300,7 +300,7 @@ module constants
integer, parameter :: NO_BIN_FOUND = -1
! Tally filter and map types
integer, parameter :: N_FILTER_TYPES = 9
integer, parameter :: N_FILTER_TYPES = 10
integer, parameter :: &
FILTER_UNIVERSE = 1, &
FILTER_MATERIAL = 2, &
@ -310,7 +310,8 @@ module constants
FILTER_MESH = 6, &
FILTER_ENERGYIN = 7, &
FILTER_ENERGYOUT = 8, &
FILTER_DISTRIBCELL = 9
FILTER_DISTRIBCELL = 9, &
FILTER_MU = 10
! Tally surface current directions
integer, parameter :: &

View file

@ -607,7 +607,8 @@ contains
! Write filter bins
if (t % filters(j) % type == FILTER_ENERGYIN .or. &
t % filters(j) % type == FILTER_ENERGYOUT) then
t % filters(j) % type == FILTER_ENERGYOUT .or. &
t % filters(j) % type == FILTER_MU) then
call su % write_data(t % filters(j) % real_bins, "bins", &
length=size(t % filters(j) % real_bins), &
group="tallies/tally " // trim(to_str(t % id)) &
@ -653,6 +654,10 @@ contains
call su % write_data("energyout", "type_name", &
group="tallies/tally " // trim(to_str(t % id)) &
// "/filter " // trim(to_str(j)))
case(FILTER_MU)
call su % write_data("mu", "type_name", &
group="tallies/tally " // trim(to_str(t % id)) &
// "/filter " // trim(to_str(j)))
end select
end do FILTER_LOOP

View file

@ -2085,6 +2085,9 @@ contains
integer :: imomstr ! Index of MOMENT_STRS & MOMENT_N_STRS
logical :: file_exists ! does tallies.xml file exist?
real(8) :: rarray3(3) ! temporary double prec. array
integer :: Nmu ! Number of angular bins
real(8) :: dmu ! Mu spacing if using automatic allocation
integer :: imu ! Loop counter for building mu filter bins
character(MAX_LINE_LEN) :: filename
character(MAX_WORD_LEN) :: word
character(MAX_WORD_LEN) :: score_name
@ -2350,8 +2353,9 @@ contains
! Determine number of bins
if (check_for_node(node_filt, "bins")) then
if (trim(temp_str) == 'energy' .or. &
trim(temp_str) == 'energyout') then
if ((trim(temp_str) == 'energy' .or. &
trim(temp_str) == 'energyout') .or. &
(trim(temp_str) == 'mu')) then
n_words = get_arraysize_double(node_filt, "bins")
else
n_words = get_arraysize_integer(node_filt, "bins")
@ -2487,6 +2491,40 @@ contains
! Set to analog estimator
t % estimator = ESTIMATOR_ANALOG
case ('mu')
! Set type of filter
t % filters(j) % type = FILTER_MU
! Set number of bins
t % filters(j) % n_bins = n_words - 1
! Allocate and store bins
allocate(t % filters(j) % real_bins(n_words))
call get_node_array(node_filt, "bins", t % filters(j) % real_bins)
! Easter egg! Allow a user to input a negative number, if it is only one
! and that will mean you subivide [-1,1] evenly with the input
! being the number of bins
if (n_words == 1) then
Nmu = abs(int(t % filters(j) % real_bins(1)))
if (Nmu > 1) then
t % filters(j) % n_bins = Nmu - 1
dmu = TWO / (real(Nmu,8) - ONE)
deallocate(t % filters(j) % real_bins)
allocate(t % filters(j) % real_bins(Nmu))
do imu = 1, Nmu
t % filters(j) % real_bins(imu) = -ONE + (imu - 1) * dmu
end do
else
call fatal_error("Must have more than one bin for mu filter &
& on tally " // trim(to_str(t % id)) // ".")
end if
end if
! Set to analog estimator
t % estimator = ESTIMATOR_ANALOG
case default
! Specified tally filter is invalid, raise error
call fatal_error("Unknown filter type '" &
@ -2656,6 +2694,7 @@ contains
! scores then strip off the n and store it as an integer to be used
! later. Then perform the select case on this modified (number
! removed) string
n_order = -1
score_name = sarray(l)
do imomstr = 1, size(MOMENT_STRS)
if (starts_with(score_name,trim(MOMENT_STRS(imomstr)))) then
@ -2701,6 +2740,23 @@ contains
end do
end if
! Check to see if the mu filter is applied and if that makes sense.
if ((.not. starts_with(score_name,'scatter')) .and. &
(.not. starts_with(score_name,'nu-scatter'))) then
if (t % find_filter(FILTER_MU) > 0) then
call fatal_error("Cannot tally " // trim(score_name) //" with a &
&change of angle (mu) filter.")
end if
! Also check to see if this is a legendre expansion or not.
! If so, we can accept this score and filter combo for p0, but not
! elsewhere.
else if (n_order > 0) then
if (t % find_filter(FILTER_MU) > 0) then
call fatal_error("Cannot tally " // trim(score_name) //" with a &
&change of angle (mu) filter unless order is 0.")
end if
end if
select case (trim(score_name))
case ('flux')
! Prohibit user from tallying flux for an individual nuclide

View file

@ -853,6 +853,17 @@ contains
write(unit_,*) ' Outgoing Energy Bins:' // trim(string)
end if
! Write any change-in-angle bins if present
j = t % find_filter(FILTER_MU)
if (j > 0) then
string = ""
do i = 1, t % filters(j) % n_bins + 1
string = trim(string) // ' ' // trim(to_str(&
t % filters(j) % real_bins(i)))
end do
write(unit_,*) ' Change-in-Angle Bins:' // trim(string)
end if
! Write nuclides bins
write(unit_,fmt='(1X,A)',advance='no') ' Nuclide Bins:'
do i = 1, t % n_nuclide_bins
@ -1734,6 +1745,7 @@ contains
filter_name(FILTER_MESH) = "Mesh"
filter_name(FILTER_ENERGYIN) = "Incoming Energy"
filter_name(FILTER_ENERGYOUT) = "Outgoing Energy"
filter_name(FILTER_MU) = "Change-in-Angle"
! Initialize names for scores
score_names(abs(SCORE_FLUX)) = "Flux"
@ -2171,7 +2183,7 @@ contains
label = "Index (" // trim(to_str(ijk(1))) // ", " // &
trim(to_str(ijk(2))) // ", " // trim(to_str(ijk(3))) // ")"
end if
case (FILTER_ENERGYIN, FILTER_ENERGYOUT)
case (FILTER_ENERGYIN, FILTER_ENERGYOUT, FILTER_MU)
E0 = t % filters(i_filter) % real_bins(bin)
E1 = t % filters(i_filter) % real_bins(bin + 1)
label = "[" // trim(to_str(E0)) // ", " // trim(to_str(E1)) // ")"

View file

@ -256,7 +256,8 @@ contains
group="tallies/tally " // trim(to_str(tally % id)) // &
"/filter " // to_str(j))
if (tally % filters(j) % type == FILTER_ENERGYIN .or. &
tally % filters(j) % type == FILTER_ENERGYOUT) then
tally % filters(j) % type == FILTER_ENERGYOUT .or. &
tally % filters(j) % type == FILTER_MU) then
call sp % write_data(tally % filters(j) % real_bins, "bins", &
group="tallies/tally " // trim(to_str(tally % id)) // &
"/filter " // to_str(j), &
@ -843,7 +844,8 @@ contains
group="tallies/tally " // trim(to_str(curr_key)) // &
"/filter " // to_str(j))
if (tally % filters(j) % type == FILTER_ENERGYIN .or. &
tally % filters(j) % type == FILTER_ENERGYOUT) then
tally % filters(j) % type == FILTER_ENERGYOUT .or. &
tally % filters(j) % type == FILTER_MU) then
call sp % read_data(tally % filters(j) % real_bins, "bins", &
group="tallies/tally " // trim(to_str(curr_key)) // &
"/filter " // to_str(j), &

View file

@ -990,7 +990,6 @@ contains
matching_bins(i) = binary_search(t % filters(i) % real_bins, &
k + 1, p % E)
end if
end select
! Check if no matching bin was found
@ -1245,6 +1244,13 @@ contains
n + 1, p % E)
end if
case (FILTER_MU)
! determine mu bin
n = t % filters(i) % n_bins
! search to find incoming energy bin
matching_bins(i) = binary_search(t % filters(i) % real_bins, &
n + 1, p % mu)
end select
! If the current filter didn't match, exit this subroutine