www.mooseframework.org
MaxQpsThread.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "MaxQpsThread.h"
16 #include "FEProblem.h"
17 
18 #include "libmesh/fe_base.h"
19 #include "libmesh/threads.h"
20 #include LIBMESH_INCLUDE_UNORDERED_SET
21 LIBMESH_DEFINE_HASH_POINTERS
22 #include "libmesh/quadrature.h"
23 
25  QuadratureType qtype,
26  Order order,
27  Order face_order)
28  : _fe_problem(fe_problem),
29  _qtype(qtype),
30  _order(order),
31  _face_order(face_order),
32  _max(0),
33  _max_shape_funcs(0)
34 {
35 }
36 
37 // Splitting Constructor
38 MaxQpsThread::MaxQpsThread(MaxQpsThread & x, Threads::split /*split*/)
40  _qtype(x._qtype),
41  _order(x._order),
43  _max(x._max),
45 {
46 }
47 
48 void
49 MaxQpsThread::operator()(const ConstElemRange & range)
50 {
51  ParallelUniqueId puid;
52  _tid = puid.id;
53 
54  // For short circuiting reinit
55  std::set<ElemType> seen_it;
56  for (const auto & elem : range)
57  {
58  // Only reinit if the element type has not previously been seen
59  if (seen_it.insert(elem->type()).second)
60  {
61  FEType fe_type(FIRST, LAGRANGE);
62  unsigned int dim = elem->dim();
63  unsigned int side = 0; // we assume that any element will have at least one side ;)
64 
65  // We cannot mess with the FE objects in Assembly, because we might need to request second
66  // derivatives
67  // later on. If we used them, we'd call reinit on them, thus making the call to request second
68  // derivatives harmful (i.e. leading to segfaults/asserts). Thus, we have to use a locally
69  // allocated object here.
70  std::unique_ptr<FEBase> fe(FEBase::build(dim, fe_type));
71 
72  // figure out the number of qps for volume
73  std::unique_ptr<QBase> qrule(QBase::build(_qtype, dim, _order));
74  fe->attach_quadrature_rule(qrule.get());
75  fe->reinit(elem);
76  if (qrule->n_points() > _max)
77  _max = qrule->n_points();
78 
79  unsigned int n_shape_funcs = fe->n_shape_functions();
80  if (n_shape_funcs > _max_shape_funcs)
81  _max_shape_funcs = n_shape_funcs;
82 
83  // figure out the number of qps for the face
84  // NOTE: user might specify higher order rule for faces, thus possibly ending up with more qps
85  // than in the volume
86  std::unique_ptr<QBase> qrule_face(QBase::build(_qtype, dim - 1, _face_order));
87  fe->attach_quadrature_rule(qrule_face.get());
88  fe->reinit(elem, side);
89  if (qrule_face->n_points() > _max)
90  _max = qrule_face->n_points();
91  }
92  }
93 }
94 
95 void
97 {
98  if (y._max > _max)
99  _max = y._max;
100 
103 }
This class determines the maximum number of Quadrature Points and Shape Functions used for a given si...
Definition: MaxQpsThread.h:33
unsigned int _max
Maximum number of qps encountered.
Definition: MaxQpsThread.h:59
unsigned int _max_shape_funcs
Maximum number of shape functions encountered.
Definition: MaxQpsThread.h:62
void join(const MaxQpsThread &y)
Definition: MaxQpsThread.C:96
void operator()(const ConstElemRange &range)
Definition: MaxQpsThread.C:49
Order _face_order
Definition: MaxQpsThread.h:54
static PetscErrorCode Vec x
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
FEProblemBase & _fe_problem
Definition: MaxQpsThread.h:50
QuadratureType _qtype
Definition: MaxQpsThread.h:52
MaxQpsThread(FEProblemBase &fe_problem, QuadratureType type, Order order, Order face_order)
Definition: MaxQpsThread.C:24
THREAD_ID _tid
Definition: MaxQpsThread.h:56