crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 14:15:42 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Merge branch 'develop' into photon-production-fix

Browse source
This commit is contained in:
amandalund 2018-09-11 08:52:23 -05:00
commit 7d163b256e
96 changed files with 4264 additions and 3912 deletions
Download patch file Download diff file Expand all files Collapse all files

17
include/openmc/capi.h
View file

@ -37,6 +37,7 @@ extern "C" {
int openmc_filter_set_type(int32_t index, const char* type);
int openmc_finalize();
int openmc_find_cell(double* xyz, int32_t* index, int32_t* instance);
int openmc_fission_bank(struct Bank** ptr, int64_t* n);
int openmc_get_cell_index(int32_t id, int32_t* index);
int openmc_get_filter_index(int32_t id, int32_t* index);
void openmc_get_filter_next_id(int32_t* id);
@ -71,7 +72,7 @@ extern "C" {
int openmc_mesh_get_params(int32_t index, double** ll, double** ur, double** width, int* n);
int openmc_mesh_set_id(int32_t index, int32_t id);
int openmc_mesh_set_dimension(int32_t index, int n, const int* dims);
int openmc_mesh_set_params(int32_t index, const double* ll, const double* ur, const double* width, int n);
int openmc_mesh_set_params(int32_t index, int n, const double* ll, const double* ur, const double* width);
int openmc_meshsurface_filter_get_mesh(int32_t index, int32_t* index_mesh);
int openmc_meshsurface_filter_set_mesh(int32_t index, int32_t index_mesh);
int openmc_next_batch(int* status);
@ -135,16 +136,11 @@ extern "C" {
extern char openmc_err_msg[256];
extern double openmc_keff;
extern double openmc_keff_std;
extern int32_t gen_per_batch;
extern int32_t n_batches;
extern int32_t n_cells;
extern int32_t n_filters;
extern int32_t n_inactive;
extern int32_t n_lattices;
extern int32_t n_materials;
extern int32_t n_meshes;
extern int n_nuclides;
extern int64_t n_particles;
extern int32_t n_plots;
extern int32_t n_realizations;
extern int32_t n_sab_tables;
@ -152,9 +148,7 @@ extern "C" {
extern int32_t n_surfaces;
extern int32_t n_tallies;
extern int32_t n_universes;
extern int openmc_run_mode;
extern bool openmc_simulation_initialized;
extern int openmc_verbosity;
// Variables that are shared by necessity (can be removed from public header
// later)
@ -164,13 +158,6 @@ extern "C" {
extern int openmc_rank;
extern int64_t openmc_work;
// Run modes
const int RUN_MODE_FIXEDSOURCE = 1;
const int RUN_MODE_EIGENVALUE = 2;
const int RUN_MODE_PLOTTING = 3;
const int RUN_MODE_PARTICLE = 4;
const int RUN_MODE_VOLUME = 5;
#ifdef __cplusplus
}
#endif

20
include/openmc/constants.h
View file

@ -11,16 +11,9 @@
namespace openmc {
// TODO: Replace with xtensor/other library?
typedef std::vector<double> double_1dvec;
typedef std::vector<std::vector<double> > double_2dvec;
typedef std::vector<std::vector<std::vector<double> > > double_3dvec;
typedef std::vector<std::vector<std::vector<std::vector<double> > > > double_4dvec;
typedef std::vector<std::vector<std::vector<std::vector<std::vector<double> > > > > double_5dvec;
typedef std::vector<std::vector<std::vector<std::vector<std::vector<std::vector<double> > > > > > double_6dvec;
typedef std::vector<int> int_1dvec;
typedef std::vector<std::vector<int> > int_2dvec;
typedef std::vector<std::vector<std::vector<int> > > int_3dvec;
using double_2dvec = std::vector<std::vector<double>>;
using double_3dvec = std::vector<std::vector<std::vector<double>>>;
using double_4dvec = std::vector<std::vector<std::vector<std::vector<double>>>>;
// ============================================================================
// VERSIONING NUMBERS
@ -431,6 +424,13 @@ enum class Interpolation {
histogram, lin_lin, lin_log, log_lin, log_log
};
// Run modes
constexpr int RUN_MODE_FIXEDSOURCE {1};
constexpr int RUN_MODE_EIGENVALUE {2};
constexpr int RUN_MODE_PLOTTING {3};
constexpr int RUN_MODE_PARTICLE {4};
constexpr int RUN_MODE_VOLUME {5};
} // namespace openmc
#endif // OPENMC_CONSTANTS_H

41
include/openmc/eigenvalue.h Normal file
View file

@ -0,0 +1,41 @@
#ifndef OPENMC_EIGENVALUE_H
#define OPENMC_EIGENVALUE_H
#include <cstdint> // for int64_t
#include <vector>
#include "xtensor/xtensor.hpp"
#include "openmc/particle.h"
namespace openmc {
//==============================================================================
// Global variables
//==============================================================================
extern std::vector<double> entropy; //!< Shannon entropy at each generation
extern xt::xtensor<double, 1> source_frac; //!< Source fraction for UFS
extern "C" int64_t n_bank;
#pragma omp threadprivate(n_bank)
//==============================================================================
// Non-member functions
//==============================================================================
//! Calculates the Shannon entropy of the fission source distribution to assess
//! source convergence
extern "C" void shannon_entropy();
//! Determines the source fraction in each UFS mesh cell and reweights the
//! source bank so that the sum of the weights is equal to n_particles. The
//! 'source_frac' variable is used later to bias the production of fission sites
extern "C" void ufs_count_sites();
//! Get UFS weight corresponding to particle's location
extern "C" double ufs_get_weight(const Particle* p);
} // namespace openmc
#endif // OPENMC_EIGENVALUE_H

87
include/openmc/hdf5_interface.h
View file

@ -14,6 +14,7 @@
#include "xtensor/xarray.hpp"
#include "openmc/position.h"
#include "openmc/error.h"
namespace openmc {
@ -46,39 +47,6 @@ hid_t file_open(const std::string& filename, char mode, bool parallel=false);
void write_string(hid_t group_id, const char* name, const std::string& buffer,
bool indep);
void
read_nd_vector(hid_t obj_id, const char* name, std::vector<double>& result,
bool must_have = false);
void
read_nd_vector(hid_t obj_id, const char* name,
std::vector<std::vector<double> >& result,
bool must_have = false);
void
read_nd_vector(hid_t obj_id, const char* name,
std::vector<std::vector<int> >& result, bool must_have = false);
void
read_nd_vector(hid_t obj_id, const char* name,
std::vector<std::vector<std::vector<double> > >& result,
bool must_have = false);
void
read_nd_vector(hid_t obj_id, const char* name,
std::vector<std::vector<std::vector<int> > >& result,
bool must_have = false);
void
read_nd_vector(hid_t obj_id, const char* name,
std::vector<std::vector<std::vector<std::vector<double> > > >& result,
bool must_have = false);
void
read_nd_vector(hid_t obj_id, const char* name,
std::vector<std::vector<std::vector<std::vector<std::vector<double> > > > >& result,
bool must_have = false);
std::vector<hsize_t> attribute_shape(hid_t obj_id, const char* name);
std::vector<std::string> dataset_names(hid_t group_id);
void ensure_exists(hid_t group_id, const char* name);
@ -236,7 +204,7 @@ read_attribute(hid_t obj_id, const char* name, std::vector<std::string>& vec)
}
//==============================================================================
// Templates/overloads for read_dataset
// Templates/overloads for read_dataset and related methods
//==============================================================================
template<typename T>
@ -294,6 +262,40 @@ void read_dataset(hid_t obj_id, const char* name, xt::xarray<T>& arr, bool indep
close_dataset(dset);
}
template <typename T, std::size_t N>
void read_dataset_as_shape(hid_t obj_id, const char* name,
xt::xtensor<T, N>& arr, bool indep=false)
{
hid_t dset = open_dataset(obj_id, name);
// Allocate new array to read data into
std::size_t size = 1;
for (const auto x : arr.shape())
size *= x;
T* buffer = new T[size];
// Read data from attribute
read_dataset(dset, nullptr, H5TypeMap<T>::type_id, buffer, indep);
// Adapt into xarray
arr = xt::adapt(buffer, size, xt::acquire_ownership(), arr.shape());
close_dataset(dset);
}
template <typename T, std::size_t N>
void read_nd_vector(hid_t obj_id, const char* name, xt::xtensor<T, N>& result,
bool must_have=false)
{
if (object_exists(obj_id, name)) {
read_dataset_as_shape(obj_id, name, result, true);
} else if (must_have) {
fatal_error(std::string("Must provide " + std::string(name) + "!"));
}
}
//==============================================================================
// Templates/overloads for write_attribute
//==============================================================================
@ -317,6 +319,14 @@ write_attribute(hid_t obj_id, const char* name, const std::array<T, N>& buffer)
write_attr(obj_id, 1, dims, name, H5TypeMap<T>::type_id, buffer.data());
}
template<typename T> inline void
write_attribute(hid_t obj_id, const char* name, const std::vector<T>& buffer)
{
hsize_t dims[] {buffer.size()};
write_attr(obj_id, 1, dims, name, H5TypeMap<T>::type_id, buffer.data());
}
//==============================================================================
// Templates/overloads for write_dataset
//==============================================================================
@ -347,6 +357,15 @@ write_dataset(hid_t obj_id, const char* name, const std::vector<T>& buffer)
write_dataset(obj_id, 1, dims, name, H5TypeMap<T>::type_id, buffer.data(), false);
}
template<typename T> inline void
write_dataset(hid_t obj_id, const char* name, const xt::xarray<T>& arr)
{
auto s = arr.shape();
std::vector<hsize_t> dims {s.cbegin(), s.cend()};
write_dataset(obj_id, dims.size(), dims.data(), name, H5TypeMap<T>::type_id,
arr.data(), false);
}
inline void
write_dataset(hid_t obj_id, const char* name, Position r)
{

131
include/openmc/mesh.h Normal file
View file

@ -0,0 +1,131 @@
//! \file mesh.h
//! \brief Mesh types used for tallies, Shannon entropy, CMFD, etc.
#ifndef OPENMC_MESH_H
#define OPENMC_MESH_H
#include <memory> // for unique_ptr
#include <vector>
#include <unordered_map>
#include "hdf5.h"
#include "pugixml.hpp"
#include "xtensor/xarray.hpp"
#include "openmc/particle.h"
#include "openmc/position.h"
namespace openmc {
//==============================================================================
//! Tessellation of n-dimensional Euclidean space by congruent squares or cubes
//==============================================================================
class RegularMesh {
public:
// Constructors
RegularMesh() = default;
RegularMesh(pugi::xml_node node);
// Methods
//! Determine which bins were crossed by a particle
//!
//! \param[in] p Particle to check
//! \param[out] bins Bins that were crossed
//! \param[out] lengths Fraction of tracklength in each bin
void bins_crossed(const Particle* p, std::vector<int>& bins,
std::vector<double>& lengths) const;
//! Determine which surface bins were crossed by a particle
//!
//! \param[in] p Particle to check
//! \param[out] bins Surface bins that were crossed
void surface_bins_crossed(const Particle* p, std::vector<int>& bins) const;
//! Get bin at a given position in space
//!
//! \param[in] r Position to get bin for
//! \return Mesh bin
int get_bin(Position r) const;
//! Get bin given mesh indices
//!
//! \param[in] Array of mesh indices
//! \return Mesh bin
int get_bin_from_indices(const int* ijk) const;
//! Get mesh indices given a position
//!
//! \param[in] r Position to get indices for
//! \param[out] ijk Array of mesh indices
//! \param[out] in_mesh Whether position is in mesh
void get_indices(Position r, int* ijk, bool* in_mesh) const;
//! Get mesh indices corresponding to a mesh bin
//!
//! \param[in] bin Mesh bin
//! \param[out] ijk Mesh indices
void get_indices_from_bin(int bin, int* ijk) const;
//! Check if a line connected by two points intersects the mesh
//!
//! \param[in] r0 Starting position
//! \param[in] r1 Ending position
//! \return Whether line connecting r0 and r1 intersects mesh
bool intersects(Position r0, Position r1) const;
//! Write mesh data to an HDF5 group
//!
//! \param[in] group HDF5 group
void to_hdf5(hid_t group) const;
//! Count number of bank sites in each mesh bin / energy bin
//!
//! \param[in] n Number of bank sites
//! \param[in] bank Array of bank sites
//! \param[in] n_energy Number of energies
//! \param[in] energies Array of energies
//! \param[out] Whether any bank sites are outside the mesh
//! \return Array indicating number of sites in each mesh/energy bin
xt::xarray<double> count_sites(int64_t n, const Bank* bank,
int n_energy, const double* energies, bool* outside) const;
int id_ {-1}; //!< User-specified ID
int n_dimension_; //!< Number of dimensions
double volume_frac_; //!< Volume fraction of each mesh element
xt::xarray<int> shape_; //!< Number of mesh elements in each dimension
xt::xarray<double> lower_left_; //!< Lower-left coordinates of mesh
xt::xarray<double> upper_right_; //!< Upper-right coordinates of mesh
xt::xarray<double> width_; //!< Width of each mesh element
private:
bool intersects_1d(Position r0, Position r1) const;
bool intersects_2d(Position r0, Position r1) const;
bool intersects_3d(Position r0, Position r1) const;
};
//==============================================================================
// Non-member functions
//==============================================================================
//! Read meshes from either settings/tallies
//! \param[in] root XML node
extern "C" void read_meshes(pugi::xml_node* root);
//! Write mesh data to an HDF5 group
//! \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

1
include/openmc/message_passing.h
View file

@ -10,6 +10,7 @@ namespace mpi {
extern int rank;
extern int n_procs;
extern bool master;
#ifdef OPENMC_MPI
extern MPI_Datatype bank;

26
include/openmc/mgxs.h
View file

@ -7,6 +7,8 @@
#include <string>
#include <vector>
#include "xtensor/xtensor.hpp"
#include "openmc/constants.h"
#include "openmc/hdf5_interface.h"
#include "openmc/xsdata.h"
@ -35,7 +37,7 @@ struct CacheData {
class Mgxs {
private:
double_1dvec kTs; // temperature in eV (k * T)
xt::xtensor<double, 1> kTs; // temperature in eV (k * T)
int scatter_format; // flag for if this is legendre, histogram, or tabular
int num_delayed_groups; // number of delayed neutron groups
int num_groups; // number of energy groups
@ -44,8 +46,8 @@ class Mgxs {
bool is_isotropic; // used to skip search for angle indices if isotropic
int n_pol;
int n_azi;
double_1dvec polar;
double_1dvec azimuthal;
std::vector<double> polar;
std::vector<double> azimuthal;
//! \brief Initializes the Mgxs object metadata
//!
@ -62,10 +64,10 @@ class Mgxs {
//! @param in_polar Polar angle grid.
//! @param in_azimuthal Azimuthal angle grid.
void
init(const std::string& in_name, double in_awr, const double_1dvec& in_kTs,
init(const std::string& in_name, double in_awr, const std::vector<double>& in_kTs,
bool in_fissionable, int in_scatter_format, int in_num_groups,
int in_num_delayed_groups, bool in_is_isotropic,
const double_1dvec& in_polar, const double_1dvec& in_azimuthal);
const std::vector<double>& in_polar, const std::vector<double>& in_azimuthal);
//! \brief Initializes the Mgxs object metadata from the HDF5 file
//!
@ -80,8 +82,8 @@ class Mgxs {
//! @param method Method of choosing nearest temperatures.
void
metadata_from_hdf5(hid_t xs_id, int in_num_groups,
int in_num_delayed_groups, const double_1dvec& temperature,
double tolerance, int_1dvec& temps_to_read, int& order_dim,
int in_num_delayed_groups, const std::vector<double>& temperature,
double tolerance, std::vector<int>& temps_to_read, int& order_dim,
int& method);
//! \brief Performs the actual act of combining the microscopic data for a
@ -93,8 +95,8 @@ class Mgxs {
//! corresponds to the temperature of interest.
//! @param this_t The temperature index of the macroscopic object.
void
combine(const std::vector<Mgxs*>& micros, const double_1dvec& scalars,
const int_1dvec& micro_ts, int this_t);
combine(const std::vector<Mgxs*>& micros, const std::vector<double>& scalars,
const std::vector<int>& micro_ts, int this_t);
//! \brief Checks to see if this and that are able to be combined
//!
@ -128,7 +130,7 @@ class Mgxs {
//! provides the number of points to use in the tabular representation.
//! @param method Method of choosing nearest temperatures.
Mgxs(hid_t xs_id, int energy_groups,
int delayed_groups, const double_1dvec& temperature, double tolerance,
int delayed_groups, const std::vector<double>& temperature, double tolerance,
int max_order, bool legendre_to_tabular,
int legendre_to_tabular_points, int& method);
@ -141,8 +143,8 @@ class Mgxs {
//! @param atom_densities Atom densities of those microscopic quantities.
//! @param tolerance Tolerance of temperature selection method.
//! @param method Method of choosing nearest temperatures.
Mgxs(const std::string& in_name, const double_1dvec& mat_kTs,
const std::vector<Mgxs*>& micros, const double_1dvec& atom_densities,
Mgxs(const std::string& in_name, const std::vector<double>& mat_kTs,
const std::vector<Mgxs*>& micros, const std::vector<double>& atom_densities,
double tolerance, int& method);
//! \brief Provides a cross section value given certain parameters

2
include/openmc/output.h
View file

@ -22,5 +22,7 @@ void header(const char* msg, int level);
extern "C" void print_overlap_check();
extern "C" void title();
} // namespace openmc
#endif // OPENMC_OUTPUT_H

2
include/openmc/particle.h
View file

@ -161,7 +161,7 @@ extern "C" {
{mark_as_lost(message.str());}
//! create a particle restart HDF5 file
void write_restart();
void write_restart() const;
};

4
include/openmc/position.h
View file

@ -1,6 +1,7 @@
#ifndef OPENMC_POSITION_H
#define OPENMC_POSITION_H
#include <cmath>
#include <vector>
namespace openmc {
@ -46,6 +47,9 @@ struct Position {
inline double dot(Position other) {
return x*other.x + y*other.y + z*other.z;
}
inline double norm() {
return std::sqrt(x*x + y*y + z*z);
}
// Data members
double x = 0.;

77
include/openmc/scattdata.h
View file

@ -6,6 +6,8 @@
#include <vector>
#include "xtensor/xtensor.hpp"
#include "openmc/constants.h"
namespace openmc {
@ -25,23 +27,25 @@ class ScattData {
protected:
//! \brief Initializes the attributes of the base class.
void
base_init(int order, const int_1dvec& in_gmin, const int_1dvec& in_gmax,
const double_2dvec& in_energy, const double_2dvec& in_mult);
base_init(int order, const xt::xtensor<int, 1>& in_gmin,
const xt::xtensor<int, 1>& in_gmax, const double_2dvec& in_energy,
const double_2dvec& in_mult);
//! \brief Combines microscopic ScattDatas into a macroscopic one.
void
base_combine(int max_order, const std::vector<ScattData*>& those_scatts,
const double_1dvec& scalars, int_1dvec& in_gmin, int_1dvec& in_gmax,
double_2dvec& sparse_mult, double_3dvec& sparse_scatter);
base_combine(size_t max_order, const std::vector<ScattData*>& those_scatts,
const std::vector<double>& scalars, xt::xtensor<int, 1>& in_gmin,
xt::xtensor<int, 1>& in_gmax, double_2dvec& sparse_mult,
double_3dvec& sparse_scatter);
public:
double_2dvec energy; // Normalized p0 matrix for sampling Eout
double_2dvec mult; // nu-scatter multiplication (nu-scatt/scatt)
double_3dvec dist; // Angular distribution
int_1dvec gmin; // minimum outgoing group
int_1dvec gmax; // maximum outgoing group
double_1dvec scattxs; // Isotropic Sigma_{s,g_{in}}
double_2dvec energy; // Normalized p0 matrix for sampling Eout
double_2dvec mult; // nu-scatter multiplication (nu-scatt/scatt)
double_3dvec dist; // Angular distribution
xt::xtensor<double, 1> gmin; // minimum outgoing group
xt::xtensor<double, 1> gmax; // maximum outgoing group
xt::xtensor<double, 1> scattxs; // Isotropic Sigma_{s,g_{in}}
//! \brief Calculates the value of normalized f(mu).
//!
@ -72,7 +76,7 @@ class ScattData {
//! @param in_mult Input sparse multiplicity matrix
//! @param coeffs Input sparse scattering matrix
virtual void
init(const int_1dvec& in_gmin, const int_1dvec& in_gmax,
init(const xt::xtensor<int, 1>& in_gmin, const xt::xtensor<int, 1>& in_gmax,
const double_2dvec& in_mult, const double_3dvec& coeffs) = 0;
//! \brief Combines the microscopic data.
@ -81,7 +85,7 @@ class ScattData {
//! @param scalars Scalars to multiply the microscopic data by.
virtual void
combine(const std::vector<ScattData*>& those_scatts,
const double_1dvec& scalars) = 0;
const std::vector<double>& scalars) = 0;
//! \brief Getter for the dimensionality of the scattering order.
//!
@ -89,7 +93,7 @@ class ScattData {
//! of points, and for Histogram this is the number of bins.
//!
//! @return The order.
virtual int
virtual size_t
get_order() = 0;
//! \brief Builds a dense scattering matrix from the constituent parts
@ -97,8 +101,8 @@ class ScattData {
//! @param max_order If Legendre this is the maximum value of "n" in "Pn"
//! requested; ignored otherwise.
//! @return The dense scattering matrix.
virtual double_3dvec
get_matrix(int max_order) = 0;
virtual xt::xtensor<double, 3>
get_matrix(size_t max_order) = 0;
//! \brief Samples the outgoing energy from the ScattData info.
//!
@ -142,12 +146,12 @@ class ScattDataLegendre: public ScattData {
public:
void
init(const int_1dvec& in_gmin, const int_1dvec& in_gmax,
init(const xt::xtensor<int, 1>& in_gmin, const xt::xtensor<int, 1>& in_gmax,
const double_2dvec& in_mult, const double_3dvec& coeffs);
void
combine(const std::vector<ScattData*>& those_scatts,
const double_1dvec& scalars);
const std::vector<double>& scalars);
//! \brief Find the maximal value of the angular distribution to use as a
// bounding box with rejection sampling.
@ -160,11 +164,11 @@ class ScattDataLegendre: public ScattData {
void
sample(int gin, int& gout, double& mu, double& wgt);
int
size_t
get_order() {return dist[0][0].size() - 1;};
double_3dvec
get_matrix(int max_order);
xt::xtensor<double, 3>
get_matrix(size_t max_order);
};
//==============================================================================
@ -176,19 +180,19 @@ class ScattDataHistogram: public ScattData {
protected:
double_1dvec mu; // Angle distribution mu bin boundaries
double dmu; // Quick storage of the spacing between the mu bin points
double_3dvec fmu; // The angular distribution histogram
xt::xtensor<double, 1> mu; // Angle distribution mu bin boundaries
double dmu; // Quick storage of the mu spacing
double_3dvec fmu; // The angular distribution histogram
public:
void
init(const int_1dvec& in_gmin, const int_1dvec& in_gmax,
init(const xt::xtensor<int, 1>& in_gmin, const xt::xtensor<int, 1>& in_gmax,
const double_2dvec& in_mult, const double_3dvec& coeffs);
void
combine(const std::vector<ScattData*>& those_scatts,
const double_1dvec& scalars);
const std::vector<double>& scalars);
double
calc_f(int gin, int gout, double mu);
@ -196,11 +200,11 @@ class ScattDataHistogram: public ScattData {
void
sample(int gin, int& gout, double& mu, double& wgt);
int
size_t
get_order() {return dist[0][0].size();};
double_3dvec
get_matrix(int max_order);
xt::xtensor<double, 3>
get_matrix(size_t max_order);
};
//==============================================================================
@ -212,9 +216,9 @@ class ScattDataTabular: public ScattData {
protected:
double_1dvec mu; // Angle distribution mu grid points
double dmu; // Quick storage of the spacing between the mu points
double_3dvec fmu; // The angular distribution function
xt::xtensor<double, 1> mu; // Angle distribution mu grid points
double dmu; // Quick storage of the mu spacing
double_3dvec fmu; // The angular distribution function
// Friend convert_legendre_to_tabular so it has access to protected
// parameters
@ -225,12 +229,12 @@ class ScattDataTabular: public ScattData {
public:
void
init(const int_1dvec& in_gmin, const int_1dvec& in_gmax,
init(const xt::xtensor<int, 1>& in_gmin, const xt::xtensor<int, 1>& in_gmax,
const double_2dvec& in_mult, const double_3dvec& coeffs);
void
combine(const std::vector<ScattData*>& those_scatts,
const double_1dvec& scalars);
const std::vector<double>& scalars);
double
calc_f(int gin, int gout, double mu);
@ -238,10 +242,11 @@ class ScattDataTabular: public ScattData {
void
sample(int gin, int& gout, double& mu, double& wgt);
int
size_t
get_order() {return dist[0][0].size();};
double_3dvec get_matrix(int max_order);
xt::xtensor<double, 3>
get_matrix(size_t max_order);
};
//==============================================================================

98
include/openmc/settings.h
View file

@ -5,6 +5,7 @@
//! \brief Settings for OpenMC
#include <array>
#include <cstdint>
#include <string>
#include "pugixml.hpp"
@ -15,39 +16,84 @@ namespace openmc {
// Global variable declarations
//==============================================================================
// Defined on Fortran side
extern "C" bool openmc_check_overlaps;
extern "C" bool openmc_particle_restart_run;
extern "C" bool openmc_photon_transport;
extern "C" bool openmc_restart_run;
extern "C" bool openmc_run_CE;
extern "C" int openmc_verbosity;
extern "C" bool openmc_write_all_tracks;
extern "C" bool openmc_write_initial_source;
namespace settings {
// Defined in .cpp
// TODO: Make strings instead of char* once Fortran is gone
extern "C" char* openmc_path_input;
extern "C" char* openmc_path_statepoint;
extern "C" char* openmc_path_sourcepoint;
extern "C" char* openmc_path_particle_restart;
extern std::string path_cross_sections;
extern std::string path_multipole;
extern std::string path_output;
// Boolean flags
extern "C" bool assume_separate; //!< assume tallies are spatially separate?
extern "C" bool check_overlaps; //!< check overlaps in geometry?
extern "C" bool cmfd_run; //!< use CMFD?
extern "C" bool confidence_intervals; //!< use confidence intervals for results?
extern "C" bool create_fission_neutrons; //!< create fission neutrons (fixed source)?
extern "C" bool entropy_on; //!< calculate Shannon entropy?
extern "C" bool legendre_to_tabular; //!< convert Legendre distributions to tabular?
extern "C" bool output_summary; //!< write summary.h5?
extern "C" bool output_tallies; //!< write tallies.out?
extern "C" bool particle_restart_run; //!< particle restart run?
extern "C" bool photon_transport; //!< photon transport turned on?
extern "C" bool reduce_tallies; //!< reduce tallies at end of batch?
extern "C" bool res_scat_on; //!< use resonance upscattering method?
extern "C" bool restart_run; //!< restart run?
extern "C" bool run_CE; //!< run with continuous-energy data?
extern "C" bool source_latest; //!< write latest source at each batch?
extern "C" bool source_separate; //!< write source to separate file?
extern "C" bool source_write; //!< write source in HDF5 files?
extern "C" bool survival_biasing; //!< use survival biasing?
extern "C" bool temperature_multipole; //!< use multipole data?
extern "C" bool trigger_on; //!< tally triggers enabled?
extern "C" bool trigger_predict; //!< predict batches for triggers?
extern "C" bool ufs_on; //!< uniform fission site method on?
extern "C" bool urr_ptables_on; //!< use unresolved resonance prob. tables?
extern "C" bool write_all_tracks; //!< write track files for every particle?
extern "C" bool write_initial_source; //!< write out initial source file?
// Paths to various files
extern std::string path_cross_sections; //!< path to cross_sections.xml
extern std::string path_input; //!< directory where main .xml files resides
extern std::string path_multipole; //!< directory containing multipole files
extern std::string path_output; //!< directory where output files are written
extern std::string path_particle_restart; //!< path to a particle restart file
extern std::string path_source;
extern std::string path_sourcepoint; //!< path to a source file
extern std::string path_statepoint; //!< path to a statepoint file
extern int temperature_method;
extern bool temperature_multipole;
extern double temperature_tolerance;
extern double temperature_default;
extern std::array<double, 2> temperature_range;
extern "C" int32_t index_entropy_mesh; //!< Index of entropy mesh in global mesh array
extern "C" int32_t index_ufs_mesh; //!< Index of UFS mesh in global mesh array
extern "C" int32_t n_batches; //!< number of (inactive+active) batches
extern "C" int32_t n_inactive; //!< number of inactive batches
extern "C" int32_t gen_per_batch; //!< number of generations per batch
extern "C" int64_t n_particles; //!< number of particles per generation
extern "C" int electron_treatment; //!< how to treat secondary electrons
extern "C" double energy_cutoff[4]; //!< Energy cutoff in [eV] for each particle type
extern "C" int legendre_to_tabular_points; //!< number of points to convert Legendres
extern "C" int max_order; //!< Maximum Legendre order for multigroup data
extern "C" int n_log_bins; //!< number of bins for logarithmic energy grid
extern "C" int n_max_batches; //!< Maximum number of batches
extern "C" int res_scat_method; //!< resonance upscattering method
extern "C" double res_scat_energy_min; //!< Min energy in [eV] for res. upscattering
extern "C" double res_scat_energy_max; //!< Max energy in [eV] for res. upscattering
extern "C" int run_mode; //!< Run mode (eigenvalue, fixed src, etc.)
extern "C" int temperature_method; //!< method for choosing temperatures
extern "C" double temperature_tolerance; //!< Tolerance in [K] on choosing temperatures
extern "C" double temperature_default; //!< Default T in [K]
extern "C" double temperature_range[2]; //!< Min/max T in [K] over which to load xs
extern "C" int trace_batch; //!< Batch to trace particle on
extern "C" int trace_gen; //!< Generation to trace particle on
extern "C" int64_t trace_particle; //!< Particle ID to enable trace on
extern "C" int trigger_batch_interval; //!< Batch interval for triggers
extern "C" int verbosity; //!< How verbose to make output
extern "C" double weight_cutoff; //!< Weight cutoff for Russian roulette
extern "C" double weight_survive; //!< Survival weight after Russian roulette
} // namespace settings
//==============================================================================
//! Read settings from XML file
//! \param[in] root XML node for <settings>
//==============================================================================
extern "C" void read_settings_xml();
extern "C" void read_settings(pugi::xml_node* root);
extern "C" void read_settings_xml_f(pugi::xml_node_struct* root_ptr);
} // namespace openmc

12
include/openmc/xml_interface.h
View file

@ -1,11 +1,14 @@
#ifndef OPENMC_XML_INTERFACE_H
#define OPENMC_XML_INTERFACE_H
#include <cstddef> // for size_t
#include <sstream> // for stringstream
#include <string>
#include <vector>
#include "pugixml.hpp"
#include "xtensor/xarray.hpp"
#include "xtensor/xadapt.hpp"
namespace openmc {
@ -38,5 +41,14 @@ std::vector<T> get_node_array(pugi::xml_node node, const char* name,
return values;
}
template <typename T>
xt::xarray<T> get_node_xarray(pugi::xml_node node, const char* name,
bool lowercase=false)
{
std::vector<T> v = get_node_array<T>(node, name, lowercase);
std::vector<std::size_t> shape = {v.size()};
return xt::adapt(v, shape);
}
} // namespace openmc
#endif // OPENMC_XML_INTERFACE_H

72
include/openmc/xsdata.h
View file

@ -7,6 +7,8 @@
#include <memory>
#include <vector>
#include "xtensor/xtensor.hpp"
#include "openmc/hdf5_interface.h"
#include "openmc/scattdata.h"
@ -22,41 +24,77 @@ class XsData {
private:
//! \brief Reads scattering data from the HDF5 file
void
scatter_from_hdf5(hid_t xsdata_grp, int n_pol, int n_azi, int energy_groups,
scatter_from_hdf5(hid_t xsdata_grp, size_t n_ang, size_t energy_groups,
int scatter_format, int final_scatter_format, int order_data,
int max_order, int legendre_to_tabular_points);
//! \brief Reads fission data from the HDF5 file
void
fission_from_hdf5(hid_t xsdata_grp, int n_pol, int n_azi, int energy_groups,
int delayed_groups, bool is_isotropic);
fission_from_hdf5(hid_t xsdata_grp, size_t n_ang, size_t energy_groups,
size_t delayed_groups, bool is_isotropic);
//! \brief Reads fission data formatted as chi and nu-fission vectors from
// the HDF5 file when beta is provided.
void
fission_vector_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups, size_t delayed_groups, bool is_isotropic);
//! \brief Reads fission data formatted as chi and nu-fission vectors from
// the HDF5 file when beta is not provided.
void
fission_vector_no_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups, size_t delayed_groups);
//! \brief Reads fission data formatted as chi and nu-fission vectors from
// the HDF5 file when no delayed data is provided.
void
fission_vector_no_delayed_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups);
//! \brief Reads fission data formatted as a nu-fission matrix from
// the HDF5 file when beta is provided.
void
fission_matrix_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups, size_t delayed_groups, bool is_isotropic);
//! \brief Reads fission data formatted as a nu-fission matrix from
// the HDF5 file when beta is not provided.
void
fission_matrix_no_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups, size_t delayed_groups);
//! \brief Reads fission data formatted as a nu-fission matrix from
// the HDF5 file when no delayed data is provided.
void
fission_matrix_no_delayed_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups);
public:
// The following quantities have the following dimensions:
// [angle][incoming group]
double_2dvec total;
double_2dvec absorption;
double_2dvec nu_fission;
double_2dvec prompt_nu_fission;
double_2dvec kappa_fission;
double_2dvec fission;
double_2dvec inverse_velocity;
xt::xtensor<double, 2> total;
xt::xtensor<double, 2> absorption;
xt::xtensor<double, 2> nu_fission;
xt::xtensor<double, 2> prompt_nu_fission;
xt::xtensor<double, 2> kappa_fission;
xt::xtensor<double, 2> fission;
xt::xtensor<double, 2> inverse_velocity;
// decay_rate has the following dimensions:
// [angle][delayed group]
double_2dvec decay_rate;
xt::xtensor<double, 2> decay_rate;
// delayed_nu_fission has the following dimensions:
// [angle][incoming group][delayed group]
double_3dvec delayed_nu_fission;
xt::xtensor<double, 3> delayed_nu_fission;
// chi_prompt has the following dimensions:
// [angle][incoming group][outgoing group]
double_3dvec chi_prompt;
xt::xtensor<double, 3> chi_prompt;
// chi_delayed has the following dimensions:
// [angle][incoming group][outgoing group][delayed group]
double_4dvec chi_delayed;
xt::xtensor<double, 4> chi_delayed;
// scatter has the following dimensions: [angle]
std::vector<std::shared_ptr<ScattData> > scatter;
std::vector<std::shared_ptr<ScattData>> scatter;
XsData() = default;
@ -68,7 +106,7 @@ class XsData {
//! @param scatter_format The scattering representation of the file.
//! @param n_pol Number of polar angles.
//! @param n_azi Number of azimuthal angles.
XsData(int num_groups, int num_delayed_groups, bool fissionable,
XsData(size_t num_groups, size_t num_delayed_groups, bool fissionable,
int scatter_format, int n_pol, int n_azi);
//! \brief Loads the XsData object from the HDF5 file
@ -101,7 +139,7 @@ class XsData {
//! @param micros Microscopic objects to combine.
//! @param scalars Scalars to multiply the microscopic data by.
void
combine(const std::vector<XsData*>& those_xs, const double_1dvec& scalars);
combine(const std::vector<XsData*>& those_xs, const std::vector<double>& scalars);
//! \brief Checks to see if this and that are able to be combined
//!