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 "KEigenvalue.h"
      22             : #include "openmc/eigenvalue.h"
      23             : #include "openmc/math_functions.h"
      24             : #include "openmc/constants.h"
      25             : 
      26             : registerMooseObject("CardinalApp", KEigenvalue);
      27             : 
      28             : InputParameters
      29        1549 : KEigenvalue::validParams()
      30             : {
      31        1549 :   InputParameters params = GeneralPostprocessor::validParams();
      32        1549 :   params += OpenMCBase::validParams();
      33        1549 :   params.addClassDescription("k eigenvalue computed by OpenMC");
      34        3098 :   params.addParam<MooseEnum>("value_type",
      35        3098 :                              getEigenvalueEnum(),
      36             :                              "Type of eigenvalue global tally to report");
      37        3098 :   params.addParam<MooseEnum>(
      38             :       "output",
      39        3098 :       getStatsOutputEnum(),
      40             :       "The value to output. Options are $k_{eff}$ (mean), the standard deviation "
      41             :       "of $k_{eff}$ (std_dev), or the relative error of $k_{eff}$ (rel_err).");
      42        1549 :   return params;
      43           0 : }
      44             : 
      45         507 : KEigenvalue::KEigenvalue(const InputParameters & parameters)
      46             :   : GeneralPostprocessor(parameters),
      47             :     OpenMCBase(this, parameters),
      48         507 :     _type(getParam<MooseEnum>("value_type").getEnum<eigenvalue::EigenvalueEnum>()),
      49        1521 :     _output(getParam<MooseEnum>("output").getEnum<statistics::OutputEnum>())
      50             : {
      51         507 :   if (openmc::settings::run_mode != openmc::RunMode::EIGENVALUE)
      52           4 :     mooseError("Eigenvalues are only computed when running OpenMC in eigenvalue mode!");
      53         503 : }
      54             : 
      55             : Real
      56         538 : KEigenvalue::getValue() const
      57             : {
      58         538 :   switch (_output)
      59             :   {
      60         428 :     case statistics::OutputEnum::Mean:
      61         428 :       return kMean();
      62             : 
      63          78 :     case statistics::OutputEnum::StDev:
      64          78 :       return KStandardDeviation();
      65             : 
      66          32 :     case statistics::OutputEnum::RelError:
      67          32 :       return kRelativeError();
      68             : 
      69           0 :     default:
      70           0 :       mooseError("Internal error: Unhandled statistics::OutputEnum enum in KEigenvalue.");
      71             :       break;
      72             :   }
      73             : }
      74             : 
      75             : Real
      76         564 : KEigenvalue::kMean() const
      77             : {
      78         564 :   int n = openmc::simulation::n_realizations;
      79             :   const auto & gt = openmc::simulation::global_tallies;
      80             : 
      81         564 :   switch (_type)
      82             :   {
      83             :     case eigenvalue::collision:
      84          48 :       return gt(openmc::GlobalTally::K_COLLISION, openmc::TallyResult::SUM) / n;
      85             :     case eigenvalue::absorption:
      86          40 :       return gt(openmc::GlobalTally::K_ABSORPTION, openmc::TallyResult::SUM) / n;
      87             :     case eigenvalue::tracklength:
      88          40 :       return gt(openmc::GlobalTally::K_TRACKLENGTH, openmc::TallyResult::SUM) / n;
      89         436 :     case eigenvalue::combined:
      90             :     {
      91         436 :       if (n <= 3)
      92           0 :         mooseError("Cannot compute combined k-effective estimate with fewer than 4 realizations!\n"
      93             :           "Please change the 'value_type' to either 'collision', 'tracklength', or 'absorption'.");
      94             : 
      95             :       double k_eff[2];
      96         436 :       openmc::openmc_get_keff(k_eff);
      97         436 :       return k_eff[0];
      98             :     }
      99           0 :     default:
     100           0 :       mooseError("Internal error: Unhandled EigenvalueEnum in KEigenvalue!");
     101             :   }
     102             : }
     103             : 
     104             : Real
     105         134 : KEigenvalue::KStandardDeviation() const
     106             : {
     107             :   const auto & gt = openmc::simulation::global_tallies;
     108         134 :   int n = openmc::simulation::n_realizations;
     109             : 
     110             :   double t_n1 = 1.0;
     111         134 :   if (openmc::settings::confidence_intervals)
     112             :   {
     113             :     double alpha = 1.0 - openmc::CONFIDENCE_LEVEL;
     114           0 :     t_n1 = openmc::t_percentile(1.0 - alpha / 2.0, n - 1);
     115             :   }
     116             : 
     117         134 :   switch (_type)
     118             :   {
     119             :     case eigenvalue::collision:
     120             :     {
     121          16 :       double mean = gt(openmc::GlobalTally::K_COLLISION, openmc::TallyResult::SUM) / n;
     122          16 :       double sum_sq = gt(openmc::GlobalTally::K_COLLISION, openmc::TallyResult::SUM_SQ);
     123          16 :       return t_n1 * stdev(mean, sum_sq, openmc::simulation::n_realizations);
     124             :     }
     125             :     case eigenvalue::absorption:
     126             :     {
     127          16 :       double mean = gt(openmc::GlobalTally::K_ABSORPTION, openmc::TallyResult::SUM) / n;
     128          16 :       double sum_sq = gt(openmc::GlobalTally::K_ABSORPTION, openmc::TallyResult::SUM_SQ);
     129          16 :       return t_n1 * stdev(mean, sum_sq, openmc::simulation::n_realizations);
     130             :     }
     131             :     case eigenvalue::tracklength:
     132             :     {
     133          16 :       double mean = gt(openmc::GlobalTally::K_TRACKLENGTH, openmc::TallyResult::SUM) / n;
     134          16 :       double sum_sq = gt(openmc::GlobalTally::K_TRACKLENGTH, openmc::TallyResult::SUM_SQ);
     135          16 :       return t_n1 * stdev(mean, sum_sq, openmc::simulation::n_realizations);
     136             :     }
     137          86 :     case eigenvalue::combined:
     138             :     {
     139          86 :       if (n <= 3)
     140           0 :         mooseError("Cannot compute combined k-effective standard deviation with fewer than 4 "
     141             :                    "realizations!");
     142             : 
     143             :       double k_eff[2];
     144          86 :       openmc::openmc_get_keff(k_eff);
     145          86 :       return k_eff[1];
     146             :     }
     147           0 :     default:
     148           0 :       mooseError("Internal error: Unhandled StandardDeviationEnum in KEigenvalue!");
     149             :   }
     150             : }
     151             : 
     152             : Real
     153          56 : KEigenvalue::kRelativeError() const
     154             : {
     155          56 :   const auto mean = kMean();
     156          56 :   return mean > 0.0 ? KStandardDeviation() / kMean() : 0.0;
     157             : }
     158             : 
     159             : #endif