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 "MaxQpsThread.h"
      11             : #include "FEProblem.h"
      12             : #include "Assembly.h"
      13             : 
      14             : #include "libmesh/fe.h"
      15             : #include "libmesh/threads.h"
      16             : #include LIBMESH_INCLUDE_UNORDERED_SET
      17             : LIBMESH_DEFINE_HASH_POINTERS
      18             : #include "libmesh/quadrature.h"
      19             : 
      20       57708 : MaxQpsThread::MaxQpsThread(FEProblemBase & fe_problem) : _fe_problem(fe_problem), _max(0) {}
      21             : 
      22             : // Splitting Constructor
      23        4977 : MaxQpsThread::MaxQpsThread(MaxQpsThread & x, Threads::split /*split*/)
      24        4977 :   : _fe_problem(x._fe_problem), _max(x._max)
      25             : {
      26        4977 : }
      27             : 
      28             : void
      29       62685 : MaxQpsThread::operator()(const libMesh::ConstElemRange & range)
      30             : {
      31       62685 :   ParallelUniqueId puid;
      32       62685 :   _tid = puid.id;
      33             : 
      34             :   // Not actually using any pre-existing data so it shouldn't matter which assembly we use
      35       62685 :   auto & assem = _fe_problem.assembly(_tid, 0);
      36             : 
      37             :   // For short circuiting reinit.  With potential block-specific qrules we
      38             :   // need to track "seen" element types by their subdomains as well.
      39       62685 :   std::set<std::pair<ElemType, SubdomainID>> seen_it;
      40             : 
      41    10387113 :   for (const auto & elem : range)
      42             :   {
      43             :     // Only reinit if the element type has not previously been seen
      44    10324428 :     if (!seen_it.insert(std::make_pair(elem->type(), elem->subdomain_id())).second)
      45    10249929 :       continue;
      46             : 
      47             :     // This ensures we can access the correct qrules if any block-specific
      48             :     // qrules have been created.
      49       74499 :     assem.setCurrentSubdomainID(elem->subdomain_id());
      50       74499 :     assem.setVolumeQRule(elem);
      51             : 
      52       74499 :     auto & qrule = assem.writeableQRule();
      53       74499 :     qrule->init(*elem);
      54       74499 :     if (qrule->n_points() > _max)
      55       61308 :       _max = qrule->n_points();
      56             : 
      57             :     // We used to check side quadrature rules too - badly (assuming
      58             :     // side 0 doesn't have a smaller rule than the others, which is
      59             :     // often untrue on prisms) and generally unnecessarily (a side
      60             :     // rule will always use fewer points than a corresponding interior
      61             :     // rule).
      62             :     //
      63             :     // We should handle the possibility of users manually specifying a
      64             :     // higher order for side than for interior integration.  Doing
      65             :     // this efficiently will need to wait for a new libMesh API,
      66             :     // though.
      67             : 
      68             :     // In initial conditions nodes are enumerated as pretend quadrature points
      69             :     // using the _qp index to access coupled variables. In order to be able to
      70             :     // use _zero (resized according to _max_qps) with _qp, we need to count nodes.
      71       74499 :     if (elem->n_nodes() > _max)
      72       17073 :       _max = elem->n_nodes();
      73             :   }
      74             : 
      75             :   // Clear the cached quadrature rules because we may add FE objects in between now and simulation
      76             :   // start and we only ensure we set all the FE quadrature rules if a quadrature rule is different
      77             :   // from the cached quadrature rule
      78       62685 :   assem.clearCachedQRules();
      79       62685 : }
      80             : 
      81             : void
      82        4977 : MaxQpsThread::join(const MaxQpsThread & y)
      83             : {
      84        4977 :   if (y._max > _max)
      85          13 :     _max = y._max;
      86        4977 : }