crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-27 21:55:41 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Move DistribcellFilter to C++

Browse source
This commit is contained in:
Sterling Harper 2018-10-11 22:27:13 -04:00
parent bf025edd58
commit f2a5e10c46
9 changed files with 121 additions and 317 deletions
Download patch file Download diff file Expand all files Collapse all files

27
include/openmc/geometry_aux.h
View file

@ -5,6 +5,7 @@
#define OPENMC_GEOMETRY_AUX_H
#include <cstdint>
#include <string>
namespace openmc {
@ -41,7 +42,7 @@ extern "C" void neighbor_lists();
//! Populate all data structures needed for distribcells.
//==============================================================================
extern "C" void prepare_distribcell(int32_t* filter_cell_list, int n);
extern "C" void prepare_distribcell();
//==============================================================================
//! Recursively search through the geometry and count cell instances.
@ -65,21 +66,6 @@ extern "C" void count_cell_instances(int32_t univ_indx);
extern "C" int
count_universe_instances(int32_t search_univ, int32_t target_univ_id);
//==============================================================================
//! Find the length necessary for a string to contain a distribcell path.
//! @param target_cell The index of the Cell in the global Cell array.
//! @param map The index of the distribcell mapping corresponding to the target
//! cell.
//! @param target_offset An instance number for a distributed cell.
//! @param root_univ The index of the root Universe in the global Universe
//! array.
//! @return The size of a character array needed to fit the distribcell path.
//==============================================================================
extern "C" int
distribcell_path_len(int32_t target_cell, int32_t map, int32_t target_offset,
int32_t root_univ);
//==============================================================================
//! Build a character array representing the path to a distribcell instance.
//! @param target_cell The index of the Cell in the global Cell array.
@ -88,13 +74,12 @@ distribcell_path_len(int32_t target_cell, int32_t map, int32_t target_offset,
//! @param target_offset An instance number for a distributed cell.
//! @param root_univ The index of the root Universe in the global Universe
//! array.
//! @param[out] path The unique traversal through the geometry tree that leads
//! to the desired instance of the target cell.
//! @return The unique traversal through the geometry tree that leads to the
//! desired instance of the target cell.
//==============================================================================
extern "C" void
distribcell_path(int32_t target_cell, int32_t map, int32_t target_offset,
int32_t root_univ, char *path);
std::string
distribcell_path(int32_t target_cell, int32_t map, int32_t target_offset);
//==============================================================================
//! Determine the maximum number of nested coordinate levels in the geometry.

93
include/openmc/tallies/tally_filter_distribcell.h Normal file
View file

@ -0,0 +1,93 @@
#ifndef OPENMC_TALLY_FILTER_DISTRIBCELL_H
#define OPENMC_TALLY_FILTER_DISTRIBCELL_H
#include <string>
#include "openmc/cell.h"
#include "openmc/error.h"
#include "openmc/geometry_aux.h" // For distribcell_path
#include "openmc/lattice.h"
#include "openmc/tallies/tally_filter.h"
namespace openmc {
class DistribcellFilter : public TallyFilter
{
public:
virtual ~DistribcellFilter() override = default;
virtual void
from_xml(pugi::xml_node node) override
{
auto cells = get_node_array<int32_t>(node, "bins");
if (cells.size() != 1) {
fatal_error("Only one cell can be specified per distribcell filter.");
}
cell_ = cells[0];
}
virtual void
initialize() override
{
auto search = cell_map.find(cell_);
if (search != cell_map.end()) {
cell_ = search->second;
n_bins_ = cells[cell_]->n_instances_;
} else {
std::stringstream err_msg;
err_msg << "Could not find cell " << cell_
<< " specified on tally filter.";
fatal_error(err_msg);
}
}
virtual void
get_all_bins(Particle* p, int estimator, TallyFilterMatch& match)
const override
{
int offset = 0;
auto distribcell_index = cells[cell_]->distribcell_index_;
for (int i = 0; i < p->n_coord; i++) {
auto& c {*cells[p->coord[i].cell]};
if (c.type_ == FILL_UNIVERSE) {
offset += c.offset_[distribcell_index];
} else if (c.type_ == FILL_LATTICE) {
auto& lat {*lattices[p->coord[i+1].lattice-1]};
int i_xyz[3] {p->coord[i+1].lattice_x,
p->coord[i+1].lattice_y,
p->coord[i+1].lattice_z};
if (lat.are_valid_indices(i_xyz)) {
offset += lat.offset(distribcell_index, i_xyz);
}
}
if (cell_ == p->coord[i].cell) {
match.bins.push_back(offset + 1);
match.weights.push_back(1);
return;
}
}
}
virtual void
to_statepoint(hid_t filter_group) const override
{
write_dataset(filter_group, "type", "distribcell");
write_dataset(filter_group, "n_bins", n_bins_);
write_dataset(filter_group, "bins", cells[cell_]->id_);
}
virtual std::string
text_label(int bin) const override
{
auto map = cells[cell_]->distribcell_index_;
auto path = distribcell_path(cell_, map, bin-1);
return "Distributed Cell " + path;
}
int32_t cell_;
};
} // namespace openmc
#endif // OPENMC_TALLY_FILTER_DISTRIBCELL_H

6
src/cell.cpp
View file

@ -803,10 +803,6 @@ extern "C" {
int32_t cell_fill(Cell* c) {return c->fill_;}
int32_t cell_n_instances(Cell* c) {return c->n_instances_;}
int cell_distribcell_index(Cell* c) {return c->distribcell_index_;}
int cell_material_size(Cell* c) {return c->material_.size();}
//TODO: off-by-one
@ -821,8 +817,6 @@ extern "C" {
double cell_sqrtkT(Cell* c, int i) {return c->sqrtkT_[i];}
int32_t cell_offset(Cell* c, int map) {return c->offset_[map];}
void extend_cells_c(int32_t n)
{
cells.reserve(cells.size() + n);

40
src/geometry_aux.cpp
View file

@ -12,6 +12,8 @@
#include "openmc/material.h"
#include "openmc/settings.h"
#include "openmc/surface.h"
#include "openmc/tallies/tally_filter.h"
#include "openmc/tallies/tally_filter_distribcell.h"
namespace openmc {
@ -181,12 +183,16 @@ neighbor_lists()
//==============================================================================
void
prepare_distribcell(int32_t* filter_cell_list, int n)
prepare_distribcell()
{
// Read the list of cells contained in distribcell filters from Fortran.
// Find all cells listed in a DistribcellFilter.
std::unordered_set<int32_t> distribcells;
for (int i = 0; i < n; i++) {
distribcells.insert(filter_cell_list[i]);
for (auto* filt : tally_filters) {
//TODO: replace this cast
auto* distrib_filt = dynamic_cast<DistribcellFilter*>(filt);
if (distrib_filt) {
distribcells.insert(distrib_filt->cell_);
}
}
// Find all cells with distributed materials or temperatures. Make sure that
@ -393,29 +399,11 @@ distribcell_path_inner(int32_t target_cell, int32_t map, int32_t target_offset,
}
}
//==============================================================================
int
distribcell_path_len(int32_t target_cell, int32_t map, int32_t target_offset,
int32_t root_univ)
std::string
distribcell_path(int32_t target_cell, int32_t map, int32_t target_offset)
{
Universe& root = *universes[root_univ];
std::string path_ {distribcell_path_inner(target_cell, map, target_offset,
root, 0)};
return path_.size() + 1;
}
//==============================================================================
void
distribcell_path(int32_t target_cell, int32_t map, int32_t target_offset,
int32_t root_univ, char* path)
{
Universe& root = *universes[root_univ];
std::string path_ {distribcell_path_inner(target_cell, map, target_offset,
root, 0)};
path_.copy(path, path_.size());
path[path_.size()] = '\0';
auto& root_univ = *universes[openmc_root_universe];
return distribcell_path_inner(target_cell, map, target_offset, root_univ, 0);
}
//==============================================================================

78
src/geometry_header.F90
View file

@ -59,20 +59,6 @@ module geometry_header
integer(C_INT32_T) :: fill
end function cell_fill_c
function cell_n_instances_c(cell_ptr) bind(C, name='cell_n_instances') &
result(n_instances)
import C_PTR, C_INT32_T
type(C_PTR), intent(in), value :: cell_ptr
integer(C_INT32_T) :: n_instances
end function cell_n_instances_c
function cell_distribcell_index_c(cell_ptr) &
bind(C, name='cell_distribcell_index') result(distribcell_index)
import C_PTR, C_INT
type(C_PTR), intent(in), value :: cell_ptr
integer(C_INT) :: distribcell_index
end function cell_distribcell_index_c
function cell_material_size_c(cell_ptr) bind(C, name='cell_material_size') &
result(n)
import C_PTR, C_INT
@ -103,14 +89,6 @@ module geometry_header
real(C_DOUBLE) :: sqrtkT
end function cell_sqrtkT_c
function cell_offset_c(cell_ptr, map) bind(C, name="cell_offset") &
result(offset)
import C_PTR, C_INT, C_INT32_T
type(C_PTR), intent(in), value :: cell_ptr
integer(C_INT), intent(in), value :: map
integer(C_INT32_T) :: offset
end function cell_offset_c
function lattice_pointer(lat_ind) bind(C) result(ptr)
import C_PTR, C_INT32_T
integer(C_INT32_T), intent(in), value :: lat_ind
@ -123,23 +101,6 @@ module geometry_header
integer(C_INT32_T) :: id
end function lattice_id_c
function lattice_are_valid_indices_c(lat_ptr, i_xyz) &
bind(C, name='lattice_are_valid_indices') result (is_valid)
import C_PTR, C_INT, C_BOOL
type(C_PTR), intent(in), value :: lat_ptr
integer(C_INT), intent(in) :: i_xyz(3)
logical(C_BOOL) :: is_valid
end function lattice_are_valid_indices_c
function lattice_offset_c(lat_ptr, map, i_xyz) &
bind(C, name='lattice_offset') result(offset)
import C_PTR, C_INT, C_INT32_T
type(C_PTR), intent(in), value :: lat_ptr
integer(C_INT), intent(in), value :: map
integer(C_INT), intent(in) :: i_xyz(3)
integer(C_INT32_T) :: offset
end function lattice_offset_c
subroutine extend_cells_c(n) bind(C)
import C_INT32_t
integer(C_INT32_T), intent(in), value :: n
@ -154,8 +115,6 @@ module geometry_header
type(C_PTR) :: ptr
contains
procedure :: id => lattice_id
procedure :: are_valid_indices => lattice_are_valid_indices
procedure :: offset => lattice_offset
end type Lattice
!===============================================================================
@ -176,13 +135,10 @@ module geometry_header
procedure :: type => cell_type
procedure :: universe => cell_universe
procedure :: fill => cell_fill
procedure :: n_instances => cell_n_instances
procedure :: distribcell_index => cell_distribcell_index
procedure :: material_size => cell_material_size
procedure :: material => cell_material
procedure :: sqrtkT_size => cell_sqrtkT_size
procedure :: sqrtkT => cell_sqrtkT
procedure :: offset => cell_offset
end type Cell
@ -208,21 +164,6 @@ contains
id = lattice_id_c(this % ptr)
end function lattice_id
function lattice_are_valid_indices(this, i_xyz) result (is_valid)
class(Lattice), intent(in) :: this
integer(C_INT), intent(in) :: i_xyz(3)
logical(C_BOOL) :: is_valid
is_valid = lattice_are_valid_indices_c(this % ptr, i_xyz)
end function lattice_are_valid_indices
function lattice_offset(this, map, i_xyz) result(offset)
class(Lattice), intent(in) :: this
integer(C_INT), intent(in) :: map
integer(C_INT), intent(in) :: i_xyz(3)
integer(C_INT32_T) :: offset
offset = lattice_offset_c(this % ptr, map, i_xyz)
end function lattice_offset
!===============================================================================
function cell_id(this) result(id)
@ -255,18 +196,6 @@ contains
fill = cell_fill_c(this % ptr)
end function cell_fill
function cell_n_instances(this) result(n_instances)
class(Cell), intent(in) :: this
integer(C_INT32_T) :: n_instances
n_instances = cell_n_instances_c(this % ptr)
end function cell_n_instances
function cell_distribcell_index(this) result(distribcell_index)
class(Cell), intent(in) :: this
integer(C_INT) :: distribcell_index
distribcell_index = cell_distribcell_index_c(this % ptr)
end function cell_distribcell_index
function cell_material_size(this) result(n)
class(Cell), intent(in) :: this
integer(C_INT) :: n
@ -293,13 +222,6 @@ contains
sqrtkT = cell_sqrtkT_c(this % ptr, i)
end function cell_sqrtkT
function cell_offset(this, map) result(offset)
class(Cell), intent(in) :: this
integer(C_INT), intent(in) :: map
integer(C_INT32_T) :: offset
offset = cell_offset_c(this % ptr, map)
end function cell_offset
!===============================================================================
! GET_TEMPERATURES returns a list of temperatures that each nuclide/S(a,b) table
! appears at in the model. Later, this list is used to determine the actual

37
src/input_xml.F90
View file

@ -59,12 +59,8 @@ module input_xml
integer(C_INT32_T), intent(in), value :: univ_indx
end subroutine count_cell_instances
subroutine prepare_distribcell_c(cell_list, n) &
bind(C, name="prepare_distribcell")
import C_INT32_T, C_INT
integer(C_INT), intent(in), value :: n
integer(C_INT32_T), intent(in) :: cell_list(n)
end subroutine prepare_distribcell_c
subroutine prepare_distribcell() bind(C)
end subroutine prepare_distribcell
subroutine read_surfaces(node_ptr) bind(C)
import C_PTR
@ -3148,33 +3144,4 @@ contains
end subroutine read_multipole_data
!===============================================================================
! PREPARE_DISTRIBCELL initializes any distribcell filters present and sets the
! offsets for distribcells
!===============================================================================
subroutine prepare_distribcell()
integer :: i, j
type(SetInt) :: cell_list ! distribcells to track
integer(C_INT32_T), allocatable :: cell_list_c(:)
! Find all cells listed in a distribcell filter.
do i = 1, n_tallies
do j = 1, size(tallies(i) % obj % filter)
select type(filt => filters(tallies(i) % obj % filter(j)) % obj)
type is (DistribcellFilter)
call cell_list % add(filt % cell)
end select
end do
end do
allocate(cell_list_c(cell_list % size()))
do i = 1, cell_list % size()
cell_list_c(i) = cell_list % get_item(i) - 1
end do
call prepare_distribcell_c(cell_list_c, cell_list % size())
end subroutine prepare_distribcell
end module input_xml

6
src/lattice.cpp
View file

@ -893,12 +893,6 @@ extern "C" {
Lattice* lattice_pointer(int lat_ind) {return lattices[lat_ind];}
int32_t lattice_id(Lattice *lat) {return lat->id_;}
bool lattice_are_valid_indices(Lattice *lat, const int i_xyz[3])
{return lat->are_valid_indices(i_xyz);}
int32_t lattice_offset(Lattice *lat, int map, const int i_xyz[3])
{return lat->offset(map, i_xyz);}
}
} // namespace openmc

3
src/tallies/tally_filter.cpp
View file

@ -6,6 +6,7 @@
#include "openmc/tallies/tally_filter_cell.h"
#include "openmc/tallies/tally_filter_cellborn.h"
#include "openmc/tallies/tally_filter_cellfrom.h"
#include "openmc/tallies/tally_filter_distribcell.h"
#include "openmc/tallies/tally_filter_mesh.h"
#include "openmc/tallies/tally_filter_meshsurface.h"
@ -85,6 +86,8 @@ extern "C" {
tally_filters.push_back(new CellbornFilter());
} else if (type_ == "cellfrom") {
tally_filters.push_back(new CellFromFilter());
} else if (type_ == "distribcell") {
tally_filters.push_back(new DistribcellFilter());
} else if (type_ == "mesh") {
tally_filters.push_back(new MeshFilter());
} else if (type_ == "meshsurface") {

148
src/tallies/tally_filter_distribcell.F90
View file

@ -1,16 +1,6 @@
module tally_filter_distribcell
use, intrinsic :: ISO_C_BINDING
use constants
use dict_header, only: EMPTY
use error
use geometry_header
use hdf5_interface
use particle_header, only: Particle
use string, only: to_str
use tally_filter_header
use xml_interface
implicit none
private
@ -20,150 +10,18 @@ module tally_filter_distribcell
! reside in.
!===============================================================================
type, public, extends(TallyFilter) :: DistribcellFilter
type, public, extends(CppTallyFilter) :: DistribcellFilter
integer :: cell
contains
procedure :: from_xml
procedure :: get_all_bins => get_all_bins_distribcell
procedure :: to_statepoint => to_statepoint_distribcell
procedure :: text_label => text_label_distribcell
procedure :: initialize => initialize_distribcell
end type DistribcellFilter
contains
subroutine from_xml(this, node)
class(DistribcellFilter), intent(inout) :: this
type(XMLNode), intent(in) :: node
integer :: n
n = node_word_count(node, "bins")
if (n /= 1) call fatal_error("Only one cell can be &
&specified per distribcell filter.")
! Store bins
call get_node_value(node, "bins", this % cell)
end subroutine from_xml
subroutine get_all_bins_distribcell(this, p, estimator, match)
class(DistribcellFilter), intent(in) :: this
type(Particle), intent(in) :: p
integer, intent(in) :: estimator
type(TallyFilterMatch), intent(inout) :: match
integer :: distribcell_index, offset, i
distribcell_index = cells(this % cell) % distribcell_index()
offset = 0
do i = 1, p % n_coord
if (cells(p % coord(i) % cell + 1) % type() == FILL_UNIVERSE) then
offset = offset + cells(p % coord(i) % cell + 1) &
% offset(distribcell_index)
elseif (cells(p % coord(i) % cell + 1) % type() == FILL_LATTICE) then
if (lattices(p % coord(i + 1) % lattice) % are_valid_indices([&
p % coord(i + 1) % lattice_x, &
p % coord(i + 1) % lattice_y, &
p % coord(i + 1) % lattice_z])) then
offset = offset + lattices(p % coord(i + 1) % lattice) &
% offset(distribcell_index, &
[p % coord(i + 1) % lattice_x, &
p % coord(i + 1) % lattice_y, &
p % coord(i + 1) % lattice_z])
end if
end if
if (this % cell == p % coord(i) % cell + 1) then
call match % bins_push_back(offset + 1)
call match % weights_push_back(ONE)
return
end if
end do
end subroutine get_all_bins_distribcell
subroutine to_statepoint_distribcell(this, filter_group)
class(DistribcellFilter), intent(in) :: this
integer(HID_T), intent(in) :: filter_group
call write_dataset(filter_group, "type", "distribcell")
call write_dataset(filter_group, "n_bins", this % n_bins)
call write_dataset(filter_group, "bins", cells(this % cell) % id())
end subroutine to_statepoint_distribcell
subroutine initialize_distribcell(this)
class(DistribcellFilter), intent(inout) :: this
integer :: id
integer :: val
! Convert id to index.
id = this % cell
val = cell_dict % get(id)
if (val /= EMPTY) then
this % cell = val
this % n_bins = cells(this % cell) % n_instances()
else
call fatal_error("Could not find cell " // trim(to_str(id)) &
&// " specified on tally filter.")
end if
call this % initialize_cpp_inner()
this % n_bins = this % n_bins_cpp()
end subroutine initialize_distribcell
function text_label_distribcell(this, bin) result(label)
class(DistribcellFilter), intent(in) :: this
integer, intent(in) :: bin
character(MAX_LINE_LEN) :: label
label = ''
call find_offset(this % cell, bin-1, label)
label = "Distributed Cell " // label
end function text_label_distribcell
!===============================================================================
! FIND_OFFSET (for distribcell) uses a given map number, a target cell ID, and
! a target offset to build a string which is the path from the base universe to
! the target cell with the given offset
!===============================================================================
subroutine find_offset(i_cell, target_offset, path)
integer, intent(in) :: i_cell ! The target cell index
integer, intent(in) :: target_offset ! Target offset
character(*), intent(inout) :: path ! Path to offset
integer :: map ! Index in maps vector
integer :: i ! Index over cells
integer(C_INT) :: path_len
character(kind=C_CHAR), allocatable, target :: path_c(:)
interface
function distribcell_path_len(target_cell, map, target_offset, root_univ)&
bind(C) result(len)
import C_INT32_T, C_INT
integer(C_INT32_T), intent(in), value :: target_cell, map, &
target_offset, root_univ
integer(C_INT) :: len
end function distribcell_path_len
subroutine distribcell_path(target_cell, map, target_offset, root_univ, &
path) bind(C)
import C_INT32_T, C_CHAR
integer(C_INT32_T), intent(in), value :: target_cell, map, &
target_offset, root_univ
character(kind=C_CHAR), intent(out) :: path(*)
end subroutine distribcell_path
end interface
! Get the distribcell index for this cell
map = cells(i_cell) % distribcell_index()
path_len = distribcell_path_len(i_cell-1, map, target_offset, &
root_universe)
allocate(path_c(path_len))
call distribcell_path(i_cell-1, map, target_offset, root_universe, &
path_c)
do i = 1, min(path_len, MAX_LINE_LEN)
if (path_c(i) == C_NULL_CHAR) exit
path(i:i) = path_c(i)
end do
end subroutine find_offset
end module tally_filter_distribcell