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 "ExtraIDIntegralVectorPostprocessor.h"
      11             : #include "MooseVariable.h"
      12             : #include "MooseMesh.h"
      13             : #include "MooseMeshUtils.h"
      14             : 
      15             : #include "libmesh/mesh_base.h"
      16             : 
      17             : registerMooseObject("MooseApp", ExtraIDIntegralVectorPostprocessor);
      18             : 
      19             : InputParameters
      20        3289 : ExtraIDIntegralVectorPostprocessor::validParams()
      21             : {
      22        3289 :   InputParameters params = ElementVariableVectorPostprocessor::validParams();
      23       13156 :   params.addParam<std::vector<MaterialPropertyName>>(
      24             :       "mat_prop", "The names of material properties that this VectorPostprocessor operates on");
      25       13156 :   params.addRequiredParam<std::vector<ExtraElementIDName>>(
      26             :       "id_name", "List of extra element ID names by which to separate integral(s).");
      27       13156 :   params.addParam<bool>("average", false, "Whether or not to compute volume average");
      28       13156 :   params.addParam<std::string>("spatial_value_name",
      29             :                                "To specify which variable or material property is to be used when "
      30             :                                "using as a spatial user object functor");
      31        6578 :   params.addClassDescription("Integrates or averages variables based on extra element IDs");
      32        3289 :   params.makeParamNotRequired("variable");
      33        3289 :   return params;
      34           0 : }
      35             : 
      36         120 : ExtraIDIntegralVectorPostprocessor::ExtraIDIntegralVectorPostprocessor(
      37         120 :     const InputParameters & parameters)
      38             :   : SpatialUserObjectFunctor<ElementVariableVectorPostprocessor>(parameters),
      39         120 :     _average(getParam<bool>("average")),
      40         480 :     _nvar(isParamValid("variable") ? coupledComponents("variable") : 0),
      41         284 :     _nprop(isParamValid("mat_prop") ? getParam<std::vector<MaterialPropertyName>>("mat_prop").size()
      42             :                                     : 0),
      43         284 :     _prop_names(isParamValid("mat_prop") ? getParam<std::vector<MaterialPropertyName>>("mat_prop")
      44             :                                          : std::vector<MaterialPropertyName>()),
      45         240 :     _extra_id(getParam<std::vector<ExtraElementIDName>>("id_name")),
      46         120 :     _n_extra_id(_extra_id.size()),
      47         240 :     _spatial_evaluation_index(std::numeric_limits<unsigned int>::max())
      48             : {
      49         120 :   if (!_nvar && !_nprop)
      50           0 :     mooseError("Neither 'variable' nor 'mat_prop' was specified.");
      51             : 
      52             :   // create map of element ids to parsed vpp ids
      53             :   _unique_vpp_ids =
      54         120 :       MooseMeshUtils::getExtraIDUniqueCombinationMap(_mesh.getMesh(), blockIDs(), _extra_id);
      55             : 
      56             :   // set up variable vector containing parsed extra ids
      57         306 :   for (unsigned int i = 0; i < _n_extra_id; ++i)
      58             :   {
      59         186 :     const auto id_type = _extra_id[i];
      60         186 :     auto & p = declareVector("Level-" + std::to_string(i) + "-" + id_type);
      61             :     // collect local map of local vpp id to extra id
      62         186 :     std::map<dof_id_type, dof_id_type> extra_ids;
      63      688796 :     for (const auto & elem : _mesh.getMesh().active_local_element_ptr_range())
      64             :     {
      65      688610 :       if (!hasBlocks(elem->subdomain_id()))
      66          40 :         continue;
      67      688570 :       auto vpp_id = _unique_vpp_ids[elem->id()];
      68      688570 :       if (extra_ids.find(vpp_id) == extra_ids.end())
      69        6467 :         extra_ids[vpp_id] = elem->get_extra_integer(getElementIDIndexByName(id_type));
      70         186 :     }
      71             :     // create global map of local vpp id to extra id
      72         186 :     comm().set_union(extra_ids);
      73         186 :     p.resize(extra_ids.size());
      74        8786 :     for (auto it : extra_ids)
      75        8600 :       p[it.first] = it.second;
      76         186 :     _extra_ids.push_back(&p);
      77         186 :   }
      78             : 
      79         120 :   std::vector<std::string> vector_names;
      80             : 
      81             :   // declare vectors containing variable integral values
      82         284 :   for (unsigned int i = 0; i < _nvar; ++i)
      83             :   {
      84         328 :     _vars.push_back(getVar("variable", i));
      85         328 :     _var_values.push_back(&coupledValue("variable", i));
      86         164 :     vector_names.push_back(_vars[i]->name());
      87         164 :     auto & p = declareVector(vector_names.back());
      88         164 :     p.resize((*_extra_ids[0]).size());
      89         164 :     _integrals.push_back(&p);
      90             :   }
      91             : 
      92             :   // declare vectors containing material property integral values
      93         164 :   for (auto & name : _prop_names)
      94             :   {
      95          44 :     _props.push_back(&getMaterialPropertyByName<Real>(name));
      96          44 :     vector_names.push_back(name);
      97          44 :     auto & p = declareVector(vector_names.back());
      98          44 :     p.resize((*_extra_ids[0]).size());
      99          44 :     _integrals.push_back(&p);
     100             :   }
     101             : 
     102         360 :   if (isParamValid("spatial_value_name"))
     103             :   {
     104          34 :     const auto & name = getParam<std::string>("spatial_value_name");
     105          17 :     _spatial_evaluation_index = vector_names.size();
     106          17 :     for (const auto i : index_range(vector_names))
     107          17 :       if (name == vector_names[i])
     108             :       {
     109          17 :         _spatial_evaluation_index = i;
     110          17 :         break;
     111             :       }
     112          17 :     if (_spatial_evaluation_index == vector_names.size())
     113           0 :       paramError("spatial_value_name",
     114             :                  "not a variable name or material property name of this vector postprocessor");
     115             :   }
     116         120 : }
     117             : 
     118             : void
     119         119 : ExtraIDIntegralVectorPostprocessor::initialize()
     120             : {
     121         326 :   for (auto & integral : _integrals)
     122         207 :     std::fill(integral->begin(), integral->end(), 0);
     123             : 
     124         119 :   if (_average)
     125          22 :     _volumes.assign((*_extra_ids[0]).size(), 0);
     126         119 : }
     127             : 
     128             : void
     129      344412 : ExtraIDIntegralVectorPostprocessor::execute()
     130             : {
     131      344412 :   if (hasBlocks(_current_elem->subdomain_id()))
     132             :   {
     133      344412 :     unsigned int i = 0;
     134      344412 :     auto ipos = _unique_vpp_ids[_current_elem->id()];
     135      860856 :     for (unsigned int ivar = 0; ivar < _nvar; ++ivar, ++i)
     136      516444 :       if (_vars[ivar]->hasBlocks(_current_elem->subdomain_id()))
     137     1033848 :         for (unsigned int qp = 0; qp < _qrule->n_points(); qp++)
     138      517404 :           (*_integrals[i])[ipos] += _JxW[qp] * _coord[qp] * (*_var_values[ivar])[qp];
     139             : 
     140      516444 :     for (unsigned int iprop = 0; iprop < _nprop; ++iprop, ++i)
     141      344064 :       for (unsigned int qp = 0; qp < _qrule->n_points(); qp++)
     142      172032 :         (*_integrals[i])[ipos] += _JxW[qp] * _coord[qp] * (*_props[iprop])[qp];
     143             : 
     144      344412 :     if (_average)
     145       43036 :       _volumes[ipos] += _current_elem_volume;
     146             :   }
     147      344412 : }
     148             : 
     149             : void
     150         107 : ExtraIDIntegralVectorPostprocessor::finalize()
     151             : {
     152         294 :   for (auto & integral : _integrals)
     153         187 :     gatherSum(*integral);
     154             : 
     155         107 :   if (_average)
     156             :   {
     157          20 :     gatherSum(_volumes);
     158             : 
     159          70 :     for (auto & integral : _integrals)
     160        2630 :       for (unsigned int i = 0; i < integral->size(); ++i)
     161        2580 :         (*integral)[i] /= _volumes[i];
     162             :   }
     163         107 : }
     164             : 
     165             : void
     166          12 : ExtraIDIntegralVectorPostprocessor::threadJoin(const UserObject & s)
     167             : {
     168          12 :   const auto & sibling = static_cast<const ExtraIDIntegralVectorPostprocessor &>(s);
     169             : 
     170          32 :   for (unsigned int i = 0; i < _integrals.size(); ++i)
     171         932 :     for (size_t j = 0; j < (*_integrals[i]).size(); ++j)
     172         912 :       (*_integrals[i])[j] += (*sibling._integrals[i])[j];
     173             : 
     174          12 :   if (_average)
     175          68 :     for (unsigned int i = 0; i < _volumes.size(); ++i)
     176          66 :       _volumes[i] += sibling._volumes[i];
     177          12 : }
     178             : 
     179             : Real
     180           0 : ExtraIDIntegralVectorPostprocessor::evaluate(const ElemArg & elem,
     181             :                                              const Moose::StateArg & libmesh_dbg_var(state)) const
     182             : {
     183             :   mooseAssert(state.state == 0, "We do not currently support evaluating at old states");
     184           0 :   return elementValue(elem.elem);
     185             : }
     186             : 
     187             : Real
     188         859 : ExtraIDIntegralVectorPostprocessor::evaluate(const ElemQpArg & elem,
     189             :                                              const Moose::StateArg & libmesh_dbg_var(state)) const
     190             : {
     191             :   mooseAssert(state.state == 0, "We do not currently support evaluating at old states");
     192         859 :   return elementValue(elem.elem);
     193             : }
     194             : 
     195             : Real
     196         859 : ExtraIDIntegralVectorPostprocessor::elementValue(const Elem * elem) const
     197             : {
     198         859 :   if (_spatial_evaluation_index == std::numeric_limits<unsigned int>::max())
     199           8 :     paramError("spatial_value_name",
     200             :                "Must set when ExtraIDIntegralVectorPostprocessor is used as a functor");
     201         855 :   const auto it = _unique_vpp_ids.find(elem->id());
     202         855 :   if (it == _unique_vpp_ids.end())
     203             :   {
     204           3 :     if (!hasBlocks(elem->subdomain_id()))
     205           6 :       mooseError(
     206             :           "Failed evaluating spatial value of element ",
     207           3 :           elem->id(),
     208             :           " in subdomain ",
     209           3 :           elem->subdomain_id(),
     210             :           ". This is likely due to the object using this "
     211             :           "ExtraIDIntegralVectorPostprocessor as a functor operates on a subdomain outside of "
     212             :           "this ExtraIDIntegralVectorPostprocessor");
     213             :     else
     214           0 :       mooseError("Failed evaluating spatial value of element ",
     215           0 :                  elem->id(),
     216             :                  ". We should not hit this error. Please contact code developers for help");
     217             :   }
     218        1704 :   return (*_integrals[_spatial_evaluation_index])[it->second];
     219             : }