crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-26 05:05:30 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Merge pull request #1120 from paulromano/namespaced-globals

Browse source 
Put global variables in a sub-namespace
This commit is contained in:
Adam Nelson 2018-11-19 12:48:39 -05:00 • committed by GitHub
commit fc4a838840
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
63 changed files with 857 additions and 695 deletions
Download patch file Download diff file Expand all files Collapse all files

4
docs/source/capi/index.rst
View file

@ -465,13 +465,15 @@ Functions
:return: Return status (negative if an error occurred)
:rtype: int
.. c:function:: int openmc_statepoint_write(const char filename[])
.. c:function:: int openmc_statepoint_write(const char filename[], const bool* write_source)
Write a statepoint file
:param filename: Name of file to create. If a null pointer is passed, a
filename is assigned automatically.
:type filename: const char[]
:param write_source: Whether to include the source bank
:type write_source: const bool*
:return: Return status (negative if an error occurs)
:rtype: int

20
docs/source/devguide/styleguide.rst
View file

@ -237,6 +237,26 @@ are written differently by convention (e.g., ``E`` for energy). Data members of
classes (but not structs) additionally have trailing underscores (e.g.,
``a_class_member_``).
The following conventions are used for variables with short names:
- ``d`` stands for "distance"
- ``E`` stands for "energy"
- ``p`` stands for "particle"
- ``r`` stands for "position"
- ``rx`` stands for "reaction"
- ``u`` stands for "direction"
- ``xs`` stands for "cross section"
All classes and non-member functions should be declared within the ``openmc``
namespace. Global variables must be declared in a namespace nested within the
``openmc`` namespace. The following sub-namespaces are in use:
- ``openmc::data``: Fundamental nuclear data (cross sections, multigroup data,
decay constants, etc.)
- ``openmc::model``: Variables related to geometry, materials, and tallies
- ``openmc::settings``: Global settings / options
- ``openmc::simulation``: Variables used only during a simulation
Accessors and mutators (get and set functions) may be named like
variables. These often correspond to actual member variables, but this is not
required. For example, ``int count()`` and ``void set_count(int count)``.

7
include/openmc/capi.h
View file

@ -135,17 +135,10 @@ extern "C" {
// Global variables
extern char openmc_err_msg[256];
extern int32_t n_cells;
extern int32_t n_lattices;
extern int32_t n_materials;
extern int n_nuclides;
extern int32_t n_plots;
extern int32_t n_realizations;
extern int32_t n_sab_tables;
extern int32_t n_sources;
extern int32_t n_surfaces;
extern int32_t n_tallies;
extern int32_t n_universes;
// Variables that are shared by necessity (can be removed from public header
// later)

17
include/openmc/cell.h
View file

@ -39,16 +39,21 @@ constexpr int32_t OP_UNION {std::numeric_limits<int32_t>::max() - 4};
// Global variables
//==============================================================================
class Cell;
class Universe;
namespace model {
extern "C" int32_t n_cells;
class Cell;
extern std::vector<Cell*> cells;
extern std::unordered_map<int32_t, int32_t> cell_map;
class Universe;
extern std::vector<Universe*> universes;
extern std::unordered_map<int32_t, int32_t> universe_map;
} // namespace model
//==============================================================================
//! A geometry primitive that fills all space and contains cells.
//==============================================================================
@ -145,13 +150,13 @@ public:
virtual ~Cell() {}
};
class CSGCell : public Cell
{
public:
CSGCell();
explicit CSGCell(pugi::xml_node cell_node);
bool
@ -163,7 +168,7 @@ public:
void to_hdf5(hid_t group_id) const;
protected:
bool contains_simple(Position r, Direction u, int32_t on_surface) const;
bool contains_complex(Position r, Direction u, int32_t on_surface) const;
@ -183,6 +188,6 @@ public:
};
#endif
} // namespace openmc
#endif // OPENMC_CELL_H

14
include/openmc/dagmc.h
View file

@ -10,12 +10,24 @@
namespace openmc {
//==============================================================================
// Global variables
//==============================================================================
namespace model {
extern moab::DagMC* DAG;
} // namespace model
//==============================================================================
// Non-member functions
//==============================================================================
extern "C" void load_dagmc_geometry();
extern "C" void free_memory_dagmc();
}
} // namespace openmc
#endif // DAGMC

4
include/openmc/eigenvalue.h
View file

@ -18,6 +18,8 @@ namespace openmc {
// Global variables
//==============================================================================
namespace simulation {
extern double keff_generation; //!< Single-generation k on each processor
extern std::array<double, 2> k_sum; //!< Used to reduce sum and sum_sq
extern std::vector<double> entropy; //!< Shannon entropy at each generation
@ -26,6 +28,8 @@ extern xt::xtensor<double, 1> source_frac; //!< Source fraction for UFS
extern "C" int64_t n_bank;
#pragma omp threadprivate(n_bank)
} // namespace simulation
//==============================================================================
// Non-member functions
//==============================================================================

6
include/openmc/geometry.h
View file

@ -13,10 +13,14 @@ namespace openmc {
// Global variables
//==============================================================================
extern "C" int openmc_root_universe;
namespace model {
extern "C" int root_universe;
extern std::vector<int64_t> overlap_check_count;
} // namespace model
//==============================================================================
//! Check for overlapping cells at a particle's position.
//==============================================================================

6
include/openmc/lattice.h
View file

@ -31,10 +31,14 @@ enum class LatticeType {
//==============================================================================
class Lattice;
extern std::vector<Lattice*> lattices;
namespace model {
extern std::vector<Lattice*> lattices;
extern std::unordered_map<int32_t, int32_t> lattice_map;
} // namespace model
//==============================================================================
//! \class Lattice
//! \brief Abstract type for ordered array of universes.

5
include/openmc/material.h
View file

@ -14,9 +14,14 @@ namespace openmc {
//==============================================================================
class Material;
namespace model {
extern std::vector<Material*> materials;
extern std::unordered_map<int32_t, int32_t> material_map;
} // namespace model
//==============================================================================
//! A substance with constituent nuclides and thermal scattering data
//==============================================================================

22
include/openmc/mesh.h
View file

@ -17,6 +17,19 @@
namespace openmc {
//==============================================================================
// Global variables
//==============================================================================
class RegularMesh;
namespace model {
extern std::vector<std::unique_ptr<RegularMesh>> meshes;
extern std::unordered_map<int32_t, int32_t> mesh_map;
} // namespace model
//==============================================================================
//! Tessellation of n-dimensional Euclidean space by congruent squares or cubes
//==============================================================================
@ -117,15 +130,6 @@ extern "C" void read_meshes(pugi::xml_node* root);
//! \param[in] group HDF5 group
extern "C" void meshes_to_hdf5(hid_t group);
//==============================================================================
// Global variables
//==============================================================================
extern std::vector<std::unique_ptr<RegularMesh>> meshes;
extern std::unordered_map<int32_t, int32_t> mesh_map;
} // namespace openmc
#endif // OPENMC_MESH_H

4
include/openmc/mgxs_interface.h
View file

@ -15,6 +15,8 @@ namespace openmc {
// Global variables
//==============================================================================
namespace data {
extern std::vector<Mgxs> nuclides_MG;
extern std::vector<Mgxs> macro_xs;
extern "C" int num_energy_groups;
@ -22,6 +24,8 @@ extern std::vector<double> energy_bins;
extern std::vector<double> energy_bin_avg;
extern std::vector<double> rev_energy_bins;
} // namespace data
//==============================================================================
// Mgxs data loading interface methods
//==============================================================================

4
include/openmc/nuclide.h
View file

@ -14,11 +14,15 @@ namespace openmc {
// Global variables
//==============================================================================
namespace data {
// Minimum/maximum transport energy for each particle type. Order corresponds to
// that of the ParticleType enum
extern std::array<double, 2> energy_min;
extern std::array<double, 2> energy_max;
} // namespace data
//===============================================================================
//! Cached microscopic cross sections for a particular nuclide at the current
//! energy

13
include/openmc/physics_mg.h
View file

@ -10,34 +10,29 @@
namespace openmc {
//TODO: Remove energy_bin_avg and material_xs parameters when they reside on
//TODO: Remove material_xs parameters when they reside on
// the C-side this should happen after materials, physics, input, and tallies
// are brought over
//! \brief samples particle behavior after a collision event.
//! \param p Particle to operate on
//! \param energy_bin_avg Average energy within each energy bin
//! \param material_xs The cross section cache for the current material
extern "C" void
collision_mg(Particle* p, const double* energy_bin_avg,
const MaterialMacroXS* material_xs);
collision_mg(Particle* p, const MaterialMacroXS* material_xs);
//! \brief samples a reaction type.
//!
//! Note that there is special logic when suvival biasing is turned on since
//! fission and disappearance are treated implicitly.
//! \param p Particle to operate on
//! \param energy_bin_avg Average energy within each energy bin
//! \param material_xs The cross section cache for the current material
void
sample_reaction(Particle* p, const double* energy_bin_avg,
const MaterialMacroXS* material_xs);
sample_reaction(Particle* p, const MaterialMacroXS* material_xs);
//! \brief Samples the scattering event
//! \param p Particle to operate on
//! \param energy_bin_avg Average energy within each energy bin
void
scatter(Particle* p, const double* energy_bin_avg);
scatter(Particle* p);
//! \brief Determines the average total, prompt and delayed neutrons produced
//! from fission and creates the appropriate bank sites.

18
include/openmc/plot.h
View file

@ -18,14 +18,14 @@ namespace openmc {
// Global variables
//===============================================================================
extern int PLOT_LEVEL_LOWEST; //!< lower bound on plot universe level
class Plot;
namespace model {
extern std::vector<Plot> plots; //!< Plot instance container
extern std::unordered_map<int, int> plot_map; //!< map of plot ids to index
extern "C" int32_t n_plots; //!< number of plots in openmc run
class Plot;
extern std::vector<Plot> plots; //!< Plot instance container
} // namespace model
//===============================================================================
// RGBColor holds color information for plotted objects
@ -124,8 +124,7 @@ void draw_mesh_lines(Plot pl, ImageData& data);
//! Write a ppm image to file using a plot object's image data
//! \param[in] plot object
//! \param[out] image data associated with the plot object
void output_ppm(Plot pl,
const ImageData& data);
void output_ppm(Plot pl, const ImageData& data);
//! Get the rgb color for a given particle position in a plot
//! \param[in] particle with position for current pixel
@ -141,7 +140,7 @@ void position_rgb(Particle p, Plot pl, RGBColor& rgb, int& id);
//! \param[out] dataset pointer to voxesl data
//! \param[out] pointer to memory space of voxel data
void voxel_init(hid_t file_id, const hsize_t* dims, hid_t* dspace,
hid_t* dset, hid_t* memspace);
hid_t* dset, hid_t* memspace);
//! Write a section of the voxel data to hdf5
//! \param[in] voxel slice
@ -149,7 +148,8 @@ void voxel_init(hid_t file_id, const hsize_t* dims, hid_t* dspace,
//! \param[out] dataset pointer to voxesl data
//! \param[out] pointer to data to write
void voxel_write_slice(int x, hid_t dspace, hid_t dset,
hid_t memspace, void* buf);
hid_t memspace, void* buf);
//! Close voxel file entities
//! \param[in] data space to close
//! \param[in] dataset to close

18
include/openmc/source.h
View file

@ -16,6 +16,18 @@
namespace openmc {
//==============================================================================
// Global variables
//==============================================================================
class SourceDistribution;
namespace model {
extern std::vector<SourceDistribution> external_sources;
} // namespace model
//==============================================================================
//! External source distribution
//==============================================================================
@ -40,12 +52,6 @@ private:
UPtrDist energy_; //!< Energy distribution
};
//==============================================================================
// Global variables
//==============================================================================
extern std::vector<SourceDistribution> external_sources;
//==============================================================================
// Functions
//==============================================================================

8
include/openmc/surface.h
View file

@ -32,13 +32,15 @@ extern "C" const int BC_PERIODIC;
// Global variables
//==============================================================================
extern "C" int32_t n_surfaces;
class Surface;
extern std::vector<Surface*> surfaces;
namespace model {
extern std::vector<Surface*> surfaces;
extern std::map<int, int> surface_map;
} // namespace model
//==============================================================================
//! Coordinates for an axis-aligned cube that bounds a geometric object.
//==============================================================================

11
include/openmc/tallies/filter.h
View file

@ -29,7 +29,7 @@ public:
// Without an explicit instantiation of vector<FilterMatch>, the Intel compiler
// will complain about the threadprivate directive on filter_matches. Note that
// this has to happen *outside* of the openmc namespace
template class std::vector<openmc::FilterMatch>;
extern template class std::vector<openmc::FilterMatch>;
namespace openmc {
@ -77,13 +77,20 @@ public:
// Global variables
//==============================================================================
extern "C" int32_t n_filters;
namespace simulation {
extern std::vector<FilterMatch> filter_matches;
#pragma omp threadprivate(filter_matches)
} // namespace simulation
namespace model {
extern "C" int32_t n_filters;
extern std::vector<std::unique_ptr<Filter>> tally_filters;
} // namespace model
//==============================================================================
extern "C" void free_memory_tally_c();

36
include/openmc/timer.h
View file

@ -5,6 +5,27 @@
namespace openmc {
//==============================================================================
// Global variables
//==============================================================================
class Timer;
namespace simulation {
extern Timer time_active;
extern Timer time_bank;
extern Timer time_bank_sample;
extern Timer time_bank_sendrecv;
extern Timer time_finalize;
extern Timer time_inactive;
extern Timer time_initialize;
extern Timer time_tallies;
extern Timer time_total;
extern Timer time_transport;
} // namespace simulation
//==============================================================================
//! Class for measuring time elapsed
//==============================================================================
@ -34,21 +55,6 @@ private:
double elapsed_ {0.0}; //!< elasped time in [s]
};
//==============================================================================
// Global variables
//==============================================================================
extern Timer time_active;
extern Timer time_bank;
extern Timer time_bank_sample;
extern Timer time_bank_sendrecv;
extern Timer time_finalize;
extern Timer time_inactive;
extern Timer time_initialize;
extern Timer time_tallies;
extern Timer time_total;
extern Timer time_transport;
//==============================================================================
// Non-member functions
//==============================================================================

2
openmc/capi/core.py
View file

@ -52,7 +52,7 @@ _dll.openmc_simulation_init.restype = c_int
_dll.openmc_simulation_init.errcheck = _error_handler
_dll.openmc_simulation_finalize.restype = c_int
_dll.openmc_simulation_finalize.errcheck = _error_handler
_dll.openmc_statepoint_write.argtypes = [POINTER(c_char_p), POINTER(c_bool)]
_dll.openmc_statepoint_write.argtypes = [c_char_p, POINTER(c_bool)]
_dll.openmc_statepoint_write.restype = c_int
_dll.openmc_statepoint_write.errcheck = _error_handler

108
src/cell.cpp
View file

@ -21,6 +21,8 @@ namespace openmc {
// Global variables
//==============================================================================
namespace model {
int32_t n_cells {0};
std::vector<Cell*> cells;
@ -29,6 +31,8 @@ std::unordered_map<int32_t, int32_t> cell_map;
std::vector<Universe*> universes;
std::unordered_map<int32_t, int32_t> universe_map;
} // namespace model
//==============================================================================
//! Convert region specification string to integer tokens.
//!
@ -197,7 +201,7 @@ Universe::to_hdf5(hid_t universes_group) const
// Write the contained cells.
if (cells_.size() > 0) {
std::vector<int32_t> cell_ids;
for (auto i_cell : cells_) cell_ids.push_back(cells[i_cell]->id_);
for (auto i_cell : cells_) cell_ids.push_back(model::cells[i_cell]->id_);
write_dataset(group, "cells", cell_ids);
}
@ -209,7 +213,7 @@ Universe::to_hdf5(hid_t universes_group) const
//==============================================================================
CSGCell::CSGCell() {} // empty constructor
CSGCell::CSGCell(pugi::xml_node cell_node)
{
if (check_for_node(cell_node, "id")) {
@ -314,7 +318,7 @@ CSGCell::CSGCell(pugi::xml_node cell_node)
// Convert user IDs to surface indices.
for (auto& r : region_) {
if (r < OP_UNION) {
r = copysign(surface_map[abs(r)] + 1, r);
r = copysign(model::surface_map[abs(r)] + 1, r);
}
}
@ -418,7 +422,7 @@ CSGCell::distance(Position r, Direction u, int32_t on_surface) const
// Calculate the distance to this surface.
// Note the off-by-one indexing
bool coincident {token == on_surface};
double d {surfaces[abs(token)-1]->distance(r, u, coincident)};
double d {model::surfaces[abs(token)-1]->distance(r, u, coincident)};
// Check if this distance is the new minimum.
if (d < min_dist) {
@ -446,7 +450,7 @@ CSGCell::to_hdf5(hid_t cell_group) const
write_string(group, "name", name_, false);
}
write_dataset(group, "universe", universes[universe_]->id_);
write_dataset(group, "universe", model::universes[universe_]->id_);
// Write the region specification.
if (!region_.empty()) {
@ -464,7 +468,7 @@ CSGCell::to_hdf5(hid_t cell_group) const
} else {
// Note the off-by-one indexing
region_spec << " "
<< copysign(surfaces[abs(token)-1]->id_, token);
<< copysign(model::surfaces[abs(token)-1]->id_, token);
}
}
write_string(group, "region", region_spec.str(), false);
@ -476,7 +480,7 @@ CSGCell::to_hdf5(hid_t cell_group) const
std::vector<int32_t> mat_ids;
for (auto i_mat : material_) {
if (i_mat != MATERIAL_VOID) {
mat_ids.push_back(materials[i_mat]->id_);
mat_ids.push_back(model::materials[i_mat]->id_);
} else {
mat_ids.push_back(MATERIAL_VOID);
}
@ -494,7 +498,7 @@ CSGCell::to_hdf5(hid_t cell_group) const
} else if (type_ == FILL_UNIVERSE) {
write_dataset(group, "fill_type", "universe");
write_dataset(group, "fill", universes[fill_]->id_);
write_dataset(group, "fill", model::universes[fill_]->id_);
if (translation_ != Position(0, 0, 0)) {
write_dataset(group, "translation", translation_);
}
@ -505,7 +509,7 @@ CSGCell::to_hdf5(hid_t cell_group) const
} else if (type_ == FILL_LATTICE) {
write_dataset(group, "fill_type", "lattice");
write_dataset(group, "lattice", lattices[fill_]->id_);
write_dataset(group, "lattice", model::lattices[fill_]->id_);
}
close_group(group);
@ -527,7 +531,7 @@ CSGCell::contains_simple(Position r, Direction u, int32_t on_surface) const
return false;
} else {
// Note the off-by-one indexing
bool sense = surfaces[abs(token)-1]->sense(r, u);
bool sense = model::surfaces[abs(token)-1]->sense(r, u);
if (sense != (token > 0)) {return false;}
}
}
@ -569,7 +573,7 @@ CSGCell::contains_complex(Position r, Direction u, int32_t on_surface) const
stack[i_stack] = false;
} else {
// Note the off-by-one indexing
bool sense = surfaces[abs(token)-1]->sense(r, u);
bool sense = model::surfaces[abs(token)-1]->sense(r, u);
stack[i_stack] = (sense == (token > 0));
}
}
@ -609,10 +613,10 @@ DAGCell::distance(Position r, Direction u, int32_t on_surface) const
} else { // indicate that particle is lost
surf_idx = -1;
}
return {dist, surf_idx};
}
bool DAGCell::contains(Position r, Direction u, int32_t on_surface) const
{
moab::ErrorCode rval;
@ -622,7 +626,7 @@ bool DAGCell::contains(Position r, Direction u, int32_t on_surface) const
double pnt[3] = {r.x, r.y, r.z};
double dir[3] = {u.x, u.y, u.z};
rval = dagmc_ptr_->point_in_volume(vol, pnt, result, dir);
MB_CHK_ERR_CONT(rval);
MB_CHK_ERR_CONT(rval);
return result;
}
@ -638,23 +642,23 @@ extern "C" void
read_cells(pugi::xml_node* node)
{
// Count the number of cells.
for (pugi::xml_node cell_node: node->children("cell")) {n_cells++;}
if (n_cells == 0) {
for (pugi::xml_node cell_node: node->children("cell")) {model::n_cells++;}
if (model::n_cells == 0) {
fatal_error("No cells found in geometry.xml!");
}
// Loop over XML cell elements and populate the array.
cells.reserve(n_cells);
for (pugi::xml_node cell_node: node->children("cell")) {
cells.push_back(new CSGCell(cell_node));
model::cells.reserve(model::n_cells);
for (pugi::xml_node cell_node : node->children("cell")) {
model::cells.push_back(new CSGCell(cell_node));
}
// Fill the cell map.
for (int i = 0; i < cells.size(); i++) {
int32_t id = cells[i]->id_;
auto search = cell_map.find(id);
if (search == cell_map.end()) {
cell_map[id] = i;
for (int i = 0; i < model::cells.size(); i++) {
int32_t id = model::cells[i]->id_;
auto search = model::cell_map.find(id);
if (search == model::cell_map.end()) {
model::cell_map[id] = i;
} else {
std::stringstream err_msg;
err_msg << "Two or more cells use the same unique ID: " << id;
@ -663,22 +667,24 @@ read_cells(pugi::xml_node* node)
}
// Populate the Universe vector and map.
for (int i = 0; i < cells.size(); i++) {
int32_t uid = cells[i]->universe_;
auto it = universe_map.find(uid);
if (it == universe_map.end()) {
universes.push_back(new Universe());
universes.back()->id_ = uid;
universes.back()->cells_.push_back(i);
universe_map[uid] = universes.size() - 1;
for (int i = 0; i < model::cells.size(); i++) {
int32_t uid = model::cells[i]->universe_;
auto it = model::universe_map.find(uid);
if (it == model::universe_map.end()) {
model::universes.push_back(new Universe());
model::universes.back()->id_ = uid;
model::universes.back()->cells_.push_back(i);
model::universe_map[uid] = model::universes.size() - 1;
} else {
universes[it->second]->cells_.push_back(i);
model::universes[it->second]->cells_.push_back(i);
}
}
universes.shrink_to_fit();
model::universes.shrink_to_fit();
// Allocate the cell overlap count if necessary.
if (settings::check_overlaps) overlap_check_count.resize(n_cells, 0);
if (settings::check_overlaps) {
model::overlap_check_count.resize(model::cells.size(), 0);
}
}
//==============================================================================
@ -688,9 +694,9 @@ read_cells(pugi::xml_node* node)
extern "C" int
openmc_cell_get_fill(int32_t index, int* type, int32_t** indices, int32_t* n)
{
if (index >= 1 && index <= cells.size()) {
if (index >= 1 && index <= model::cells.size()) {
//TODO: off-by-one
Cell& c {*cells[index - 1]};
Cell& c {*model::cells[index - 1]};
*type = c.type_;
if (c.type_ == FILL_MATERIAL) {
*indices = c.material_.data();
@ -710,9 +716,9 @@ extern "C" int
openmc_cell_set_fill(int32_t index, int type, int32_t n,
const int32_t* indices)
{
if (index >= 1 && index <= cells.size()) {
if (index >= 1 && index <= model::cells.size()) {
//TODO: off-by-one
Cell& c {*cells[index - 1]};
Cell& c {*model::cells[index - 1]};
if (type == FILL_MATERIAL) {
c.type_ = FILL_MATERIAL;
c.material_.clear();
@ -720,7 +726,7 @@ openmc_cell_set_fill(int32_t index, int type, int32_t n,
int i_mat = indices[i];
if (i_mat == MATERIAL_VOID) {
c.material_.push_back(MATERIAL_VOID);
} else if (i_mat >= 1 && i_mat <= materials.size()) {
} else if (i_mat >= 1 && i_mat <= model::materials.size()) {
//TODO: off-by-one
c.material_.push_back(i_mat - 1);
} else {
@ -745,9 +751,9 @@ openmc_cell_set_fill(int32_t index, int type, int32_t n,
extern "C" int
openmc_cell_set_temperature(int32_t index, double T, const int32_t* instance)
{
if (index >= 1 && index <= cells.size()) {
if (index >= 1 && index <= model::cells.size()) {
//TODO: off-by-one
Cell& c {*cells[index - 1]};
Cell& c {*model::cells[index - 1]};
if (instance) {
if (*instance >= 0 && *instance < c.sqrtkT_.size()) {
@ -775,7 +781,7 @@ openmc_cell_set_temperature(int32_t index, double T, const int32_t* instance)
//==============================================================================
extern "C" {
Cell* cell_pointer(int32_t cell_ind) {return cells[cell_ind];}
Cell* cell_pointer(int32_t cell_ind) {return model::cells[cell_ind];}
int32_t cell_id(Cell* c) {return c->id_;}
@ -785,9 +791,9 @@ extern "C" {
c->id_ = id;
// Find the index of this cell and update the cell map.
for (int i = 0; i < cells.size(); i++) {
if (cells[i] == c) {
cell_map[id] = i;
for (int i = 0; i < model::cells.size(); i++) {
if (model::cells[i] == c) {
model::cell_map[id] = i;
break;
}
}
@ -830,17 +836,17 @@ extern "C" {
void extend_cells_c(int32_t n)
{
cells.reserve(cells.size() + n);
model::cells.reserve(model::cells.size() + n);
for (int32_t i = 0; i < n; i++) {
cells.push_back(new CSGCell());
model::cells.push_back(new CSGCell());
}
n_cells = cells.size();
model::n_cells = model::cells.size();
}
int32_t universe_id(int i_univ) {return universes[i_univ]->id_;}
int32_t universe_id(int i_univ) {return model::universes[i_univ]->id_;}
void universes_to_hdf5(hid_t universes_group)
{for (Universe* u : universes) u->to_hdf5(universes_group);}
{for (Universe* u : model::universes) u->to_hdf5(universes_group);}
}

2
src/cmfd_execute.cpp
View file

@ -24,7 +24,7 @@ cmfd_populate_sourcecounts(int n_energy, const double* energies,
openmc_source_bank(&source_bank, &n);
// Get source counts in each mesh bin / energy bin
auto& m = meshes.at(settings::index_cmfd_mesh);
auto& m = model::meshes.at(settings::index_cmfd_mesh);
xt::xarray<double> counts = m->count_sites(simulation::work, source_bank, n_energy, energies, outside);
// Copy data from the xarray into the source counts array

79
src/dagmc.cpp
View file

@ -1,5 +1,6 @@
#include "openmc/dagmc.h"
#include "openmc/cell.h"
#include "openmc/error.h"
#include "openmc/string_utils.h"
#include "openmc/settings.h"
@ -13,20 +14,24 @@
namespace openmc {
namespace model {
moab::DagMC* DAG;
} // namespace model
void load_dagmc_geometry()
{
if (!DAG) {
DAG = new moab::DagMC();
if (!model::DAG) {
model::DAG = new moab::DagMC();
}
int32_t dagmc_univ_id = 0; // universe is always 0 for DAGMC
moab::ErrorCode rval = DAG->load_file("dagmc.h5m");
moab::ErrorCode rval = model::DAG->load_file("dagmc.h5m");
MB_CHK_ERR_CONT(rval);
rval = DAG->init_OBBTree();
rval = model::DAG->init_OBBTree();
MB_CHK_ERR_CONT(rval);
std::vector<std::string> prop_keywords;
@ -34,48 +39,45 @@ void load_dagmc_geometry()
prop_keywords.push_back("boundary");
std::map<std::string, std::string> ph;
DAG->parse_properties(prop_keywords, ph, ":");
model::DAG->parse_properties(prop_keywords, ph, ":");
MB_CHK_ERR_CONT(rval);
// initialize cell objects
n_cells = DAG->num_entities(3);
model::n_cells = model::DAG->num_entities(3);
// Allocate the cell overlap count if necessary.
if (settings::check_overlaps) overlap_check_count.resize(n_cells, 0);
for (int i = 0; i < n_cells; i++) {
moab::EntityHandle vol_handle = DAG->entity_by_index(3, i+1);
for (int i = 0; i < model::n_cells; i++) {
moab::EntityHandle vol_handle = model::DAG->entity_by_index(3, i+1);
// set cell ids using global IDs
DAGCell* c = new DAGCell();
c->id_ = DAG->id_by_index(3, i+1);
c->dagmc_ptr_ = DAG;
c->id_ = model::DAG->id_by_index(3, i+1);
c->dagmc_ptr_ = model::DAG;
c->universe_ = dagmc_univ_id; // set to zero for now
c->fill_ = C_NONE; // no fill, single universe
cells.push_back(c);
cell_map[c->id_] = i;
model::cells.push_back(c);
model::cell_map[c->id_] = i;
// Populate the Universe vector and dict
auto it = universe_map.find(dagmc_univ_id);
if (it == universe_map.end()) {
universes.push_back(new Universe());
universes.back()-> id_ = dagmc_univ_id;
universes.back()->cells_.push_back(i);
universe_map[dagmc_univ_id] = universes.size() - 1;
auto it = model::universe_map.find(dagmc_univ_id);
if (it == model::universe_map.end()) {
model::universes.push_back(new Universe());
model::universes.back()-> id_ = dagmc_univ_id;
model::universes.back()->cells_.push_back(i);
model::universe_map[dagmc_univ_id] = model::universes.size() - 1;
} else {
universes[it->second]->cells_.push_back(i);
model::universes[it->second]->cells_.push_back(i);
}
if (DAG->is_implicit_complement(vol_handle)) {
if (model::DAG->is_implicit_complement(vol_handle)) {
// assuming implicit complement is void for now
c->material_.push_back(MATERIAL_VOID);
continue;
}
if (DAG->has_prop(vol_handle, "mat")){
if (model::DAG->has_prop(vol_handle, "mat")){
std::string mat_value;
rval = DAG->prop_value(vol_handle, "mat", mat_value);
rval = model::DAG->prop_value(vol_handle, "mat", mat_value);
MB_CHK_ERR_CONT(rval);
to_lower(mat_value);
@ -91,21 +93,26 @@ void load_dagmc_geometry()
}
}
// Allocate the cell overlap count if necessary.
if (settings::check_overlaps) {
model::overlap_check_count.resize(model::cells.size(), 0);
}
// initialize surface objects
n_surfaces = DAG->num_entities(2);
surfaces.resize(n_surfaces);
int n_surfaces = model::DAG->num_entities(2);
model::surfaces.resize(n_surfaces);
for (int i = 0; i < n_surfaces; i++) {
moab::EntityHandle surf_handle = DAG->entity_by_index(2, i+1);
moab::EntityHandle surf_handle = model::DAG->entity_by_index(2, i+1);
// set cell ids using global IDs
DAGSurface* s = new DAGSurface();
s->id_ = DAG->id_by_index(2, i+1);
s->dagmc_ptr_ = DAG;
s->id_ = model::DAG->id_by_index(2, i+1);
s->dagmc_ptr_ = model::DAG;
if (DAG->has_prop(surf_handle, "boundary")) {
if (model::DAG->has_prop(surf_handle, "boundary")) {
std::string bc_value;
rval = DAG->prop_value(surf_handle, "boundary", bc_value);
rval = model::DAG->prop_value(surf_handle, "boundary", bc_value);
MB_CHK_ERR_CONT(rval);
to_lower(bc_value);
@ -128,8 +135,8 @@ void load_dagmc_geometry()
}
// add to global array and map
surfaces[i] = s;
surface_map[s->id_] = s->id_;
model::surfaces[i] = s;
model::surface_map[s->id_] = s->id_;
}
return;
@ -137,7 +144,7 @@ void load_dagmc_geometry()
void free_memory_dagmc()
{
delete DAG;
delete model::DAG;
}
}

78
src/eigenvalue.cpp
View file

@ -30,11 +30,15 @@ namespace openmc {
// Global variables
//==============================================================================
namespace simulation {
double keff_generation;
std::array<double, 2> k_sum;
std::vector<double> entropy;
xt::xtensor<double, 1> source_frac;
} // namespace simulation
//==============================================================================
// Non-member functions
//==============================================================================
@ -44,15 +48,15 @@ void calculate_generation_keff()
auto gt = global_tallies();
// Get keff for this generation by subtracting off the starting value
keff_generation = gt(K_TRACKLENGTH, RESULT_VALUE) - keff_generation;
simulation::keff_generation = gt(K_TRACKLENGTH, RESULT_VALUE) - simulation::keff_generation;
double keff_reduced;
#ifdef OPENMC_MPI
// Combine values across all processors
MPI_Allreduce(&keff_generation, &keff_reduced, 1, MPI_DOUBLE,
MPI_Allreduce(&simulation::keff_generation, &keff_reduced, 1, MPI_DOUBLE,
MPI_SUM, mpi::intracomm);
#else
keff_reduced = keff_generation;
keff_reduced = simulation::keff_generation;
#endif
// Normalize single batch estimate of k
@ -63,7 +67,7 @@ void calculate_generation_keff()
void synchronize_bank()
{
time_bank.start();
simulation::time_bank.start();
// Get pointers to source/fission bank
Bank* source_bank;
@ -83,21 +87,21 @@ void synchronize_bank()
#ifdef OPENMC_MPI
int64_t start = 0;
MPI_Exscan(&n_bank, &start, 1, MPI_INT64_T, MPI_SUM, mpi::intracomm);
MPI_Exscan(&simulation::n_bank, &start, 1, MPI_INT64_T, MPI_SUM, mpi::intracomm);
// While we would expect the value of start on rank 0 to be 0, the MPI
// standard says that the receive buffer on rank 0 is undefined and not
// significant
if (mpi::rank == 0) start = 0;
int64_t finish = start + n_bank;
int64_t finish = start + simulation::n_bank;
int64_t total = finish;
MPI_Bcast(&total, 1, MPI_INT64_T, mpi::n_procs - 1, mpi::intracomm);
#else
int64_t start = 0;
int64_t finish = n_bank;
int64_t total = n_bank;
int64_t finish = simulation::n_bank;
int64_t total = simulation::n_bank;
#endif
// If there are not that many particles per generation, it's possible that no
@ -105,7 +109,7 @@ void synchronize_bank()
// extra logic to treat this circumstance, we really want to ensure the user
// runs enough particles to avoid this in the first place.
if (n_bank == 0) {
if (simulation::n_bank == 0) {
fatal_error("No fission sites banked on MPI rank " + std::to_string(mpi::rank));
}
@ -127,7 +131,7 @@ void synchronize_bank()
}
double p_sample = static_cast<double>(sites_needed) / total;
time_bank_sample.start();
simulation::time_bank_sample.start();
// ==========================================================================
// SAMPLE N_PARTICLES FROM FISSION BANK AND PLACE IN TEMP_SITES
@ -136,7 +140,7 @@ void synchronize_bank()
int64_t index_temp = 0;
std::vector<Bank> temp_sites(3*simulation::work);
for (int64_t i = 0; i < n_bank; ++i) {
for (int64_t i = 0; i < simulation::n_bank; ++i) {
// If there are less than n_particles particles banked, automatically add
// int(n_particles/total) sites to temp_sites. For example, if you need
// 1000 and 300 were banked, this would add 3 source sites per banked site
@ -191,7 +195,7 @@ void synchronize_bank()
// fission bank
sites_needed = settings::n_particles - finish;
for (int i = 0; i < sites_needed; ++i) {
temp_sites[index_temp] = fission_bank[n_bank - sites_needed + i];
temp_sites[index_temp] = fission_bank[simulation::n_bank - sites_needed + i];
++index_temp;
}
}
@ -200,8 +204,8 @@ void synchronize_bank()
finish = simulation::work_index[mpi::rank + 1];
}
time_bank_sample.stop();
time_bank_sendrecv.start();
simulation::time_bank_sample.stop();
simulation::time_bank_sendrecv.start();
#ifdef OPENMC_MPI
// ==========================================================================
@ -299,8 +303,8 @@ void synchronize_bank()
std::copy(temp_sites.data(), temp_sites.data() + settings::n_particles, source_bank);
#endif
time_bank_sendrecv.stop();
time_bank.stop();
simulation::time_bank_sendrecv.stop();
simulation::time_bank.stop();
}
void calculate_average_keff()
@ -320,11 +324,11 @@ void calculate_average_keff()
simulation::keff = simulation::k_generation[i];
} else {
// Sample mean of keff
k_sum[0] += simulation::k_generation[i];
k_sum[1] += std::pow(simulation::k_generation[i], 2);
simulation::k_sum[0] += simulation::k_generation[i];
simulation::k_sum[1] += std::pow(simulation::k_generation[i], 2);
// Determine mean
simulation::keff = k_sum[0] / n;
simulation::keff = simulation::k_sum[0] / n;
if (n > 1) {
double t_value;
@ -337,7 +341,7 @@ void calculate_average_keff()
}
// Standard deviation of the sample mean of k
simulation::keff_std = t_value * std::sqrt((k_sum[1]/n -
simulation::keff_std = t_value * std::sqrt((simulation::k_sum[1]/n -
std::pow(simulation::keff, 2)) / (n - 1));
}
}
@ -493,7 +497,7 @@ int openmc_get_keff(double* k_combined)
void shannon_entropy()
{
// Get pointer to entropy mesh
auto& m = meshes[settings::index_entropy_mesh];
auto& m = model::meshes[settings::index_entropy_mesh];
// Get pointer to fission bank
Bank* fission_bank;
@ -503,7 +507,7 @@ void shannon_entropy()
// Get source weight in each mesh bin
bool sites_outside;
xt::xtensor<double, 1> p = m->count_sites(
n_bank, fission_bank, 0, nullptr, &sites_outside);
simulation::n_bank, fission_bank, 0, nullptr, &sites_outside);
// display warning message if there were sites outside entropy box
if (sites_outside) {
@ -523,20 +527,20 @@ void shannon_entropy()
}
// Add value to vector
entropy.push_back(H);
simulation::entropy.push_back(H);
}
}
void ufs_count_sites()
{
auto &m = meshes[settings::index_ufs_mesh];
auto &m = model::meshes[settings::index_ufs_mesh];
if (simulation::current_batch == 1 && simulation::current_gen == 1) {
// On the first generation, just assume that the source is already evenly
// distributed so that effectively the production of fission sites is not
// biased
auto s = xt::view(source_frac, xt::all());
auto s = xt::view(simulation::source_frac, xt::all());
s = m->volume_frac_;
} else {
@ -547,7 +551,7 @@ void ufs_count_sites()
// count number of source sites in each ufs mesh cell
bool sites_outside;
source_frac = m->count_sites(simulation::work, source_bank, 0, nullptr,
simulation::source_frac = m->count_sites(simulation::work, source_bank, 0, nullptr,
&sites_outside);
// Check for sites outside of the mesh
@ -558,12 +562,12 @@ void ufs_count_sites()
#ifdef OPENMC_MPI
// Send source fraction to all processors
int n_bins = xt::prod(m->shape_)();
MPI_Bcast(source_frac.data(), n_bins, MPI_DOUBLE, 0, mpi::intracomm);
MPI_Bcast(simulation::source_frac.data(), n_bins, MPI_DOUBLE, 0, mpi::intracomm);
#endif
// Normalize to total weight to get fraction of source in each cell
double total = xt::sum(source_frac)();
source_frac /= total;
double total = xt::sum(simulation::source_frac)();
simulation::source_frac /= total;
// Since the total starting weight is not equal to n_particles, we need to
// renormalize the weight of the source sites
@ -575,7 +579,7 @@ void ufs_count_sites()
double ufs_get_weight(const Particle* p)
{
auto& m = meshes[settings::index_ufs_mesh];
auto& m = model::meshes[settings::index_ufs_mesh];
// Determine indices on ufs mesh for current location
// TODO: off by one
@ -585,8 +589,8 @@ double ufs_get_weight(const Particle* p)
fatal_error("Source site outside UFS mesh!");
}
if (source_frac(mesh_bin) != 0.0) {
return m->volume_frac_ / source_frac(mesh_bin);
if (simulation::source_frac(mesh_bin) != 0.0) {
return m->volume_frac_ / simulation::source_frac(mesh_bin);
} else {
return 1.0;
}
@ -598,7 +602,7 @@ extern "C" void write_eigenvalue_hdf5(hid_t group)
write_dataset(group, "generations_per_batch", settings::gen_per_batch);
write_dataset(group, "k_generation", simulation::k_generation);
if (settings::entropy_on) {
write_dataset(group, "entropy", entropy);
write_dataset(group, "entropy", simulation::entropy);
}
write_dataset(group, "k_col_abs", simulation::k_col_abs);
write_dataset(group, "k_col_tra", simulation::k_col_tra);
@ -615,7 +619,7 @@ extern "C" void read_eigenvalue_hdf5(hid_t group)
simulation::k_generation.resize(n);
read_dataset(group, "k_generation", simulation::k_generation);
if (settings::entropy_on) {
read_dataset(group, "entropy", entropy);
read_dataset(group, "entropy", simulation::entropy);
}
read_dataset(group, "k_col_abs", simulation::k_col_abs);
read_dataset(group, "k_col_tra", simulation::k_col_tra);
@ -628,14 +632,14 @@ extern "C" void read_eigenvalue_hdf5(hid_t group)
extern "C" double entropy_c(int i)
{
return entropy.at(i - 1);
return simulation::entropy.at(i - 1);
}
extern "C" void entropy_clear()
{
entropy.clear();
simulation::entropy.clear();
}
extern "C" void k_sum_reset() { k_sum.fill(0.0); }
extern "C" void k_sum_reset() { simulation::k_sum.fill(0.0); }
} // namespace openmc

8
src/finalize.cpp
View file

@ -75,10 +75,10 @@ int openmc_finalize()
simulation::satisfy_triggers = false;
simulation::total_gen = 0;
energy_max = {INFTY, INFTY};
energy_min = {0.0, 0.0};
data::energy_max = {INFTY, INFTY};
data::energy_min = {0.0, 0.0};
n_tallies = 0;
openmc_root_universe = -1;
model::root_universe = -1;
openmc_set_seed(DEFAULT_SEED);
// Deallocate arrays
@ -113,7 +113,7 @@ int openmc_reset()
simulation::k_col_abs = 0.0;
simulation::k_col_tra = 0.0;
simulation::k_abs_tra = 0.0;
k_sum = {0.0, 0.0};
simulation::k_sum = {0.0, 0.0};
// Reset timers
reset_timers();

53
src/geometry.cpp
View file

@ -14,8 +14,21 @@
namespace openmc {
//==============================================================================
// Global variables
//==============================================================================
namespace model {
int root_universe {-1};
std::vector<int64_t> overlap_check_count;
} // namespace model
//==============================================================================
// Non-member functions
//==============================================================================
extern "C" bool
@ -24,21 +37,21 @@ check_cell_overlap(Particle* p) {
// Loop through each coordinate level
for (int j = 0; j < n_coord; j++) {
Universe& univ = *universes[p->coord[j].universe];
Universe& univ = *model::universes[p->coord[j].universe];
int n = univ.cells_.size();
// Loop through each cell on this level
for (auto index_cell : univ.cells_) {
Cell& c = *cells[index_cell];
Cell& c = *model::cells[index_cell];
if (c.contains(p->coord[j].xyz, p->coord[j].uvw, p->surface)) {
if (index_cell != p->coord[j].cell) {
std::stringstream err_msg;
err_msg << "Overlapping cells detected: " << c.id_ << ", "
<< cells[p->coord[j].cell]->id_ << " on universe "
<< model::cells[p->coord[j].cell]->id_ << " on universe "
<< univ.id_;
fatal_error(err_msg);
}
++overlap_check_count[index_cell];
++model::overlap_check_count[index_cell];
}
}
}
@ -57,8 +70,8 @@ find_cell(Particle* p, int search_surf) {
// Determine universe (if not yet set, use root universe)
int i_universe = p->coord[p->n_coord-1].universe;
if (i_universe == C_NONE) {
p->coord[p->n_coord-1].universe = openmc_root_universe;
i_universe = openmc_root_universe;
p->coord[p->n_coord-1].universe = model::root_universe;
i_universe = model::root_universe;
}
// If a surface was indicated, only search cells from the neighbor list of
@ -66,12 +79,12 @@ find_cell(Particle* p, int search_surf) {
// the positive or negative side of the surface should be searched.
const std::vector<int>* search_cells;
if (search_surf > 0) {
search_cells = &surfaces[search_surf-1]->neighbor_pos_;
search_cells = &model::surfaces[search_surf-1]->neighbor_pos_;
} else if (search_surf < 0) {
search_cells = &surfaces[-search_surf-1]->neighbor_neg_;
search_cells = &model::surfaces[-search_surf-1]->neighbor_neg_;
} else {
// No surface was indicated, search all cells in the universe.
search_cells = &universes[i_universe]->cells_;
search_cells = &model::universes[i_universe]->cells_;
}
// Find which cell of this universe the particle is in.
@ -81,17 +94,17 @@ find_cell(Particle* p, int search_surf) {
i_cell = (*search_cells)[i];
// Make sure the search cell is in the same universe.
if (cells[i_cell]->universe_ != i_universe) continue;
if (model::cells[i_cell]->universe_ != i_universe) continue;
Position r {p->coord[p->n_coord-1].xyz};
Direction u {p->coord[p->n_coord-1].uvw};
int32_t surf = p->surface;
if (cells[i_cell]->contains(r, u, surf)) {
if (model::cells[i_cell]->contains(r, u, surf)) {
p->coord[p->n_coord-1].cell = i_cell;
if (settings::verbosity >= 10 || simulation::trace) {
std::stringstream msg;
msg << " Entering cell " << cells[i_cell]->id_;
msg << " Entering cell " << model::cells[i_cell]->id_;
write_message(msg, 1);
}
found = true;
@ -100,7 +113,7 @@ find_cell(Particle* p, int search_surf) {
}
if (found) {
Cell& c {*cells[i_cell]};
Cell& c {*model::cells[i_cell]};
if (c.type_ == FILL_MATERIAL) {
//=======================================================================
//! Found a material cell which means this is the lowest coord level.
@ -109,11 +122,11 @@ find_cell(Particle* p, int search_surf) {
if (c.material_.size() > 1 || c.sqrtkT_.size() > 1) {
int offset = 0;
for (int i = 0; i < p->n_coord; i++) {
Cell& c_i {*cells[p->coord[i].cell]};
Cell& c_i {*model::cells[p->coord[i].cell]};
if (c_i.type_ == FILL_UNIVERSE) {
offset += c_i.offset_[c.distribcell_index_];
} else if (c_i.type_ == FILL_LATTICE) {
Lattice& lat {*lattices[p->coord[i+1].lattice-1]};
Lattice& lat {*model::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};
@ -198,7 +211,7 @@ find_cell(Particle* p, int search_surf) {
//========================================================================
//! Found a lower lattice, update this coord level then search the next.
Lattice& lat {*lattices[c.fill_]};
Lattice& lat {*model::lattices[c.fill_]};
// Determine lattice indices.
Position r {p->coord[p->n_coord-1].xyz};
@ -251,7 +264,7 @@ find_cell(Particle* p, int search_surf) {
extern "C" void
cross_lattice(Particle* p, int lattice_translation[3])
{
Lattice& lat {*lattices[p->coord[p->n_coord-1].lattice-1]};
Lattice& lat {*model::lattices[p->coord[p->n_coord-1].lattice-1]};
if (settings::verbosity >= 10 || simulation::trace) {
std::stringstream msg;
@ -326,7 +339,7 @@ distance_to_boundary(Particle* p, double* dist, int* surface_crossed,
for (int i = 0; i < p->n_coord; i++) {
Position r {p->coord[i].xyz};
Direction u {p->coord[i].uvw};
Cell& c {*cells[p->coord[i].cell]};
Cell& c {*model::cells[p->coord[i].cell]};
// Find the oncoming surface in this cell and the distance to it.
auto surface_distance = c.distance(r, u, p->surface);
@ -335,7 +348,7 @@ distance_to_boundary(Particle* p, double* dist, int* surface_crossed,
// Find the distance to the next lattice tile crossing.
if (p->coord[i].lattice != F90_NONE) {
Lattice& lat {*lattices[p->coord[i].lattice-1]};
Lattice& lat {*model::lattices[p->coord[i].lattice-1]};
std::array<int, 3> i_xyz {p->coord[i].lattice_x, p->coord[i].lattice_y,
p->coord[i].lattice_z};
//TODO: refactor so both lattice use the same position argument (which
@ -377,7 +390,7 @@ distance_to_boundary(Particle* p, double* dist, int* surface_crossed,
*surface_crossed = level_surf_cross;
} else {
Position r_hit = r + d_surf * u;
Surface& surf {*surfaces[std::abs(level_surf_cross)-1]};
Surface& surf {*model::surfaces[std::abs(level_surf_cross)-1]};
Direction norm = surf.normal(r_hit);
if (u.dot(norm) > 0) {
*surface_crossed = std::abs(level_surf_cross);

124
src/geometry_aux.cpp
View file

@ -24,15 +24,15 @@ void
adjust_indices()
{
// Adjust material/fill idices.
for (Cell* c : cells) {
for (Cell* c : model::cells) {
if (c->fill_ != C_NONE) {
int32_t id = c->fill_;
auto search_univ = universe_map.find(id);
auto search_lat = lattice_map.find(id);
if (search_univ != universe_map.end()) {
auto search_univ = model::universe_map.find(id);
auto search_lat = model::lattice_map.find(id);
if (search_univ != model::universe_map.end()) {
c->type_ = FILL_UNIVERSE;
c->fill_ = search_univ->second;
} else if (search_lat != lattice_map.end()) {
} else if (search_lat != model::lattice_map.end()) {
c->type_ = FILL_LATTICE;
c->fill_ = search_lat->second;
} else {
@ -46,8 +46,8 @@ adjust_indices()
for (auto it = c->material_.begin(); it != c->material_.end(); it++) {
int32_t mid = *it;
if (mid != MATERIAL_VOID) {
auto search = material_map.find(mid);
if (search != material_map.end()) {
auto search = model::material_map.find(mid);
if (search != model::material_map.end()) {
*it = search->second;
} else {
std::stringstream err_msg;
@ -61,9 +61,9 @@ adjust_indices()
}
// Change cell.universe values from IDs to indices.
for (Cell* c : cells) {
auto search = universe_map.find(c->universe_);
if (search != universe_map.end()) {
for (Cell* c : model::cells) {
auto search = model::universe_map.find(c->universe_);
if (search != model::universe_map.end()) {
c->universe_ = search->second;
} else {
std::stringstream err_msg;
@ -74,7 +74,7 @@ adjust_indices()
}
// Change all lattice universe values from IDs to indices.
for (Lattice* l : lattices) {
for (Lattice* l : model::lattices) {
l->adjust_indices();
}
}
@ -84,7 +84,7 @@ adjust_indices()
void
assign_temperatures()
{
for (Cell* c : cells) {
for (Cell* c : model::cells) {
// Ignore non-material cells and cells with defined temperature.
if (c->material_.size() == 0) continue;
if (c->sqrtkT_.size() > 0) continue;
@ -96,9 +96,9 @@ assign_temperatures()
c->sqrtkT_.push_back(0);
} else {
if (materials[i_mat]->temperature_ >= 0) {
if (model::materials[i_mat]->temperature_ >= 0) {
// This material has a default temperature; use that value.
auto T = materials[i_mat]->temperature_;
auto T = model::materials[i_mat]->temperature_;
c->sqrtkT_.push_back(std::sqrt(K_BOLTZMANN * T));
} else {
// Use the global default temperature.
@ -117,12 +117,12 @@ find_root_universe()
{
// Find all the universes listed as a cell fill.
std::unordered_set<int32_t> fill_univ_ids;
for (Cell* c : cells) {
for (Cell* c : model::cells) {
fill_univ_ids.insert(c->fill_);
}
// Find all the universes contained in a lattice.
for (Lattice* lat : lattices) {
for (Lattice* lat : model::lattices) {
for (auto it = lat->begin(); it != lat->end(); ++it) {
fill_univ_ids.insert(*it);
}
@ -134,8 +134,8 @@ find_root_universe()
// Figure out which universe is not in the set. This is the root universe.
bool root_found {false};
int32_t root_univ;
for (int32_t i = 0; i < universes.size(); i++) {
auto search = fill_univ_ids.find(universes[i]->id_);
for (int32_t i = 0; i < model::universes.size(); i++) {
auto search = fill_univ_ids.find(model::universes[i]->id_);
if (search == fill_univ_ids.end()) {
if (root_found) {
fatal_error("Two or more universes are not used as fill universes, so "
@ -160,21 +160,21 @@ neighbor_lists()
{
write_message("Building neighboring cells lists for each surface...", 6);
for (int i = 0; i < cells.size(); i++) {
for (auto token : cells[i]->region_) {
for (int i = 0; i < model::cells.size(); i++) {
for (auto token : model::cells[i]->region_) {
// Skip operator tokens.
if (std::abs(token) >= OP_UNION) continue;
// This token is a surface index. Add the cell to the surface's list.
if (token > 0) {
surfaces[std::abs(token)-1]->neighbor_pos_.push_back(i);
model::surfaces[std::abs(token)-1]->neighbor_pos_.push_back(i);
} else {
surfaces[std::abs(token)-1]->neighbor_neg_.push_back(i);
model::surfaces[std::abs(token)-1]->neighbor_neg_.push_back(i);
}
}
}
for (Surface* surf : surfaces) {
for (Surface* surf : model::surfaces) {
surf->neighbor_pos_.shrink_to_fit();
surf->neighbor_neg_.shrink_to_fit();
}
@ -187,7 +187,7 @@ prepare_distribcell()
{
// Find all cells listed in a DistribcellFilter.
std::unordered_set<int32_t> distribcells;
for (auto& filt : tally_filters) {
for (auto& filt : model::tally_filters) {
auto* distrib_filt = dynamic_cast<DistribcellFilter*>(filt.get());
if (distrib_filt) {
distribcells.insert(distrib_filt->cell_);
@ -196,8 +196,8 @@ prepare_distribcell()
// Find all cells with distributed materials or temperatures. Make sure that
// the number of materials/temperatures matches the number of cell instances.
for (int i = 0; i < cells.size(); i++) {
Cell& c {*cells[i]};
for (int i = 0; i < model::cells.size(); i++) {
Cell& c {*model::cells[i]};
if (c.material_.size() > 1) {
if (c.material_.size() != c.n_instances_) {
@ -228,10 +228,10 @@ prepare_distribcell()
// unique distribcell array index.
int distribcell_index = 0;
std::vector<int32_t> target_univ_ids;
for (Universe* u : universes) {
for (Universe* u : model::universes) {
for (auto cell_indx : u->cells_) {
if (distribcells.find(cell_indx) != distribcells.end()) {
cells[cell_indx]->distribcell_index_ = distribcell_index;
model::cells[cell_indx]->distribcell_index_ = distribcell_index;
target_univ_ids.push_back(u->id_);
++distribcell_index;
}
@ -240,22 +240,22 @@ prepare_distribcell()
// Allocate the cell and lattice offset tables.
int n_maps = target_univ_ids.size();
for (Cell* c : cells) {
for (Cell* c : model::cells) {
if (c->type_ != FILL_MATERIAL) {
c->offset_.resize(n_maps, C_NONE);
}
}
for (Lattice* lat : lattices) {
for (Lattice* lat : model::lattices) {
lat->allocate_offset_table(n_maps);
}
// Fill the cell and lattice offset tables.
for (int map = 0; map < target_univ_ids.size(); map++) {
auto target_univ_id = target_univ_ids[map];
for (Universe* univ : universes) {
for (Universe* univ : model::universes) {
int32_t offset {0}; // TODO: is this a bug? It matches F90 implementation.
for (int32_t cell_indx : univ->cells_) {
Cell& c = *cells[cell_indx];
Cell& c = *model::cells[cell_indx];
if (c.type_ == FILL_UNIVERSE) {
c.offset_[map] = offset;
@ -263,7 +263,7 @@ prepare_distribcell()
offset += count_universe_instances(search_univ, target_univ_id);
} else if (c.type_ == FILL_LATTICE) {
Lattice& lat = *lattices[c.fill_];
Lattice& lat = *model::lattices[c.fill_];
offset = lat.fill_offset_table(offset, target_univ_id, map);
}
}
@ -276,8 +276,8 @@ prepare_distribcell()
void
count_cell_instances(int32_t univ_indx)
{
for (int32_t cell_indx : universes[univ_indx]->cells_) {
Cell& c = *cells[cell_indx];
for (int32_t cell_indx : model::universes[univ_indx]->cells_) {
Cell& c = *model::cells[cell_indx];
++c.n_instances_;
if (c.type_ == FILL_UNIVERSE) {
@ -286,7 +286,7 @@ count_cell_instances(int32_t univ_indx)
} else if (c.type_ == FILL_LATTICE) {
// This cell contains a lattice. Recurse into the lattice universes.
Lattice& lat = *lattices[c.fill_];
Lattice& lat = *model::lattices[c.fill_];
for (auto it = lat.begin(); it != lat.end(); ++it) {
count_cell_instances(*it);
}
@ -300,20 +300,20 @@ int
count_universe_instances(int32_t search_univ, int32_t target_univ_id)
{
// If this is the target, it can't contain itself.
if (universes[search_univ]->id_ == target_univ_id) {
if (model::universes[search_univ]->id_ == target_univ_id) {
return 1;
}
int count {0};
for (int32_t cell_indx : universes[search_univ]->cells_) {
Cell& c = *cells[cell_indx];
for (int32_t cell_indx : model::universes[search_univ]->cells_) {
Cell& c = *model::cells[cell_indx];
if (c.type_ == FILL_UNIVERSE) {
int32_t next_univ = c.fill_;
count += count_universe_instances(next_univ, target_univ_id);
} else if (c.type_ == FILL_LATTICE) {
Lattice& lat = *lattices[c.fill_];
Lattice& lat = *model::lattices[c.fill_];
for (auto it = lat.begin(); it != lat.end(); ++it) {
int32_t next_univ = *it;
count += count_universe_instances(next_univ, target_univ_id);
@ -338,7 +338,7 @@ distribcell_path_inner(int32_t target_cell, int32_t map, int32_t target_offset,
// write to the path and return.
for (int32_t cell_indx : search_univ.cells_) {
if ((cell_indx == target_cell) && (offset == target_offset)) {
Cell& c = *cells[cell_indx];
Cell& c = *model::cells[cell_indx];
path << "c" << c.id_;
return path.str();
}
@ -350,7 +350,7 @@ distribcell_path_inner(int32_t target_cell, int32_t map, int32_t target_offset,
std::vector<std::int32_t>::const_reverse_iterator cell_it
{search_univ.cells_.crbegin()};
for (; cell_it != search_univ.cells_.crend(); ++cell_it) {
Cell& c = *cells[*cell_it];
Cell& c = *model::cells[*cell_it];
// Material cells don't contain other cells so ignore them.
if (c.type_ != FILL_MATERIAL) {
@ -358,7 +358,7 @@ distribcell_path_inner(int32_t target_cell, int32_t map, int32_t target_offset,
if (c.type_ == FILL_UNIVERSE) {
temp_offset = offset + c.offset_[map];
} else {
Lattice& lat = *lattices[c.fill_];
Lattice& lat = *model::lattices[c.fill_];
int32_t indx = lat.universes_.size()*map + lat.begin().indx_;
temp_offset = offset + lat.offsets_[indx];
}
@ -370,18 +370,18 @@ distribcell_path_inner(int32_t target_cell, int32_t map, int32_t target_offset,
}
// Add the cell to the path string.
Cell& c = *cells[*cell_it];
Cell& c = *model::cells[*cell_it];
path << "c" << c.id_ << "->";
if (c.type_ == FILL_UNIVERSE) {
// Recurse into the fill cell.
offset += c.offset_[map];
path << distribcell_path_inner(target_cell, map, target_offset,
*universes[c.fill_], offset);
*model::universes[c.fill_], offset);
return path.str();
} else {
// Recurse into the lattice cell.
Lattice& lat = *lattices[c.fill_];
Lattice& lat = *model::lattices[c.fill_];
path << "l" << lat.id_;
for (ReverseLatticeIter it = lat.rbegin(); it != lat.rend(); ++it) {
int32_t indx = lat.universes_.size()*map + it.indx_;
@ -390,7 +390,7 @@ distribcell_path_inner(int32_t target_cell, int32_t map, int32_t target_offset,
offset = temp_offset;
path << "(" << lat.index_to_string(it.indx_) << ")->";
path << distribcell_path_inner(target_cell, map, target_offset,
*universes[*it], offset);
*model::universes[*it], offset);
return path.str();
}
}
@ -401,7 +401,7 @@ distribcell_path_inner(int32_t target_cell, int32_t map, int32_t target_offset,
std::string
distribcell_path(int32_t target_cell, int32_t map, int32_t target_offset)
{
auto& root_univ = *universes[openmc_root_universe];
auto& root_univ = *model::universes[model::root_universe];
return distribcell_path_inner(target_cell, map, target_offset, root_univ, 0);
}
@ -412,13 +412,13 @@ maximum_levels(int32_t univ)
{
int levels_below {0};
for (int32_t cell_indx : universes[univ]->cells_) {
Cell& c = *cells[cell_indx];
for (int32_t cell_indx : model::universes[univ]->cells_) {
Cell& c = *model::cells[cell_indx];
if (c.type_ == FILL_UNIVERSE) {
int32_t next_univ = c.fill_;
levels_below = std::max(levels_below, maximum_levels(next_univ));
} else if (c.type_ == FILL_LATTICE) {
Lattice& lat = *lattices[c.fill_];
Lattice& lat = *model::lattices[c.fill_];
for (auto it = lat.begin(); it != lat.end(); ++it) {
int32_t next_univ = *it;
levels_below = std::max(levels_below, maximum_levels(next_univ));
@ -435,20 +435,20 @@ maximum_levels(int32_t univ)
void
free_memory_geometry_c()
{
for (Cell* c : cells) {delete c;}
cells.clear();
cell_map.clear();
n_cells = 0;
for (Cell* c : model::cells) {delete c;}
model::cells.clear();
model::cell_map.clear();
model::n_cells = 0;
for (Universe* u : universes) {delete u;}
universes.clear();
universe_map.clear();
for (Universe* u : model::universes) {delete u;}
model::universes.clear();
model::universe_map.clear();
for (Lattice* lat : lattices) {delete lat;}
lattices.clear();
lattice_map.clear();
for (Lattice* lat : model::lattices) {delete lat;}
model::lattices.clear();
model::lattice_map.clear();
overlap_check_count.clear();
model::overlap_check_count.clear();
}
} // namespace openmc

2
src/geometry_header.F90
View file

@ -143,7 +143,7 @@ module geometry_header
end type Cell
! array index of the root universe
integer(C_INT), bind(C, name='openmc_root_universe') :: root_universe = -1
integer(C_INT), bind(C) :: root_universe
integer(C_INT32_T), bind(C) :: n_cells ! # of cells
integer(C_INT32_T), bind(C) :: n_universes ! # of universes

8
src/initialize.cpp
View file

@ -54,8 +54,8 @@ int openmc_init(int argc, char* argv[], const void* intracomm)
if (err) return err;
// Start total and initialization timer
time_total.start();
time_initialize.start();
simulation::time_total.start();
simulation::time_initialize.start();
#ifdef _OPENMP
// If OMP_SCHEDULE is not set, default to a static schedule
@ -79,7 +79,7 @@ int openmc_init(int argc, char* argv[], const void* intracomm)
if (settings::particle_restart_run) settings::run_mode = RUN_MODE_PARTICLE;
// Stop initialization timer
time_initialize.stop();
simulation::time_initialize.stop();
return 0;
}
@ -210,7 +210,7 @@ parse_command_line(int argc, char* argv[])
}
omp_set_num_threads(simulation::n_threads);
#else
if (openmc_master)
if (mpi::master)
warning("Ignoring number of threads specified on command line.");
#endif

8
src/input_xml.F90
View file

@ -347,8 +347,8 @@ contains
call read_surfaces(root % ptr)
! Allocate surfaces array
allocate(surfaces(n_surfaces))
do i = 1, n_surfaces
allocate(surfaces(surfaces_size()))
do i = 1, size(surfaces)
surfaces(i) % ptr = surface_pointer(i - 1);
if (surfaces(i) % bc() /= BC_TRANSMIT) boundary_exists = .true.
@ -505,9 +505,9 @@ contains
integer :: i
! Allocate surfaces array
allocate(surfaces(n_surfaces))
allocate(surfaces(surfaces_size()))
do i = 1, n_surfaces
do i = 1, size(surfaces)
surfaces(i) % ptr = surface_pointer(i - 1);
end do

37
src/lattice.cpp
View file

@ -18,10 +18,13 @@ namespace openmc {
// Global variables
//==============================================================================
std::vector<Lattice*> lattices;
namespace model {
std::vector<Lattice*> lattices;
std::unordered_map<int32_t, int32_t> lattice_map;
}
//==============================================================================
// Lattice implementation
//==============================================================================
@ -65,8 +68,8 @@ Lattice::adjust_indices()
// Adjust the indices for the universes array.
for (LatticeIter it = begin(); it != end(); ++it) {
int uid = *it;
auto search = universe_map.find(uid);
if (search != universe_map.end()) {
auto search = model::universe_map.find(uid);
if (search != model::universe_map.end()) {
*it = search->second;
} else {
std::stringstream err_msg;
@ -78,8 +81,8 @@ Lattice::adjust_indices()
// Adjust the index for the outer universe.
if (outer_ != NO_OUTER_UNIVERSE) {
auto search = universe_map.find(outer_);
if (search != universe_map.end()) {
auto search = model::universe_map.find(outer_);
if (search != model::universe_map.end()) {
outer_ = search->second;
} else {
std::stringstream err_msg;
@ -118,7 +121,7 @@ Lattice::to_hdf5(hid_t lattices_group) const
}
if (outer_ != NO_OUTER_UNIVERSE) {
int32_t outer_id = universes[outer_]->id_;
int32_t outer_id = model::universes[outer_]->id_;
write_dataset(lat_group, "outer", outer_id);
} else {
write_dataset(lat_group, "outer", outer_);
@ -370,7 +373,7 @@ RectLattice::to_hdf5_inner(hid_t lat_group) const
for (int j = 0; j < nx; j++) {
int indx1 = nx*ny*m + nx*k + j;
int indx2 = nx*ny*m + nx*(ny-k-1) + j;
out[indx2] = universes[universes_[indx1]]->id_;
out[indx2] = model::universes[universes_[indx1]]->id_;
}
}
}
@ -387,7 +390,7 @@ RectLattice::to_hdf5_inner(hid_t lat_group) const
for (int j = 0; j < nx; j++) {
int indx1 = nx*k + j;
int indx2 = nx*(ny-k-1) + j;
out[indx2] = universes[universes_[indx1]]->id_;
out[indx2] = model::universes[universes_[indx1]]->id_;
}
}
@ -847,7 +850,7 @@ HexLattice::to_hdf5_inner(hid_t lat_group) const
// This array position is never used; put a -1 to indicate this.
out[indx] = -1;
} else {
out[indx] = universes[universes_[indx]]->id_;
out[indx] = model::universes[universes_[indx]]->id_;
}
}
}
@ -865,18 +868,18 @@ extern "C" void
read_lattices(pugi::xml_node *node)
{
for (pugi::xml_node lat_node : node->children("lattice")) {
lattices.push_back(new RectLattice(lat_node));
model::lattices.push_back(new RectLattice(lat_node));
}
for (pugi::xml_node lat_node : node->children("hex_lattice")) {
lattices.push_back(new HexLattice(lat_node));
model::lattices.push_back(new HexLattice(lat_node));
}
// Fill the lattice map.
for (int i_lat = 0; i_lat < lattices.size(); i_lat++) {
int id = lattices[i_lat]->id_;
auto in_map = lattice_map.find(id);
if (in_map == lattice_map.end()) {
lattice_map[id] = i_lat;
for (int i_lat = 0; i_lat < model::lattices.size(); i_lat++) {
int id = model::lattices[i_lat]->id_;
auto in_map = model::lattice_map.find(id);
if (in_map == model::lattice_map.end()) {
model::lattice_map[id] = i_lat;
} else {
std::stringstream err_msg;
err_msg << "Two or more lattices use the same unique ID: " << id;
@ -890,7 +893,7 @@ read_lattices(pugi::xml_node *node)
//==============================================================================
extern "C" {
Lattice* lattice_pointer(int lat_ind) {return lattices[lat_ind];}
Lattice* lattice_pointer(int lat_ind) {return model::lattices[lat_ind];}
int32_t lattice_id(Lattice *lat) {return lat->id_;}
}

24
src/main.cpp
View file

@ -1,13 +1,15 @@
#ifdef OPENMC_MPI
#include "mpi.h"
#include <mpi.h>
#endif
#include "openmc/capi.h"
#include "openmc/constants.h"
#include "openmc/error.h"
#include "openmc/message_passing.h"
#include "openmc/settings.h"
int main(int argc, char* argv[]) {
using namespace openmc;
int err;
// Initialize run -- when run with MPI, pass communicator
@ -21,30 +23,30 @@ int main(int argc, char* argv[]) {
// This happens for the -h and -v flags
return 0;
} else if (err) {
openmc::fatal_error(openmc_err_msg);
fatal_error(openmc_err_msg);
}
// start problem based on mode
switch (openmc::settings::run_mode) {
case openmc::RUN_MODE_FIXEDSOURCE:
case openmc::RUN_MODE_EIGENVALUE:
switch (settings::run_mode) {
case RUN_MODE_FIXEDSOURCE:
case RUN_MODE_EIGENVALUE:
err = openmc_run();
break;
case openmc::RUN_MODE_PLOTTING:
case RUN_MODE_PLOTTING:
err = openmc_plot_geometry();
break;
case openmc::RUN_MODE_PARTICLE:
if (openmc_master) err = openmc_particle_restart();
case RUN_MODE_PARTICLE:
if (mpi::master) err = openmc_particle_restart();
break;
case openmc::RUN_MODE_VOLUME:
case RUN_MODE_VOLUME:
err = openmc_calculate_volumes();
break;
}
if (err) openmc::fatal_error(openmc_err_msg);
if (err) fatal_error(openmc_err_msg);
// Finalize and free up memory
err = openmc_finalize();
if (err) openmc::fatal_error(openmc_err_msg);
if (err) fatal_error(openmc_err_msg);
// If MPI is in use and enabled, terminate it
#ifdef OPENMC_MPI

40
src/material.cpp
View file

@ -13,9 +13,13 @@ namespace openmc {
// Global variables
//==============================================================================
namespace model {
std::vector<Material*> materials;
std::unordered_map<int32_t, int32_t> material_map;
} // namespace model
//==============================================================================
// Material implementation
//==============================================================================
@ -46,16 +50,16 @@ read_materials(pugi::xml_node* node)
{
// Loop over XML material elements and populate the array.
for (pugi::xml_node material_node : node->children("material")) {
materials.push_back(new Material(material_node));
model::materials.push_back(new Material(material_node));
}
materials.shrink_to_fit();
model::materials.shrink_to_fit();
// Populate the material map.
for (int i = 0; i < materials.size(); i++) {
int32_t mid = materials[i]->id_;
auto search = material_map.find(mid);
if (search == material_map.end()) {
material_map[mid] = i;
for (int i = 0; i < model::materials.size(); i++) {
int32_t mid = model::materials[i]->id_;
auto search = model::material_map.find(mid);
if (search == model::material_map.end()) {
model::material_map[mid] = i;
} else {
std::stringstream err_msg;
err_msg << "Two or more materials use the same unique ID: " << mid;
@ -71,8 +75,8 @@ read_materials(pugi::xml_node* node)
extern "C" int
openmc_material_get_volume(int32_t index, double* volume)
{
if (index >= 1 && index <= materials.size()) {
Material* m = materials[index - 1];
if (index >= 1 && index <= model::materials.size()) {
Material* m = model::materials[index - 1];
if (m->volume_ >= 0.0) {
*volume = m->volume_;
return 0;
@ -91,8 +95,8 @@ openmc_material_get_volume(int32_t index, double* volume)
extern "C" int
openmc_material_set_volume(int32_t index, double volume)
{
if (index >= 1 && index <= materials.size()) {
Material* m = materials[index - 1];
if (index >= 1 && index <= model::materials.size()) {
Material* m = model::materials[index - 1];
if (volume >= 0.0) {
m->volume_ = volume;
return 0;
@ -111,7 +115,7 @@ openmc_material_set_volume(int32_t index, double volume)
//==============================================================================
extern "C" {
Material* material_pointer(int32_t indx) {return materials[indx];}
Material* material_pointer(int32_t indx) {return model::materials[indx];}
int32_t material_id(Material* mat) {return mat->id_;}
@ -119,7 +123,7 @@ extern "C" {
{
mat->id_ = id;
//TODO: off-by-one
material_map[id] = index - 1;
model::material_map[id] = index - 1;
}
bool material_fissionable(Material* mat) {return mat->fissionable;}
@ -131,17 +135,17 @@ extern "C" {
void extend_materials_c(int32_t n)
{
materials.reserve(materials.size() + n);
model::materials.reserve(model::materials.size() + n);
for (int32_t i = 0; i < n; i++) {
materials.push_back(new Material());
model::materials.push_back(new Material());
}
}
void free_memory_material_c()
{
for (Material *mat : materials) {delete mat;}
materials.clear();
material_map.clear();
for (Material *mat : model::materials) {delete mat;}
model::materials.clear();
model::material_map.clear();
}
}

61
src/mesh.cpp
View file

@ -29,10 +29,13 @@ namespace openmc {
// Global variables
//==============================================================================
std::vector<std::unique_ptr<RegularMesh>> meshes;
namespace model {
std::vector<std::unique_ptr<RegularMesh>> meshes;
std::unordered_map<int32_t, int32_t> mesh_map;
} // namespace model
//==============================================================================
// RegularMesh implementation
//==============================================================================
@ -44,7 +47,7 @@ RegularMesh::RegularMesh(pugi::xml_node node)
id_ = std::stoi(get_node_value(node, "id"));
// Check to make sure 'id' hasn't been used
if (mesh_map.find(id_) != mesh_map.end()) {
if (model::mesh_map.find(id_) != model::mesh_map.end()) {
fatal_error("Two or more meshes use the same unique ID: " +
std::to_string(id_));
}
@ -724,11 +727,11 @@ xt::xarray<double> RegularMesh::count_sites(int64_t n, const Bank* bank,
extern "C" int
openmc_extend_meshes(int32_t n, int32_t* index_start, int32_t* index_end)
{
if (index_start) *index_start = meshes.size();
if (index_start) *index_start = model::meshes.size();
for (int i = 0; i < n; ++i) {
meshes.emplace_back(new RegularMesh{});
model::meshes.emplace_back(new RegularMesh{});
}
if (index_end) *index_end = meshes.size() - 1;
if (index_end) *index_end = model::meshes.size() - 1;
return 0;
}
@ -737,8 +740,8 @@ openmc_extend_meshes(int32_t n, int32_t* index_start, int32_t* index_end)
extern "C" int
openmc_get_mesh_index(int32_t id, int32_t* index)
{
auto pair = mesh_map.find(id);
if (pair == mesh_map.end()) {
auto pair = model::mesh_map.find(id);
if (pair == model::mesh_map.end()) {
set_errmsg("No mesh exists with ID=" + std::to_string(id) + ".");
return OPENMC_E_INVALID_ID;
}
@ -750,11 +753,11 @@ openmc_get_mesh_index(int32_t id, int32_t* index)
extern "C" int
openmc_mesh_get_id(int32_t index, int32_t* id)
{
if (index < 0 || index >= meshes.size()) {
if (index < 0 || index >= model::meshes.size()) {
set_errmsg("Index in meshes array is out of bounds.");
return OPENMC_E_OUT_OF_BOUNDS;
}
*id = meshes[index]->id_;
*id = model::meshes[index]->id_;
return 0;
}
@ -762,12 +765,12 @@ openmc_mesh_get_id(int32_t index, int32_t* id)
extern "C" int
openmc_mesh_set_id(int32_t index, int32_t id)
{
if (index < 0 || index >= meshes.size()) {
if (index < 0 || index >= model::meshes.size()) {
set_errmsg("Index in meshes array is out of bounds.");
return OPENMC_E_OUT_OF_BOUNDS;
}
meshes[index]->id_ = id;
mesh_map[id] = index;
model::meshes[index]->id_ = id;
model::mesh_map[id] = index;
return 0;
}
@ -775,12 +778,12 @@ openmc_mesh_set_id(int32_t index, int32_t id)
extern "C" int
openmc_mesh_get_dimension(int32_t index, int** dims, int* n)
{
if (index < 0 || index >= meshes.size()) {
if (index < 0 || index >= model::meshes.size()) {
set_errmsg("Index in meshes array is out of bounds.");
return OPENMC_E_OUT_OF_BOUNDS;
}
*dims = meshes[index]->shape_.data();
*n = meshes[index]->n_dimension_;
*dims = model::meshes[index]->shape_.data();
*n = model::meshes[index]->n_dimension_;
return 0;
}
@ -788,14 +791,14 @@ openmc_mesh_get_dimension(int32_t index, int** dims, int* n)
extern "C" int
openmc_mesh_set_dimension(int32_t index, int n, const int* dims)
{
if (index < 0 || index >= meshes.size()) {
if (index < 0 || index >= model::meshes.size()) {
set_errmsg("Index in meshes array is out of bounds.");
return OPENMC_E_OUT_OF_BOUNDS;
}
// Copy dimension
std::vector<std::size_t> shape = {static_cast<std::size_t>(n)};
auto& m = meshes[index];
auto& m = model::meshes[index];
m->shape_ = xt::adapt(dims, n, xt::no_ownership(), shape);
m->n_dimension_ = m->shape_.size();
@ -806,12 +809,12 @@ openmc_mesh_set_dimension(int32_t index, int n, const int* dims)
extern "C" int
openmc_mesh_get_params(int32_t index, double** ll, double** ur, double** width, int* n)
{
if (index < 0 || index >= meshes.size()) {
if (index < 0 || index >= model::meshes.size()) {
set_errmsg("Index in meshes array is out of bounds.");
return OPENMC_E_OUT_OF_BOUNDS;
}
auto& m = meshes[index];
auto& m = model::meshes[index];
if (m->lower_left_.dimension() == 0) {
set_errmsg("Mesh parameters have not been set.");
return OPENMC_E_ALLOCATE;
@ -829,12 +832,12 @@ extern "C" int
openmc_mesh_set_params(int32_t index, int n, const double* ll, const double* ur,
const double* width)
{
if (index < 0 || index >= meshes.size()) {
if (index < 0 || index >= model::meshes.size()) {
set_errmsg("Index in meshes array is out of bounds.");
return OPENMC_E_OUT_OF_BOUNDS;
}
auto& m = meshes[index];
auto& m = model::meshes[index];
std::vector<std::size_t> shape = {static_cast<std::size_t>(n)};
if (ll && ur) {
m->lower_left_ = xt::adapt(ll, n, xt::no_ownership(), shape);
@ -864,10 +867,10 @@ void read_meshes(pugi::xml_node* root)
{
for (auto node : root->children("mesh")) {
// Read mesh and add to vector
meshes.emplace_back(new RegularMesh{node});
model::meshes.emplace_back(new RegularMesh{node});
// Map ID to position in vector
mesh_map[meshes.back()->id_] = meshes.size() - 1;
model::mesh_map[model::meshes.back()->id_] = model::meshes.size() - 1;
}
}
@ -875,13 +878,13 @@ void meshes_to_hdf5(hid_t group)
{
// Write number of meshes
hid_t meshes_group = create_group(group, "meshes");
int32_t n_meshes = meshes.size();
int32_t n_meshes = model::meshes.size();
write_attribute(meshes_group, "n_meshes", n_meshes);
if (n_meshes > 0) {
// Write IDs of meshes
std::vector<int> ids;
for (const auto& m : meshes) {
for (const auto& m : model::meshes) {
m->to_hdf5(meshes_group);
ids.push_back(m->id_);
}
@ -896,9 +899,9 @@ void meshes_to_hdf5(hid_t group)
//==============================================================================
extern "C" {
int n_meshes() { return meshes.size(); }
int n_meshes() { return model::meshes.size(); }
RegularMesh* mesh_ptr(int i) { return meshes.at(i).get(); }
RegularMesh* mesh_ptr(int i) { return model::meshes.at(i).get(); }
int32_t mesh_id(RegularMesh* m) { return m->id_; }
@ -936,8 +939,8 @@ extern "C" {
void free_memory_mesh()
{
meshes.clear();
mesh_map.clear();
model::meshes.clear();
model::mesh_map.clear();
}
}

8
src/mgxs.cpp
View file

@ -22,10 +22,18 @@
namespace openmc {
//==============================================================================
// Global variables
//==============================================================================
namespace data {
// Storage for the MGXS data
std::vector<Mgxs> nuclides_MG;
std::vector<Mgxs> macro_xs;
} // namespace data
//==============================================================================
// Mgxs base-class methods
//==============================================================================

42
src/mgxs_interface.cpp
View file

@ -13,10 +13,14 @@ namespace openmc {
// Global variable definitions
//==============================================================================
namespace data {
std::vector<double> energy_bins;
std::vector<double> energy_bin_avg;
std::vector<double> rev_energy_bins;
} // namesapce data
//==============================================================================
// Mgxs data loading interface methods
//==============================================================================
@ -44,7 +48,7 @@ add_mgxs_c(hid_t file_id, const char* name, int energy_groups,
Mgxs mg(xs_grp, energy_groups, delayed_groups, temperature, tolerance,
max_order, legendre_to_tabular, legendre_to_tabular_points, method);
nuclides_MG.push_back(mg);
data::nuclides_MG.push_back(mg);
close_group(xs_grp);
}
@ -55,7 +59,7 @@ query_fissionable_c(int n_nuclides, const int i_nuclides[])
{
bool result = false;
for (int n = 0; n < n_nuclides; n++) {
if (nuclides_MG[i_nuclides[n] - 1].fissionable) result = true;
if (data::nuclides_MG[i_nuclides[n] - 1].fissionable) result = true;
}
return result;
}
@ -79,16 +83,16 @@ create_macro_xs_c(const char* mat_name, int n_nuclides, const int i_nuclides[],
// material
std::vector<Mgxs*> mgxs_ptr(n_nuclides);
for (int n = 0; n < n_nuclides; n++) {
mgxs_ptr[n] = &nuclides_MG[i_nuclides[n] - 1];
mgxs_ptr[n] = &data::nuclides_MG[i_nuclides[n] - 1];
}
Mgxs macro(mat_name, temperature, mgxs_ptr, atom_densities_vec,
tolerance, method);
macro_xs.emplace_back(macro);
data::macro_xs.emplace_back(macro);
} else {
// Preserve the ordering of materials by including a blank entry
Mgxs macro;
macro_xs.emplace_back(macro);
data::macro_xs.emplace_back(macro);
}
}
@ -97,18 +101,18 @@ create_macro_xs_c(const char* mat_name, int n_nuclides, const int i_nuclides[],
void read_mg_cross_sections_header_c(hid_t file_id)
{
ensure_exists(file_id, "energy_groups", true);
read_attribute(file_id, "energy_groups", num_energy_groups);
read_attribute(file_id, "energy_groups", data::num_energy_groups);
ensure_exists(file_id, "group structure", true);
read_attribute(file_id, "group structure", rev_energy_bins);
read_attribute(file_id, "group structure", data::rev_energy_bins);
// Reverse energy bins
std::copy(rev_energy_bins.crbegin(), rev_energy_bins.crend(),
std::back_inserter(energy_bins));
std::copy(data::rev_energy_bins.crbegin(), data::rev_energy_bins.crend(),
std::back_inserter(data::energy_bins));
// Create average energies
for (int i = 0; i < energy_bins.size() - 1; ++i) {
energy_bin_avg.push_back(0.5*(energy_bins[i] + energy_bins[i+1]));
for (int i = 0; i < data::energy_bins.size() - 1; ++i) {
data::energy_bin_avg.push_back(0.5*(data::energy_bins[i] + data::energy_bins[i+1]));
}
// Add entries into libraries for MG data
@ -134,7 +138,7 @@ void
calculate_xs_c(int i_mat, int gin, double sqrtkT, const double uvw[3],
double& total_xs, double& abs_xs, double& nu_fiss_xs)
{
macro_xs[i_mat - 1].calculate_xs(gin - 1, sqrtkT, uvw, total_xs, abs_xs,
data::macro_xs[i_mat - 1].calculate_xs(gin - 1, sqrtkT, uvw, total_xs, abs_xs,
nu_fiss_xs);
}
@ -159,7 +163,7 @@ get_nuclide_xs_c(int index, int xstype, int gin, int* gout, double* mu, int* dg)
} else {
dg_c_p = dg;
}
return nuclides_MG[index - 1].get_xs(xstype, gin - 1, gout_c_p, mu, dg_c_p);
return data::nuclides_MG[index - 1].get_xs(xstype, gin - 1, gout_c_p, mu, dg_c_p);
}
//==============================================================================
@ -183,7 +187,7 @@ get_macro_xs_c(int index, int xstype, int gin, int* gout, double* mu, int* dg)
} else {
dg_c_p = dg;
}
return macro_xs[index - 1].get_xs(xstype, gin - 1, gout_c_p, mu, dg_c_p);
return data::macro_xs[index - 1].get_xs(xstype, gin - 1, gout_c_p, mu, dg_c_p);
}
//==============================================================================
@ -192,7 +196,7 @@ void
set_nuclide_angle_index_c(int index, const double uvw[3])
{
// Update the values
nuclides_MG[index - 1].set_angle_index(uvw);
data::nuclides_MG[index - 1].set_angle_index(uvw);
}
//==============================================================================
@ -201,7 +205,7 @@ void
set_macro_angle_index_c(int index, const double uvw[3])
{
// Update the values
macro_xs[index - 1].set_angle_index(uvw);
data::macro_xs[index - 1].set_angle_index(uvw);
}
//==============================================================================
@ -210,7 +214,7 @@ void
set_nuclide_temperature_index_c(int index, double sqrtkT)
{
// Update the values
nuclides_MG[index - 1].set_temperature_index(sqrtkT);
data::nuclides_MG[index - 1].set_temperature_index(sqrtkT);
}
//==============================================================================
@ -225,7 +229,7 @@ get_name_c(int index, int name_len, char* name)
std::strcpy(name, str.c_str());
// Now get the data and copy to the C-string
str = nuclides_MG[index - 1].name;
str = data::nuclides_MG[index - 1].name;
std::strcpy(name, str.c_str());
// Finally, remove the null terminator
@ -237,7 +241,7 @@ get_name_c(int index, int name_len, char* name)
double
get_awr_c(int index)
{
return nuclides_MG[index - 1].awr;
return data::nuclides_MG[index - 1].awr;
}
} // namespace openmc

8
src/nuclide.cpp
View file

@ -6,9 +6,13 @@ namespace openmc {
// Global variables
//==============================================================================
namespace data {
std::array<double, 2> energy_min {0.0, 0.0};
std::array<double, 2> energy_max {INFTY, INFTY};
} // namespace data
//==============================================================================
// Fortran compatibility functions
//==============================================================================
@ -16,8 +20,8 @@ std::array<double, 2> energy_max {INFTY, INFTY};
extern "C" void
set_particle_energy_bounds(int particle, double E_min, double E_max)
{
energy_min[particle - 1] = E_min;
energy_max[particle - 1] = E_max;
data::energy_min[particle - 1] = E_min;
data::energy_max[particle - 1] = E_max;
}
} // namespace openmc

32
src/output.cpp
View file

@ -51,7 +51,7 @@ std::string time_stamp()
<< ":" << now->tm_min << ":" << now->tm_sec;
return ts.str();
}
//==============================================================================
//===============================================================================
@ -62,7 +62,7 @@ void print_plot() {
header("PLOTTING SUMMARY", 5);
for (auto pl : plots) {
for (auto pl : model::plots) {
// Plot id
std::cout << "Plot ID: " << pl.id_ << "\n";
// Plot filename
@ -74,7 +74,7 @@ void print_plot() {
if (PlotType::slice == pl.type_) {
std::cout << "Plot Type: Slice" << "\n";
} else if (PlotType::voxel == pl.type_) {
std::cout << "Plot Type: Voxel" << "\n";
std::cout << "Plot Type: Voxel" << "\n";
}
// Plot parameters
@ -98,9 +98,9 @@ void print_plot() {
if (PlotColorBy::cells == pl.color_by_) {
std::cout << "Coloring: Cells" << "\n";
} else if (PlotColorBy::mats == pl.color_by_) {
std::cout << "Coloring: Materials" << "\n";
std::cout << "Coloring: Materials" << "\n";
}
if (PlotType::slice == pl.type_) {
switch(pl.basis_) {
case PlotBasis::xy:
@ -122,29 +122,29 @@ void print_plot() {
}
std::cout << "\n";
}
}
void
print_overlap_check() {
#ifdef OPENMC_MPI
std::vector<int64_t> temp(overlap_check_count);
int err = MPI_Reduce(temp.data(), overlap_check_count.data(),
overlap_check_count.size(), MPI_INT64_T, MPI_SUM, 0,
std::vector<int64_t> temp(model::overlap_check_count);
int err = MPI_Reduce(temp.data(), model::overlap_check_count.data(),
model::overlap_check_count.size(), MPI_INT64_T, MPI_SUM, 0,
mpi::intracomm);
#endif
if (openmc_master) {
if (mpi::master) {
header("cell overlap check summary", 1);
std::cout << " Cell ID No. Overlap Checks\n";
std::vector<int32_t> sparse_cell_ids;
for (int i = 0; i < n_cells; i++) {
std::cout << " " << std::setw(8) << cells[i]->id_ << std::setw(17)
<< overlap_check_count[i] << "\n";
if (overlap_check_count[i] < 10) {
sparse_cell_ids.push_back(cells[i]->id_);
for (int i = 0; i < model::cells.size(); i++) {
std::cout << " " << std::setw(8) << model::cells[i]->id_ << std::setw(17)
<< model::overlap_check_count[i] << "\n";
if (model::overlap_check_count[i] < 10) {
sparse_cell_ids.push_back(model::cells[i]->id_);
}
}

2
src/particle.cpp
View file

@ -113,7 +113,7 @@ Particle::from_source(const Bank* src)
} else {
g = static_cast<int>(src->E);
last_g = static_cast<int>(src->E);
E = energy_bin_avg[g - 1];
E = data::energy_bin_avg[g - 1];
}
last_E = E;
}

22
src/physics_mg.cpp
View file

@ -21,14 +21,13 @@
namespace openmc {
void
collision_mg(Particle* p, const double* energy_bin_avg,
const MaterialMacroXS* material_xs)
collision_mg(Particle* p, const MaterialMacroXS* material_xs)
{
// Add to the collision counter for the particle
p->n_collision++;
// Sample the reaction type
sample_reaction(p, energy_bin_avg, material_xs);
sample_reaction(p, material_xs);
// Display information about collision
if ((settings::verbosity >= 10) || (simulation::trace)) {
@ -39,21 +38,20 @@ collision_mg(Particle* p, const double* energy_bin_avg,
}
void
sample_reaction(Particle* p, const double* energy_bin_avg,
const MaterialMacroXS* material_xs)
sample_reaction(Particle* p, const MaterialMacroXS* material_xs)
{
// Create fission bank sites. Note that while a fission reaction is sampled,
// it never actually "happens", i.e. the weight of the particle does not
// change when sampling fission sites. The following block handles all
// absorption (including fission)
if (materials[p->material - 1]->fissionable) {
if (model::materials[p->material - 1]->fissionable) {
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
Bank* result_bank;
int64_t result_bank_size;
// Get pointer to fission bank from Fortran side
openmc_fission_bank(&result_bank, &result_bank_size);
create_fission_sites(p, result_bank, &n_bank, result_bank_size,
create_fission_sites(p, result_bank, &simulation::n_bank, result_bank_size,
material_xs);
} else if ((settings::run_mode == RUN_MODE_FIXEDSOURCE) &&
(settings::create_fission_neutrons)) {
@ -72,7 +70,7 @@ sample_reaction(Particle* p, const double* energy_bin_avg,
if (!p->alive) return;
// Sample a scattering event to determine the energy of the exiting neutron
scatter(p, energy_bin_avg);
scatter(p);
// Play Russian roulette if survival biasing is turned on
if (settings::survival_biasing) {
@ -82,14 +80,14 @@ sample_reaction(Particle* p, const double* energy_bin_avg,
}
void
scatter(Particle* p, const double* energy_bin_avg)
scatter(Particle* p)
{
// Adjust indices for Fortran to C++ indexing
// TODO: Remove when no longer needed
int gin = p->last_g - 1;
int gout = p->g - 1;
int i_mat = p->material - 1;
macro_xs[i_mat].sample_scatter(gin, gout, p->mu, p->wgt);
data::macro_xs[i_mat].sample_scatter(gin, gout, p->mu, p->wgt);
// Adjust return value for fortran indexing
// TODO: Remove when no longer needed
@ -99,7 +97,7 @@ scatter(Particle* p, const double* energy_bin_avg)
rotate_angle_c(p->coord[0].uvw, p->mu, nullptr);
// Update energy value for downstream compatability (in tallying)
p->E = energy_bin_avg[gout];
p->E = data::energy_bin_avg[gout];
// Set event component
p->event = EVENT_SCATTER;
@ -180,7 +178,7 @@ create_fission_sites(Particle* p, Bank* bank_array, int64_t* size_bank,
// the energy in the fission bank
int dg;
int gout;
macro_xs[p->material - 1].sample_fission_energy(p->g - 1, dg, gout);
data::macro_xs[p->material - 1].sample_fission_energy(p->g - 1, dg, gout);
bank_array[i].E = static_cast<double>(gout + 1);
bank_array[i].delayed_group = dg + 1;

75
src/plot.cpp
View file

@ -9,6 +9,7 @@
#include "openmc/geometry.h"
#include "openmc/cell.h"
#include "openmc/material.h"
#include "openmc/message_passing.h"
#include "openmc/string_utils.h"
#include "openmc/mesh.h"
#include "openmc/output.h"
@ -19,19 +20,23 @@
namespace openmc {
//==============================================================================
// Constants
//==============================================================================
const RGBColor WHITE {255, 255, 255};
constexpr int PLOT_LEVEL_LOWEST {-1}; //!< lower bound on plot universe level
//==============================================================================
// Global variables
//==============================================================================
int PLOT_LEVEL_LOWEST = -1;
std::unordered_map<int, int> plot_map;
int n_plots;
namespace model {
std::vector<Plot> plots;
std::unordered_map<int, int> plot_map;
const RGBColor WHITE = {255, 255, 255};
} // namespace model
//==============================================================================
// RUN_PLOT controls the logic for making one or many plots
@ -42,7 +47,7 @@ int openmc_plot_geometry()
{
int err;
for (auto pl : plots) {
for (auto pl : model::plots) {
std::stringstream ss;
ss << "Processing plot " << pl.id_ << ": "
<< pl.path_plot_ << "...";
@ -63,11 +68,10 @@ int openmc_plot_geometry()
void
read_plots(pugi::xml_node* plots_node)
{
n_plots = 0;
for (auto node : plots_node->children("plot")) {
Plot pl(node);
plots.push_back(pl);
plot_map[pl.id_] = n_plots++;
model::plots.push_back(pl);
model::plot_map[pl.id_] = model::plots.size() - 1;
}
}
@ -122,7 +126,7 @@ void create_ppm(Plot pl)
p.initialize();
std::copy(xyz, xyz+3, p.coord[0].xyz);
std::copy(dir, dir+3, p.coord[0].uvw);
p.coord[0].universe = openmc_root_universe;
p.coord[0].universe = model::root_universe;
#pragma omp for
for (int y = 0; y < height; y++) {
@ -155,7 +159,7 @@ Plot::set_id(pugi::xml_node plot_node)
}
// Check to make sure 'id' hasn't been used
if (plot_map.find(id_) != plot_map.end()) {
if (model::plot_map.find(id_) != model::plot_map.end()) {
std::stringstream err_msg;
err_msg << "Two or more plots use the same unique ID: " << id_;
fatal_error(err_msg.str());
@ -243,7 +247,7 @@ Plot::set_bg_color(pugi::xml_node plot_node)
if (check_for_node(plot_node, "background")) {
std::vector<int> bg_rgb = get_node_array<int>(plot_node, "background");
if (PlotType::voxel == type_) {
if (openmc_master) {
if (mpi::master) {
std::stringstream err_msg;
err_msg << "Background color ignored in voxel plot "
<< id_;
@ -361,30 +365,27 @@ Plot::set_default_colors(pugi::xml_node plot_node)
}
if ("cell" == pl_color_by) {
color_by_ = PlotColorBy::cells;
colors_.resize(n_cells);
for (int i = 0; i < n_cells; i++) {
colors_[i] = random_color();
}
colors_.resize(model::cells.size());
} else if("material" == pl_color_by) {
color_by_ = PlotColorBy::mats;
colors_.resize(n_materials);
for (int i = 0; i < materials.size(); i++) {
colors_[i] = random_color();
}
colors_.resize(model::materials.size());
} else {
std::stringstream err_msg;
err_msg << "Unsupported plot color type '" << pl_color_by
<< "' in plot " << id_;
fatal_error(err_msg);
}
for (auto& c : colors_) {
c = random_color();
}
}
void
Plot::set_user_colors(pugi::xml_node plot_node)
{
if (!plot_node.select_nodes("color").empty() && PlotType::voxel == type_) {
if (openmc_master) {
if (mpi::master) {
std::stringstream err_msg;
err_msg << "Color specifications ignored in voxel plot "
<< id_;
@ -412,8 +413,8 @@ Plot::set_user_colors(pugi::xml_node plot_node)
}
// Add RGB
if (PlotColorBy::cells == color_by_) {
if (cell_map.find(col_id) != cell_map.end()) {
col_id = cell_map[col_id];
if (model::cell_map.find(col_id) != model::cell_map.end()) {
col_id = model::cell_map[col_id];
colors_[col_id] = user_rgb;
} else {
std::stringstream err_msg;
@ -422,8 +423,8 @@ Plot::set_user_colors(pugi::xml_node plot_node)
fatal_error(err_msg);
}
} else if (PlotColorBy::mats == color_by_) {
if (material_map.find(col_id) != material_map.end()) {
col_id = material_map[col_id];
if (model::material_map.find(col_id) != model::material_map.end()) {
col_id = model::material_map[col_id];
colors_[col_id] = user_rgb;
} else {
std::stringstream err_msg;
@ -552,7 +553,7 @@ Plot::set_mask(pugi::xml_node plot_node)
if (!mask_nodes.empty()) {
if (PlotType::voxel == type_) {
if (openmc_master) {
if (mpi::master) {
std::stringstream wrn_msg;
wrn_msg << "Mask ignored in voxel plot " << id_;
warning(wrn_msg);
@ -575,8 +576,8 @@ Plot::set_mask(pugi::xml_node plot_node)
// in the cell and material arrays
for (auto& col_id : iarray) {
if (PlotColorBy::cells == color_by_) {
if (cell_map.find(col_id) != cell_map.end()) {
col_id = cell_map[col_id];
if (model::cell_map.find(col_id) != model::cell_map.end()) {
col_id = model::cell_map[col_id];
}
else {
std::stringstream err_msg;
@ -585,8 +586,8 @@ Plot::set_mask(pugi::xml_node plot_node)
fatal_error(err_msg);
}
} else if (PlotColorBy::mats == color_by_) {
if (material_map.find(col_id) != material_map.end()) {
col_id = material_map[col_id];
if (model::material_map.find(col_id) != model::material_map.end()) {
col_id = model::material_map[col_id];
}
else {
std::stringstream err_msg;
@ -660,7 +661,7 @@ void position_rgb(Particle p, Plot pl, RGBColor& rgb, int& id)
} else {
if (PlotColorBy::mats == pl.color_by_) {
// Assign color based on material
Cell* c = cells[p.coord[j].cell];
Cell* c = model::cells[p.coord[j].cell];
if (c->type_ == FILL_UNIVERSE) {
// If we stopped on a middle universe level, treat as if not found
rgb = pl.not_found_;
@ -671,12 +672,12 @@ void position_rgb(Particle p, Plot pl, RGBColor& rgb, int& id)
id = -1;
} else {
rgb = pl.colors_[p.material - 1];
id = materials[p.material - 1]->id_;
id = model::materials[p.material - 1]->id_;
}
} else if (PlotColorBy::cells == pl.color_by_) {
// Assign color based on cell
rgb = pl.colors_[p.coord[j].cell];
id = cells[p.coord[j].cell]->id_;
id = model::cells[p.coord[j].cell]->id_;
}
} // endif found_cell
}
@ -759,7 +760,7 @@ void draw_mesh_lines(Plot pl, ImageData& data)
width[1] = xyz_ur_plot[1] - xyz_ll_plot[1];
width[2] = xyz_ur_plot[2] - xyz_ll_plot[2];
auto& m = meshes[pl.index_meshlines_mesh_];
auto& m = model::meshes[pl.index_meshlines_mesh_];
int ijk_ll[3], ijk_ur[3];
bool in_mesh;
@ -851,7 +852,7 @@ void create_voxel(Plot pl)
p.initialize();
std::copy(ll.begin(), ll.begin()+ll.size(), p.coord[0].xyz);
std::copy(dir, dir+3, p.coord[0].uvw);
p.coord[0].universe = openmc_root_universe;
p.coord[0].universe = model::root_universe;
// Open binary plot file for writing
std::ofstream of;

39
src/settings.cpp
View file

@ -17,6 +17,7 @@
#include "openmc/error.h"
#include "openmc/file_utils.h"
#include "openmc/mesh.h"
#include "openmc/message_passing.h"
#include "openmc/output.h"
#include "openmc/random_lcg.h"
#include "openmc/simulation.h"
@ -150,7 +151,7 @@ void get_run_parameters(pugi::xml_node node_base)
// Preallocate space for keff and entropy by generation
int m = settings::n_max_batches * settings::gen_per_batch;
simulation::k_generation.reserve(m);
entropy.reserve(m);
simulation::entropy.reserve(m);
// Get the trigger information for keff
if (check_for_node(node_base, "keff_trigger")) {
@ -231,7 +232,7 @@ void read_settings_xml()
// To this point, we haven't displayed any output since we didn't know what
// the verbosity is. Now that we checked for it, show the title if necessary
if (openmc_master) {
if (mpi::master) {
if (verbosity >= 2) title();
}
write_message("Reading settings XML file...", 5);
@ -398,7 +399,7 @@ void read_settings_xml()
omp_set_num_threads(simulation::n_threads);
}
#else
if (openmc_master) warning("OpenMC was not compiled with OpenMP support; "
if (mpi::master) warning("OpenMC was not compiled with OpenMP support; "
"ignoring number of threads.");
#endif
}
@ -415,17 +416,17 @@ void read_settings_xml()
// Get point to list of <source> elements and make sure there is at least one
for (pugi::xml_node node : root.children("source")) {
external_sources.emplace_back(node);
model::external_sources.emplace_back(node);
}
// If no source specified, default to isotropic point source at origin with Watt spectrum
if (external_sources.empty()) {
if (model::external_sources.empty()) {
SourceDistribution source {
UPtrSpace{new SpatialPoint({0.0, 0.0, 0.0})},
UPtrAngle{new Isotropic()},
UPtrDist{new Watt(0.988, 2.249e-6)}
};
external_sources.push_back(std::move(source));
model::external_sources.push_back(std::move(source));
}
// Check if we want to write out source
@ -506,12 +507,12 @@ void read_settings_xml()
// Shannon Entropy mesh
if (check_for_node(root, "entropy_mesh")) {
int temp = std::stoi(get_node_value(root, "entropy_mesh"));
if (mesh_map.find(temp) == mesh_map.end()) {
if (model::mesh_map.find(temp) == model::mesh_map.end()) {
std::stringstream msg;
msg << "Mesh " << temp << " specified for Shannon entropy does not exist.";
fatal_error(msg);
}
index_entropy_mesh = mesh_map.at(temp);
index_entropy_mesh = model::mesh_map.at(temp);
} else if (check_for_node(root, "entropy")) {
warning("Specifying a Shannon entropy mesh via the <entropy> element "
@ -520,18 +521,18 @@ void read_settings_xml()
// Read entropy mesh from <entropy>
auto node_entropy = root.child("entropy");
meshes.emplace_back(new RegularMesh{node_entropy});
model::meshes.emplace_back(new RegularMesh{node_entropy});
// Set entropy mesh index
index_entropy_mesh = meshes.size() - 1;
index_entropy_mesh = model::meshes.size() - 1;
// Assign ID and set mapping
meshes.back()->id_ = 10000;
mesh_map[10000] = index_entropy_mesh;
model::meshes.back()->id_ = 10000;
model::mesh_map[10000] = index_entropy_mesh;
}
if (index_entropy_mesh >= 0) {
auto& m = *meshes[index_entropy_mesh];
auto& m = *model::meshes[index_entropy_mesh];
if (m.shape_.dimension() == 0) {
// If the user did not specify how many mesh cells are to be used in
// each direction, we automatically determine an appropriate number of
@ -551,13 +552,13 @@ void read_settings_xml()
// Uniform fission source weighting mesh
if (check_for_node(root, "ufs_mesh")) {
auto temp = std::stoi(get_node_value(root, "ufs_mesh"));
if (mesh_map.find(temp) == mesh_map.end()) {
if (model::mesh_map.find(temp) == model::mesh_map.end()) {
std::stringstream msg;
msg << "Mesh " << temp << " specified for uniform fission site method "
"does not exist.";
fatal_error(msg);
}
index_ufs_mesh = mesh_map.at(temp);
index_ufs_mesh = model::mesh_map.at(temp);
} else if (check_for_node(root, "uniform_fs")) {
warning("Specifying a UFS mesh via the <uniform_fs> element "
@ -566,14 +567,14 @@ void read_settings_xml()
// Read entropy mesh from <entropy>
auto node_ufs = root.child("uniform_fs");
meshes.emplace_back(new RegularMesh{node_ufs});
model::meshes.emplace_back(new RegularMesh{node_ufs});
// Set entropy mesh index
index_ufs_mesh = meshes.size() - 1;
index_ufs_mesh = model::meshes.size() - 1;
// Assign ID and set mapping
meshes.back()->id_ = 10001;
mesh_map[10001] = index_entropy_mesh;
model::meshes.back()->id_ = 10001;
model::mesh_map[10001] = index_entropy_mesh;
}
if (index_ufs_mesh >= 0) {

43
src/simulation.cpp
View file

@ -25,10 +25,10 @@ extern "C" bool cmfd_on;
extern "C" void accumulate_tallies();
extern "C" void allocate_banks();
extern "C" void allocate_tally_results();
extern "C" void check_triggers();
extern "C" void cmfd_init_batch();
extern "C" void cmfd_tally_init();
extern "C" void configure_tallies();
extern "C" void execute_cmfd();
extern "C" void init_tally_routines();
extern "C" void join_bank_from_threads();
@ -81,7 +81,7 @@ int openmc_simulation_init()
// Allocate array for matching filter bins
#pragma omp parallel
{
filter_matches.resize(n_filters);
simulation::filter_matches.resize(model::tally_filters.size());
}
// Set up tally procedure pointers
@ -92,7 +92,7 @@ int openmc_simulation_init()
allocate_banks();
// Allocate tally results arrays if they're not allocated yet
configure_tallies();
allocate_tally_results();
// Activate the CMFD tallies
cmfd_tally_init();
@ -104,8 +104,9 @@ int openmc_simulation_init()
// will potentially populate k_generation and entropy)
simulation::current_batch = 0;
simulation::k_generation.clear();
entropy.clear();
simulation::entropy.clear();
simulation::need_depletion_rx = false;
openmc_reset();
// If this is a restart run, load the state point data and binary source
// file
@ -139,12 +140,12 @@ int openmc_simulation_finalize()
if (!simulation::initialized) return 0;
// Stop active batch timer and start finalization timer
time_active.stop();
time_finalize.start();
simulation::time_active.stop();
simulation::time_finalize.start();
#pragma omp parallel
{
filter_matches.clear();
simulation::filter_matches.clear();
}
// Deallocate Fortran variables, set tallies to inactive
@ -166,8 +167,8 @@ int openmc_simulation_finalize()
}
// Stop timers and show timing statistics
time_finalize.stop();
time_total.stop();
simulation::time_finalize.stop();
simulation::time_total.stop();
if (mpi::master) {
if (settings::verbosity >= 6) print_runtime();
if (settings::verbosity >= 4) print_results();
@ -200,7 +201,7 @@ int openmc_next_batch(int* status)
initialize_generation();
// Start timer for transport
time_transport.start();
simulation::time_transport.start();
// ====================================================================
// LOOP OVER PARTICLES
@ -218,7 +219,7 @@ int openmc_next_batch(int* status)
}
// Accumulate time for transport
time_transport.stop();
simulation::time_transport.stop();
finalize_generation();
}
@ -305,10 +306,10 @@ void initialize_batch()
// Manage active/inactive timers and activate tallies if necessary.
if (first_inactive) {
time_inactive.start();
simulation::time_inactive.start();
} else if (first_active) {
time_inactive.stop();
time_active.start();
simulation::time_inactive.stop();
simulation::time_active.start();
for (int i = 1; i <= n_tallies; ++i) {
// TODO: change one-based index
openmc_tally_set_active(i, true);
@ -327,9 +328,9 @@ void initialize_batch()
void finalize_batch()
{
// Reduce tallies onto master process and accumulate
time_tallies.start();
simulation::time_tallies.start();
accumulate_tallies();
time_tallies.stop();
simulation::time_tallies.stop();
// Reset global tally results
if (simulation::current_batch <= settings::n_inactive) {
@ -384,13 +385,13 @@ void initialize_generation()
{
if (settings::run_mode == RUN_MODE_EIGENVALUE) {
// Reset number of fission bank sites
n_bank = 0;
simulation::n_bank = 0;
// Count source sites if using uniform fission source weighting
if (settings::ufs_on) ufs_count_sites();
// Store current value of tracklength k
keff_generation = global_tallies()(K_TRACKLENGTH, RESULT_VALUE);
simulation::keff_generation = global_tallies()(K_TRACKLENGTH, RESULT_VALUE);
}
}
@ -564,4 +565,10 @@ extern "C" void k_generation_clear() { simulation::k_generation.clear(); }
extern "C" void k_generation_reserve(int i) { simulation::k_generation.reserve(i); }
extern "C" int64_t work_index(int rank) { return simulation::work_index[rank]; }
// This function was moved here to get around a bug on macOS whereby an invalid
// pointer is returned for the threadprivate filter_matches
extern "C" FilterMatch* filter_match_pointer(int indx) {
return &simulation::filter_matches[indx];
}
} // namespace openmc

34
src/source.cpp
View file

@ -27,8 +27,12 @@ namespace openmc {
// Global variables
//==============================================================================
namespace model {
std::vector<SourceDistribution> external_sources;
}
//==============================================================================
// SourceDistribution implementation
//==============================================================================
@ -167,9 +171,9 @@ Bank SourceDistribution::sample() const
if (space_box) {
if (space_box->only_fissionable()) {
// Determine material
auto c = cells[cell_index - 1];
auto c = model::cells[cell_index - 1];
int32_t mat_index = c->material_[instance];
auto m = materials[mat_index];
auto m = model::materials[mat_index];
if (mat_index == MATERIAL_VOID) {
found = false;
@ -207,10 +211,10 @@ Bank SourceDistribution::sample() const
auto energy_ptr = dynamic_cast<Discrete*>(energy_.get());
if (energy_ptr) {
auto energies = xt::adapt(energy_ptr->x());
if (xt::any(energies > energy_max[p-1])) {
if (xt::any(energies > data::energy_max[p-1])) {
fatal_error("Source energy above range of energies of at least "
"one cross section table");
} else if (xt::any(energies < energy_min[p-1])) {
} else if (xt::any(energies < data::energy_min[p-1])) {
fatal_error("Source energy below range of energies of at least "
"one cross section table");
}
@ -221,7 +225,7 @@ Bank SourceDistribution::sample() const
site.E = energy_->sample();
// Resample if energy falls outside minimum or maximum particle energy
if (site.E < energy_max[p-1] && site.E > energy_min[p-1]) break;
if (site.E < data::energy_max[p-1] && site.E > data::energy_min[p-1]) break;
}
// Set delayed group
@ -299,28 +303,28 @@ Bank sample_external_source()
// Determine total source strength
double total_strength = 0.0;
for (auto& s : external_sources)
for (auto& s : model::external_sources)
total_strength += s.strength();
// Sample from among multiple source distributions
int i = 0;
if (external_sources.size() > 1) {
if (model::external_sources.size() > 1) {
double xi = prn()*total_strength;
double c = 0.0;
for (; i < external_sources.size(); ++i) {
c += external_sources[i].strength();
for (; i < model::external_sources.size(); ++i) {
c += model::external_sources[i].strength();
if (xi < c) break;
}
}
// Sample source site from i-th source distribution
Bank site {external_sources[i].sample()};
Bank site {model::external_sources[i].sample()};
// If running in MG, convert site % E to group
if (!settings::run_CE) {
site.E = lower_bound_index(rev_energy_bins.begin(), rev_energy_bins.end(),
site.E);
site.E = num_energy_groups - site.E;
site.E = lower_bound_index(data::rev_energy_bins.begin(),
data::rev_energy_bins.end(), site.E);
site.E = data::num_energy_groups - site.E;
}
// Set the random number generator back to the tracking stream.
@ -335,13 +339,13 @@ Bank sample_external_source()
extern "C" void free_memory_source()
{
external_sources.clear();
model::external_sources.clear();
}
extern "C" double total_source_strength()
{
double strength = 0.0;
for (const auto& s : external_sources) {
for (const auto& s : model::external_sources) {
strength += s.strength();
}
return strength;

4
src/state_point.F90
View file

@ -55,7 +55,7 @@ contains
!===============================================================================
function openmc_statepoint_write(filename, write_source) result(err) bind(C)
type(C_PTR), intent(in), optional :: filename
type(C_PTR), value :: filename
logical(C_BOOL), intent(in), optional :: write_source
integer(C_INT) :: err
@ -91,7 +91,7 @@ contains
err = 0
! Set the filename
if (present(filename)) then
if (c_associated(filename)) then
call c_f_pointer(filename, string, [MAX_FILE_LEN])
filename_ = to_f_string(string)
else

8
src/state_point.cpp
View file

@ -113,7 +113,7 @@ write_source_bank(hid_t group_id, Bank* source_bank)
#else
if (openmc_master) {
if (mpi::master) {
// Create dataset big enough to hold all source sites
hsize_t dims[] {static_cast<hsize_t>(settings::n_particles)};
hid_t dspace = H5Screate_simple(1, dims, nullptr);
@ -320,11 +320,11 @@ void restart_set_keff()
{
if (simulation::restart_batch > settings::n_inactive) {
for (int i = settings::n_inactive; i < simulation::restart_batch; ++i) {
k_sum[0] += simulation::k_generation[i];
k_sum[1] += std::pow(simulation::k_generation[i], 2);
simulation::k_sum[0] += simulation::k_generation[i];
simulation::k_sum[1] += std::pow(simulation::k_generation[i], 2);
}
int n = settings::gen_per_batch*n_realizations;
simulation::keff = k_sum[0] / n;
simulation::keff = simulation::k_sum[0] / n;
} else {
simulation::keff = simulation::k_generation.back();
}

18
src/summary.cpp
View file

@ -10,7 +10,7 @@ extern "C" void
write_geometry(hid_t file_id) {
auto geom_group = create_group(file_id, "geometry");
#ifdef DAGMC
if (settings::dagmc) {
write_attribute(geom_group, "dagmc", 1);
@ -18,25 +18,25 @@ write_geometry(hid_t file_id) {
}
#endif
write_attribute(geom_group, "n_cells", cells.size());
write_attribute(geom_group, "n_surfaces", surfaces.size());
write_attribute(geom_group, "n_universes", universes.size());
write_attribute(geom_group, "n_lattices", lattices.size());
write_attribute(geom_group, "n_cells", model::cells.size());
write_attribute(geom_group, "n_surfaces", model::surfaces.size());
write_attribute(geom_group, "n_universes", model::universes.size());
write_attribute(geom_group, "n_lattices", model::lattices.size());
auto cells_group = create_group(geom_group, "cells");
for (Cell* c : cells) c->to_hdf5(cells_group);
for (Cell* c : model::cells) c->to_hdf5(cells_group);
close_group(cells_group);
auto surfaces_group = create_group(geom_group, "surfaces");
for (Surface* surf : surfaces) surf->to_hdf5(surfaces_group);
for (Surface* surf : model::surfaces) surf->to_hdf5(surfaces_group);
close_group(surfaces_group);
auto universes_group = create_group(geom_group, "universes");
for (Universe* u : universes) u->to_hdf5(universes_group);
for (Universe* u : model::universes) u->to_hdf5(universes_group);
close_group(universes_group);
auto lattices_group = create_group(geom_group, "lattices");
for (Lattice* lat : lattices) lat->to_hdf5(lattices_group);
for (Lattice* lat : model::lattices) lat->to_hdf5(lattices_group);
close_group(lattices_group);
close_group(geom_group);

73
src/surface.cpp
View file

@ -25,12 +25,13 @@ extern "C" const int BC_PERIODIC {3};
// Global variables
//==============================================================================
int32_t n_surfaces;
namespace model {
std::vector<Surface*> surfaces;
std::map<int, int> surface_map;
} // namespace model
//==============================================================================
// Helper functions for reading the "coeffs" node of an XML surface element
//==============================================================================
@ -1063,13 +1064,14 @@ extern "C" void
read_surfaces(pugi::xml_node* node)
{
// Count the number of surfaces.
for (pugi::xml_node surf_node: node->children("surface")) {n_surfaces++;}
int n_surfaces = 0;
for (pugi::xml_node surf_node : node->children("surface")) {n_surfaces++;}
if (n_surfaces == 0) {
fatal_error("No surfaces found in geometry.xml!");
}
// Loop over XML surface elements and populate the array.
surfaces.reserve(n_surfaces);
model::surfaces.reserve(n_surfaces);
{
pugi::xml_node surf_node;
int i_surf;
@ -1078,40 +1080,40 @@ read_surfaces(pugi::xml_node* node)
std::string surf_type = get_node_value(surf_node, "type", true, true);
if (surf_type == "x-plane") {
surfaces.push_back(new SurfaceXPlane(surf_node));
model::surfaces.push_back(new SurfaceXPlane(surf_node));
} else if (surf_type == "y-plane") {
surfaces.push_back(new SurfaceYPlane(surf_node));
model::surfaces.push_back(new SurfaceYPlane(surf_node));
} else if (surf_type == "z-plane") {
surfaces.push_back(new SurfaceZPlane(surf_node));
model::surfaces.push_back(new SurfaceZPlane(surf_node));
} else if (surf_type == "plane") {
surfaces.push_back(new SurfacePlane(surf_node));
model::surfaces.push_back(new SurfacePlane(surf_node));
} else if (surf_type == "x-cylinder") {
surfaces.push_back(new SurfaceXCylinder(surf_node));
model::surfaces.push_back(new SurfaceXCylinder(surf_node));
} else if (surf_type == "y-cylinder") {
surfaces.push_back(new SurfaceYCylinder(surf_node));
model::surfaces.push_back(new SurfaceYCylinder(surf_node));
} else if (surf_type == "z-cylinder") {
surfaces.push_back(new SurfaceZCylinder(surf_node));
model::surfaces.push_back(new SurfaceZCylinder(surf_node));
} else if (surf_type == "sphere") {
surfaces.push_back(new SurfaceSphere(surf_node));
model::surfaces.push_back(new SurfaceSphere(surf_node));
} else if (surf_type == "x-cone") {
surfaces.push_back(new SurfaceXCone(surf_node));
model::surfaces.push_back(new SurfaceXCone(surf_node));
} else if (surf_type == "y-cone") {
surfaces.push_back(new SurfaceYCone(surf_node));
model::surfaces.push_back(new SurfaceYCone(surf_node));
} else if (surf_type == "z-cone") {
surfaces.push_back(new SurfaceZCone(surf_node));
model::surfaces.push_back(new SurfaceZCone(surf_node));
} else if (surf_type == "quadric") {
surfaces.push_back(new SurfaceQuadric(surf_node));
model::surfaces.push_back(new SurfaceQuadric(surf_node));
} else {
std::stringstream err_msg;
@ -1122,11 +1124,11 @@ read_surfaces(pugi::xml_node* node)
}
// Fill the surface map.
for (int i_surf = 0; i_surf < n_surfaces; i_surf++) {
int id = surfaces[i_surf]->id_;
auto in_map = surface_map.find(id);
if (in_map == surface_map.end()) {
surface_map[id] = i_surf;
for (int i_surf = 0; i_surf < model::surfaces.size(); i_surf++) {
int id = model::surfaces[i_surf]->id_;
auto in_map = model::surface_map.find(id);
if (in_map == model::surface_map.end()) {
model::surface_map[id] = i_surf;
} else {
std::stringstream err_msg;
err_msg << "Two or more surfaces use the same unique ID: " << id;
@ -1138,10 +1140,10 @@ read_surfaces(pugi::xml_node* node)
double xmin {INFTY}, xmax {-INFTY}, ymin {INFTY}, ymax {-INFTY},
zmin {INFTY}, zmax {-INFTY};
int i_xmin, i_xmax, i_ymin, i_ymax, i_zmin, i_zmax;
for (int i_surf = 0; i_surf < n_surfaces; i_surf++) {
if (surfaces[i_surf]->bc_ == BC_PERIODIC) {
for (int i_surf = 0; i_surf < model::surfaces.size(); i_surf++) {
if (model::surfaces[i_surf]->bc_ == BC_PERIODIC) {
// Downcast to the PeriodicSurface type.
Surface* surf_base = surfaces[i_surf];
Surface* surf_base = model::surfaces[i_surf];
PeriodicSurface* surf = dynamic_cast<PeriodicSurface*>(surf_base);
// Make sure this surface inherits from PeriodicSurface.
@ -1183,10 +1185,10 @@ read_surfaces(pugi::xml_node* node)
}
// Set i_periodic for periodic BC surfaces.
for (int i_surf = 0; i_surf < n_surfaces; i_surf++) {
if (surfaces[i_surf]->bc_ == BC_PERIODIC) {
for (int i_surf = 0; i_surf < model::surfaces.size(); i_surf++) {
if (model::surfaces[i_surf]->bc_ == BC_PERIODIC) {
// Downcast to the PeriodicSurface type.
Surface* surf_base = surfaces[i_surf];
Surface* surf_base = model::surfaces[i_surf];
PeriodicSurface* surf = dynamic_cast<PeriodicSurface*>(surf_base);
// Also try downcasting to the SurfacePlane type (which must be handled
@ -1216,7 +1218,7 @@ read_surfaces(pugi::xml_node* node)
}
} else {
// Convert the surface id to an index.
surf->i_periodic_ = surface_map[surf->i_periodic_];
surf->i_periodic_ = model::surface_map[surf->i_periodic_];
}
} else {
// This is a SurfacePlane. We won't try to find it's partner if the
@ -1228,12 +1230,12 @@ read_surfaces(pugi::xml_node* node)
fatal_error(err_msg);
} else {
// Convert the surface id to an index.
surf->i_periodic_ = surface_map[surf->i_periodic_];
surf->i_periodic_ = model::surface_map[surf->i_periodic_];
}
}
// Make sure the opposite surface is also periodic.
if (surfaces[surf->i_periodic_]->bc_ != BC_PERIODIC) {
if (model::surfaces[surf->i_periodic_]->bc_ != BC_PERIODIC) {
std::stringstream err_msg;
err_msg << "Could not find matching surface for periodic boundary "
"condition on surface " << surf->id_;
@ -1248,7 +1250,7 @@ read_surfaces(pugi::xml_node* node)
//==============================================================================
extern "C" {
Surface* surface_pointer(int surf_ind) {return surfaces[surf_ind];}
Surface* surface_pointer(int surf_ind) {return model::surfaces[surf_ind];}
int surface_id(Surface* surf) {return surf->id_;}
@ -1288,11 +1290,12 @@ extern "C" {
void free_memory_surfaces_c()
{
for (Surface* surf : surfaces) {delete surf;}
surfaces.clear();
n_surfaces = 0;
surface_map.clear();
for (Surface* surf : model::surfaces) {delete surf;}
model::surfaces.clear();
model::surface_map.clear();
}
int surfaces_size() { return model::surfaces.size(); }
}
} // namespace openmc

9
src/surface_header.F90
View file

@ -78,10 +78,15 @@ module surface_header
end type Surface
integer(C_INT32_T), bind(C) :: n_surfaces ! # of surfaces
type(Surface), allocatable, target :: surfaces(:)
interface
function surfaces_size() result(sz) bind(C)
import C_INT
integer(C_INT) :: sz
end function
end interface
contains
pure function surface_id(this) result(id)

60
src/tallies/filter.cpp
View file

@ -26,6 +26,8 @@
#include "openmc/tallies/filter_universe.h"
#include "openmc/tallies/filter_zernike.h"
// explicit template instantiation definition
template class std::vector<openmc::FilterMatch>;
namespace openmc {
@ -33,8 +35,13 @@ namespace openmc {
// Global variables
//==============================================================================
namespace simulation {
std::vector<FilterMatch> filter_matches;
} // namespace simulation
namespace model {
std::vector<std::unique_ptr<Filter>> tally_filters;
} // namespace model
//==============================================================================
// Non-member functions
@ -45,10 +52,10 @@ free_memory_tally_c()
{
#pragma omp parallel
{
filter_matches.clear();
simulation::filter_matches.clear();
}
tally_filters.clear();
model::tally_filters.clear();
}
//==============================================================================
@ -56,8 +63,7 @@ free_memory_tally_c()
//==============================================================================
extern "C" {
FilterMatch* filter_match_pointer(int indx)
{return &filter_matches[indx];}
// filter_match_point moved to simulation.cpp
void
filter_match_bins_push_back(FilterMatch* match, int val)
@ -96,53 +102,53 @@ extern "C" {
{
std::string type_ {type};
if (type_ == "azimuthal") {
tally_filters.push_back(std::make_unique<AzimuthalFilter>());
model::tally_filters.push_back(std::make_unique<AzimuthalFilter>());
} else if (type_ == "cell") {
tally_filters.push_back(std::make_unique<CellFilter>());
model::tally_filters.push_back(std::make_unique<CellFilter>());
} else if (type_ == "cellborn") {
tally_filters.push_back(std::make_unique<CellbornFilter>());
model::tally_filters.push_back(std::make_unique<CellbornFilter>());
} else if (type_ == "cellfrom") {
tally_filters.push_back(std::make_unique<CellFromFilter>());
model::tally_filters.push_back(std::make_unique<CellFromFilter>());
} else if (type_ == "distribcell") {
tally_filters.push_back(std::make_unique<DistribcellFilter>());
model::tally_filters.push_back(std::make_unique<DistribcellFilter>());
} else if (type_ == "delayedgroup") {
tally_filters.push_back(std::make_unique<DelayedGroupFilter>());
model::tally_filters.push_back(std::make_unique<DelayedGroupFilter>());
} else if (type_ == "energyfunction") {
tally_filters.push_back(std::make_unique<EnergyFunctionFilter>());
model::tally_filters.push_back(std::make_unique<EnergyFunctionFilter>());
} else if (type_ == "energy") {
tally_filters.push_back(std::make_unique<EnergyFilter>());
model::tally_filters.push_back(std::make_unique<EnergyFilter>());
} else if (type_ == "energyout") {
tally_filters.push_back(std::make_unique<EnergyoutFilter>());
model::tally_filters.push_back(std::make_unique<EnergyoutFilter>());
} else if (type_ == "legendre") {
tally_filters.push_back(std::make_unique<LegendreFilter>());
model::tally_filters.push_back(std::make_unique<LegendreFilter>());
} else if (type_ == "material") {
tally_filters.push_back(std::make_unique<MaterialFilter>());
model::tally_filters.push_back(std::make_unique<MaterialFilter>());
} else if (type_ == "mesh") {
tally_filters.push_back(std::make_unique<MeshFilter>());
model::tally_filters.push_back(std::make_unique<MeshFilter>());
} else if (type_ == "meshsurface") {
tally_filters.push_back(std::make_unique<MeshSurfaceFilter>());
model::tally_filters.push_back(std::make_unique<MeshSurfaceFilter>());
} else if (type_ == "mu") {
tally_filters.push_back(std::make_unique<MuFilter>());
model::tally_filters.push_back(std::make_unique<MuFilter>());
} else if (type_ == "particle") {
tally_filters.push_back(std::make_unique<ParticleFilter>());
model::tally_filters.push_back(std::make_unique<ParticleFilter>());
} else if (type_ == "polar") {
tally_filters.push_back(std::make_unique<PolarFilter>());
model::tally_filters.push_back(std::make_unique<PolarFilter>());
} else if (type_ == "surface") {
tally_filters.push_back(std::make_unique<SurfaceFilter>());
model::tally_filters.push_back(std::make_unique<SurfaceFilter>());
} else if (type_ == "spatiallegendre") {
tally_filters.push_back(std::make_unique<SpatialLegendreFilter>());
model::tally_filters.push_back(std::make_unique<SpatialLegendreFilter>());
} else if (type_ == "sphericalharmonics") {
tally_filters.push_back(std::make_unique<SphericalHarmonicsFilter>());
model::tally_filters.push_back(std::make_unique<SphericalHarmonicsFilter>());
} else if (type_ == "universe") {
tally_filters.push_back(std::make_unique<UniverseFilter>());
model::tally_filters.push_back(std::make_unique<UniverseFilter>());
} else if (type_ == "zernike") {
tally_filters.push_back(std::make_unique<ZernikeFilter>());
model::tally_filters.push_back(std::make_unique<ZernikeFilter>());
} else if (type_ == "zernikeradial") {
tally_filters.push_back(std::make_unique<ZernikeRadialFilter>());
model::tally_filters.push_back(std::make_unique<ZernikeRadialFilter>());
} else {
return nullptr;
}
return tally_filters.back().get();
return model::tally_filters.back().get();
}
void filter_from_xml(Filter* filt, pugi::xml_node* node)

8
src/tallies/filter_cell.cpp
View file

@ -21,8 +21,8 @@ CellFilter::initialize()
{
// Convert cell IDs to indices of the global array.
for (auto& c : cells_) {
auto search = cell_map.find(c);
if (search != cell_map.end()) {
auto search = model::cell_map.find(c);
if (search != model::cell_map.end()) {
c = search->second;
} else {
std::stringstream err_msg;
@ -56,7 +56,7 @@ CellFilter::to_statepoint(hid_t filter_group) const
{
Filter::to_statepoint(filter_group);
std::vector<int32_t> cell_ids;
for (auto c : cells_) cell_ids.push_back(cells[c]->id_);
for (auto c : cells_) cell_ids.push_back(model::cells[c]->id_);
write_dataset(filter_group, "bins", cell_ids);
}
@ -64,7 +64,7 @@ std::string
CellFilter::text_label(int bin) const
{
//TODO: off-by-one
return "Cell " + std::to_string(cells[cells_[bin-1]]->id_);
return "Cell " + std::to_string(model::cells[cells_[bin-1]]->id_);
}
//==============================================================================

2
src/tallies/filter_cellborn.cpp
View file

@ -20,7 +20,7 @@ std::string
CellbornFilter::text_label(int bin) const
{
//TODO: off-by-one
return "Birth Cell " + std::to_string(cells[cells_[bin-1]]->id_);
return "Birth Cell " + std::to_string(model::cells[cells_[bin-1]]->id_);
}
} // namespace openmc

2
src/tallies/filter_cellfrom.cpp
View file

@ -22,7 +22,7 @@ std::string
CellFromFilter::text_label(int bin) const
{
//TODO: off-by-one
return "Cell from " + std::to_string(cells[cells_[bin-1]]->id_);
return "Cell from " + std::to_string(model::cells[cells_[bin-1]]->id_);
}
} // namespace openmc

16
src/tallies/filter_distribcell.cpp
View file

@ -22,10 +22,10 @@ void
DistribcellFilter::initialize()
{
// Convert the cell ID to an index of the global array.
auto search = cell_map.find(cell_);
if (search != cell_map.end()) {
auto search = model::cell_map.find(cell_);
if (search != model::cell_map.end()) {
cell_ = search->second;
n_bins_ = cells[cell_]->n_instances_;
n_bins_ = model::cells[cell_]->n_instances_;
} else {
std::stringstream err_msg;
err_msg << "Could not find cell " << cell_
@ -39,13 +39,13 @@ DistribcellFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
int offset = 0;
auto distribcell_index = cells[cell_]->distribcell_index_;
auto distribcell_index = model::cells[cell_]->distribcell_index_;
for (int i = 0; i < p->n_coord; i++) {
auto& c {*cells[p->coord[i].cell]};
auto& c {*model::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]};
auto& lat {*model::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};
@ -66,13 +66,13 @@ void
DistribcellFilter::to_statepoint(hid_t filter_group) const
{
Filter::to_statepoint(filter_group);
write_dataset(filter_group, "bins", cells[cell_]->id_);
write_dataset(filter_group, "bins", model::cells[cell_]->id_);
}
std::string
DistribcellFilter::text_label(int bin) const
{
auto map = cells[cell_]->distribcell_index_;
auto map = model::cells[cell_]->distribcell_index_;
//TODO: off-by-one
auto path = distribcell_path(cell_, map, bin-1);
return "Distributed Cell " + path;

10
src/tallies/filter_energy.cpp
View file

@ -25,10 +25,10 @@ EnergyFilter::from_xml(pugi::xml_node node)
// (after flipping for the different ordering of the library and tallying
// systems).
if (!settings::run_CE) {
if (n_bins_ == num_energy_groups) {
if (n_bins_ == data::num_energy_groups) {
matches_transport_groups_ = true;
for (auto i = 0; i < n_bins_ + 1; i++) {
if (rev_energy_bins[i] != bins_[i]) {
if (data::rev_energy_bins[i] != bins_[i]) {
matches_transport_groups_ = false;
break;
}
@ -44,10 +44,10 @@ const
if (p->g != F90_NONE && matches_transport_groups_) {
if (estimator == ESTIMATOR_TRACKLENGTH) {
//TODO: off-by-one
match.bins_.push_back(num_energy_groups - p->g + 1);
match.bins_.push_back(data::num_energy_groups - p->g + 1);
} else {
//TODO: off-by-one
match.bins_.push_back(num_energy_groups - p->last_g + 1);
match.bins_.push_back(data::num_energy_groups - p->last_g + 1);
}
match.weights_.push_back(1.0);
@ -90,7 +90,7 @@ EnergyoutFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
if (p->g != F90_NONE && matches_transport_groups_) {
match.bins_.push_back(num_energy_groups - p->g + 1);
match.bins_.push_back(data::num_energy_groups - p->g + 1);
match.weights_.push_back(1.0);
} else {

8
src/tallies/filter_material.cpp
View file

@ -21,8 +21,8 @@ MaterialFilter::initialize()
{
// Convert material IDs to indices of the global array.
for (auto& m : materials_) {
auto search = material_map.find(m);
if (search != material_map.end()) {
auto search = model::material_map.find(m);
if (search != model::material_map.end()) {
m = search->second;
} else {
std::stringstream err_msg;
@ -55,7 +55,7 @@ MaterialFilter::to_statepoint(hid_t filter_group) const
{
Filter::to_statepoint(filter_group);
std::vector<int32_t> material_ids;
for (auto c : materials_) material_ids.push_back(materials[c]->id_);
for (auto c : materials_) material_ids.push_back(model::materials[c]->id_);
write_dataset(filter_group, "bins", material_ids);
}
@ -63,7 +63,7 @@ std::string
MaterialFilter::text_label(int bin) const
{
//TODO: off-by-one
return "Material " + std::to_string(materials[materials_[bin-1]]->id_);
return "Material " + std::to_string(model::materials[materials_[bin-1]]->id_);
}
//==============================================================================

16
src/tallies/filter_mesh.cpp
View file

@ -20,8 +20,8 @@ MeshFilter::from_xml(pugi::xml_node node)
}
auto id = bins_[0];
auto search = mesh_map.find(id);
if (search != mesh_map.end()) {
auto search = model::mesh_map.find(id);
if (search != model::mesh_map.end()) {
set_mesh(search->second);
} else{
std::stringstream err_msg;
@ -35,14 +35,14 @@ MeshFilter::get_all_bins(const Particle* p, int estimator, FilterMatch& match)
const
{
if (estimator != ESTIMATOR_TRACKLENGTH) {
auto bin = meshes[mesh_]->get_bin(p->coord[0].xyz);
auto bin = model::meshes[mesh_]->get_bin(p->coord[0].xyz);
if (bin >= 0) {
//TODO: off-by-one
match.bins_.push_back(bin);
match.weights_.push_back(1.0);
}
} else {
meshes[mesh_]->bins_crossed(p, match.bins_, match.weights_);
model::meshes[mesh_]->bins_crossed(p, match.bins_, match.weights_);
}
}
@ -50,13 +50,13 @@ void
MeshFilter::to_statepoint(hid_t filter_group) const
{
Filter::to_statepoint(filter_group);
write_dataset(filter_group, "bins", meshes[mesh_]->id_);
write_dataset(filter_group, "bins", model::meshes[mesh_]->id_);
}
std::string
MeshFilter::text_label(int bin) const
{
auto& mesh = *meshes[mesh_];
auto& mesh = *model::meshes[mesh_];
int n_dim = mesh.n_dimension_;
int ijk[n_dim];
@ -76,7 +76,7 @@ MeshFilter::set_mesh(int32_t mesh)
{
mesh_ = mesh;
n_bins_ = 1;
for (auto dim : meshes[mesh_]->shape_) n_bins_ *= dim;
for (auto dim : model::meshes[mesh_]->shape_) n_bins_ *= dim;
}
//==============================================================================
@ -123,7 +123,7 @@ openmc_mesh_filter_set_mesh(int32_t index, int32_t index_mesh)
}
// Check the mesh index.
if (index_mesh < 0 || index_mesh >= meshes.size()) {
if (index_mesh < 0 || index_mesh >= model::meshes.size()) {
set_errmsg("Index in 'meshes' array is out of bounds.");
return OPENMC_E_OUT_OF_BOUNDS;
}

8
src/tallies/filter_meshsurface.cpp
View file

@ -11,14 +11,14 @@ void
MeshSurfaceFilter::get_all_bins(const Particle* p, int estimator,
FilterMatch& match) const
{
meshes[mesh_]->surface_bins_crossed(p, match.bins_);
model::meshes[mesh_]->surface_bins_crossed(p, match.bins_);
for (auto b : match.bins_) match.weights_.push_back(1.0);
}
std::string
MeshSurfaceFilter::text_label(int bin) const
{
auto& mesh = *meshes[mesh_];
auto& mesh = *model::meshes[mesh_];
int n_dim = mesh.n_dimension_;
// Get flattend mesh index and surface index.
@ -76,8 +76,8 @@ void
MeshSurfaceFilter::set_mesh(int32_t mesh)
{
mesh_ = mesh;
n_bins_ = 4 * meshes[mesh_]->n_dimension_;
for (auto dim : meshes[mesh_]->shape_) n_bins_ *= dim;
n_bins_ = 4 * model::meshes[mesh_]->n_dimension_;
for (auto dim : model::meshes[mesh_]->shape_) n_bins_ *= dim;
}
//==============================================================================

8
src/tallies/filter_surface.cpp
View file

@ -20,8 +20,8 @@ SurfaceFilter::initialize()
{
// Convert surface IDs to indices of the global array.
for (auto& s : surfaces_) {
auto search = surface_map.find(s);
if (search != surface_map.end()) {
auto search = model::surface_map.find(s);
if (search != model::surface_map.end()) {
s = search->second;
} else {
std::stringstream err_msg;
@ -58,7 +58,7 @@ SurfaceFilter::to_statepoint(hid_t filter_group) const
{
Filter::to_statepoint(filter_group);
std::vector<int32_t> surface_ids;
for (auto c : surfaces_) surface_ids.push_back(surfaces[c]->id_);
for (auto c : surfaces_) surface_ids.push_back(model::surfaces[c]->id_);
write_dataset(filter_group, "bins", surface_ids);
}
@ -66,7 +66,7 @@ std::string
SurfaceFilter::text_label(int bin) const
{
//TODO: off-by-one
return "Surface " + std::to_string(surfaces[surfaces_[bin-1]]->id_);
return "Surface " + std::to_string(model::surfaces[surfaces_[bin-1]]->id_);
}
} // namespace openmc

8
src/tallies/filter_universe.cpp
View file

@ -20,8 +20,8 @@ UniverseFilter::initialize()
{
// Convert universe IDs to indices of the global array.
for (auto& u : universes_) {
auto search = universe_map.find(u);
if (search != universe_map.end()) {
auto search = model::universe_map.find(u);
if (search != model::universe_map.end()) {
u = search->second;
} else {
std::stringstream err_msg;
@ -56,7 +56,7 @@ UniverseFilter::to_statepoint(hid_t filter_group) const
{
Filter::to_statepoint(filter_group);
std::vector<int32_t> universe_ids;
for (auto u : universes_) universe_ids.push_back(universes[u]->id_);
for (auto u : universes_) universe_ids.push_back(model::universes[u]->id_);
write_dataset(filter_group, "bins", universe_ids);
}
@ -64,7 +64,7 @@ std::string
UniverseFilter::text_label(int bin) const
{
//TODO: off-by-one
return "Universe " + std::to_string(universes[universes_[bin-1]]->id_);
return "Universe " + std::to_string(model::universes[universes_[bin-1]]->id_);
}
} // namespace openmc

10
src/tallies/tally.F90
View file

@ -3792,11 +3792,6 @@ contains
! Accumulate on master only unless run is not reduced then do it on all
if (master .or. (.not. reduce_tallies)) then
! Accumulate results for each tally
do i = 1, active_tallies % size()
call tallies(active_tallies % data(i)) % obj % accumulate()
end do
if (run_mode == MODE_EIGENVALUE) then
if (current_batch > n_inactive) then
! Accumulate products of different estimators of k
@ -3820,6 +3815,11 @@ contains
end do
end if
! Accumulate results for each tally
do i = 1, active_tallies % size()
call tallies(active_tallies % data(i)) % obj % accumulate()
end do
end subroutine accumulate_tallies
!===============================================================================

49
src/tallies/tally_header.F90
View file

@ -6,7 +6,7 @@ module tally_header
use error
use dict_header, only: DictIntInt
use hdf5_interface, only: HID_T, HSIZE_T
use message_passing, only: n_procs
use message_passing, only: n_procs, master
use nuclide_header, only: nuclide_dict
use settings, only: reduce_tallies, run_mode
use stl_vector, only: VectorInt
@ -17,7 +17,7 @@ module tally_header
implicit none
private
public :: configure_tallies
public :: allocate_tally_results
public :: free_memory_tally
public :: openmc_extend_tallies
public :: openmc_get_tally_index
@ -182,25 +182,27 @@ contains
this % n_realizations = this % n_realizations + n_procs
end if
! Calculate total source strength for normalization
if (run_mode == MODE_FIXEDSOURCE) then
total_source = total_source_strength()
else
total_source = ONE
end if
if (master .or. (.not. reduce_tallies)) then
! Calculate total source strength for normalization
if (run_mode == MODE_FIXEDSOURCE) then
total_source = total_source_strength()
else
total_source = ONE
end if
! Accumulate each result
do j = 1, size(this % results, 3)
do i = 1, size(this % results, 2)
val = this % results(RESULT_VALUE, i, j)/total_weight * total_source
this % results(RESULT_VALUE, i, j) = ZERO
! Accumulate each result
do j = 1, size(this % results, 3)
do i = 1, size(this % results, 2)
val = this % results(RESULT_VALUE, i, j)/total_weight * total_source
this % results(RESULT_VALUE, i, j) = ZERO
this % results(RESULT_SUM, i, j) = &
this % results(RESULT_SUM, i, j) + val
this % results(RESULT_SUM_SQ, i, j) = &
this % results(RESULT_SUM_SQ, i, j) + val*val
this % results(RESULT_SUM, i, j) = &
this % results(RESULT_SUM, i, j) + val
this % results(RESULT_SUM_SQ, i, j) = &
this % results(RESULT_SUM_SQ, i, j) + val*val
end do
end do
end do
end if
end subroutine tally_accumulate
subroutine tally_write_results_hdf5(this, group_id)
@ -275,9 +277,6 @@ contains
allocate(this % results(3, this % total_score_bins, this % n_filter_bins))
end if
! Initialize results array to zero
this % results(:,:,:) = ZERO
end subroutine tally_allocate_results
@ -392,21 +391,21 @@ contains
! tallies.xml file.
!===============================================================================
subroutine configure_tallies() bind(C)
subroutine allocate_tally_results() bind(C)
integer :: i
! Allocate and initialize global tallies
! Allocate global tallies
if (.not. allocated(global_tallies)) then
allocate(global_tallies(3, N_GLOBAL_TALLIES))
end if
global_tallies(:,:) = ZERO
! Allocate results arrays for tallies
do i = 1, n_tallies
call tallies(i) % obj % allocate_results()
end do
end subroutine configure_tallies
end subroutine allocate_tally_results
!===============================================================================
! FREE_MEMORY_TALLY deallocates global arrays defined in this module

74
src/timer.cpp
View file

@ -2,6 +2,25 @@
namespace openmc {
//==============================================================================
// Global variables
//==============================================================================
namespace simulation {
Timer time_active;
Timer time_bank;
Timer time_bank_sample;
Timer time_bank_sendrecv;
Timer time_finalize;
Timer time_inactive;
Timer time_initialize;
Timer time_tallies;
Timer time_total;
Timer time_transport;
} // namespace simulation
//==============================================================================
// Timer implementation
//==============================================================================
@ -34,35 +53,20 @@ double Timer::elapsed()
}
}
//==============================================================================
// Global variables
//==============================================================================
Timer time_active;
Timer time_bank;
Timer time_bank_sample;
Timer time_bank_sendrecv;
Timer time_finalize;
Timer time_inactive;
Timer time_initialize;
Timer time_tallies;
Timer time_total;
Timer time_transport;
//==============================================================================
// Fortran compatibility
//==============================================================================
extern "C" double time_active_elapsed() { return time_active.elapsed(); }
extern "C" double time_bank_elapsed() { return time_bank.elapsed(); }
extern "C" double time_bank_sample_elapsed() { return time_bank_sample.elapsed(); }
extern "C" double time_bank_sendrecv_elapsed() { return time_bank_sendrecv.elapsed(); }
extern "C" double time_finalize_elapsed() { return time_finalize.elapsed(); }
extern "C" double time_inactive_elapsed() { return time_inactive.elapsed(); }
extern "C" double time_initialize_elapsed() { return time_initialize.elapsed(); }
extern "C" double time_tallies_elapsed() { return time_tallies.elapsed(); }
extern "C" double time_total_elapsed() { return time_total.elapsed(); }
extern "C" double time_transport_elapsed() { return time_transport.elapsed(); }
extern "C" double time_active_elapsed() { return simulation::time_active.elapsed(); }
extern "C" double time_bank_elapsed() { return simulation::time_bank.elapsed(); }
extern "C" double time_bank_sample_elapsed() { return simulation::time_bank_sample.elapsed(); }
extern "C" double time_bank_sendrecv_elapsed() { return simulation::time_bank_sendrecv.elapsed(); }
extern "C" double time_finalize_elapsed() { return simulation::time_finalize.elapsed(); }
extern "C" double time_inactive_elapsed() { return simulation::time_inactive.elapsed(); }
extern "C" double time_initialize_elapsed() { return simulation::time_initialize.elapsed(); }
extern "C" double time_tallies_elapsed() { return simulation::time_tallies.elapsed(); }
extern "C" double time_total_elapsed() { return simulation::time_total.elapsed(); }
extern "C" double time_transport_elapsed() { return simulation::time_transport.elapsed(); }
//==============================================================================
// Non-member functions
@ -70,16 +74,16 @@ extern "C" double time_transport_elapsed() { return time_transport.elapsed(); }
void reset_timers()
{
time_active.reset();
time_bank.reset();
time_bank_sample.reset();
time_bank_sendrecv.reset();
time_finalize.reset();
time_inactive.reset();
time_initialize.reset();
time_tallies.reset();
time_total.reset();
time_transport.reset();
simulation::time_active.reset();
simulation::time_bank.reset();
simulation::time_bank_sample.reset();
simulation::time_bank_sendrecv.reset();
simulation::time_finalize.reset();
simulation::time_inactive.reset();
simulation::time_initialize.reset();
simulation::time_tallies.reset();
simulation::time_total.reset();
simulation::time_transport.reset();
}
} // namespace openmc

5
src/tracking.F90
View file

@ -33,10 +33,9 @@ module tracking
implicit none
interface
subroutine collision_mg(p, energy_bin_avg, material_xs) bind(C)
subroutine collision_mg(p, material_xs) bind(C)
import Particle, C_DOUBLE, MaterialMacroXS
type(Particle), intent(inout) :: p
real(C_DOUBLE), intent(in) :: energy_bin_avg(*)
type(MaterialMacroXS), intent(in) :: material_xs
end subroutine collision_mg
@ -235,7 +234,7 @@ contains
if (run_CE) then
call collision(p)
else
call collision_mg(p, energy_bin_avg, material_xs)
call collision_mg(p, material_xs)
end if
! Score collision estimator tallies -- this is done after a collision