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 "IndependentGaussianMH.h"
      11             : #include "Normal.h"
      12             : #include "TruncatedNormal.h"
      13             : 
      14             : registerMooseObject("StochasticToolsApp", IndependentGaussianMH);
      15             : 
      16             : InputParameters
      17          46 : IndependentGaussianMH::validParams()
      18             : {
      19          46 :   InputParameters params = PMCMCBase::validParams();
      20          46 :   params.addClassDescription("Perform M-H MCMC sampling with independent Gaussian propoposals.");
      21          92 :   params.addRequiredParam<ReporterName>("seed_inputs",
      22             :                                         "Reporter with seed inputs values for the next proposals.");
      23          92 :   params.addRequiredParam<std::vector<Real>>("std_prop",
      24             :                                              "Standard deviations for making the next proposal.");
      25          46 :   return params;
      26           0 : }
      27             : 
      28          26 : IndependentGaussianMH::IndependentGaussianMH(const InputParameters & parameters)
      29             :   : PMCMCBase(parameters),
      30          26 :     _seed_inputs(getReporterValue<std::vector<Real>>("seed_inputs")),
      31          78 :     _std_prop(getParam<std::vector<Real>>("std_prop"))
      32             : {
      33             :   // Error check for sizes of proposal stds
      34          26 :   if (_std_prop.size() != _priors.size())
      35           2 :     paramError("std_prop",
      36             :                "The number of proposal stds, initial values, and priors should be the same.");
      37          24 : }
      38             : 
      39             : void
      40         144 : IndependentGaussianMH::proposeSamples()
      41             : {
      42         144 :   std::vector<Real> old_sample = (_t_step + 1 > decisionStep()) ? _seed_inputs : _initial_values;
      43         864 :   for (unsigned int j = 0; j < _num_parallel_proposals; ++j)
      44        2160 :     for (unsigned int i = 0; i < _priors.size(); ++i)
      45             :     {
      46        1440 :       if (_lower_bound)
      47           0 :         _new_samples[j][i] = TruncatedNormal::quantile(
      48           0 :             random(), old_sample[i], _std_prop[i], (*_lower_bound)[i], (*_upper_bound)[i]);
      49             :       else
      50        1440 :         _new_samples[j][i] = Normal::quantile(random(), old_sample[i], _std_prop[i]);
      51             :     }
      52         144 : }