Line data    Source code

       1             : /********************************************************************/
       2             : /*                  SOFTWARE COPYRIGHT NOTIFICATION                 */
       3             : /*                             Cardinal                             */
       4             : /*                                                                  */
       5             : /*                  (c) 2021 UChicago Argonne, LLC                  */
       6             : /*                        ALL RIGHTS RESERVED                       */
       7             : /*                                                                  */
       8             : /*                 Prepared by UChicago Argonne, LLC                */
       9             : /*               Under Contract No. DE-AC02-06CH11357               */
      10             : /*                With the U. S. Department of Energy               */
      11             : /*                                                                  */
      12             : /*             Prepared by Battelle Energy Alliance, LLC            */
      13             : /*               Under Contract No. DE-AC07-05ID14517               */
      14             : /*                With the U. S. Department of Energy               */
      15             : /*                                                                  */
      16             : /*                 See LICENSE for full restrictions                */
      17             : /********************************************************************/
      18             : 
      19             : #ifdef ENABLE_OPENMC_COUPLING
      20             : 
      21             : #include "ElementOpticalDepthIndicator.h"
      22             : 
      23             : #include "CardinalEnums.h"
      24             : #include "TallyBase.h"
      25             : 
      26             : registerMooseObject("CardinalApp", ElementOpticalDepthIndicator);
      27             : 
      28             : InputParameters
      29          58 : ElementOpticalDepthIndicator::validParams()
      30             : {
      31          58 :   auto params = OpenMCIndicator::validParams();
      32          58 :   params.addClassDescription(
      33             :       "A class which returns the estimate of a given element's optical depth under the assumption "
      34             :       "that "
      35             :       "elements with a large optical depth experience large solution gradients.");
      36         116 :   params.addRequiredParam<MooseEnum>("rxn_rate",
      37         116 :                                      getSingleTallyScoreEnum(),
      38             :                                      "The reaction rate to use for computing the optical depth.");
      39         116 :   params.addParam<MooseEnum>(
      40             :       "h_type",
      41         174 :       MooseEnum("min max cube_root", "max"),
      42             :       "The estimate for the length of the element used to compute the optical depth. Options are "
      43             :       "the "
      44             :       "minimum vertex separation (min), the maximum vertex separation (max), and the cube root of "
      45             :       "the element volume (cube_root).");
      46         116 :   params.addParam<bool>(
      47             :       "invert",
      48         116 :       false,
      49             :       "Whether the optical depth is computed as the optical depth (false) or the inverse of the "
      50             :       "optical depth (true).");
      51             : 
      52          58 :   return params;
      53           0 : }
      54             : 
      55          26 : ElementOpticalDepthIndicator::ElementOpticalDepthIndicator(const InputParameters & parameters)
      56             :   : OpenMCIndicator(parameters),
      57          26 :     _h_type(getParam<MooseEnum>("h_type").getEnum<HType>()),
      58          78 :     _invert(getParam<bool>("invert"))
      59             : {
      60          52 :   std::string score = getParam<MooseEnum>("rxn_rate");
      61             :   std::replace(score.begin(), score.end(), '_', '-');
      62             : 
      63             :   // Error check to make sure the score is a reaction rate score and to make sure one of the
      64             :   // [Tallies] has added the score and a flux score.
      65          26 :   if (!_openmc_problem->isReactionRateScore(score))
      66           2 :     paramError(
      67             :         "rxn_rate",
      68           2 :         "At present the ElementOpticalDepthIndicator only works with reaction rate scores. " +
      69           2 :             std::string(getParam<MooseEnum>("rxn_rate")) + " is not a valid reaction rate score.");
      70             : 
      71          24 :   if (!_openmc_problem->hasScore(score))
      72           2 :     paramError("rxn_rate",
      73           2 :                "The problem does not contain any score named " +
      74           2 :                    std::string(getParam<MooseEnum>("rxn_rate")) +
      75             :                    "! Please "
      76             :                    "ensure that one of your [Tallies] is scoring the requested reaction rate.");
      77             : 
      78          22 :   if (!_openmc_problem->hasScore("flux"))
      79           2 :     mooseError("In order to use an ElementOpticalDepthIndicator one of your [Tallies] must add a "
      80             :                "flux score.");
      81             : 
      82             :   // Check to ensure the reaction rate / flux variables are CONSTANT MONOMIALS.
      83             :   bool const_mon = true;
      84          60 :   for (const auto v : _openmc_problem->getTallyScoreVariables(score, _tid, "", true))
      85             :     const_mon &= v->feType() == FEType(CONSTANT, MONOMIAL);
      86          60 :   for (const auto v : _openmc_problem->getTallyScoreVariables("flux", _tid, "", true))
      87             :     const_mon &= v->feType() == FEType(CONSTANT, MONOMIAL);
      88             : 
      89          20 :   if (!const_mon)
      90           0 :     paramError("rxn_rate",
      91             :                "ElementOpticalDepthIndicator only supports CONSTANT MONOMIAL field variables. "
      92             :                "Please ensure your [Tallies] are adding CONSTANT MONOMIAL field variables.");
      93             : 
      94             :   // Grab the reaction rate / flux variables from the [Tallies].
      95          20 :   _rxn_rates = _openmc_problem->getTallyScoreVariableValues(score, _tid, "", true);
      96          40 :   _scalar_fluxes = _openmc_problem->getTallyScoreVariableValues("flux", _tid, "", true);
      97          20 : }
      98             : 
      99             : void
     100      281984 : ElementOpticalDepthIndicator::computeIndicator()
     101             : {
     102             :   Real rxn_rate = 0.0;
     103             :   Real scalar_flux = 0.0;
     104             : 
     105      563968 :   for (const auto & var : _rxn_rates)
     106      281984 :     rxn_rate += (*var)[0];
     107             : 
     108      563968 :   for (const auto & var : _scalar_fluxes)
     109      281984 :     scalar_flux += (*var)[0];
     110             : 
     111      281984 :   auto od = scalar_flux < libMesh::TOLERANCE ? 0.0 : rxn_rate / scalar_flux;
     112             : 
     113      281984 :   switch (_h_type)
     114             :   {
     115        1024 :     case HType::Min:
     116        1024 :       od /= _current_elem->hmin();
     117        1024 :       break;
     118      279424 :     case HType::Max:
     119      279424 :       od /= _current_elem->hmax();
     120      279424 :       break;
     121        1536 :     case HType::CubeRoot:
     122        1536 :       od /= std::cbrt(_current_elem->volume());
     123        1536 :       break;
     124           0 :     default:
     125           0 :       mooseError("Internal error: unhandled HType enum state in ElementOpticalDepthIndicator.");
     126             :       break;
     127             :   }
     128             : 
     129      281984 :   if (_invert && od > libMesh::TOLERANCE)
     130        1536 :     _field_var.setNodalValue(1.0 / od);
     131             :   else
     132      280448 :     _field_var.setNodalValue(od);
     133      281984 : }
     134             : 
     135             : #endif