Skip to content

Add chi cross sections to MGXS machinery #3740

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
yardasol wants to merge 1 commit into
base: develop
Choose a base branch
from
+83 −10
Open

Add chi cross sections to MGXS machinery #3740

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions include/openmc/constants.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -285,6 +285,7 @@ enum class MgxsType {
PROMPT_NU_FISSION,
DELAYED_NU_FISSION,
NU_FISSION,
CHI,
CHI_PROMPT,
CHI_DELAYED
};
Expand Down
3 changes: 3 additions & 0 deletions include/openmc/xsdata.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -83,6 +83,9 @@ class XsData {
// delayed_nu_fission has the following dimensions:
// [angle][delayed group][incoming group]
xt::xtensor<double, 3> delayed_nu_fission;
// chi has the following dimensions:
// [angle][incoming group][outgoing group]
xt::xtensor<double, 3> chi;
// chi_prompt has the following dimensions:
// [angle][incoming group][outgoing group]
xt::xtensor<double, 3> chi_prompt;
Expand Down
38 changes: 38 additions & 0 deletions src/mgxs.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -434,31 +434,39 @@ double Mgxs::get_xs(MgxsType xstype, int gin, const int* gout, const double* mu,
{
XsData* xs_t = &xs[t];
double val;
std::string mgxs_type;
switch (xstype) {
case MgxsType::TOTAL:
val = xs_t->total(a, gin);
mgxs_type = "total";
break;
case MgxsType::NU_FISSION:
val = fissionable ? xs_t->nu_fission(a, gin) : 0.;
mgxs_type = "nu-fission";
break;
case MgxsType::ABSORPTION:
val = xs_t->absorption(a, gin);
mgxs_type = "absorption";
;
break;
case MgxsType::FISSION:
val = fissionable ? xs_t->fission(a, gin) : 0.;
mgxs_type = "fission";
break;
case MgxsType::KAPPA_FISSION:
val = fissionable ? xs_t->kappa_fission(a, gin) : 0.;
mgxs_type = "kappa-fission";
break;
case MgxsType::NU_SCATTER:
case MgxsType::SCATTER:
case MgxsType::NU_SCATTER_FMU:
case MgxsType::SCATTER_FMU:
val = xs_t->scatter[a]->get_xs(xstype, gin, gout, mu);
mgxs_type = "scatter_data";
break;
case MgxsType::PROMPT_NU_FISSION:
val = fissionable ? xs_t->prompt_nu_fission(a, gin) : 0.;
mgxs_type = "prompt-nu-fission";
break;
case MgxsType::DELAYED_NU_FISSION:
if (fissionable) {
Expand All @@ -473,6 +481,23 @@ double Mgxs::get_xs(MgxsType xstype, int gin, const int* gout, const double* mu,
} else {
val = 0.;
}
mgxs_type = "delayed-nu-fission";
break;
case MgxsType::CHI:
if (fissionable) {
if (gout != nullptr) {
val = xs_t->chi(a, gin, *gout);
} else {
// provide an outgoing group-wise sum
val = 0.;
for (int g = 0; g < xs_t->chi.shape()[2]; g++) {
val += xs_t->chi(a, gin, g);
}
}
} else {
val = 0.;
}
mgxs_type = "chi";
break;
case MgxsType::CHI_PROMPT:
if (fissionable) {
Expand All @@ -488,6 +513,7 @@ double Mgxs::get_xs(MgxsType xstype, int gin, const int* gout, const double* mu,
} else {
val = 0.;
}
mgxs_type = "chi-prompt";
break;
case MgxsType::CHI_DELAYED:
if (fissionable) {
Expand Down Expand Up @@ -515,20 +541,32 @@ double Mgxs::get_xs(MgxsType xstype, int gin, const int* gout, const double* mu,
} else {
val = 0.;
}
mgxs_type = "chi-delayed";
break;
case MgxsType::INVERSE_VELOCITY:
val = xs_t->inverse_velocity(a, gin);
mgxs_type = "inverse-velocity";
break;
case MgxsType::DECAY_RATE:
if (dg != nullptr) {
val = xs_t->decay_rate(a, *dg);
} else {
val = xs_t->decay_rate(a, 0);
}
mgxs_type = "decay-rate";
break;
default:
val = 0.;
}
// TODO: need to add more robust handling of zero-values cross sections that
// produce NANs on normalization
if (val != val) {
warning(
fmt::format("The {} cross section for this material has a nan present "
"(possibly due to a division by zero). Setting to zero...",
mgxs_type));
val = 0;
}
return val;
}

Expand Down
3 changes: 1 addition & 2 deletions src/random_ray/flat_source_domain.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1156,8 +1156,7 @@ void FlatSourceDomain::flatten_xs()
m.get_xs(MgxsType::FISSION, g_out, NULL, NULL, NULL, t, a);
sigma_f_.push_back(sigma_f);

double chi =
m.get_xs(MgxsType::CHI_PROMPT, g_out, &g_out, NULL, NULL, t, a);
double chi = m.get_xs(MgxsType::CHI, g_out, &g_out, NULL, NULL, t, a);
if (!std::isfinite(chi)) {
// MGXS interface may return NaN in some cases, such as when material
// is fissionable but has very small sigma_f.
Expand Down
48 changes: 40 additions & 8 deletions src/xsdata.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -56,8 +56,9 @@ XsData::XsData(bool fissionable, AngleDistributionType scatter_format,
// allocate delayed_nu_fission; [temperature][angle][delay group][in group]
delayed_nu_fission = xt::zeros<double>(shape);

// chi_prompt; [temperature][angle][in group][out group]
// chi and chi_prompt; [temperature][angle][in group][out group]
shape = {n_ang, n_g_, n_g_};
chi = xt::zeros<double>(shape);
chi_prompt = xt::zeros<double>(shape);

// chi_delayed; [temperature][angle][delay group][in group][out group]
Expand Down Expand Up @@ -133,8 +134,13 @@ void XsData::fission_vector_beta_from_hdf5(
// Normalize chi by summing over the outgoing groups for each incoming angle
temp_chi /= xt::view(xt::sum(temp_chi, {1}), xt::all(), xt::newaxis());

// Store chi in chi
chi = xt::view(temp_chi, xt::all(), xt::newaxis(), xt::all());

// Now every incoming group in prompt_chi and delayed_chi is the normalized
// chi we just made
// TODO: This is incorrect. This makes chi_prompt and chi_delayed identical
// to chi.
chi_prompt = xt::view(temp_chi, xt::all(), xt::newaxis(), xt::all());
chi_delayed =
xt::view(temp_chi, xt::all(), xt::newaxis(), xt::newaxis(), xt::all());
Expand Down Expand Up @@ -176,20 +182,31 @@ void XsData::fission_vector_beta_from_hdf5(

void XsData::fission_vector_no_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang)
{
// Data is provided separately as prompt + delayed nu-fission and chi
// If chi is included in the dataset, we should store it!
if (object_exists(xsdata_grp, "chi")) {
xt::xtensor<double, 2> temp_chi({n_ang, n_g_}, 0.);
read_nd_vector(xsdata_grp, "chi", temp_chi, true);
// Normalize chi by summing over the outgoing groups for each incoming angle
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
chi = xt::view(temp_chi, xt::all(), xt::newaxis(), xt::all());
}

// Data is provided separately as prompt + delayed nu-fission and chi
// Get chi-prompt
xt::xtensor<double, 2> temp_chi_p({n_ang, n_g_}, 0.);
read_nd_vector(xsdata_grp, "chi-prompt", temp_chi_p, true);

// Normalize chi by summing over the outgoing groups for each incoming angle
// Normalize chi-prompt by summing over the outgoing groups for each incoming
// angle
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, n_dg_, n_g_}, 0.);
read_nd_vector(xsdata_grp, "chi-delayed", temp_chi_d, true);

// Normalize chi by summing over the outgoing groups for each incoming angle
// Normalize chi-delayed by summing over the outgoing groups for each incoming
// angle
temp_chi_d /=
xt::view(xt::sum(temp_chi_d, {2}), xt::all(), xt::all(), xt::newaxis());

Expand Down Expand Up @@ -217,6 +234,7 @@ void XsData::fission_vector_no_delayed_from_hdf5(hid_t xsdata_grp, size_t n_ang)
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
chi = xt::view(temp_chi, xt::all(), xt::newaxis(), xt::all());
chi_prompt = xt::view(temp_chi, xt::all(), xt::newaxis(), xt::all());

// Get nu-fission directly
Expand Down Expand Up @@ -268,6 +286,10 @@ void XsData::fission_matrix_beta_from_hdf5(
xt::view(temp_beta, xt::all(), xt::all(), xt::newaxis(), xt::newaxis()) *
xt::view(temp_matrix, xt::all(), xt::newaxis(), xt::all(), xt::all());

// Store chi
chi =
chi_prompt * (1. - temp_beta_sum) + xt::sum(temp_beta * chi_delayed, {1});

} else if (beta_ndims == ndim_target + 1) {
xt::xtensor<double, 3> temp_beta({n_ang, n_dg_, n_g_}, 0.);
read_nd_vector(xsdata_grp, "beta", temp_beta, true);
Expand All @@ -293,14 +315,20 @@ void XsData::fission_matrix_beta_from_hdf5(
chi_delayed =
xt::view(temp_beta, xt::all(), xt::all(), xt::all(), xt::newaxis()) *
xt::view(temp_matrix, xt::all(), xt::newaxis(), xt::all(), xt::all());

// Store chi
chi =
chi_prompt * (1. - temp_beta_sum) + xt::sum(temp_beta * chi_delayed, {1});
}

// Normalize both chis
// Normalize chis
chi_prompt /=
xt::view(xt::sum(chi_prompt, {2}), xt::all(), xt::all(), xt::newaxis());

chi_delayed /= xt::view(
xt::sum(chi_delayed, {3}), xt::all(), xt::all(), xt::all(), xt::newaxis());

chi /= xt::view(xt::sum(chi, {2}), xt::all(), xt::all(), xt::newaxis());
}

void XsData::fission_matrix_no_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang)
Expand All @@ -326,8 +354,8 @@ void XsData::fission_matrix_no_beta_from_hdf5(hid_t xsdata_grp, size_t n_ang)
// delayed_nu_fission is the sum over outgoing groups
delayed_nu_fission = xt::sum(temp_matrix_d, {3});

// chi_prompt is this matrix but normalized over outgoing groups, which we
// have already stored in prompt_nu_fission
// chi_delayee is this matrix but normalized over outgoing groups, which we
// have already stored in delayed_nu_fission
chi_delayed = temp_matrix_d / xt::view(delayed_nu_fission, xt::all(),
xt::all(), xt::all(), xt::newaxis());
}
Expand All @@ -346,6 +374,8 @@ void XsData::fission_matrix_no_delayed_from_hdf5(hid_t xsdata_grp, size_t n_ang)

// chi_prompt is this matrix but normalized over outgoing groups, which we
// have already stored in prompt_nu_fission
chi = temp_matrix /
xt::view(prompt_nu_fission, xt::all(), xt::all(), xt::newaxis());
chi_prompt = temp_matrix /
xt::view(prompt_nu_fission, xt::all(), xt::all(), xt::newaxis());
}
Expand Down Expand Up @@ -525,6 +555,7 @@ void XsData::combine(
kappa_fission += scalar * that->kappa_fission;
fission += scalar * that->fission;
delayed_nu_fission += scalar * that->delayed_nu_fission;
chi += scalar * that->chi;
chi_prompt += scalar *
xt::view(xt::sum(that->prompt_nu_fission, {1}), xt::all(),
xt::newaxis(), xt::newaxis()) *
Expand All @@ -537,8 +568,9 @@ void XsData::combine(
decay_rate += scalar * that->decay_rate;
}

// Ensure the chi_prompt and chi_delayed are normalized to 1 for each
// Ensure chi, chi_prompt and chi_delayed are normalized to 1 for each
// azimuthal angle and delayed group (for chi_delayed)
chi /= xt::view(xt::sum(chi, {2}), xt::all(), xt::all(), xt::newaxis());
chi_prompt /=
xt::view(xt::sum(chi_prompt, {2}), xt::all(), xt::all(), xt::newaxis());
chi_delayed /= xt::view(
Expand Down
Loading