From 7ae1b024ed044c3e098cdba17d4e5f536e1dd779 Mon Sep 17 00:00:00 2001 From: liangjg Date: Mon, 18 Mar 2019 19:20:51 -0400 Subject: [PATCH] do heating tallies for both neutron and photon --- include/openmc/constants.h | 3 ++- openmc/capi/tally.py | 2 +- src/output.cpp | 1 + src/photon.cpp | 10 ++++++--- src/reaction.cpp | 2 ++ src/tallies/tally.cpp | 5 +++-- src/tallies/tally_scoring.cpp | 41 +++++++++++++++++++++++++++++++++++ 7 files changed, 57 insertions(+), 7 deletions(-) diff --git a/include/openmc/constants.h b/include/openmc/constants.h index d9307abe2..dd9ecf211 100644 --- a/include/openmc/constants.h +++ b/include/openmc/constants.h @@ -230,7 +230,7 @@ constexpr int N_XD {204}; constexpr int N_XT {205}; constexpr int N_X3HE {206}; constexpr int N_XA {207}; -constexpr int HEATING {301}; +constexpr int NEUTRON_HEATING {301}; constexpr int DAMAGE_ENERGY {444}; constexpr int COHERENT {502}; constexpr int INCOHERENT {504}; @@ -366,6 +366,7 @@ constexpr int SCORE_INVERSE_VELOCITY {-13}; // flux-weighted inverse velocity constexpr int SCORE_FISS_Q_PROMPT {-14}; // prompt fission Q-value constexpr int SCORE_FISS_Q_RECOV {-15}; // recoverable fission Q-value constexpr int SCORE_DECAY_RATE {-16}; // delayed neutron precursor decay rate +constexpr int SCORE_HEATING {-17}; // nuclear heating (neutron or photon) // Tally map bin finding constexpr int NO_BIN_FOUND {-1}; diff --git a/openmc/capi/tally.py b/openmc/capi/tally.py index eec3a2b66..9529a31f2 100644 --- a/openmc/capi/tally.py +++ b/openmc/capi/tally.py @@ -89,7 +89,7 @@ _SCORES = { -5: 'absorption', -6: 'fission', -7: 'nu-fission', -8: 'kappa-fission', -9: 'current', -10: 'events', -11: 'delayed-nu-fission', -12: 'prompt-nu-fission', -13: 'inverse-velocity', -14: 'fission-q-prompt', - -15: 'fission-q-recoverable', -16: 'decay-rate' + -15: 'fission-q-recoverable', -16: 'decay-rate', -17: 'heating' } _ESTIMATORS = { 1: 'analog', 2: 'tracklength', 3: 'collision' diff --git a/src/output.cpp b/src/output.cpp index df0976d94..0c7ebcbda 100644 --- a/src/output.cpp +++ b/src/output.cpp @@ -618,6 +618,7 @@ const std::unordered_map score_names = { {SCORE_FISS_Q_PROMPT, "Prompt fission power"}, {SCORE_FISS_Q_RECOV, "Recoverable fission power"}, {SCORE_CURRENT, "Current"}, + {SCORE_HEATING, "Heating"}, }; //! Create an ASCII output file showing all tally results. diff --git a/src/photon.cpp b/src/photon.cpp index 86bf4390a..20000c4a9 100644 --- a/src/photon.cpp +++ b/src/photon.cpp @@ -97,9 +97,13 @@ PhotonInteraction::PhotonInteraction(hid_t group, int i_element) close_group(rgroup); // Read heating - rgroup = open_group(group, "heating"); - read_dataset(rgroup, "xs", heating_); - close_group(rgroup); + if (object_exists(group, "heating")) { + rgroup = open_group(group, "heating"); + read_dataset(rgroup, "xs", heating_); + close_group(rgroup); + } else { + heating_ = xt::zeros_like(energy_); + } // Read subshell photoionization cross section and atomic relaxation data rgroup = open_group(group, "subshells"); diff --git a/src/reaction.cpp b/src/reaction.cpp index d2405d899..222822230 100644 --- a/src/reaction.cpp +++ b/src/reaction.cpp @@ -121,6 +121,8 @@ std::string reaction_name(int mt) return "fission-q-prompt"; } else if (mt == SCORE_FISS_Q_RECOV) { return "fission-q-recoverable"; + } else if (mt == SCORE_HEATING) { + return "heating"; // Normal ENDF-based reactions } else if (mt == TOTAL_XS) { diff --git a/src/tallies/tally.cpp b/src/tallies/tally.cpp index 682db3540..d26725db2 100644 --- a/src/tallies/tally.cpp +++ b/src/tallies/tally.cpp @@ -109,6 +109,9 @@ score_str_to_int(std::string score_str) if (score_str == "fission-q-recoverable") return SCORE_FISS_Q_RECOV; + if (score_str == "heating") + return SCORE_HEATING; + if (score_str == "current") return SCORE_CURRENT; @@ -207,8 +210,6 @@ score_str_to_int(std::string score_str) return N_X3HE; if (score_str == "(n,Xa)" || score_str == "He4-production") return N_XA; - if (score_str == "heating") - return HEATING; if (score_str == "damage-energy") return DAMAGE_ENERGY; diff --git a/src/tallies/tally_scoring.cpp b/src/tallies/tally_scoring.cpp index 7e235fa45..d574ed109 100644 --- a/src/tallies/tally_scoring.cpp +++ b/src/tallies/tally_scoring.cpp @@ -7,6 +7,7 @@ #include "openmc/material.h" #include "openmc/mgxs_interface.h" #include "openmc/nuclide.h" +#include "openmc/photon.h" #include "openmc/reaction_product.h" #include "openmc/search.h" #include "openmc/settings.h" @@ -1148,6 +1149,46 @@ score_general_ce(Particle* p, int i_tally, int start_index, break; + case SCORE_HEATING: + if (p->type_ == Particle::Type::neutron) { + score_bin = NEUTRON_HEATING; + // No break here, continue to run the default neutron MT scoring + } else if (p->type_ == Particle::Type::photon) { + if (tally.estimator_ != ESTIMATOR_ANALOG) { + // Calculate photon heating cross section on-the-fly + score = 0.; + if (i_nuclide >= 0) { + // Find the element corresponding to the nuclide + auto name = data::nuclides[i_nuclide]->name_; + int pos = name.find_first_of("0123456789"); + std::string element = name.substr(0, pos); + int i_element = data::element_map[element]; + auto& heating {data::elements[i_element].heating_}; + auto i_grid = simulation::micro_photon_xs[i_element].index_grid; + auto f = simulation::micro_photon_xs[i_element].interp_factor; + score = std::exp(heating(i_grid) + f * (heating(i_grid+1) - + heating(i_grid))) * atom_density * flux; + } else { + if (p->material_ != MATERIAL_VOID) { + const Material& material {*model::materials[p->material_]}; + for (auto i = 0; i < material.nuclide_.size(); ++i) { + auto i_element = material.element_[i]; + auto atom_density = material.atom_density_(i); + auto& heating {data::elements[i_element].heating_}; + auto i_grid = simulation::micro_photon_xs[i_element].index_grid; + auto f = simulation::micro_photon_xs[i_element].interp_factor; + score += std::exp(heating(i_grid) + f * (heating(i_grid+1) - + heating(i_grid))) * atom_density * flux; + } + } + } + } + break; + } else { + // Nuclear heating of any other particles not implmented + break; + } + default: if (tally.estimator_ == ESTIMATOR_ANALOG) { // Any other score is assumed to be a MT number. Thus, we just need