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 "QuadraturePointsPositions.h"
      11             : #include "libmesh/quadrature_gauss.h"
      12             : #include "SetupQuadratureAction.h"
      13             : 
      14             : using namespace libMesh;
      15             : 
      16             : registerMooseObject("MooseApp", QuadraturePointsPositions);
      17             : 
      18             : InputParameters
      19       14317 : QuadraturePointsPositions::validParams()
      20             : {
      21       14317 :   InputParameters params = Positions::validParams();
      22       14317 :   params.addClassDescription("Positions of element quadrature points.");
      23       14317 :   params += BlockRestrictable::validParams();
      24             : 
      25       42951 :   params.addParam<MooseEnum>("quadrature_type",
      26       28634 :                              SetupQuadratureAction::getQuadratureTypesEnum(),
      27             :                              "Type of the quadrature rule");
      28       42951 :   params.addParam<MooseEnum>("quadrature_order",
      29       28634 :                              SetupQuadratureAction::getQuadratureOrderEnum(),
      30             :                              "Order of the volumetric quadrature. If unspecified, defaults to the "
      31             :                              "local element default order. This is not the problem default, which "
      32             :                              "is based on the order of the variables in the system.");
      33             : 
      34             :   // Element centroids could be sorted by XYZ or by id. Default to not sorting
      35       14317 :   params.set<bool>("auto_sort") = false;
      36             :   // Gathered locally, should be broadcast on every process
      37       14317 :   params.set<bool>("auto_broadcast") = true;
      38             : 
      39       14317 :   return params;
      40           0 : }
      41             : 
      42          26 : QuadraturePointsPositions::QuadraturePointsPositions(const InputParameters & parameters)
      43             :   : Positions(parameters),
      44             :     BlockRestrictable(this),
      45          52 :     _mesh(_fe_problem.mesh()),
      46          26 :     _q_type(Moose::stringToEnum<QuadratureType>(getParam<MooseEnum>("quadrature_type"))),
      47          52 :     _q_order(Moose::stringToEnum<Order>(getParam<MooseEnum>("quadrature_order")))
      48             : {
      49             :   // Mesh is ready at construction
      50          26 :   initialize();
      51             :   // Trigger synchronization as the initialization is distributed
      52          26 :   finalize();
      53          26 : }
      54             : 
      55             : void
      56          26 : QuadraturePointsPositions::initialize()
      57             : {
      58          26 :   clearPositions();
      59             : 
      60             :   // By default, initialize should be called on meshChanged()
      61             :   // Gathering of positions is local, reporter system makes sure to make it global
      62          26 :   if (blockRestricted())
      63             :   {
      64          13 :     _positions_2d.resize(numBlocks());
      65          13 :     unsigned int b_index = 0;
      66          26 :     for (const auto & sub_id : blockIDs())
      67             :     {
      68          33 :       for (const auto & elem : _mesh.getMesh().active_local_subdomain_elements_ptr_range(sub_id))
      69             :       {
      70             :         // Get a quadrature going of the requested type and order
      71          10 :         const FEFamily mapping_family = FEMap::map_fe_type(*elem);
      72          10 :         const FEType fe_type(elem->default_order(), mapping_family);
      73          10 :         std::unique_ptr<FEBase> fe = FEBase::build(elem->dim(), fe_type);
      74             :         const auto q_order =
      75          10 :             (_q_order == INVALID_ORDER) ? fe_type.default_quadrature_order() : _q_order;
      76          10 :         auto qrule = QBase::build(_q_type, elem->dim(), q_order);
      77          10 :         fe->attach_quadrature_rule(qrule.get());
      78          10 :         const auto & q_points = fe->get_xyz();
      79          10 :         fe->reinit(elem);
      80             : 
      81          90 :         for (const auto & q : q_points)
      82             :         {
      83          80 :           _positions.emplace_back(q);
      84          80 :           _positions_2d[b_index].emplace_back(q);
      85             :         }
      86          23 :       }
      87          13 :       b_index += 1;
      88             :     }
      89             :   }
      90             :   else
      91             :   {
      92          13 :     _positions.resize(_mesh.getMesh().n_active_local_elem());
      93          53 :     for (const auto & elem : _mesh.getMesh().active_local_element_ptr_range())
      94             :     {
      95             :       // Get a quadrature going on the element
      96          20 :       const FEFamily mapping_family = FEMap::map_fe_type(*elem);
      97          20 :       const FEType fe_type(elem->default_order(), mapping_family);
      98          20 :       std::unique_ptr<FEBase> fe = FEBase::build(elem->dim(), fe_type);
      99             :       const auto q_order =
     100          20 :           (_q_order == INVALID_ORDER) ? fe_type.default_quadrature_order() : _q_order;
     101          20 :       auto qrule = QBase::build(_q_type, elem->dim(), q_order);
     102          20 :       fe->attach_quadrature_rule(qrule.get());
     103          20 :       const auto & q_points = fe->get_xyz();
     104          20 :       fe->reinit(elem);
     105             : 
     106         180 :       for (const auto & q : q_points)
     107         160 :         _positions.emplace_back(q);
     108          33 :     }
     109             :   }
     110          26 :   _initialized = true;
     111          26 : }