Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://mooseframework.inl.gov
       3             : //*
       4             : //* All rights reserved, see COPYRIGHT for full restrictions
       5             : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
       6             : //*
       7             : //* Licensed under LGPL 2.1, please see LICENSE for details
       8             : //* https://www.gnu.org/licenses/lgpl-2.1.html
       9             : 
      10             : #include "SobolStatistics.h"
      11             : #include "SobolSampler.h"
      12             : #include "SobolCalculators.h"
      13             : #include "BootstrapCalculators.h"
      14             : 
      15             : // MOOSE includes
      16             : #include "MooseVariable.h"
      17             : #include "ThreadedElementLoopBase.h"
      18             : #include "ThreadedNodeLoop.h"
      19             : 
      20             : #include "libmesh/quadrature.h"
      21             : 
      22             : #include <numeric>
      23             : 
      24             : registerMooseObjectDeprecated("StochasticToolsApp", SobolStatistics, "11/03/2021 12:00");
      25             : 
      26             : InputParameters
      27          68 : SobolStatistics::validParams()
      28             : {
      29          68 :   InputParameters params = GeneralVectorPostprocessor::validParams();
      30          68 :   params.addClassDescription(
      31             :       "Compute SOBOL statistics values of a given VectorPostprocessor objects and vectors.");
      32         136 :   params.addParam<SamplerName>("sampler", "SobolSampler object.");
      33         136 :   params.addParam<VectorPostprocessorName>(
      34             :       "results", "StochasticResults object containing data to use for calculation.");
      35             : 
      36          68 :   params.addParam<std::vector<Real>>(
      37             :       "ci_levels",
      38          68 :       std::vector<Real>(),
      39             :       "A vector of confidence levels to consider, values must be in (0, 1).");
      40         136 :   params.addParam<unsigned int>(
      41             :       "ci_replicates",
      42         136 :       10000,
      43             :       "The number of replicates to use when computing confidence level intervals.");
      44         136 :   params.addParam<unsigned int>("ci_seed",
      45         136 :                                 1,
      46             :                                 "The random number generator seed used for creating replicates "
      47             :                                 "while computing confidence level intervals.");
      48          68 :   return params;
      49           0 : }
      50             : 
      51          34 : SobolStatistics::SobolStatistics(const InputParameters & parameters)
      52             :   : GeneralVectorPostprocessor(parameters),
      53          34 :     _sobol_sampler(getSampler<SobolSampler>("sampler")),
      54          68 :     _ci_levels(getParam<std::vector<Real>>("ci_levels")),
      55          68 :     _ci_replicates(getParam<unsigned int>("ci_replicates")),
      56         102 :     _ci_seed(getParam<unsigned int>("ci_seed"))
      57             : {
      58          34 : }
      59             : 
      60             : void
      61          34 : SobolStatistics::initialSetup()
      62             : {
      63          34 :   const VectorPostprocessorName & vpp_name = getParam<VectorPostprocessorName>("results");
      64          34 :   const VectorPostprocessor & vpp_object = _fe_problem.getVectorPostprocessorObjectByName(vpp_name);
      65          34 :   const std::set<std::string> & vpp_vectors = vpp_object.getVectorNames();
      66          34 :   _sobol_ci_vectors.resize(_ci_levels.size());
      67          68 :   for (const auto & vec_name : vpp_vectors)
      68             :   {
      69             :     ReporterName rname(vpp_name, vec_name);
      70          34 :     _result_vectors.push_back(std::make_pair(
      71          34 :         &getReporterValueByName<VectorPostprocessorValue>(rname), vpp_object.isDistributed()));
      72          34 :     _sobol_stat_vectors.push_back(&declareVector(vpp_name + "_" + vec_name));
      73             : 
      74         187 :     for (const auto & l : index_range(_ci_levels))
      75             :     {
      76         153 :       std::stringstream vname; /// Vectors computed by this object
      77         306 :       vname << vpp_name << "_" << vec_name << "_" << _ci_levels[l] * 100 << "CI";
      78         153 :       _sobol_ci_vectors[l].push_back(&declareVector(vname.str()));
      79         153 :     }
      80             :   }
      81          34 : }
      82             : 
      83             : void
      84          34 : SobolStatistics::execute()
      85             : {
      86          68 :   TIME_SECTION("execute", 3, "Executing Sobol Statistics");
      87             : 
      88             :   StochasticTools::SobolCalculator<std::vector<Real>, Real> calc(
      89          34 :       *this, "SOBOL", _sobol_sampler.resample());
      90          34 :   auto boot_calc = _ci_levels.empty()
      91          34 :                        ? nullptr
      92          68 :                        : makeBootstrapCalculator(MooseEnum("percentile", "percentile"),
      93             :                                                  *this,
      94             :                                                  _ci_levels,
      95          17 :                                                  _ci_replicates,
      96          17 :                                                  _ci_seed,
      97          17 :                                                  calc);
      98          68 :   for (std::size_t i = 0; i < _result_vectors.size(); ++i)
      99             :   {
     100          34 :     const std::size_t ncol = _sobol_sampler.resample() ? 2 * _sobol_sampler.getNumberOfCols() + 2
     101           0 :                                                        : _sobol_sampler.getNumberOfCols() + 2;
     102             :     const std::vector<std::vector<Real>> data =
     103          34 :         StochasticTools::reshapeVector(*(_result_vectors[i].first), ncol, /*row_major =*/true);
     104          34 :     (*_sobol_stat_vectors[i]) = calc.compute(data, _result_vectors[i].second);
     105             : 
     106          34 :     if (boot_calc)
     107             :     {
     108          17 :       std::vector<std::vector<Real>> sobol_ci = boot_calc->compute(data, _result_vectors[i].second);
     109         116 :       for (const auto & l : index_range(sobol_ci))
     110          99 :         (*_sobol_ci_vectors[l][i]) = std::move(sobol_ci[l]);
     111          17 :     }
     112          34 :   }
     113          68 : }