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 "NestedMonteCarloSampler.h"
      11             : #include "Distribution.h"
      12             : 
      13             : registerMooseObjectAliased("StochasticToolsApp", NestedMonteCarloSampler, "NestedMonteCarlo");
      14             : 
      15             : InputParameters
      16         100 : NestedMonteCarloSampler::validParams()
      17             : {
      18         100 :   InputParameters params = Sampler::validParams();
      19         100 :   params.addClassDescription("Monte Carlo sampler for nested loops of parameters.");
      20         200 :   params.addRequiredParam<std::vector<dof_id_type>>(
      21             :       "num_rows",
      22             :       "The number of rows for each loop of parameters. The first number represents the outermost "
      23             :       "loop.");
      24         200 :   params.addRequiredParam<std::vector<std::vector<DistributionName>>>(
      25             :       "distributions",
      26             :       "Sets of distribution names to be sampled. Each set defines the parameters for the nested "
      27             :       "loop, with the first set being the outermost loop.");
      28         100 :   return params;
      29           0 : }
      30             : 
      31          60 : NestedMonteCarloSampler::NestedMonteCarloSampler(const InputParameters & parameters)
      32          60 :   : Sampler(parameters)
      33             : {
      34             :   // Grab inputs and make sure size is consistent
      35          60 :   const auto & dnames = getParam<std::vector<std::vector<DistributionName>>>("distributions");
      36         120 :   const auto & nrows = getParam<std::vector<dof_id_type>>("num_rows");
      37          60 :   if (dnames.size() != nrows.size())
      38           4 :     paramError("distributions",
      39             :                "There must be a set of distributions for each loop defined by 'num_rows'.");
      40             : 
      41             :   // Gather distribution pointers and fill in loop index
      42             :   const std::size_t nloop = dnames.size();
      43         224 :   for (const auto & n : make_range(nloop))
      44         504 :     for (const auto & name : dnames[n])
      45             :     {
      46         336 :       _distributions.push_back(&getDistributionByName(name));
      47         336 :       _loop_index.push_back(n);
      48             :     }
      49             : 
      50             :   // Compute what row indices need to recompute which columns
      51          56 :   _loop_mod.resize(nloop);
      52          56 :   std::partial_sum(
      53          56 :       nrows.rbegin(), nrows.rend(), _loop_mod.rbegin(), std::multiplies<dof_id_type>());
      54             :   _loop_mod.erase(_loop_mod.begin());
      55          56 :   _loop_mod.push_back(1);
      56             : 
      57             :   // Allocate row storage
      58          56 :   _row_data.resize(_distributions.size());
      59             : 
      60          56 :   setNumberOfRows(std::accumulate(nrows.begin(), nrows.end(), 1, std::multiplies<dof_id_type>()));
      61          56 :   setNumberOfCols(_distributions.size());
      62          56 : }
      63             : 
      64             : void
      65          40 : NestedMonteCarloSampler::sampleSetUp(const SampleMode mode)
      66             : {
      67          40 :   if (mode == Sampler::SampleMode::GLOBAL || getNumberOfRows() == 0)
      68           0 :     return;
      69             : 
      70             :   dof_id_type curr_row = 0;
      71          70 :   for (const auto & mod : _loop_mod)
      72             :   {
      73          70 :     if (getLocalRowBegin() % mod == 0)
      74             :       break;
      75             : 
      76          30 :     const dof_id_type target_row = std::floor(getLocalRowBegin() / mod) * mod;
      77          30 :     advanceGenerators((target_row - curr_row) * getNumberOfCols());
      78         210 :     for (const auto & j : make_range(getNumberOfCols()))
      79         180 :       computeSample(target_row, j);
      80          30 :     curr_row = target_row + 1;
      81             :   }
      82             :   restoreGeneratorState();
      83             : }
      84             : 
      85             : Real
      86       12180 : NestedMonteCarloSampler::computeSample(dof_id_type row_index, dof_id_type col_index)
      87             : {
      88       12180 :   const Real rn = getRand();
      89             :   const auto & loop = _loop_index[col_index];
      90       12180 :   if (row_index % _loop_mod[loop] == 0)
      91        4750 :     _row_data[col_index] = _distributions[col_index]->quantile(rn);
      92       12180 :   return _row_data[col_index];
      93             : }