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 "SpatialAverageBase.h"
      11             : 
      12             : // MOOSE includes
      13             : #include "MooseVariableFE.h"
      14             : 
      15             : #include "libmesh/quadrature.h"
      16             : 
      17             : InputParameters
      18       28580 : SpatialAverageBase::validParams()
      19             : {
      20       28580 :   InputParameters params = ElementVectorPostprocessor::validParams();
      21       28580 :   params.addParam<unsigned int>("bin_number", 50, "Number of histogram bins");
      22       28580 :   params.addCoupledVar("variable", "Variables to average radially");
      23       28580 :   params.addRequiredParam<Real>("radius", "Radius to average out to");
      24       28580 :   params.addParam<Point>("origin", Point(), "Origin of the cylinder");
      25       85740 :   params.addParam<Real>(
      26       57160 :       "empty_bin_value", 0.0, "Value to assign to bins into which no datapoints fall");
      27       28580 :   return params;
      28           0 : }
      29             : 
      30          26 : SpatialAverageBase::SpatialAverageBase(const InputParameters & parameters)
      31             :   : ElementVectorPostprocessor(parameters),
      32          26 :     _nbins(getParam<unsigned int>("bin_number")),
      33          26 :     _radius(getParam<Real>("radius")),
      34          26 :     _origin(getParam<Point>("origin")),
      35          26 :     _deltaR(_radius / _nbins),
      36          26 :     _nvals(coupledComponents("variable")),
      37          26 :     _values(coupledValues("variable")),
      38          26 :     _empty_bin_value(getParam<Real>("empty_bin_value")),
      39          26 :     _bin_center(declareVector("radius")),
      40          26 :     _counts(_nbins),
      41          52 :     _average(_nvals)
      42             : {
      43          26 :   if (coupledComponents("variable") != 1)
      44           0 :     mooseError("SpatialAverageBase works on exactly one coupled variable");
      45             : 
      46             :   // Note: We associate the local variable "i" with nbins and "j" with nvals throughout.
      47             : 
      48             :   // couple variables initialize vectors
      49          52 :   for (MooseIndex(_average) j = 0; j < _nvals; ++j)
      50          26 :     _average[j] = &declareVector(coupledName("variable", j));
      51             : 
      52             :   // initialize the bin center value vector
      53          26 :   _bin_center.resize(_nbins);
      54         286 :   for (MooseIndex(_counts) i = 0; i < _nbins; ++i)
      55         260 :     _bin_center[i] = (i + 0.5) * _deltaR;
      56          26 : }
      57             : 
      58             : void
      59          48 : SpatialAverageBase::initialize()
      60             : {
      61             :   // reset the histogram
      62          96 :   for (auto vec_ptr : _average)
      63          48 :     vec_ptr->assign(_nbins, 0.0);
      64             : 
      65             :   // reset bin counts
      66          48 :   _counts.assign(_nbins, 0.0);
      67          48 : }
      68             : 
      69             : void
      70       17600 : SpatialAverageBase::execute()
      71             : {
      72             :   // loop over quadrature points
      73      158400 :   for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
      74             :   {
      75             :     // compute target bin
      76      140800 :     auto bin = computeDistance() / _deltaR;
      77             : 
      78             :     // add the volume contributed by the current quadrature point
      79      140800 :     if (bin >= 0 && bin < static_cast<int>(_nbins))
      80             :     {
      81      152320 :       for (MooseIndex(_nvals) j = 0; j < _nvals; ++j)
      82       76160 :         (*_average[j])[bin] += (*_values[j])[_qp];
      83             : 
      84       76160 :       _counts[bin]++;
      85             :     }
      86             :   }
      87       17600 : }
      88             : 
      89             : void
      90          44 : SpatialAverageBase::finalize()
      91             : {
      92          44 :   gatherSum(_counts);
      93             : 
      94          88 :   for (MooseIndex(_average) j = 0; j < _nvals; ++j)
      95             :   {
      96          44 :     gatherSum(*_average[j]);
      97             : 
      98         484 :     for (MooseIndex(_counts) i = 0; i < _nbins; ++i)
      99         440 :       (*_average[j])[i] =
     100         440 :           _counts[i] > 0 ? (*_average[j])[i] / static_cast<Real>(_counts[i]) : _empty_bin_value;
     101             :   }
     102          44 : }
     103             : 
     104             : void
     105           4 : SpatialAverageBase::threadJoin(const UserObject & y)
     106             : {
     107           4 :   const auto & uo = static_cast<const SpatialAverageBase &>(y);
     108             : 
     109          44 :   for (MooseIndex(_counts) i = 0; i < _nbins; ++i)
     110             :   {
     111          40 :     _counts[i] += uo._counts[i];
     112             : 
     113          80 :     for (MooseIndex(_average) j = 0; j < _nvals; ++j)
     114          40 :       (*_average[j])[i] += (*uo._average[j])[i];
     115             :   }
     116           4 : }