From ea0a91b0e6c88eacdaca7978f6a944804ec998cb Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Sun, 14 Oct 2018 13:28:03 -0400 Subject: [PATCH] Move SurfaceFilter to C++ --- include/openmc/tallies/tally_filter_surface.h | 85 +++++++++++++++++ src/input_xml.F90 | 5 - src/surface_header.F90 | 5 - src/tallies/tally_filter.cpp | 3 + src/tallies/tally_filter_surface.F90 | 95 +------------------ 5 files changed, 89 insertions(+), 104 deletions(-) create mode 100644 include/openmc/tallies/tally_filter_surface.h diff --git a/include/openmc/tallies/tally_filter_surface.h b/include/openmc/tallies/tally_filter_surface.h new file mode 100644 index 0000000000..a80fcaf129 --- /dev/null +++ b/include/openmc/tallies/tally_filter_surface.h @@ -0,0 +1,85 @@ +#ifndef OPENMC_TALLY_FILTER_SURFACE_H +#define OPENMC_TALLY_FILTER_SURFACE_H + +#include +#include +#include +#include + +#include "openmc/error.h" +#include "openmc/surface.h" +#include "openmc/tallies/tally_filter.h" + + +namespace openmc { + +class SurfaceFilter : public TallyFilter +{ +public: + virtual std::string type() const override {return "surface";} + + virtual ~SurfaceFilter() override = default; + + virtual void + from_xml(pugi::xml_node node) override + { + surfaces_ = get_node_array(node, "bins"); + n_bins_ = surfaces_.size(); + } + + virtual void + initialize() override + { + for (auto& s : surfaces_) { + auto search = surface_map.find(s); + if (search != surface_map.end()) { + s = search->second; + } else { + std::stringstream err_msg; + err_msg << "Could not find surface " << s + << " specified on tally filter."; + fatal_error(err_msg); + } + } + + for (int i = 0; i < surfaces_.size(); i++) { + map_[surfaces_[i]] = i; + } + } + + virtual void + get_all_bins(Particle* p, int estimator, TallyFilterMatch& match) + const override + { + auto search = map_.find(std::abs(p->surface)-1); + if (search != map_.end()) { + match.bins.push_back(search->second + 1); + if (p->surface < 0) { + match.weights.push_back(-1); + } else { + match.weights.push_back(1); + } + } + } + + virtual void + to_statepoint(hid_t filter_group) const override + { + TallyFilter::to_statepoint(filter_group); + std::vector surface_ids; + for (auto c : surfaces_) surface_ids.push_back(surfaces[c]->id_); + write_dataset(filter_group, "bins", surface_ids); + } + + virtual std::string + text_label(int bin) const override + { + return "Surface " + std::to_string(surfaces[surfaces_[bin-1]]->id_); + } + + std::vector surfaces_; + std::unordered_map map_; +}; + +} // namespace openmc +#endif // OPENMC_TALLY_FILTER_SURFACE_H diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 8a3682c055..d4fd5dfa45 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -468,9 +468,6 @@ contains surfaces(i) % ptr = surface_pointer(i - 1); if (surfaces(i) % bc() /= BC_TRANSMIT) boundary_exists = .true. - - ! Add surface to dictionary - call surface_dict % set(surfaces(i) % id(), i) end do ! Check to make sure a boundary condition was applied to at least one @@ -628,8 +625,6 @@ contains do i = 1, n_surfaces surfaces(i) % ptr = surface_pointer(i - 1); - ! Add surface to dictionary - call surface_dict % set(surfaces(i) % id(), i) end do end subroutine allocate_surfaces diff --git a/src/surface_header.F90 b/src/surface_header.F90 index 9d7e114613..a2ef9c0563 100644 --- a/src/surface_header.F90 +++ b/src/surface_header.F90 @@ -2,7 +2,6 @@ module surface_header use, intrinsic :: ISO_C_BINDING - use dict_header, only: DictIntInt use hdf5_interface implicit none @@ -83,9 +82,6 @@ module surface_header type(Surface), allocatable, target :: surfaces(:) - ! Dictionary that maps user IDs to indices in 'surfaces' - type(DictIntInt) :: surface_dict - contains pure function surface_id(this) result(id) @@ -128,7 +124,6 @@ contains subroutine free_memory_surfaces() if (allocated(surfaces)) deallocate(surfaces) - call surface_dict % clear() call free_memory_surfaces_c() end subroutine free_memory_surfaces diff --git a/src/tallies/tally_filter.cpp b/src/tallies/tally_filter.cpp index c5d5600a1c..d66505d470 100644 --- a/src/tallies/tally_filter.cpp +++ b/src/tallies/tally_filter.cpp @@ -13,6 +13,7 @@ #include "openmc/tallies/tally_filter_meshsurface.h" #include "openmc/tallies/tally_filter_mu.h" #include "openmc/tallies/tally_filter_polar.h" +#include "openmc/tallies/tally_filter_surface.h" #include "openmc/tallies/tally_filter_universe.h" @@ -105,6 +106,8 @@ extern "C" { tally_filters.push_back(new MuFilter()); } else if (type_ == "polar") { tally_filters.push_back(new PolarFilter()); + } else if (type_ == "surface") { + tally_filters.push_back(new SurfaceFilter()); } else if (type_ == "universe") { tally_filters.push_back(new UniverseFilter()); } else { diff --git a/src/tallies/tally_filter_surface.F90 b/src/tallies/tally_filter_surface.F90 index d545d46721..1dcb56022b 100644 --- a/src/tallies/tally_filter_surface.F90 +++ b/src/tallies/tally_filter_surface.F90 @@ -1,109 +1,16 @@ module tally_filter_surface - use, intrinsic :: ISO_C_BINDING - - use constants, only: ONE, MAX_LINE_LEN - use dict_header, only: EMPTY - use error, only: fatal_error - use hdf5_interface - use surface_header - use particle_header, only: Particle - use string, only: to_str use tally_filter_header - use xml_interface implicit none - private !=============================================================================== ! SURFACEFILTER specifies which surface particles are crossing !=============================================================================== - type, public, extends(TallyFilter) :: SurfaceFilter - integer, allocatable :: surfaces(:) - + type, extends(CppTallyFilter) :: SurfaceFilter ! True if this filter is used for surface currents logical :: current = .false. - contains - procedure :: from_xml - procedure :: get_all_bins => get_all_bins_surface - procedure :: to_statepoint => to_statepoint_surface - procedure :: text_label => text_label_surface - procedure :: initialize => initialize_surface end type SurfaceFilter -contains - - subroutine from_xml(this, node) - class(SurfaceFilter), intent(inout) :: this - type(XMLNode), intent(in) :: node - - integer :: n - - n = node_word_count(node, "bins") - - ! Allocate and store bins - this % n_bins = n - allocate(this % surfaces(n)) - call get_node_array(node, "bins", this % surfaces) - end subroutine from_xml - - subroutine get_all_bins_surface(this, p, estimator, match) - class(SurfaceFilter), intent(in) :: this - type(Particle), intent(in) :: p - integer, intent(in) :: estimator - type(TallyFilterMatch), intent(inout) :: match - - integer :: i - - do i = 1, this % n_bins - if (abs(p % surface) == this % surfaces(i)) then - call match % bins_push_back(i) - if (p % surface < 0) then - call match % weights_push_back(-ONE) - else - call match % weights_push_back(ONE) - end if - exit - end if - end do - - end subroutine get_all_bins_surface - - subroutine to_statepoint_surface(this, filter_group) - class(SurfaceFilter), intent(in) :: this - integer(HID_T), intent(in) :: filter_group - - call write_dataset(filter_group, "type", "surface") - call write_dataset(filter_group, "n_bins", this % n_bins) - call write_dataset(filter_group, "bins", this % surfaces) - end subroutine to_statepoint_surface - - subroutine initialize_surface(this) - class(SurfaceFilter), intent(inout) :: this - - integer :: i, id - integer :: val - - ! Convert ids to indices. - do i = 1, this % n_bins - id = this % surfaces(i) - val = surface_dict % get(id) - if (val /= EMPTY) then - this % surfaces(i) = val - else - call fatal_error("Could not find surface " // trim(to_str(id)) & - &// " specified on tally filter.") - end if - end do - end subroutine initialize_surface - - function text_label_surface(this, bin) result(label) - class(SurfaceFilter), intent(in) :: this - integer, intent(in) :: bin - character(MAX_LINE_LEN) :: label - - label = "Surface " // to_str(surfaces(this % surfaces(bin)) % id()) - end function text_label_surface - end module tally_filter_surface