Remove redundant surface tally subroutine. Fix broken tests

This commit is contained in:
Paul Romano 2018-03-15 22:30:07 -05:00
parent c51091a1d7
commit 206166cc41
12 changed files with 63 additions and 166 deletions

View file

@ -55,6 +55,9 @@ _dll.openmc_material_filter_set_bins.errcheck = _error_handler
_dll.openmc_mesh_filter_set_mesh.argtypes = [c_int32, c_int32]
_dll.openmc_mesh_filter_set_mesh.restype = c_int
_dll.openmc_mesh_filter_set_mesh.errcheck = _error_handler
_dll.openmc_meshsurface_filter_set_mesh.argtypes = [c_int32, c_int32]
_dll.openmc_meshsurface_filter_set_mesh.restype = c_int
_dll.openmc_meshsurface_filter_set_mesh.errcheck = _error_handler
class Filter(_FortranObjectWithID):
@ -188,6 +191,10 @@ class MeshFilter(Filter):
filter_type = 'mesh'
class MeshSurfaceFilter(Filter):
filter_type = 'meshsurface'
class MuFilter(Filter):
filter_type = 'mu'
@ -216,6 +223,7 @@ _FILTER_TYPE_MAP = {
'energyfunction': EnergyFunctionFilter,
'material': MaterialFilter,
'mesh': MeshFilter,
'meshsurface': MeshSurfaceFilter,
'mu': MuFilter,
'polar': PolarFilter,
'surface': SurfaceFilter,

View file

@ -276,8 +276,9 @@ contains
! Clear active tally lists
call active_analog_tallies % clear()
call active_tracklength_tallies % clear()
call active_current_tallies % clear()
call active_meshsurf_tallies % clear()
call active_collision_tallies % clear()
call active_surface_tallies % clear()
call active_tallies % clear()
! Reset timers

View file

@ -35,14 +35,14 @@ element tallies {
element filter {
(element id { xsd:int } | attribute id { xsd:int }) &
(
( (element type { ( "cell" | "cellfrom" | "cellborn" | "material" |
"universe" | "surface" | "distribcell" | "mesh" | "energy" |
"energyout" | "mu" | "polar" | "azimuthal" | "delayedgroup" |
"energyfunction") } |
attribute type { ( "cell" | "cellfrom" | "cellborn" | "material" |
"universe" | "surface" | "distribcell" | "mesh" | "energy" |
"energyout" | "mu" | "polar" | "azimuthal" | "delayedgroup" |
"energyfunction") }) &
( (element type { ( "cell" | "cellfrom" | "cellborn" | "material" |
"universe" | "surface" | "distribcell" | "mesh" | "energy" |
"energyout" | "mu" | "polar" | "azimuthal" | "delayedgroup" |
"energyfunction" | "meshsurface") } |
attribute type { ( "cell" | "cellfrom" | "cellborn" | "material" |
"universe" | "surface" | "distribcell" | "mesh" | "energy" |
"energyout" | "mu" | "polar" | "azimuthal" | "delayedgroup" |
"energyfunction" | "meshsurface") }) &
(element bins { list { xsd:double+ } } |
attribute bins { list { xsd:double+ } })
) |

View file

@ -182,6 +182,7 @@
<value>azimuthal</value>
<value>delayedgroup</value>
<value>energyfunction</value>
<value>meshsurface</value>
</choice>
</element>
<attribute name="type">
@ -201,6 +202,7 @@
<value>azimuthal</value>
<value>delayedgroup</value>
<value>energyfunction</value>
<value>meshsurface</value>
</choice>
</attribute>
</choice>

View file

@ -3064,9 +3064,9 @@ contains
! tally total or partial currents between two cells
!===============================================================================
subroutine score_surface_tally(p)
type(Particle), intent(in) :: p
subroutine score_surface_tally(p, tally_vec)
type(Particle), intent(in) :: p
type(VectorInt), intent(in) :: tally_vec
integer :: i
integer :: i_tally
@ -3086,9 +3086,9 @@ contains
! No collision, so no weight change when survival biasing
flux = p % wgt
TALLY_LOOP: do i = 1, active_surface_tallies % size()
TALLY_LOOP: do i = 1, tally_vec % size()
! Get index of tally and pointer to tally
i_tally = active_surface_tallies % data(i)
i_tally = tally_vec % data(i)
associate (t => tallies(i_tally) % obj)
! Find all valid bins in each filter if they have not already been found
@ -3188,134 +3188,6 @@ contains
end subroutine score_surface_tally
!===============================================================================
! SCORE_SURFACE_CURRENT tallies surface crossings in a mesh tally by manually
! determining which mesh surfaces were crossed
!===============================================================================
subroutine score_surface_current(p)
type(Particle), intent(in) :: p
integer :: i
integer :: i_tally
integer :: i_filt
integer :: i_bin
integer :: q ! loop index for scoring bins
integer :: k ! working index for expand and score
integer :: score_bin ! scoring bin, e.g. SCORE_FLUX
integer :: score_index ! scoring bin index
integer :: j ! loop index for scoring bins
integer :: filter_index ! single index for single bin
real(8) :: flux ! collision estimate of flux
real(8) :: filter_weight ! combined weight of all filters
real(8) :: score ! analog tally score
logical :: finished ! found all valid bin combinations
! No collision, so no weight change when survival biasing
flux = p % wgt
TALLY_LOOP: do i = 1, active_current_tallies % size()
! Get index of tally and pointer to tally
i_tally = active_current_tallies % data(i)
associate (t => tallies(i_tally) % obj)
! Find all valid bins in each filter if they have not already been found
! for a previous tally.
do j = 1, size(t % filter)
i_filt = t % filter(j)
if (.not. filter_matches(i_filt) % bins_present) then
call filter_matches(i_filt) % bins % clear()
call filter_matches(i_filt) % weights % clear()
call filters(i_filt) % obj % get_all_bins(p, t % estimator, &
filter_matches(i_filt))
filter_matches(i_filt) % bins_present = .true.
end if
! If there are no valid bins for this filter, then there is nothing to
! score and we can move on to the next tally.
if (filter_matches(i_filt) % bins % size() == 0) cycle TALLY_LOOP
! Set the index of the bin used in the first filter combination
filter_matches(i_filt) % i_bin = 1
end do
! ========================================================================
! Loop until we've covered all valid bins on each of the filters.
FILTER_LOOP: do
! Reset scoring index and weight
filter_index = 1
filter_weight = ONE
! Determine scoring index and weight for this filter combination
do j = 1, size(t % filter)
i_filt = t % filter(j)
i_bin = filter_matches(i_filt) % i_bin
filter_index = filter_index + (filter_matches(i_filt) % bins % &
data(i_bin) - 1) * t % stride(j)
filter_weight = filter_weight * filter_matches(i_filt) % weights % &
data(i_bin)
end do
! Determine score
score = flux * filter_weight
! Currently only one score type
k = 0
SCORE_LOOP: do q = 1, t % n_user_score_bins
k = k + 1
! determine what type of score bin
score_bin = t % score_bins(q)
! determine scoring bin index, no offset from nuclide bins
score_index = q
call expand_and_score(p, t, score_index, filter_index, score_bin, &
score, k)
end do SCORE_LOOP
! ======================================================================
! Filter logic
! Increment the filter bins, starting with the last filter to find the
! next valid bin combination
finished = .true.
do j = size(t % filter), 1, -1
i_filt = t % filter(j)
if (filter_matches(i_filt) % i_bin < filter_matches(i_filt) % &
bins % size()) then
filter_matches(i_filt) % i_bin = filter_matches(i_filt) % i_bin + 1
finished = .false.
exit
else
filter_matches(i_filt) % i_bin = 1
end if
end do
! Once we have finished all valid bins for each of the filters, exit
! the loop.
if (finished) exit FILTER_LOOP
end do FILTER_LOOP
end associate
! If the user has specified that we can assume all tallies are spatially
! separate, this implies that once a tally has been scored to, we needn't
! check the others. This cuts down on overhead when there are many
! tallies specified
if (assume_separate) exit TALLY_LOOP
end do TALLY_LOOP
! Reset filter matches flag
filter_matches(:) % bins_present = .false.
end subroutine score_surface_current
!===============================================================================
! APPLY_DERIVATIVE_TO_SCORE multiply the given score by its relative derivative
!===============================================================================
@ -4219,7 +4091,7 @@ contains
call active_collision_tallies % clear()
call active_tracklength_tallies % clear()
call active_surface_tallies % clear()
call active_current_tallies % clear()
call active_meshsurf_tallies % clear()
do i = 1, n_tallies
associate (t => tallies(i) % obj)
@ -4237,7 +4109,7 @@ contains
call active_collision_tallies % push_back(i)
end if
elseif (t % type == TALLY_MESH_SURFACE) then
call active_current_tallies % push_back(i)
call active_meshsurf_tallies % push_back(i)
elseif (t % type == TALLY_SURFACE) then
call active_surface_tallies % push_back(i)
end if

View file

@ -144,7 +144,7 @@ module tally_header
! Active tally lists
type(VectorInt), public :: active_analog_tallies
type(VectorInt), public :: active_tracklength_tallies
type(VectorInt), public :: active_current_tallies
type(VectorInt), public :: active_meshsurf_tallies
type(VectorInt), public :: active_collision_tallies
type(VectorInt), public :: active_tallies
type(VectorInt), public :: active_surface_tallies
@ -418,7 +418,7 @@ contains
! Deallocate tally node lists
call active_analog_tallies % clear()
call active_tracklength_tallies % clear()
call active_current_tallies % clear()
call active_meshsurf_tallies % clear()
call active_collision_tallies % clear()
call active_surface_tallies % clear()
call active_tallies % clear()

View file

@ -19,9 +19,9 @@ module tracking
use surface_header
use tally_header
use tally, only: score_analog_tally, score_tracklength_tally, &
score_collision_tally, score_surface_current, &
score_track_derivative, score_surface_tally, &
score_collision_derivative, zero_flux_derivs
score_collision_tally, score_surface_tally, &
score_track_derivative, zero_flux_derivs, &
score_collision_derivative
use track_output, only: initialize_particle_track, write_particle_track, &
add_particle_track, finalize_particle_track
@ -185,7 +185,8 @@ contains
p % event = EVENT_SURFACE
end if
! Score cell to cell partial currents
if(active_surface_tallies % size() > 0) call score_surface_tally(p)
if(active_surface_tallies % size() > 0) &
call score_surface_tally(p, active_surface_tallies)
else
! ====================================================================
! PARTICLE HAS COLLISION
@ -200,7 +201,8 @@ contains
! since the direction of the particle will change and we need to use the
! pre-collision direction to figure out what mesh surfaces were crossed
if (active_current_tallies % size() > 0) call score_surface_current(p)
if (active_meshsurf_tallies % size() > 0) &
call score_surface_tally(p, active_meshsurf_tallies)
! Clear surface component
p % surface = NONE
@ -317,12 +319,12 @@ contains
! forward slightly so that if the mesh boundary is on the surface, it is
! still processed
if (active_current_tallies % size() > 0) then
if (active_meshsurf_tallies % size() > 0) then
! TODO: Find a better solution to score surface currents than
! physically moving the particle forward slightly
p % coord(1) % xyz = p % coord(1) % xyz + TINY_BIT * p % coord(1) % uvw
call score_surface_current(p)
call score_surface_tally(p, active_meshsurf_tallies)
end if
! Score to global leakage tally
@ -350,10 +352,10 @@ contains
! particle to change -- artificially move the particle slightly back in
! case the surface crossing is coincident with a mesh boundary
if (active_current_tallies % size() > 0) then
if (active_meshsurf_tallies % size() > 0) then
xyz = p % coord(1) % xyz
p % coord(1) % xyz = p % coord(1) % xyz - TINY_BIT * p % coord(1) % uvw
call score_surface_current(p)
call score_surface_tally(p, active_meshsurf_tallies)
p % coord(1) % xyz = xyz
end if
@ -407,10 +409,10 @@ contains
! Score surface currents since reflection causes the direction of the
! particle to change -- artificially move the particle slightly back in
! case the surface crossing is coincident with a mesh boundary
if (active_current_tallies % size() > 0) then
if (active_meshsurf_tallies % size() > 0) then
xyz = p % coord(1) % xyz
p % coord(1) % xyz = p % coord(1) % xyz - TINY_BIT * p % coord(1) % uvw
call score_surface_current(p)
call score_surface_tally(p, active_meshsurf_tallies)
p % coord(1) % xyz = xyz
end if

View file

@ -327,18 +327,27 @@
<filter id="1" type="mesh">
<bins>1</bins>
</filter>
<filter id="4" type="meshsurface">
<bins>1</bins>
</filter>
<filter id="2" type="mesh">
<bins>2</bins>
</filter>
<filter id="5" type="meshsurface">
<bins>2</bins>
</filter>
<filter id="3" type="mesh">
<bins>3</bins>
</filter>
<filter id="6" type="meshsurface">
<bins>3</bins>
</filter>
<tally id="1" name="tally 1">
<filters>1</filters>
<scores>total</scores>
</tally>
<tally id="2" name="tally 2">
<filters>1</filters>
<filters>4</filters>
<scores>current</scores>
</tally>
<tally id="3" name="tally 3">
@ -346,7 +355,7 @@
<scores>total</scores>
</tally>
<tally id="4" name="tally 4">
<filters>2</filters>
<filters>5</filters>
<scores>current</scores>
</tally>
<tally id="5" name="tally 5">
@ -354,7 +363,7 @@
<scores>total</scores>
</tally>
<tally id="6" name="tally 6">
<filters>3</filters>
<filters>6</filters>
<scores>current</scores>
</tally>
</tallies>

View file

@ -1 +1 @@
804d161cb8eae506d3247a533d122f44a01d3cedd566b3c65c71a0b51326dd9b97f8bbf45af7304b500476f3f854d91b10ccad94122d15f23641b05b835fada6
46950c046648faaa5ff3cb7b4fdd03667ae3c6da96a7ed8121761291de452cf88481f53e967ed52407d77d90109ae24839967a22ae216531a0e1491f5051ea72

View file

@ -30,6 +30,9 @@ class FilterMeshTestHarness(HashedPyAPITestHarness):
mesh_1d_filter = openmc.MeshFilter(mesh_1d)
mesh_2d_filter = openmc.MeshFilter(mesh_2d)
mesh_3d_filter = openmc.MeshFilter(mesh_3d)
meshsurf_1d_filter = openmc.MeshSurfaceFilter(mesh_1d)
meshsurf_2d_filter = openmc.MeshSurfaceFilter(mesh_2d)
meshsurf_3d_filter = openmc.MeshSurfaceFilter(mesh_3d)
# Initialized the tallies
tally = openmc.Tally(name='tally 1')
@ -38,7 +41,7 @@ class FilterMeshTestHarness(HashedPyAPITestHarness):
self._model.tallies.append(tally)
tally = openmc.Tally(name='tally 2')
tally.filters = [mesh_1d_filter]
tally.filters = [meshsurf_1d_filter]
tally.scores = ['current']
self._model.tallies.append(tally)
@ -48,7 +51,7 @@ class FilterMeshTestHarness(HashedPyAPITestHarness):
self._model.tallies.append(tally)
tally = openmc.Tally(name='tally 4')
tally.filters = [mesh_2d_filter]
tally.filters = [meshsurf_2d_filter]
tally.scores = ['current']
self._model.tallies.append(tally)
@ -58,7 +61,7 @@ class FilterMeshTestHarness(HashedPyAPITestHarness):
self._model.tallies.append(tally)
tally = openmc.Tally(name='tally 6')
tally.filters = [mesh_3d_filter]
tally.filters = [meshsurf_3d_filter]
tally.scores = ['current']
self._model.tallies.append(tally)

View file

@ -1 +1 @@
4675d5101f4f829369c39cb33d654430836b934ab07c165777ba6e214bdf3a698b8082f4f9bb9e78f1f0e495b30ea02cf9b3d14622c59915d818d678a1e5b7b1
138b312cdaa822c9b62f757b3259522004b679c4ed289a92798b29a6442c26d12c53256635be273f13e3703816ff50a1b9f52d79770eade01e482384ac0d389f

View file

@ -9,7 +9,7 @@
</mesh>
<filter id="1">
<type>mesh</type>
<type>meshsurface</type>
<bins>1</bins>
</filter>