Add energy deposition score

The energy-deposition score is computed using three helper
functions, previously introduced in #1329 to score neutron
heating and recoverable fission energy release. This score
also relies on the modified fission energy release dataset
introduced in 2cecbd5d0f
and stored at Nuclide::modified_fission_q_.

Given the helper functions, the implementation is straightforward,
as weighting, biasing, and other tally factors are included.
The score is total_heating - fission_heating + modified_fission_q.
The fission heating is removed as its contribution to the total
heating is also included in the modified_fission_q part.

The score has been added as -18 and is included in:
- constants.h
- openmc/capi/tally.py
- output.cpp
- reaction.cpp
- tally.cpp

Analysis demonstrates that this tally is able to obtain
energy deposited in fissile material and non-fissile materials.
No regression tests have been added at this point.

include/openmc/constants.h

@@ -371,6 +371,7 @@ 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)
+constexpr int SCORE_ENERGY_DEPOSITION  {-18}; // energy deposition
 
 // Tally map bin finding
 constexpr int NO_BIN_FOUND {-1};

openmc/capi/tally.py

@@ -89,7 +89,8 @@ _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', -17: 'heating'
+    -15: 'fission-q-recoverable', -16: 'decay-rate', -17: 'heating',
+    -18: "energy-deposition",
 }
 _ESTIMATORS = {
     1: 'analog', 2: 'tracklength', 3: 'collision'

src/output.cpp

@@ -621,6 +621,7 @@ const std::unordered_map<int, const char*> score_names = {
   {SCORE_FISS_Q_RECOV,       "Recoverable fission power"},
   {SCORE_CURRENT,            "Current"},
   {SCORE_HEATING,            "Heating"},
+  {SCORE_ENERGY_DEPOSITION,  "Energy Deposition"},
 };
 
 //! Create an ASCII output file showing all tally results.

src/reaction.cpp

@@ -121,6 +121,8 @@ std::string reaction_name(int mt)
     return  "fission-q-recoverable";
   } else if (mt == SCORE_HEATING) {
     return  "heating";
+  } else if (mt == SCORE_ENERGY_DEPOSITION) {
+    return "energy-deposition";
 
   // Normal ENDF-based reactions
   } else if (mt == TOTAL_XS) {

src/tallies/tally.cpp

@@ -113,6 +113,9 @@ score_str_to_int(std::string score_str)
   if (score_str == "heating")
     return SCORE_HEATING;
 
+  if (score_str == "energy-deposition")
+    return SCORE_ENERGY_DEPOSITION;
+
   if (score_str == "current")
     return SCORE_CURRENT;
 

src/tallies/tally_scoring.cpp

@@ -183,6 +183,10 @@ double get_nuc_fission_q(const Nuclide& nuc, const Particle* p, int score_bin)
     if (nuc.fission_q_recov_) {
       return (*nuc.fission_q_recov_)(p->E_last_);
     }
+  } else if (score_bin == SCORE_ENERGY_DEPOSITION) {
+    if (nuc.modified_fission_q_) {
+      return (*nuc.modified_fission_q_)(p->E_last_);
+    }
   }
   return 0.0;
 }
@@ -1209,7 +1213,6 @@ score_general_ce(Particle* p, int i_tally, int start_index,
       }
       break;
 
-
     case SCORE_HEATING:
       score = 0.;
       if (p->type_ == Particle::Type::neutron) {
@@ -1263,6 +1266,16 @@ score_general_ce(Particle* p, int i_tally, int start_index,
       }
       break;
 
+    case SCORE_ENERGY_DEPOSITION:
+      // Need to obtain three items: total heating, fission heating
+      // and the modified fission q recoverable
+      // Fission heating is removed because this energy is included in 
+      // the modified q recoverable data
+      score = score_neutron_heating( p, tally, flux, NEUTRON_HEATING, i_nuclide, atom_density)
+        - score_neutron_heating(p, tally, flux, FISSION_HEATING, i_nuclide, atom_density)
+        + score_fission_q(p, SCORE_ENERGY_DEPOSITION, tally, flux, i_nuclide, atom_density);
+      break;
+
     default:
       if (tally.estimator_ == ESTIMATOR_ANALOG) {
         // Any other score is assumed to be a MT number. Thus, we just need
@@ -1316,7 +1329,7 @@ score_general_ce(Particle* p, int i_tally, int start_index,
       }
     }
 
-    // Add derivative information on score for differnetial tallies.
+    // Add derivative information on score for differential tallies.
     if (tally.deriv_ != C_NONE)
       apply_derivative_to_score(p, i_tally, i_nuclide, atom_density, score_bin,
         score);