From 39bd9b7b69f3e895ee4ef4afd208c7dc759106be Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Fri, 25 Jan 2019 15:36:52 -0500 Subject: [PATCH] Move indices like tally % filter_energyin to C++ --- include/openmc/tallies/tally.h | 9 ++++ src/input_xml.F90 | 6 +-- src/tallies/tally.F90 | 26 +++++----- src/tallies/tally.cpp | 44 ++++++++++++++++ src/tallies/tally_header.F90 | 92 ++++++++++++++++++++++++++-------- src/tallies/trigger.F90 | 6 +-- 6 files changed, 142 insertions(+), 41 deletions(-) diff --git a/include/openmc/tallies/tally.h b/include/openmc/tallies/tally.h index 0f68afe3f..5f4deb34c 100644 --- a/include/openmc/tallies/tally.h +++ b/include/openmc/tallies/tally.h @@ -30,6 +30,15 @@ public: int32_t n_filter_bins() const {return n_filter_bins_;} + // We need to have quick access to some filters. The following gives indices + // for various filters that could be in the tally or C_NONE if they are not + // present. + int energyin_filter_ {C_NONE}; + int energyout_filter_ {C_NONE}; + int delayedgroup_filter_ {C_NONE}; + int surface_filter_ {C_NONE}; + int mesh_filter_ {C_NONE}; + private: std::vector filters_; //< Filter indices in global filters array diff --git a/src/input_xml.F90 b/src/input_xml.F90 index a5dae6c82..304ed9f03 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -1167,10 +1167,10 @@ contains deallocate(temp_filter) ! Check for the presence of certain filter types - has_energyout = (t % energyout_filter > 0) - has_delayedgroup = (t % delayedgroup_filter > 0) + has_energyout = (t % energyout_filter() > 0) + has_delayedgroup = (t % delayedgroup_filter() > 0) has_legendre = .false. - has_surface = (t % surface_filter > 0) + has_surface = (t % surface_filter() > 0) has_cell = .false. has_cellfrom = .false. has_meshsurface = .false. diff --git a/src/tallies/tally.F90 b/src/tallies/tally.F90 index 02dacb251..0232ed602 100644 --- a/src/tallies/tally.F90 +++ b/src/tallies/tally.F90 @@ -296,7 +296,7 @@ contains if (t % estimator == ESTIMATOR_ANALOG) then if (survival_biasing .or. p % fission) then - if (t % energyout_filter > 0) then + if (t % energyout_filter() > 0) then ! Normally, we only need to make contributions to one scoring ! bin. However, in the case of fission, since multiple fission ! neutrons were emitted with different energies, multiple @@ -342,7 +342,7 @@ contains if (t % estimator == ESTIMATOR_ANALOG) then if (survival_biasing .or. p % fission) then - if (t % energyout_filter > 0) then + if (t % energyout_filter() > 0) then ! Normally, we only need to make contributions to one scoring ! bin. However, in the case of fission, since multiple fission ! neutrons were emitted with different energies, multiple @@ -406,11 +406,11 @@ contains if (material_xs % absorption == ZERO) cycle SCORE_LOOP ! Set the delayedgroup filter index and the number of delayed group bins - dg_filter = t % delayedgroup_filter + dg_filter = t % delayedgroup_filter() if (t % estimator == ESTIMATOR_ANALOG) then if (survival_biasing .or. p % fission) then - if (t % energyout_filter > 0) then + if (t % energyout_filter() > 0) then ! Normally, we only need to make contributions to one scoring ! bin. However, in the case of fission, since multiple fission ! neutrons were emitted with different energies, multiple @@ -600,7 +600,7 @@ contains if (material_xs % absorption == ZERO) cycle SCORE_LOOP ! Set the delayedgroup filter index - dg_filter = t % delayedgroup_filter + dg_filter = t % delayedgroup_filter() if (t % estimator == ESTIMATOR_ANALOG) then if (survival_biasing) then @@ -1509,7 +1509,7 @@ contains if (t % estimator == ESTIMATOR_ANALOG) then if (survival_biasing .or. p % fission) then - if (t % energyout_filter > 0) then + if (t % energyout_filter() > 0) then ! Normally, we only need to make contributions to one scoring ! bin. However, in the case of fission, since multiple fission ! neutrons were emitted with different energies, multiple @@ -1563,7 +1563,7 @@ contains if (t % estimator == ESTIMATOR_ANALOG) then if (survival_biasing .or. p % fission) then - if (t % energyout_filter > 0) then + if (t % energyout_filter() > 0) then ! Normally, we only need to make contributions to one scoring ! bin. However, in the case of fission, since multiple fission ! neutrons were emitted with different energies, multiple @@ -1617,11 +1617,11 @@ contains case (SCORE_DELAYED_NU_FISSION) ! Set the delayedgroup filter index and the number of delayed group bins - dg_filter = t % delayedgroup_filter + dg_filter = t % delayedgroup_filter() if (t % estimator == ESTIMATOR_ANALOG) then if (survival_biasing .or. p % fission) then - if (t % energyout_filter > 0) then + if (t % energyout_filter() > 0) then ! Normally, we only need to make contributions to one scoring ! bin. However, in the case of fission, since multiple fission ! neutrons were emitted with different energies, multiple @@ -1760,7 +1760,7 @@ contains case (SCORE_DECAY_RATE) ! Set the delayedgroup filter index and the number of delayed group bins - dg_filter = t % delayedgroup_filter + dg_filter = t % delayedgroup_filter() if (t % estimator == ESTIMATOR_ANALOG) then if (survival_biasing) then @@ -2350,7 +2350,7 @@ contains integer :: g_out ! energy group of fission bank site ! save original outgoing energy bin and score index - i = t % filter(t % energyout_filter) + i = t % filter(t % energyout_filter()) i_bin = filter_matches(i) % i_bin bin_energyout = filter_matches(i) % bins_data(i_bin) @@ -2429,7 +2429,7 @@ contains else if (score_bin == SCORE_DELAYED_NU_FISSION .and. g /= 0) then ! Get the index of delayed group filter - j = t % delayedgroup_filter + j = t % delayedgroup_filter() ! if the delayed group filter is present, tally to corresponding ! delayed group bin if it exists @@ -2519,7 +2519,7 @@ contains integer :: filter_index ! index for matching filter bin combination ! save original delayed group bin - i_filt = t % filter(t % delayedgroup_filter) + i_filt = t % filter(t % delayedgroup_filter()) i_bin = filter_matches(i_filt) % i_bin bin_original = filter_matches(i_filt) % bins_data(i_bin) call filter_matches(i_filt) % bins_set_data(i_bin, d_bin) diff --git a/src/tallies/tally.cpp b/src/tallies/tally.cpp index 1d816bdac..ada469ee5 100644 --- a/src/tallies/tally.cpp +++ b/src/tallies/tally.cpp @@ -4,6 +4,10 @@ #include "openmc/constants.h" #include "openmc/error.h" #include "openmc/tallies/filter.h" +#include "openmc/tallies/filter_energy.h" +#include "openmc/tallies/filter_delayedgroup.h" +#include "openmc/tallies/filter_surface.h" +#include "openmc/tallies/filter_mesh.h" #include "openmc/message_passing.h" #include "xtensor/xadapt.hpp" @@ -41,6 +45,31 @@ Tally::set_filters(const int32_t filter_indices[], int n) // Copy in the given filter indices. filters_.assign(filter_indices, filter_indices + n); + for (int i = 0; i < n; ++i) { + auto i_filt = filters_[i]; + //TODO: off-by-one + if (i_filt < 1 || i_filt > model::tally_filters.size()) + throw std::out_of_range("Index in tally filter array out of bounds."); + + //TODO: off-by-one + const auto* filt = model::tally_filters[i_filt-1].get(); + + //TODO: off-by-one on each index + if (dynamic_cast(filt)) { + if (dynamic_cast(filt)) { + energyout_filter_ = i + 1; + } else { + energyin_filter_ = i + 1; + } + } else if (dynamic_cast(filt)) { + delayedgroup_filter_ = i + 1; + } else if (dynamic_cast(filt)) { + surface_filter_ = i + 1; + } else if (dynamic_cast(filt)) { + mesh_filter_ = i + 1; + } + } + // Set the strides. Filters are traversed in reverse so that the last filter // has the shortest stride in memory and the first filter has the longest // stride. @@ -216,6 +245,21 @@ extern "C" { int32_t tally_get_n_filter_bins_c(Tally* tally) {return tally->n_filter_bins();} + + int tally_get_energyin_filter_c(Tally* tally) + {return tally->energyin_filter_;} + + int tally_get_energyout_filter_c(Tally* tally) + {return tally->energyout_filter_;} + + int tally_get_delayedgroup_filter_c(Tally* tally) + {return tally->delayedgroup_filter_;} + + int tally_get_surface_filter_c(Tally* tally) + {return tally->surface_filter_;} + + int tally_get_mesh_filter_c(Tally* tally) + {return tally->mesh_filter_;} } } // namespace openmc diff --git a/src/tallies/tally_header.F90 b/src/tallies/tally_header.F90 index 173e4959b..f1ae558f5 100644 --- a/src/tallies/tally_header.F90 +++ b/src/tallies/tally_header.F90 @@ -59,15 +59,6 @@ module tally_header logical :: active = .false. logical :: depletion_rx = .false. ! has depletion reactions, e.g. (n,2n) - ! We need to have quick access to some filters. The following gives indices - ! for various filters that could be in the tally or -1 if they are not - ! present. - integer :: energyin_filter = -1 - integer :: energyout_filter = -1 - integer :: delayedgroup_filter = -1 - integer :: surface_filter = -1 - integer :: mesh_filter = -1 - ! Individual nuclides to tally integer :: n_nuclide_bins = 0 integer, allocatable :: nuclide_bins(:) @@ -105,6 +96,11 @@ module tally_header procedure :: filter => tally_get_filter procedure :: stride => tally_get_stride procedure :: n_filter_bins => tally_get_n_filter_bins + procedure :: energyin_filter => tally_get_energyin_filter + procedure :: energyout_filter => tally_get_energyout_filter + procedure :: delayedgroup_filter => tally_get_delayedgroup_filter + procedure :: surface_filter => tally_get_surface_filter + procedure :: mesh_filter => tally_get_mesh_filter end type TallyObject type, public :: TallyContainer @@ -312,22 +308,9 @@ contains ! Set the filter index in the tally find_filter array select type (filt => filters(k) % obj) - type is (SurfaceFilter) - this % surface_filter = i - - type is (MeshFilter) - this % mesh_filter = i - - type is (EnergyFilter) - this % energyin_filter = i - type is (EnergyoutFilter) - this % energyout_filter = i this % estimator = ESTIMATOR_ANALOG - type is (DelayedGroupFilter) - this % delayedgroup_filter = i - type is (LegendreFilter) this % estimator = ESTIMATOR_ANALOG @@ -410,6 +393,71 @@ contains n_filter_bins = tally_get_n_filter_bins_c(this % ptr) end function + function tally_get_energyin_filter(this) result(filt) + class(TallyObject) :: this + integer(C_INT) :: filt + interface + function tally_get_energyin_filter_c(tally) result(filt) bind(C) + import C_PTR, C_INT + type(C_PTR), value :: tally + integer(C_INT) :: filt + end function + end interface + filt = tally_get_energyin_filter_c(this % ptr) + end function + + function tally_get_energyout_filter(this) result(filt) + class(TallyObject) :: this + integer(C_INT) :: filt + interface + function tally_get_energyout_filter_c(tally) result(filt) bind(C) + import C_PTR, C_INT + type(C_PTR), value :: tally + integer(C_INT) :: filt + end function + end interface + filt = tally_get_energyout_filter_c(this % ptr) + end function + + function tally_get_delayedgroup_filter(this) result(filt) + class(TallyObject) :: this + integer(C_INT) :: filt + interface + function tally_get_delayedgroup_filter_c(tally) result(filt) bind(C) + import C_PTR, C_INT + type(C_PTR), value :: tally + integer(C_INT) :: filt + end function + end interface + filt = tally_get_delayedgroup_filter_c(this % ptr) + end function + + function tally_get_surface_filter(this) result(filt) + class(TallyObject) :: this + integer(C_INT) :: filt + interface + function tally_get_surface_filter_c(tally) result(filt) bind(C) + import C_PTR, C_INT + type(C_PTR), value :: tally + integer(C_INT) :: filt + end function + end interface + filt = tally_get_surface_filter_c(this % ptr) + end function + + function tally_get_mesh_filter(this) result(filt) + class(TallyObject) :: this + integer(C_INT) :: filt + interface + function tally_get_mesh_filter_c(tally) result(filt) bind(C) + import C_PTR, C_INT + type(C_PTR), value :: tally + integer(C_INT) :: filt + end function + end interface + filt = tally_get_mesh_filter_c(this % ptr) + end function + !=============================================================================== ! CONFIGURE_TALLIES initializes several data structures related to tallies. This ! is called after the basic tally data has already been read from the diff --git a/src/tallies/trigger.F90 b/src/tallies/trigger.F90 index d792b7e98..31281a64f 100644 --- a/src/tallies/trigger.F90 +++ b/src/tallies/trigger.F90 @@ -260,8 +260,8 @@ contains type(RegularMesh) :: m ! surface current mesh ! Get pointer to mesh - i_filter_mesh = t % mesh_filter - i_filter_surf = t % surface_filter + i_filter_mesh = t % mesh_filter() + i_filter_surf = t % surface_filter() select type(filt => filters(i_filter_mesh) % obj) type is (MeshFilter) m = meshes(filt % mesh()) @@ -274,7 +274,7 @@ contains end do ! determine how many energyin bins there are - i_filter_ein = t % energyin_filter + i_filter_ein = t % energyin_filter() if (i_filter_ein > 0) then print_ebin = .true. n = filters(t % filter(i_filter_ein)) % obj % n_bins