added specialized fission-getting methods

This commit is contained in:
Adam G Nelson 2018-09-02 18:59:28 -04:00
parent 35def7aac2
commit d6f94f7ce8
11 changed files with 678 additions and 368 deletions

View file

@ -159,7 +159,7 @@ class Mgxs {
//! @param dg delayed group index; use nullptr if irrelevant.
//! @return Requested cross section value.
double
get_xs(int xstype, int gin, int* gout, double* mu, int* dg);
get_xs(const int xstype, const int gin, int* gout, double* mu, int* dg);
//! \brief Samples the fission neutron energy and if prompt or delayed.
//!
@ -167,7 +167,7 @@ class Mgxs {
//! @param dg Sampled delayed group index.
//! @param gout Sampled outgoing energy group.
void
sample_fission_energy(int gin, int& dg, int& gout);
sample_fission_energy(const int gin, int& dg, int& gout);
//! \brief Samples the outgoing energy and angle from a scatter event.
//!
@ -176,7 +176,7 @@ class Mgxs {
//! @param mu Sampled cosine of the change-in-angle.
//! @param wgt Weight of the particle to be adjusted.
void
sample_scatter(int gin, int& gout, double& mu, double& wgt);
sample_scatter(const int gin, int& gout, double& mu, double& wgt);
//! \brief Calculates cross section quantities needed for tracking.
//!
@ -187,14 +187,14 @@ class Mgxs {
//! @param abs_xs Resultant absorption cross section.
//! @param nu_fiss_xs Resultant nu-fission cross section.
void
calculate_xs(int gin, double sqrtkT, const double uvw[3],
calculate_xs(const int gin, const double sqrtkT, const double uvw[3],
double& total_xs, double& abs_xs, double& nu_fiss_xs);
//! \brief Sets the temperature index in cache given a temperature
//!
//! @param sqrtkT Temperature of the material.
void
set_temperature_index(double sqrtkT);
set_temperature_index(const double sqrtkT);
//! \brief Sets the angle index in cache given a direction
//!

View file

@ -33,6 +33,42 @@ class XsData {
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);
//! \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);
//! \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:

View file

@ -183,7 +183,7 @@ contains
! Assign temperatures to cells that don't have temperatures already assigned
call assign_temperatures()
! Determine desired txemperatures for each nuclide and S(a,b) table
! Determine desired temperatures for each nuclide and S(a,b) table
call get_temperatures(nuc_temps, sab_temps)
! Check to make sure there are not too many nested coordinate levels in the

View file

@ -423,7 +423,7 @@ Mgxs::combine(const std::vector<Mgxs*>& micros, const std::vector<double>& scala
//==============================================================================
double
Mgxs::get_xs(int xstype, int gin, int* gout, double* mu, int* dg)
Mgxs::get_xs(const int xstype, const int gin, int* gout, double* mu, int* dg)
{
// This method assumes that the temperature and angle indices are set
#ifdef _OPENMP
@ -535,7 +535,7 @@ Mgxs::get_xs(int xstype, int gin, int* gout, double* mu, int* dg)
//==============================================================================
void
Mgxs::sample_fission_energy(int gin, int& dg, int& gout)
Mgxs::sample_fission_energy(const int gin, int& dg, int& gout)
{
// This method assumes that the temperature and angle indices are set
#ifdef _OPENMP
@ -599,7 +599,7 @@ Mgxs::sample_fission_energy(int gin, int& dg, int& gout)
//==============================================================================
void
Mgxs::sample_scatter(int gin, int& gout, double& mu, double& wgt)
Mgxs::sample_scatter(const int gin, int& gout, double& mu, double& wgt)
{
// This method assumes that the temperature and angle indices are set
// Sample the data
@ -614,7 +614,7 @@ Mgxs::sample_scatter(int gin, int& gout, double& mu, double& wgt)
//==============================================================================
void
Mgxs::calculate_xs(int gin, double sqrtkT, const double uvw[3],
Mgxs::calculate_xs(const int gin, const double sqrtkT, const double uvw[3],
double& total_xs, double& abs_xs, double& nu_fiss_xs)
{
// Set our indices
@ -649,7 +649,7 @@ Mgxs::equiv(const Mgxs& that)
//==============================================================================
void
Mgxs::set_temperature_index(double sqrtkT)
Mgxs::set_temperature_index(const double sqrtkT)
{
// See if we need to find the new index
#ifdef _OPENMP

View file

@ -121,6 +121,241 @@ XsData::from_hdf5(hid_t xsdata_grp, bool fissionable, int scatter_format,
//==============================================================================
void
XsData::fission_vector_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups, size_t delayed_groups)
{
// Data is provided as nu-fission and chi with a beta for delayed info
// Get chi
xt::xtensor<double, 2> temp_chi({n_ang, energy_groups}, 0.);
read_nd_vector(xsdata_grp, "chi", temp_chi, true);
// Normalize chi by summing over the outgoing groups for each incoming angle
temp_chi = temp_chi / xt::view(xt::sum(temp_chi, {1}), xt::all(), xt::newaxis());
// Now every incoming group in prompt_chi and delayed_chi is the normalized
// chi we just made
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, gin, gout) = temp_chi(a, gout);
for (size_t dg = 0; dg < delayed_groups; dg++) {
chi_delayed(a, gin, gout, dg) = temp_chi(a, gout);
}
}
}
}
// Get nu-fission
xt::xtensor<double, 2> temp_nufiss({n_ang, energy_groups}, 0.);
read_nd_vector(xsdata_grp, "nu-fission", temp_nufiss, true);
// Get beta
xt::xtensor<double, 2> temp_beta({n_ang, delayed_groups}, 0.);
read_nd_vector(xsdata_grp, "beta", temp_beta, true);
// Set prompt_nu_fission = (1. - beta_total)*nu_fission
prompt_nu_fission = temp_nufiss * (1. - xt::sum(temp_beta, {1}));
// Set delayed_nu_fission as beta * nu_fission
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
delayed_nu_fission(a, gin, dg) = temp_beta(a, dg) * temp_nufiss(a, gin);
}
}
}
}
void
XsData::fission_vector_no_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups, size_t delayed_groups)
{
// Data is provided separately as prompt + delayed nu-fission and chi
// Get chi-prompt
xt::xtensor<double, 2> temp_chi_p({n_ang, energy_groups}, 0.);
read_nd_vector(xsdata_grp, "chi-prompt", temp_chi_p, true);
// Normalize chi by summing over the outgoing groups for each incoming angle
temp_chi_p = temp_chi_p / xt::view(xt::sum(temp_chi_p, {1}), xt::all(), xt::newaxis());
// Get chi-delayed
xt::xtensor<double, 3> temp_chi_d({n_ang, energy_groups, delayed_groups}, 0.);
read_nd_vector(xsdata_grp, "chi-delayed", temp_chi_d, true);
// Normalize chi by summing over the outgoing groups for each incoming angle
temp_chi_d = temp_chi_d / xt::view(xt::sum(temp_chi_d, {1}), xt::all(),
xt::newaxis(), xt::all());
// Now assign the prompt and delayed chis by replicating for each incoming group
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, gin, gout) = temp_chi_p(a, gout);
for (size_t dg = 0; dg < delayed_groups; dg++) {
chi_delayed(a, gin, gout, dg) = temp_chi_d(a, gout, dg);
}
}
}
}
// Get prompt and delayed nu-fission directly
read_nd_vector(xsdata_grp, "prompt-nu-fission", prompt_nu_fission, true);
read_nd_vector(xsdata_grp, "delayed-nu-fission", delayed_nu_fission, true);
}
void
XsData::fission_vector_no_delayed_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups)
{
// No beta is provided and there is no prompt/delay distinction.
// Therefore, the code only considers the data as prompt.
// Get chi
xt::xtensor<double, 2> temp_chi({n_ang, energy_groups}, 0.);
read_nd_vector(xsdata_grp, "chi", temp_chi, true);
// Normalize chi by summing over the outgoing groups for each incoming angle
temp_chi = temp_chi / xt::view(xt::sum(temp_chi, {1}), xt::all(), xt::newaxis());
// Now every incoming group in self.chi is the normalized chi we just made
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, gin, gout) = temp_chi(a, gout);
}
}
}
// Get nu-fission directly
if (object_exists(xsdata_grp, "prompt-nu-fission")) {
read_nd_vector(xsdata_grp, "prompt-nu-fission", prompt_nu_fission, true);
} else {
read_nd_vector(xsdata_grp, "nu-fission", prompt_nu_fission, true);
}
}
//==============================================================================
void
XsData::fission_matrix_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups, size_t delayed_groups)
{
// Data is provided as nu-fission and chi with a beta for delayed info
// Get nu-fission matrix
xt::xtensor<double, 3> temp_matrix({n_ang, energy_groups, energy_groups});
read_nd_vector(xsdata_grp, "nu-fission", temp_matrix, true);
// Get beta
xt::xtensor<double, 3> temp_beta({n_ang, energy_groups, delayed_groups}, 0.);
read_nd_vector(xsdata_grp, "beta", temp_beta, true);
// prompt_nu_fission is the sum of this matrix over outgoing groups and
// multiplied by (1 - beta_tot)
prompt_nu_fission = xt::sum(temp_matrix, {2}) * (1. - xt::sum(temp_beta, {2}));
// delayed_nu_fission is the sum of this matrix over outgoing groups and
// multiplied by beta
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
double out_sum = 0.;
for (size_t gout = 0; gout < energy_groups; gout++) {
out_sum += temp_matrix(a, gin, gout);
}
for (size_t dg = 0; dg < delayed_groups; dg++) {
delayed_nu_fission(a, gin, dg) = temp_beta(a, dg) * out_sum;
}
}
}
// Store chi-prompt
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
double beta = 0.;
for (size_t dg = 0; dg < delayed_groups; dg++) {
beta += temp_beta(a, gin, dg);
}
chi_prompt(a, gin, gout) = (1.0 - beta) * temp_matrix(a, gin, gout);
}
}
}
// Store chi-delayed
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
chi_delayed(a, gin, gout, dg) = temp_beta(a, gin, dg) * temp_matrix(a, gin, gout);
}
}
}
}
//Normalize both
chi_prompt = chi_prompt / xt::view(xt::sum(chi_prompt, {2}), xt::all(), xt::newaxis());
chi_delayed = chi_delayed / xt::view(xt::sum(chi_delayed, {2}), xt::all(),
xt::newaxis(), xt::all());
}
void
XsData::fission_matrix_no_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups, size_t delayed_groups)
{
// Data is provided separately as prompt + delayed nu-fission and chi
// Get the prompt nu-fission matrix
xt::xtensor<double, 3> temp_matrix_p({n_ang, energy_groups, energy_groups});
read_nd_vector(xsdata_grp, "prompt-nu-fission", temp_matrix_p, true);
// prompt_nu_fission is the sum over outgoing groups
prompt_nu_fission = xt::sum(temp_matrix_p, {2});
// chi_prompt is this matrix but normalized over outgoing groups, which we
// have already stored in prompt_nu_fission
chi_prompt = temp_matrix_p / prompt_nu_fission;
// Get the delayed nu-fission matrix
xt::xtensor<double, 4> temp_matrix_d({n_ang, energy_groups, energy_groups,
delayed_groups});
read_nd_vector(xsdata_grp, "delayed-nu-fission", temp_matrix_d, true);
// delayed_nu_fission is the sum over outgoing groups
delayed_nu_fission = xt::sum(temp_matrix_d, {2});
// chi_prompt is this matrix but normalized over outgoing groups, which we
// have already stored in prompt_nu_fission
chi_delayed = temp_matrix_d / delayed_nu_fission;
}
void
XsData::fission_matrix_no_delayed_from_hdf5(hid_t xsdata_grp, size_t n_ang,
size_t energy_groups)
{
// No beta is provided and there is no prompt/delay distinction.
// Therefore, the code only considers the data as prompt.
// Get nu-fission matrix
xt::xtensor<double, 3> temp_matrix({n_ang, energy_groups, energy_groups});
if (object_exists(xsdata_grp, "prompt-nu-fission")) {
read_nd_vector(xsdata_grp, "prompt-nu-fission", temp_matrix, true);
} else {
read_nd_vector(xsdata_grp, "nu-fission", temp_matrix, true);
}
// prompt_nu_fission is the sum over outgoing groups
prompt_nu_fission = xt::sum(temp_matrix, {2});
// chi_prompt is this matrix but normalized over outgoing groups, which we
// have already stored in prompt_nu_fission
chi_prompt = temp_matrix / prompt_nu_fission;
}
//==============================================================================
void
XsData::fission_from_hdf5(hid_t xsdata_grp, size_t n_ang, size_t energy_groups,
size_t delayed_groups, bool is_isotropic)
@ -129,369 +364,31 @@ XsData::fission_from_hdf5(hid_t xsdata_grp, size_t n_ang, size_t energy_groups,
read_nd_vector(xsdata_grp, "fission", fission);
read_nd_vector(xsdata_grp, "kappa-fission", kappa_fission);
// Set/get beta
xt::xtensor<double, 3> temp_beta({n_ang, energy_groups, delayed_groups}, 0.);
if (object_exists(xsdata_grp, "beta")) {
hid_t xsdata = open_dataset(xsdata_grp, "beta");
size_t ndims = dataset_ndims(xsdata);
// raise ndims to make the isotropic ndims the same as angular
if (is_isotropic) ndims += 2;
if (ndims == 3) {
// Beta is input as [delayed group]
xt::xtensor<double, 1> temp_arr({n_ang * delayed_groups}, 0.);
read_nd_vector(xsdata_grp, "beta", temp_arr);
// Broadcast to all incoming groups
size_t temp_idx = 0;
for (size_t a = 0; a < n_ang; a++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
// Set the first group index and copy the rest
temp_beta(a, 0, dg) = temp_arr[temp_idx++];
for (size_t gin = 1; gin < energy_groups; gin++) {
temp_beta(a, gin, dg) = temp_beta(a, 0, dg);
}
}
}
} else if (ndims == 4) {
// Beta is input as [in group][delayed group]
read_nd_vector(xsdata_grp, "beta", temp_beta);
} else {
fatal_error("beta must be provided as a 3D or 4D array!");
}
}
// If chi is provided, set chi-prompt and chi-delayed
// Get the data; the strategy for doing so depends on if the data is provided
// as a nu-fission matrix or a set of chi and nu-fission vectors
if (object_exists(xsdata_grp, "chi")) {
xt::xtensor<double, 2> temp_arr ({n_ang, energy_groups}, 0.);
read_nd_vector(xsdata_grp, "chi", temp_arr);
for (size_t a = 0; a < n_ang; a++) {
// First set the first group
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, 0, gout) = temp_arr(a, gout);
}
// Now normalize this data
double chi_sum = 0.;
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_sum += chi_prompt(a, 0, gout);
}
if (chi_sum <= 0.) {
fatal_error("Encountered chi for a group that is <= 0!");
}
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, 0, gout) /= chi_sum;
}
// And extend to the remaining incoming groups
for (size_t gin = 1; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, gin, gout) = chi_prompt(a, 0, gout);
}
}
// Finally set chi-delayed equal to chi-prompt
// Set chi-delayed to chi-prompt
for(size_t gin = 0; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
chi_delayed(a, gin, gout, dg) = chi_prompt(a, gin, gout);
}
}
}
}
}
// If nu-fission is provided, set prompt- and delayed-nu-fission;
// if nu-fission is a matrix, set chi-prompt and chi-delayed.
if (object_exists(xsdata_grp, "nu-fission")) {
hid_t xsdata = open_dataset(xsdata_grp, "nu-fission");
size_t ndims = dataset_ndims(xsdata);
// raise ndims to make the isotropic ndims the same as angular
if (is_isotropic) ndims += 2;
if (ndims == 3) {
// nu-fission is a 3-d array
read_nd_vector(xsdata_grp, "nu-fission", prompt_nu_fission);
// set delayed-nu-fission and correct prompt-nu-fission with beta
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
delayed_nu_fission(a, gin, dg) =
temp_beta(a, gin, dg) * prompt_nu_fission(a, gin);
}
// Correct the prompt-nu-fission using the delayed neutron fraction
if (delayed_groups > 0) {
double beta_sum = 0.;
for (size_t gin = 0; gin < energy_groups; gin++) {
beta_sum += temp_beta(a, gin);
}
prompt_nu_fission(a, gin) *= (1. - beta_sum);
}
}
}
} else if (ndims == 4) {
// nu-fission is a matrix
read_nd_vector(xsdata_grp, "nu-fission", chi_prompt);
// Normalize the chi info so the CDF is 1.
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
double chi_sum = 0.;
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_sum += chi_prompt(a, gin, gout);
}
// Set the vector nu-fission from the matrix nu-fission
prompt_nu_fission(a, gin) = chi_sum;
if (chi_sum >= 0.) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, gin, gout) /= chi_sum;
}
} else {
fatal_error("Encountered chi for a group that is <= 0!");
}
}
// set all of chi-delayed to chi-prompt
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
chi_delayed(a, gin, gout, dg) = chi_prompt(a, gin, gout);
}
}
}
// Set the delayed-nu-fission and correct prompt-nu-fission with beta
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
delayed_nu_fission(a, gin, dg) = temp_beta(a, gin, dg) *
prompt_nu_fission(a, gin);
}
// Correct prompt-nu-fission using the delayed neutron fraction
if (delayed_groups > 0) {
double beta_sum = 0.;
for (size_t dg = 0; dg < delayed_groups; dg++) {
beta_sum += temp_beta(a, gin, dg);
}
prompt_nu_fission(a, gin) *= (1. - beta_sum);
}
}
}
if (delayed_groups == 0) {
fission_vector_no_delayed_from_hdf5(xsdata_grp, n_ang, energy_groups);
} else {
fatal_error("nu-fission must be provided as a 3D or 4D array!");
}
close_dataset(xsdata);
}
// If chi-prompt is provided, set chi-prompt
if (object_exists(xsdata_grp, "chi-prompt")) {
xt::xtensor<double, 2> temp_arr({n_ang, energy_groups}, 0.);
read_nd_vector(xsdata_grp, "chi-prompt", temp_arr);
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, gin, gout) = temp_arr(a, gout);
}
// Normalize chi so its CDF goes to 1
double chi_sum = 0.;
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_sum += chi_prompt(a, gin, gout);
}
if (chi_sum >= 0.) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, gin, gout) /= chi_sum;
}
} else {
fatal_error("Encountered chi-prompt for a group that is <= 0.!");
}
if (object_exists(xsdata_grp, "beta")) {
fission_vector_beta_from_hdf5(xsdata_grp, n_ang, energy_groups,
delayed_groups);
} else {
fission_vector_no_beta_from_hdf5(xsdata_grp, n_ang, energy_groups,
delayed_groups);
}
}
}
// If chi-delayed is provided, set chi-delayed
if (object_exists(xsdata_grp, "chi-delayed")) {
hid_t xsdata = open_dataset(xsdata_grp, "chi-delayed");
size_t ndims = dataset_ndims(xsdata);
// raise ndims to make the isotropic ndims the same as angular
if (is_isotropic) ndims += 2;
close_dataset(xsdata);
if (ndims == 3) {
// chi-delayed is a [in group] vector
xt::xtensor<double, 2> temp_arr({n_ang, energy_groups}, 0.);
read_nd_vector(xsdata_grp, "chi-delayed", temp_arr);
for (size_t a = 0; a < n_ang; a++) {
// normalize the chi CDF to 1
double chi_sum = 0.;
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_sum += temp_arr(a, gout);
}
if (chi_sum <= 0.) {
fatal_error("Encountered chi-delayed for a group that is <= 0!");
}
// set chi-delayed
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t gout = 0; gout < energy_groups; gout++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
chi_delayed(a, gin, gout, dg) = temp_arr(a, gout) / chi_sum;
}
}
}
}
} else if (ndims == 4) {
// chi_delayed is a matrix
read_nd_vector(xsdata_grp, "chi-delayed", chi_delayed);
// Normalize the chi info so the CDF is 1.
for (size_t a = 0; a < n_ang; a++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
double chi_sum = 0.;
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_sum += chi_delayed(a, gin, gout, dg);
}
if (chi_sum > 0.) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_delayed(a, gin, gout, dg) /= chi_sum;
}
} else {
fatal_error("Encountered chi-delayed for a group that is <= 0!");
}
}
}
}
} else {
if (delayed_groups == 0) {
fission_matrix_no_delayed_from_hdf5(xsdata_grp, n_ang, energy_groups);
} else {
fatal_error("chi-delayed must be provided as a 3D or 4D array!");
}
}
// Get prompt-nu-fission, if present
if (object_exists(xsdata_grp, "prompt-nu-fission")) {
hid_t xsdata = open_dataset(xsdata_grp, "prompt-nu-fission");
size_t ndims = dataset_ndims(xsdata);
// raise ndims to make the isotropic ndims the same as angular
if (is_isotropic) ndims += 2;
close_dataset(xsdata);
if (ndims == 3) {
// prompt-nu-fission is a [in group] vector
read_nd_vector(xsdata_grp, "prompt-nu-fission", prompt_nu_fission);
} else if (ndims == 4) {
// prompt nu fission is a matrix,
// so set prompt_nu_fiss & chi_prompt
xt::xtensor<double, 3> temp_arr({n_ang, energy_groups, energy_groups}, 0.);
read_nd_vector(xsdata_grp, "prompt-nu-fission", temp_arr);
// The prompt_nu_fission vector from the matrix form
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
double prompt_sum = 0.;
for (size_t gout = 0; gout < energy_groups; gout++) {
prompt_sum += temp_arr(a, gin, gout);
}
prompt_nu_fission(a, gin) = prompt_sum;
}
// The chi_prompt data is just the normalized fission matrix
for (size_t gin= 0; gin < energy_groups; gin++) {
if (prompt_nu_fission(a, gin) > 0.) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_prompt(a, gin, gout) =
temp_arr(a, gin, gout) / prompt_nu_fission(a, gin);
}
} else {
fatal_error("Encountered chi-prompt for a group that is <= 0!");
}
}
if (object_exists(xsdata_grp, "beta")) {
fission_matrix_beta_from_hdf5(xsdata_grp, n_ang, energy_groups,
delayed_groups);
} else {
fission_matrix_no_beta_from_hdf5(xsdata_grp, n_ang, energy_groups,
delayed_groups);
}
} else {
fatal_error("prompt-nu-fission must be provided as a 3D or 4D array!");
}
}
// Get delayed-nu-fission, if present
if (object_exists(xsdata_grp, "delayed-nu-fission")) {
hid_t xsdata = open_dataset(xsdata_grp, "delayed-nu-fission");
size_t ndims = dataset_ndims(xsdata);
close_dataset(xsdata);
// raise ndims to make the isotropic ndims the same as angular
if (is_isotropic) ndims += 2;
if (ndims == 3) {
// delayed-nu-fission is an [in group] vector
if (temp_beta(0, 0, 0) == 0.) {
fatal_error("cannot set delayed-nu-fission with a 1D array if "
"beta is not provided");
}
xt::xtensor<double, 2> temp_arr({n_ang, energy_groups}, 0.);
read_nd_vector(xsdata_grp, "delayed-nu-fission", temp_arr);
for (size_t a = 0; a < n_ang; a++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
// Set delayed-nu-fission using beta
delayed_nu_fission(a, gin, dg) =
temp_beta(a, gin, dg) * temp_arr(a, gin);
}
}
}
} else if (ndims == 4) {
read_nd_vector(xsdata_grp, "delayed-nu-fission",
delayed_nu_fission);
} else if (ndims == 5) {
// This will contain delayed-nu-fission and chi-delayed data
xt::xtensor<double, 4> temp_arr({n_ang, energy_groups, energy_groups,
delayed_groups}, 0.);
read_nd_vector(xsdata_grp, "delayed-nu-fission", temp_arr);
// Set the 3D delayed-nu-fission matrix and 4D chi-delayed matrix
// from the 4D delayed-nu-fission matrix
for (size_t a = 0; a < n_ang; a++) {
for (size_t dg = 0; dg < delayed_groups; dg++) {
for (size_t gin = 0; gin < energy_groups; gin++) {
double gout_sum = 0.;
for (size_t gout = 0; gout < energy_groups; gout++) {
gout_sum += temp_arr(a, gin, gout, dg);
chi_delayed(a, gin, gout, dg) = temp_arr(a, gin, gout, dg);
}
delayed_nu_fission(a, gin, dg) = gout_sum;
// Normalize chi-delayed
if (gout_sum > 0.) {
for (size_t gout = 0; gout < energy_groups; gout++) {
chi_delayed(a, gin, gout, dg) /= gout_sum;
}
} else {
fatal_error("Encountered chi-delayed for a group that is <= 0!");
}
}
}
}
} else {
fatal_error("prompt-nu-fission must be provided as a 3D, 4D, or 5D "
"array!");
}
}

View file

@ -0,0 +1,33 @@
<?xml version='1.0' encoding='utf-8'?>
<geometry>
<cell id="1" material="1" region="1 -2" universe="0" />
<surface boundary="reflective" coeffs="0.0" id="1" type="x-plane" />
<surface boundary="vacuum" coeffs="929.45" id="2" type="x-plane" />
</geometry>
<?xml version='1.0' encoding='utf-8'?>
<materials>
<cross_sections>2g.h5</cross_sections>
<material id="1" name="mat_1">
<density units="macro" value="1.0" />
<macroscopic name="mat_1" />
</material>
</materials>
<?xml version='1.0' encoding='utf-8'?>
<settings>
<run_mode>eigenvalue</run_mode>
<particles>1000</particles>
<batches>10</batches>
<inactive>5</inactive>
<source strength="1.0">
<space type="box">
<parameters>-929.45 -1e+50 -1e+50 929.45 1e+50 1e+50</parameters>
</space>
</source>
<output>
<summary>false</summary>
</output>
<energy_mode>multi-group</energy_mode>
<tabular_legendre>
<enable>false</enable>
</tabular_legendre>
</settings>

View file

@ -0,0 +1,2 @@
k-combined:
9.949396E-01 2.047218E-02

View file

@ -0,0 +1,167 @@
import os
import numpy as np
import openmc
from tests.testing_harness import PyAPITestHarness
def create_library():
# Instantiate the energy group data and file object
groups = openmc.mgxs.EnergyGroups(group_edges=[0.0, 0.625, 20.0e6])
mg_cross_sections_file = openmc.MGXSLibrary(groups)
# Make the base, isotropic data
nu = [2.50, 2.50]
fiss = np.array([0.002817, 0.097])
capture = [0.008708, 0.02518]
absorption = np.add(capture, fiss)
scatter = np.array(
[[[0.31980, 0.06694], [0.004555, -0.0003972]],
[[0.00000, 0.00000], [0.424100, 0.05439000]]])
total = [0.33588, 0.54628]
chi = [1., 0.]
mat_1 = openmc.XSdata('mat_1', groups)
mat_1.order = 1
mat_1.set_nu_fission(np.multiply(nu, fiss))
mat_1.set_absorption(absorption)
mat_1.set_scatter_matrix(scatter)
mat_1.set_total(total)
mat_1.set_chi(chi)
mg_cross_sections_file.add_xsdata(mat_1)
# Make a version of mat-1 which has a tabular representation of the
# scattering vice Legendre with 33 points
mat_2 = mat_1.convert_scatter_format('tabular', 33)
mat_2.name = 'mat_2'
mg_cross_sections_file.add_xsdata(mat_2)
# Make a version of mat-1 which has a histogram representation of the
# scattering vice Legendre with 33 bins
mat_3 = mat_1.convert_scatter_format('histogram', 33)
mat_3.name = 'mat_3'
mg_cross_sections_file.add_xsdata(mat_3)
# Make a version which uses a fission matrix vice chi & nu-fission
mat_4 = openmc.XSdata('mat_4', groups)
mat_4.order = 1
mat_4.set_nu_fission(np.outer(np.multiply(nu, fiss), chi))
mat_4.set_absorption(absorption)
mat_4.set_scatter_matrix(scatter)
mat_4.set_total(total)
mg_cross_sections_file.add_xsdata(mat_4)
# Make an angle-dependent version of mat_1 with 2 polar and 2 azim. angles
mat_5 = mat_1.convert_representation('angle', 2, 2)
mat_5.name = 'mat_5'
mg_cross_sections_file.add_xsdata(mat_5)
# Make a copy of mat_1 for testing microscopic cross sections
mat_6 = openmc.XSdata('mat_6', groups)
mat_6.order = 1
mat_6.set_nu_fission(np.multiply(nu, fiss))
mat_6.set_absorption(absorption)
mat_6.set_scatter_matrix(scatter)
mat_6.set_total(total)
mat_6.set_chi(chi)
mg_cross_sections_file.add_xsdata(mat_6)
# Write the file
mg_cross_sections_file.export_to_hdf5('2g.h5')
def create_model():
create_library()
# # Make Materials
materials_file = openmc.Materials()
mat_names = ['base leg', 'base tab', 'base hist', 'base matrix', 'base ang']
macros = []
mats = []
for i in range(len(mat_names)):
macros.append(openmc.Macroscopic('mat_' + str(i + 1)))
mats.append(openmc.Material(name=mat_names[i]))
mats[-1].set_density('macro', 1.0)
mats[-1].add_macroscopic(macros[-1])
# Add in the microscopic data
mats.append(openmc.Material(name='micro'))
mats[-1].set_density("sum")
mats[-1].add_nuclide("mat_1", 0.5)
mats[-1].add_nuclide("mat_6", 0.5)
materials_file += mats
materials_file.cross_sections = '2g.h5'
# # Make Geometry
rad_outer = 929.45
# Set a cell boundary to exist for every material above (exclude the 0)
rads = np.linspace(0., rad_outer, len(mats) + 1, endpoint=True)[1:]
# Instantiate Universe
root = openmc.Universe(universe_id=0, name='root universe')
cells = []
surfs = []
surfs.append(openmc.XPlane(x0=0., boundary_type='reflective'))
for r, rad in enumerate(rads):
if r == len(rads) - 1:
surfs.append(openmc.XPlane(x0=rad, boundary_type='vacuum'))
else:
surfs.append(openmc.XPlane(x0=rad))
# Instantiate Cells
cells = []
for c in range(len(surfs) - 1):
cells.append(openmc.Cell())
cells[-1].region = (+surfs[c] & -surfs[c + 1])
cells[-1].fill = mats[c]
# Register Cells with Universe
root.add_cells(cells)
# Instantiate a Geometry, register the root Universe, and export to XML
geometry_file = openmc.Geometry(root)
# # Make Settings
# Instantiate a Settings object, set all runtime parameters
settings_file = openmc.Settings()
settings_file.energy_mode = "multi-group"
settings_file.tabular_legendre = {'enable': False}
settings_file.batches = 10
settings_file.inactive = 5
settings_file.particles = 1000
# Build source distribution
INF = 1000.
bounds = [0., -INF, -INF, rads[0], INF, INF]
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:])
settings_file.source = openmc.source.Source(space=uniform_dist)
settings_file.output = {'summary': False}
model = openmc.model.Model()
model.geometry = geometry_file
model.materials = materials_file
model.settings = settings_file
return model
class MGXSTestHarness(PyAPITestHarness):
def _cleanup(self):
super()._cleanup()
f = '2g .h5'
if os.path.exists(f):
os.remove(f)
def test_mg_benchmark():
model = create_model()
harness = PyAPITestHarness('statepoint.10.h5', model)
harness.main()

View file

@ -0,0 +1,175 @@
import os
import numpy as np
import openmc
from tests.testing_harness import PyAPITestHarness
def create_library():
# Instantiate the energy group data and file object
groups = openmc.mgxs.EnergyGroups(group_edges=[0.0, 0.625, 20.0e6])
n_dg = 2
mg_cross_sections_file = openmc.MGXSLibrary(groups)
mg_cross_sections_file.num_delayed_groups = n_dg
beta = np.array([0.003, 0.003])
one_m_beta = 1. - np.sum(beta)
nu = [2.50, 2.50]
fiss = np.array([0.002817, 0.097])
capture = [0.008708, 0.02518]
absorption = np.add(capture, fiss)
scatter = np.array(
[[[0.31980, 0.06694], [0.004555, -0.0003972]],
[[0.00000, 0.00000], [0.424100, 0.05439000]]])
total = [0.33588, 0.54628]
chi = [1., 0.]
# Make the base data that uses chi & nu-fission vectors with a beta
mat_1 = openmc.XSdata('mat_1', groups)
mat_1.order = 1
mat_1.num_delayed_groups = 2
mat_1.set_beta(beta)
mat_1.set_nu_fission(np.multiply(nu, fiss))
mat_1.set_absorption(absorption)
mat_1.set_scatter_matrix(scatter)
mat_1.set_total(total)
mat_1.set_chi(chi)
mg_cross_sections_file.add_xsdata(mat_1)
# Make a version that uses prompt and delayed version of nufiss and chi
mat_2 = openmc.XSdata('mat_2', groups)
mat_2.order = 1
mat_2.num_delayed_groups = 2
mat_2.set_prompt_nu_fission(one_m_beta * np.multiply(nu, fiss))
delay_nu_fiss = np.zeros((n_dg, groups.num_groups))
for dg in range(n_dg):
for g in range(groups.num_groups):
delay_nu_fiss[dg, g] = beta[dg] * nu[g] * fiss[g]
mat_2.set_delayed_nu_fission(delay_nu_fiss)
mat_2.set_absorption(absorption)
mat_2.set_scatter_matrix(scatter)
mat_2.set_total(total)
mat_2.set_chi_prompt(chi)
mat_2.set_chi_delayed(np.stack([chi] * n_dg))
mg_cross_sections_file.add_xsdata(mat_2)
# Make a version that uses a nu-fission matrix with a beta
mat_3 = openmc.XSdata('mat_3', groups)
mat_3.order = 1
mat_3.num_delayed_groups = 2
mat_3.set_beta(beta)
mat_3.set_nu_fission(np.outer(np.multiply(nu, fiss), chi))
mat_3.set_absorption(absorption)
mat_3.set_scatter_matrix(scatter)
mat_3.set_total(total)
mg_cross_sections_file.add_xsdata(mat_3)
# Make a version that uses prompt and delayed version of the nufiss matrix
mat_4 = openmc.XSdata('mat_4', groups)
mat_4.order = 1
mat_4.num_delayed_groups = 2
mat_4.set_prompt_nu_fission(one_m_beta * np.outer(np.multiply(nu, fiss), chi))
delay_nu_fiss = np.zeros((n_dg, groups.num_groups, groups.num_groups))
for dg in range(n_dg):
for g in range(groups.num_groups):
for go in range(groups.num_groups):
delay_nu_fiss[dg, g, go] = beta[dg] * nu[g] * fiss[g] * chi[go]
mat_4.set_delayed_nu_fission(delay_nu_fiss)
mat_4.set_absorption(absorption)
mat_4.set_scatter_matrix(scatter)
mat_4.set_total(total)
mg_cross_sections_file.add_xsdata(mat_4)
# Write the file
mg_cross_sections_file.export_to_hdf5('2g.h5')
def create_model():
create_library()
# # Make Materials
materials_file = openmc.Materials()
mat_names = ['vec beta', 'vec no beta', 'matrix beta', 'matrix no beta']
macros = []
mats = []
for i in range(len(mat_names)):
macros.append(openmc.Macroscopic('mat_' + str(i + 1)))
mats.append(openmc.Material(name=mat_names[i]))
mats[-1].set_density('macro', 1.0)
mats[-1].add_macroscopic(macros[-1])
materials_file += mats
materials_file.cross_sections = '2g.h5'
# # Make Geometry
rad_outer = 929.45
# Set a cell boundary to exist for every material above (exclude the 0)
rads = np.linspace(0., rad_outer, len(mats) + 1, endpoint=True)[1:]
# Instantiate Universe
root = openmc.Universe(universe_id=0, name='root universe')
cells = []
surfs = []
surfs.append(openmc.XPlane(x0=0., boundary_type='reflective'))
for r, rad in enumerate(rads):
if r == len(rads) - 1:
surfs.append(openmc.XPlane(x0=rad, boundary_type='vacuum'))
else:
surfs.append(openmc.XPlane(x0=rad))
# Instantiate Cells
cells = []
for c in range(len(surfs) - 1):
cells.append(openmc.Cell())
cells[-1].region = (+surfs[c] & -surfs[c + 1])
cells[-1].fill = mats[c]
# Register Cells with Universe
root.add_cells(cells)
# Instantiate a Geometry, register the root Universe, and export to XML
geometry_file = openmc.Geometry(root)
# # Make Settings
# Instantiate a Settings object, set all runtime parameters
settings_file = openmc.Settings()
settings_file.energy_mode = "multi-group"
settings_file.tabular_legendre = {'enable': False}
settings_file.batches = 10
settings_file.inactive = 5
settings_file.particles = 1000
# Build source distribution
INF = 1000.
bounds = [0., -INF, -INF, rads[0], INF, INF]
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:])
settings_file.source = openmc.source.Source(space=uniform_dist)
settings_file.output = {'summary': False}
model = openmc.model.Model()
model.geometry = geometry_file
model.materials = materials_file
model.settings = settings_file
return model
class MGXSTestHarness(PyAPITestHarness):
def _cleanup(self):
super()._cleanup()
f = '2g .h5'
if os.path.exists(f):
os.remove(f)
def test_mg_benchmark():
model = create_model()
harness = PyAPITestHarness('statepoint.10.h5', model)
harness.main()