libMesh
fe_abstract.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 // libmesh includes
19 #include "libmesh/fe.h"
20 #include "libmesh/libmesh_logging.h"
21 #include "libmesh/enum_elem_type.h"
22 #include "libmesh/boundary_info.h"
23 #include "libmesh/mesh_base.h"
24 #include "libmesh/dense_matrix.h"
25 #include "libmesh/dense_vector.h"
26 #include "libmesh/dof_map.h"
27 #include "libmesh/elem.h"
28 #include "libmesh/fe_interface.h"
29 #include "libmesh/numeric_vector.h"
30 #include "libmesh/periodic_boundaries.h"
31 #include "libmesh/periodic_boundary.h"
32 #include "libmesh/quadrature.h"
33 #include "libmesh/quadrature_gauss.h"
34 #include "libmesh/remote_elem.h"
35 #include "libmesh/tensor_value.h"
36 #include "libmesh/threads.h"
37 #include "libmesh/enum_to_string.h"
38 
39 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
40 #include "libmesh/inf_fe.h"
41 #include "libmesh/fe_interface_macros.h"
42 #endif
43 
44 namespace libMesh
45 {
46 
47 FEAbstract::FEAbstract(const unsigned int d,
48  const FEType & fet) :
49  _fe_map( FEMap::build(fet) ),
50  dim(d),
51  calculations_started(false),
52  calculate_dual(false),
53  calculate_default_dual_coeff(true),
54  calculate_nothing(false),
55  calculate_map(false),
56  calculate_phi(false),
57  calculate_dphi(false),
58  calculate_d2phi(false),
59  calculate_curl_phi(false),
60  calculate_div_phi(false),
61  calculate_dphiref(false),
62  fe_type(fet),
63  _elem_type(INVALID_ELEM),
64  _elem(nullptr),
65  _elem_p_level(0),
66  _p_level(0),
67  qrule(nullptr),
68  shapes_on_quadrature(false),
69  _n_total_qp(0),
70  _add_p_level_in_reinit(true)
71 {
72 }
73 
74 
75 FEAbstract::~FEAbstract() = default;
76 
77 
78 std::unique_ptr<FEAbstract> FEAbstract::build(const unsigned int dim,
79  const FEType & fet)
80 {
81  switch (dim)
82  {
83  // 0D
84  case 0:
85  {
86  switch (fet.family)
87  {
88  case CLOUGH:
89  return std::make_unique<FE<0,CLOUGH>>(fet);
90 
91  case HERMITE:
92  return std::make_unique<FE<0,HERMITE>>(fet);
93 
94  case LAGRANGE:
95  return std::make_unique<FE<0,LAGRANGE>>(fet);
96 
97  case LAGRANGE_VEC:
98  return std::make_unique<FE<0,LAGRANGE_VEC>>(fet);
99 
100  case L2_LAGRANGE:
101  return std::make_unique<FE<0,L2_LAGRANGE>>(fet);
102 
103  case L2_LAGRANGE_VEC:
104  return std::make_unique<FE<0,L2_LAGRANGE_VEC>>(fet);
105 
106  case HIERARCHIC_VEC:
107  return std::make_unique<FE<0,HIERARCHIC_VEC>>(fet);
108 
109  case HIERARCHIC:
110  return std::make_unique<FE<0,HIERARCHIC>>(fet);
111 
112  case L2_HIERARCHIC:
113  return std::make_unique<FE<0,L2_HIERARCHIC>>(fet);
114 
115  case L2_HIERARCHIC_VEC:
116  return std::make_unique<FE<0,L2_HIERARCHIC_VEC>>(fet);
117 
118  case SIDE_HIERARCHIC:
119  return std::make_unique<FE<0,SIDE_HIERARCHIC>>(fet);
120 
121  case MONOMIAL:
122  return std::make_unique<FE<0,MONOMIAL>>(fet);
123 
124  case MONOMIAL_VEC:
125  return std::make_unique<FE<0,MONOMIAL_VEC>>(fet);
126 
127 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
128  case SZABAB:
129  return std::make_unique<FE<0,SZABAB>>(fet);
130 
131  case BERNSTEIN:
132  return std::make_unique<FE<0,BERNSTEIN>>(fet);
133 
134  case RATIONAL_BERNSTEIN:
135  return std::make_unique<FE<0,RATIONAL_BERNSTEIN>>(fet);
136 #endif
137 
138  case XYZ:
139  return std::make_unique<FEXYZ<0>>(fet);
140 
141  case SCALAR:
142  return std::make_unique<FEScalar<0>>(fet);
143 
144  case NEDELEC_ONE:
145  return std::make_unique<FENedelecOne<0>>(fet);
146 
147  case RAVIART_THOMAS:
148  return std::make_unique<FERaviartThomas<0>>(fet);
149 
150  case L2_RAVIART_THOMAS:
151  return std::make_unique<FEL2RaviartThomas<0>>(fet);
152 
153  case SUBDIVISION:
154  return std::make_unique<FESubdivision>(fet);
155 
156  default:
157  libmesh_error_msg("ERROR: Bad FEType.family= " << Utility::enum_to_string(fet.family));
158  }
159  }
160  // 1D
161  case 1:
162  {
163  switch (fet.family)
164  {
165  case CLOUGH:
166  return std::make_unique<FE<1,CLOUGH>>(fet);
167 
168  case HERMITE:
169  return std::make_unique<FE<1,HERMITE>>(fet);
170 
171  case LAGRANGE:
172  return std::make_unique<FE<1,LAGRANGE>>(fet);
173 
174  case LAGRANGE_VEC:
175  return std::make_unique<FE<1,LAGRANGE_VEC>>(fet);
176 
177  case L2_LAGRANGE:
178  return std::make_unique<FE<1,L2_LAGRANGE>>(fet);
179 
180  case L2_LAGRANGE_VEC:
181  return std::make_unique<FE<1,L2_LAGRANGE_VEC>>(fet);
182 
183  case HIERARCHIC_VEC:
184  return std::make_unique<FE<1,HIERARCHIC_VEC>>(fet);
185 
186  case HIERARCHIC:
187  return std::make_unique<FE<1,HIERARCHIC>>(fet);
188 
189  case L2_HIERARCHIC:
190  return std::make_unique<FE<1,L2_HIERARCHIC>>(fet);
191 
192  case L2_HIERARCHIC_VEC:
193  return std::make_unique<FE<1,L2_HIERARCHIC_VEC>>(fet);
194 
195  case SIDE_HIERARCHIC:
196  return std::make_unique<FE<1,SIDE_HIERARCHIC>>(fet);
197 
198  case MONOMIAL:
199  return std::make_unique<FE<1,MONOMIAL>>(fet);
200 
201  case MONOMIAL_VEC:
202  return std::make_unique<FE<1,MONOMIAL_VEC>>(fet);
203 
204 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
205  case SZABAB:
206  return std::make_unique<FE<1,SZABAB>>(fet);
207 
208  case BERNSTEIN:
209  return std::make_unique<FE<1,BERNSTEIN>>(fet);
210 
211  case RATIONAL_BERNSTEIN:
212  return std::make_unique<FE<1,RATIONAL_BERNSTEIN>>(fet);
213 #endif
214 
215  case XYZ:
216  return std::make_unique<FEXYZ<1>>(fet);
217 
218  case SCALAR:
219  return std::make_unique<FEScalar<1>>(fet);
220 
221  case NEDELEC_ONE:
222  return std::make_unique<FENedelecOne<1>>(fet);
223 
224  case RAVIART_THOMAS:
225  return std::make_unique<FERaviartThomas<1>>(fet);
226 
227  case L2_RAVIART_THOMAS:
228  return std::make_unique<FEL2RaviartThomas<1>>(fet);
229 
230  case SUBDIVISION:
231  return std::make_unique<FESubdivision>(fet);
232 
233  default:
234  libmesh_error_msg("ERROR: Bad FEType.family= " << Utility::enum_to_string(fet.family));
235  }
236  }
237 
238 
239  // 2D
240  case 2:
241  {
242  switch (fet.family)
243  {
244  case CLOUGH:
245  return std::make_unique<FE<2,CLOUGH>>(fet);
246 
247  case HERMITE:
248  return std::make_unique<FE<2,HERMITE>>(fet);
249 
250  case LAGRANGE:
251  return std::make_unique<FE<2,LAGRANGE>>(fet);
252 
253  case LAGRANGE_VEC:
254  return std::make_unique<FE<2,LAGRANGE_VEC>>(fet);
255 
256  case L2_LAGRANGE:
257  return std::make_unique<FE<2,L2_LAGRANGE>>(fet);
258 
259  case L2_LAGRANGE_VEC:
260  return std::make_unique<FE<2,L2_LAGRANGE_VEC>>(fet);
261 
262  case HIERARCHIC_VEC:
263  return std::make_unique<FE<2,HIERARCHIC_VEC>>(fet);
264 
265  case HIERARCHIC:
266  return std::make_unique<FE<2,HIERARCHIC>>(fet);
267 
268  case L2_HIERARCHIC:
269  return std::make_unique<FE<2,L2_HIERARCHIC>>(fet);
270 
271  case L2_HIERARCHIC_VEC:
272  return std::make_unique<FE<2,L2_HIERARCHIC_VEC>>(fet);
273 
274  case SIDE_HIERARCHIC:
275  return std::make_unique<FE<2,SIDE_HIERARCHIC>>(fet);
276 
277  case MONOMIAL:
278  return std::make_unique<FE<2,MONOMIAL>>(fet);
279 
280  case MONOMIAL_VEC:
281  return std::make_unique<FE<2,MONOMIAL_VEC>>(fet);
282 
283 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
284  case SZABAB:
285  return std::make_unique<FE<2,SZABAB>>(fet);
286 
287  case BERNSTEIN:
288  return std::make_unique<FE<2,BERNSTEIN>>(fet);
289 
290  case RATIONAL_BERNSTEIN:
291  return std::make_unique<FE<2,RATIONAL_BERNSTEIN>>(fet);
292 #endif
293 
294  case XYZ:
295  return std::make_unique<FEXYZ<2>>(fet);
296 
297  case SCALAR:
298  return std::make_unique<FEScalar<2>>(fet);
299 
300  case NEDELEC_ONE:
301  return std::make_unique<FENedelecOne<2>>(fet);
302 
303  case RAVIART_THOMAS:
304  return std::make_unique<FERaviartThomas<2>>(fet);
305 
306  case L2_RAVIART_THOMAS:
307  return std::make_unique<FEL2RaviartThomas<2>>(fet);
308 
309  case SUBDIVISION:
310  return std::make_unique<FESubdivision>(fet);
311 
312  default:
313  libmesh_error_msg("ERROR: Bad FEType.family= " << Utility::enum_to_string(fet.family));
314  }
315  }
316 
317 
318  // 3D
319  case 3:
320  {
321  switch (fet.family)
322  {
323  case CLOUGH:
324  libmesh_error_msg("ERROR: Clough-Tocher elements currently only support 1D and 2D");
325 
326  case HERMITE:
327  return std::make_unique<FE<3,HERMITE>>(fet);
328 
329  case LAGRANGE:
330  return std::make_unique<FE<3,LAGRANGE>>(fet);
331 
332  case LAGRANGE_VEC:
333  return std::make_unique<FE<3,LAGRANGE_VEC>>(fet);
334 
335  case L2_LAGRANGE:
336  return std::make_unique<FE<3,L2_LAGRANGE>>(fet);
337 
338  case L2_LAGRANGE_VEC:
339  return std::make_unique<FE<3,L2_LAGRANGE_VEC>>(fet);
340 
341  case HIERARCHIC_VEC:
342  return std::make_unique<FE<3,HIERARCHIC_VEC>>(fet);
343 
344  case HIERARCHIC:
345  return std::make_unique<FE<3,HIERARCHIC>>(fet);
346 
347  case L2_HIERARCHIC:
348  return std::make_unique<FE<3,L2_HIERARCHIC>>(fet);
349 
350  case L2_HIERARCHIC_VEC:
351  return std::make_unique<FE<3,L2_HIERARCHIC_VEC>>(fet);
352 
353  case SIDE_HIERARCHIC:
354  return std::make_unique<FE<3,SIDE_HIERARCHIC>>(fet);
355 
356  case MONOMIAL:
357  return std::make_unique<FE<3,MONOMIAL>>(fet);
358 
359  case MONOMIAL_VEC:
360  return std::make_unique<FE<3,MONOMIAL_VEC>>(fet);
361 
362 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
363  case SZABAB:
364  return std::make_unique<FE<3,SZABAB>>(fet);
365 
366  case BERNSTEIN:
367  return std::make_unique<FE<3,BERNSTEIN>>(fet);
368 
369  case RATIONAL_BERNSTEIN:
370  return std::make_unique<FE<3,RATIONAL_BERNSTEIN>>(fet);
371 #endif
372 
373  case XYZ:
374  return std::make_unique<FEXYZ<3>>(fet);
375 
376  case SCALAR:
377  return std::make_unique<FEScalar<3>>(fet);
378 
379  case NEDELEC_ONE:
380  return std::make_unique<FENedelecOne<3>>(fet);
381 
382  case RAVIART_THOMAS:
383  return std::make_unique<FERaviartThomas<3>>(fet);
384 
385  case L2_RAVIART_THOMAS:
386  return std::make_unique<FEL2RaviartThomas<3>>(fet);
387 
388  default:
389  libmesh_error_msg("ERROR: Bad FEType.family= " << Utility::enum_to_string(fet.family));
390  }
391  }
392 
393  default:
394  libmesh_error_msg("Invalid dimension dim = " << dim);
395  }
396 }
397 
398 
399 
400 void FEAbstract::get_refspace_nodes(const ElemType itemType, std::vector<Point> & nodes)
401 {
402  const unsigned int n_nodes = Elem::type_to_n_nodes_map[itemType];
403  if (n_nodes == invalid_uint)
404  libmesh_error_msg("Number of nodes is not well-defined for " <<
405  Utility::enum_to_string(itemType));
406 
407  nodes.resize(n_nodes);
408  switch(itemType)
409  {
410  case NODEELEM:
411  {
412  nodes[0] = Point (0.,0.,0.);
413  return;
414  }
415  case EDGE3:
416  {
417  nodes[2] = Point (0.,0.,0.);
418  libmesh_fallthrough();
419  }
420  case EDGE2:
421  {
422  nodes[0] = Point (-1.,0.,0.);
423  nodes[1] = Point (1.,0.,0.);
424  return;
425  }
426  case EDGE4: // not nested with EDGE3
427  {
428  nodes[0] = Point (-1.,0.,0.);
429  nodes[1] = Point (1.,0.,0.);
430  nodes[2] = Point (-1./3.,0.,0.);
431  nodes[3] - Point (1./3.,0.,0.);
432  return;
433  }
434  case TRI7:
435  {
436  nodes[6] = Point (1./3.,1./3.,0.);
437  libmesh_fallthrough();
438  }
439  case TRI6:
440  {
441  nodes[3] = Point (.5,0.,0.);
442  nodes[4] = Point (.5,.5,0.);
443  nodes[5] = Point (0.,.5,0.);
444  libmesh_fallthrough();
445  }
446  case TRI3:
447  case TRISHELL3:
448  {
449  nodes[0] = Point (0.,0.,0.);
450  nodes[1] = Point (1.,0.,0.);
451  nodes[2] = Point (0.,1.,0.);
452  return;
453  }
454  case QUAD9:
455  case QUADSHELL9:
456  {
457  nodes[8] = Point (0.,0.,0.);
458  libmesh_fallthrough();
459  }
460  case QUAD8:
461  case QUADSHELL8:
462  {
463  nodes[4] = Point (0.,-1.,0.);
464  nodes[5] = Point (1.,0.,0.);
465  nodes[6] = Point (0.,1.,0.);
466  nodes[7] = Point (-1.,0.,0.);
467  libmesh_fallthrough();
468  }
469  case QUAD4:
470  case QUADSHELL4:
471  {
472  nodes[0] = Point (-1.,-1.,0.);
473  nodes[1] = Point (1.,-1.,0.);
474  nodes[2] = Point (1.,1.,0.);
475  nodes[3] = Point (-1.,1.,0.);
476  return;
477  }
478  case TET14:
479  {
480  nodes[10] = Point (1/Real(3),1/Real(3),0.);
481  nodes[11] = Point (1/Real(3),0.,1/Real(3));
482  nodes[12] = Point (1/Real(3),1/Real(3),1/Real(3));
483  nodes[13] = Point (0.,1/Real(3),1/Real(3));
484  libmesh_fallthrough();
485  }
486  case TET10:
487  {
488  nodes[4] = Point (.5,0.,0.);
489  nodes[5] = Point (.5,.5,0.);
490  nodes[6] = Point (0.,.5,0.);
491  nodes[7] = Point (0.,0.,.5);
492  nodes[8] = Point (.5,0.,.5);
493  nodes[9] = Point (0.,.5,.5);
494  libmesh_fallthrough();
495  }
496  case TET4:
497  {
498  nodes[0] = Point (0.,0.,0.);
499  nodes[1] = Point (1.,0.,0.);
500  nodes[2] = Point (0.,1.,0.);
501  nodes[3] = Point (0.,0.,1.);
502  return;
503  }
504  case HEX27:
505  {
506  nodes[20] = Point (0.,0.,-1.);
507  nodes[21] = Point (0.,-1.,0.);
508  nodes[22] = Point (1.,0.,0.);
509  nodes[23] = Point (0.,1.,0.);
510  nodes[24] = Point (-1.,0.,0.);
511  nodes[25] = Point (0.,0.,1.);
512  nodes[26] = Point (0.,0.,0.);
513  libmesh_fallthrough();
514  }
515  case HEX20:
516  {
517  nodes[8] = Point (0.,-1.,-1.);
518  nodes[9] = Point (1.,0.,-1.);
519  nodes[10] = Point (0.,1.,-1.);
520  nodes[11] = Point (-1.,0.,-1.);
521  nodes[12] = Point (-1.,-1.,0.);
522  nodes[13] = Point (1.,-1.,0.);
523  nodes[14] = Point (1.,1.,0.);
524  nodes[15] = Point (-1.,1.,0.);
525  nodes[16] = Point (0.,-1.,1.);
526  nodes[17] = Point (1.,0.,1.);
527  nodes[18] = Point (0.,1.,1.);
528  nodes[19] = Point (-1.,0.,1.);
529  libmesh_fallthrough();
530  }
531  case HEX8:
532  {
533  nodes[0] = Point (-1.,-1.,-1.);
534  nodes[1] = Point (1.,-1.,-1.);
535  nodes[2] = Point (1.,1.,-1.);
536  nodes[3] = Point (-1.,1.,-1.);
537  nodes[4] = Point (-1.,-1.,1.);
538  nodes[5] = Point (1.,-1.,1.);
539  nodes[6] = Point (1.,1.,1.);
540  nodes[7] = Point (-1.,1.,1.);
541  return;
542  }
543  case PRISM21:
544  {
545  nodes[20] = Point (1/Real(3),1/Real(3),0);
546  libmesh_fallthrough();
547  }
548  case PRISM20:
549  {
550  nodes[18] = Point (1/Real(3),1/Real(3),-1);
551  nodes[19] = Point (1/Real(3),1/Real(3),1);
552  libmesh_fallthrough();
553  }
554  case PRISM18:
555  {
556  nodes[15] = Point (.5,0.,0.);
557  nodes[16] = Point (.5,.5,0.);
558  nodes[17] = Point (0.,.5,0.);
559  libmesh_fallthrough();
560  }
561  case PRISM15:
562  {
563  nodes[6] = Point (.5,0.,-1.);
564  nodes[7] = Point (.5,.5,-1.);
565  nodes[8] = Point (0.,.5,-1.);
566  nodes[9] = Point (0.,0.,0.);
567  nodes[10] = Point (1.,0.,0.);
568  nodes[11] = Point (0.,1.,0.);
569  nodes[12] = Point (.5,0.,1.);
570  nodes[13] = Point (.5,.5,1.);
571  nodes[14] = Point (0.,.5,1.);
572  libmesh_fallthrough();
573  }
574  case PRISM6:
575  {
576  nodes[0] = Point (0.,0.,-1.);
577  nodes[1] = Point (1.,0.,-1.);
578  nodes[2] = Point (0.,1.,-1.);
579  nodes[3] = Point (0.,0.,1.);
580  nodes[4] = Point (1.,0.,1.);
581  nodes[5] = Point (0.,1.,1.);
582  return;
583  }
584  case PYRAMID18:
585  {
586  // triangle centers
587  nodes[14] = Point (-2/Real(3),0.,1/Real(3));
588  nodes[15] = Point (0.,2/Real(3),1/Real(3));
589  nodes[16] = Point (2/Real(3),0.,1/Real(3));
590  nodes[17] = Point (0.,-2/Real(3),1/Real(3));
591 
592  libmesh_fallthrough();
593  }
594  case PYRAMID14:
595  {
596  // base center
597  nodes[13] = Point (0.,0.,0.);
598 
599  libmesh_fallthrough();
600  }
601  case PYRAMID13:
602  {
603  // base midedge
604  nodes[5] = Point (0.,-1.,0.);
605  nodes[6] = Point (1.,0.,0.);
606  nodes[7] = Point (0.,1.,0.);
607  nodes[8] = Point (-1,0.,0.);
608 
609  // lateral midedge
610  nodes[9] = Point (-.5,-.5,.5);
611  nodes[10] = Point (.5,-.5,.5);
612  nodes[11] = Point (.5,.5,.5);
613  nodes[12] = Point (-.5,.5,.5);
614 
615  libmesh_fallthrough();
616  }
617  case PYRAMID5:
618  {
619  // base corners
620  nodes[0] = Point (-1.,-1.,0.);
621  nodes[1] = Point (1.,-1.,0.);
622  nodes[2] = Point (1.,1.,0.);
623  nodes[3] = Point (-1.,1.,0.);
624  // apex
625  nodes[4] = Point (0.,0.,1.);
626  return;
627  }
628 
629  default:
630  libmesh_error_msg("ERROR: Unknown element type " << Utility::enum_to_string(itemType));
631  }
632 }
633 
634 
635 
636 #ifdef LIBMESH_ENABLE_DEPRECATED
637 bool FEAbstract::on_reference_element(const Point & p, const ElemType t, const Real eps)
638 {
639  // Use Elem::on_reference_element() instead
640  libmesh_deprecated();
641 
642  libmesh_assert_greater_equal (eps, 0.);
643 
644  const Real xi = p(0);
645 #if LIBMESH_DIM > 1
646  const Real eta = p(1);
647 #else
648  const Real eta = 0.;
649 #endif
650 #if LIBMESH_DIM > 2
651  const Real zeta = p(2);
652 #else
653  const Real zeta = 0.;
654 #endif
655 
656  switch (t)
657  {
658  case NODEELEM:
659  {
660  return (!xi && !eta && !zeta);
661  }
662  case EDGE2:
663  case EDGE3:
664  case EDGE4:
665  {
666  // The reference 1D element is [-1,1].
667  if ((xi >= -1.-eps) &&
668  (xi <= 1.+eps))
669  return true;
670 
671  return false;
672  }
673 
674 
675  case TRI3:
676  case TRISHELL3:
677  case TRI6:
678  case TRI7:
679  {
680  // The reference triangle is isosceles
681  // and is bound by xi=0, eta=0, and xi+eta=1.
682  if ((xi >= 0.-eps) &&
683  (eta >= 0.-eps) &&
684  ((xi + eta) <= 1.+eps))
685  return true;
686 
687  return false;
688  }
689 
690 
691  case QUAD4:
692  case QUADSHELL4:
693  case QUAD8:
694  case QUADSHELL8:
695  case QUAD9:
696  case QUADSHELL9:
697  {
698  // The reference quadrilateral element is [-1,1]^2.
699  if ((xi >= -1.-eps) &&
700  (xi <= 1.+eps) &&
701  (eta >= -1.-eps) &&
702  (eta <= 1.+eps))
703  return true;
704 
705  return false;
706  }
707 
708 
709  case TET4:
710  case TET10:
711  case TET14:
712  {
713  // The reference tetrahedral is isosceles
714  // and is bound by xi=0, eta=0, zeta=0,
715  // and xi+eta+zeta=1.
716  if ((xi >= 0.-eps) &&
717  (eta >= 0.-eps) &&
718  (zeta >= 0.-eps) &&
719  ((xi + eta + zeta) <= 1.+eps))
720  return true;
721 
722  return false;
723  }
724 
725 
726  case HEX8:
727  case HEX20:
728  case HEX27:
729  {
730  /*
731  if ((xi >= -1.) &&
732  (xi <= 1.) &&
733  (eta >= -1.) &&
734  (eta <= 1.) &&
735  (zeta >= -1.) &&
736  (zeta <= 1.))
737  return true;
738  */
739 
740  // The reference hexahedral element is [-1,1]^3.
741  if ((xi >= -1.-eps) &&
742  (xi <= 1.+eps) &&
743  (eta >= -1.-eps) &&
744  (eta <= 1.+eps) &&
745  (zeta >= -1.-eps) &&
746  (zeta <= 1.+eps))
747  {
748  // libMesh::out << "Strange Point:\n";
749  // p.print();
750  return true;
751  }
752 
753  return false;
754  }
755 
756  case PRISM6:
757  case PRISM15:
758  case PRISM18:
759  case PRISM20:
760  case PRISM21:
761  {
762  // Figure this one out...
763  // inside the reference triangle with zeta in [-1,1]
764  if ((xi >= 0.-eps) &&
765  (eta >= 0.-eps) &&
766  (zeta >= -1.-eps) &&
767  (zeta <= 1.+eps) &&
768  ((xi + eta) <= 1.+eps))
769  return true;
770 
771  return false;
772  }
773 
774 
775  case PYRAMID5:
776  case PYRAMID13:
777  case PYRAMID14:
778  case PYRAMID18:
779  {
780  // Check that the point is on the same side of all the faces
781  // by testing whether:
782  //
783  // n_i.(x - x_i) <= 0
784  //
785  // for each i, where:
786  // n_i is the outward normal of face i,
787  // x_i is a point on face i.
788  if ((-eta - 1. + zeta <= 0.+eps) &&
789  ( xi - 1. + zeta <= 0.+eps) &&
790  ( eta - 1. + zeta <= 0.+eps) &&
791  ( -xi - 1. + zeta <= 0.+eps) &&
792  ( zeta >= 0.-eps))
793  return true;
794 
795  return false;
796  }
797 
798 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
799  case INFHEX8:
800  case INFHEX16:
801  case INFHEX18:
802  {
803  // The reference infhex8 is a [-1,1]^3.
804  if ((xi >= -1.-eps) &&
805  (xi <= 1.+eps) &&
806  (eta >= -1.-eps) &&
807  (eta <= 1.+eps) &&
808  (zeta >= -1.-eps) &&
809  (zeta <= 1.+eps))
810  {
811  return true;
812  }
813  return false;
814  }
815 
816  case INFPRISM6:
817  case INFPRISM12:
818  {
819  // inside the reference triangle with zeta in [-1,1]
820  if ((xi >= 0.-eps) &&
821  (eta >= 0.-eps) &&
822  (zeta >= -1.-eps) &&
823  (zeta <= 1.+eps) &&
824  ((xi + eta) <= 1.+eps))
825  {
826  return true;
827  }
828 
829  return false;
830  }
831 #endif
832 
833  default:
834  libmesh_error_msg("ERROR: Unknown element type " << Utility::enum_to_string(t));
835  }
836 
837  // If we get here then the point is _not_ in the
838  // reference element. Better return false.
839 
840  return false;
841 }
842 #endif // LIBMESH_ENABLE_DEPRECATED
843 
844 
845 
846 void FEAbstract::print_JxW(std::ostream & os) const
847 {
848  this->_fe_map->print_JxW(os);
849 }
850 
851 
852 
853 void FEAbstract::print_xyz(std::ostream & os) const
854 {
855  this->_fe_map->print_xyz(os);
856 }
857 
858 
859 void FEAbstract::print_info(std::ostream & os) const
860 {
861  os << "phi[i][j]: Shape function i at quadrature pt. j" << std::endl;
862  this->print_phi(os);
863 
864  os << "dphi[i][j]: Shape function i's gradient at quadrature pt. j" << std::endl;
865  this->print_dphi(os);
866 
867  os << "XYZ locations of the quadrature pts." << std::endl;
868  this->print_xyz(os);
869 
870  os << "Values of JxW at the quadrature pts." << std::endl;
871  this->print_JxW(os);
872 }
873 
874 
875 std::ostream & operator << (std::ostream & os, const FEAbstract & fe)
876 {
877  fe.print_info(os);
878  return os;
879 }
880 
881 
882 
883 #ifdef LIBMESH_ENABLE_AMR
884 
885 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
887  const Elem * elem)
888 {
889  libmesh_assert(elem);
890 
891  const unsigned int Dim = elem->dim();
892 
893  // Only constrain elements in 2,3D.
894  if (Dim == 1)
895  return;
896 
897  // Only constrain active and ancestor elements
898  if (elem->subactive())
899  return;
900 
901 
902 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
903  if (elem->infinite())
904  {
905  const FEType fe_t(elem->default_order(), FEMap::map_fe_type(*elem));
906 
907  // expand the infinite_compute_constraint in its template-arguments.
908  switch(Dim)
909  {
910  case 2:
911  {
912  inf_fe_family_mapping_switch(2, inf_compute_node_constraints (constraints, elem) , ,; break;);
913  break;
914  }
915  case 3:
916  {
917  inf_fe_family_mapping_switch(3, inf_compute_node_constraints (constraints, elem) , ,; break;);
918  break;
919  }
920  default:
921  libmesh_error_msg("Invalid dim = " << Dim);
922  }
923  return;
924  }
925 
926 #endif
927  const FEFamily mapping_family = FEMap::map_fe_type(*elem);
928  const FEType fe_type(elem->default_side_order(), mapping_family);
929 
930  // Pull objects out of the loop to reduce heap operations
931  std::vector<const Node *> my_nodes, parent_nodes;
932  std::unique_ptr<const Elem> my_side, parent_side;
933 
934  // Look at the element faces. Check to see if we need to
935  // build constraints.
936  for (auto s : elem->side_index_range())
937  if (elem->neighbor_ptr(s) != nullptr &&
938  elem->neighbor_ptr(s) != remote_elem)
939  if (elem->neighbor_ptr(s)->level() < elem->level()) // constrain dofs shared between
940  { // this element and ones coarser
941  // than this element.
942  // Get pointers to the elements of interest and its parent.
943  const Elem * parent = elem->parent();
944 
945  // This can't happen... Only level-0 elements have nullptr
946  // parents, and no level-0 elements can be at a higher
947  // level than their neighbors!
948  libmesh_assert(parent);
949 
950  elem->build_side_ptr(my_side, s);
951  parent->build_side_ptr(parent_side, s);
952 
953  const unsigned int n_side_nodes = my_side->n_nodes();
954 
955  my_nodes.clear();
956  my_nodes.reserve (n_side_nodes);
957  parent_nodes.clear();
958  parent_nodes.reserve (n_side_nodes);
959 
960  for (unsigned int n=0; n != n_side_nodes; ++n)
961  my_nodes.push_back(my_side->node_ptr(n));
962 
963  for (unsigned int n=0; n != n_side_nodes; ++n)
964  parent_nodes.push_back(parent_side->node_ptr(n));
965 
966  for (unsigned int my_side_n=0;
967  my_side_n < n_side_nodes;
968  my_side_n++)
969  {
970  // We can have an FE type that supports an order
971  // partially, such that sides do not support the same
972  // order. E.g. we say that a LAGRANGE PRISM21 supports
973  // "third" order to distinguish its shape functions from
974  // a PRISM18, but the QUAD9 sides will still only
975  // support second order.
976  FEType side_fe_type = fe_type;
977  const int side_max_order =
978  FEInterface::max_order(fe_type, my_side->type());
979 
980  if ((int)fe_type.order > side_max_order)
981  side_fe_type.order = side_max_order;
982 
983  // Do not use the p_level(), if any, that is inherited by the side.
984  libmesh_assert_less
985  (my_side_n,
986  FEInterface::n_dofs(side_fe_type, /*extra_order=*/0,
987  my_side.get()));
988 
989  const Node * my_node = my_nodes[my_side_n];
990 
991  // The support point of the DOF
992  const Point & support_point = *my_node;
993 
994  // Figure out where my node lies on their reference element.
995  const Point mapped_point = FEMap::inverse_map(Dim-1,
996  parent_side.get(),
997  support_point);
998 
999  // Compute the parent's side shape function values.
1000  for (unsigned int their_side_n=0;
1001  their_side_n < n_side_nodes;
1002  their_side_n++)
1003  {
1004  // Do not use the p_level(), if any, that is inherited by the side.
1005  libmesh_assert_less
1006  (their_side_n,
1007  FEInterface::n_dofs(side_fe_type,
1008  /*extra_order=*/0,
1009  parent_side.get()));
1010 
1011  const Node * their_node = parent_nodes[their_side_n];
1012  libmesh_assert(their_node);
1013 
1014  // Do not use the p_level(), if any, that is inherited by the side.
1015  const Real their_value = FEInterface::shape(side_fe_type,
1016  /*extra_order=*/0,
1017  parent_side.get(),
1018  their_side_n,
1019  mapped_point);
1020 
1021  const Real their_mag = std::abs(their_value);
1022 #ifdef DEBUG
1023  // Protect for the case u_i ~= u_j,
1024  // in which case i better equal j.
1025  if (their_mag > 0.999)
1026  {
1027  libmesh_assert_equal_to (my_node, their_node);
1028  libmesh_assert_less (std::abs(their_value - 1.), 0.001);
1029  }
1030  else
1031 #endif
1032  // To make nodal constraints useful for constructing
1033  // sparsity patterns faster, we need to get EVERY
1034  // POSSIBLE constraint coupling identified, even if
1035  // there is no coupling in the isoparametric
1036  // Lagrange case.
1037  if (their_mag < 1.e-5)
1038  {
1039  // since we may be running this method concurrently
1040  // on multiple threads we need to acquire a lock
1041  // before modifying the shared constraint_row object.
1042  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1043 
1044  // A reference to the constraint row.
1045  NodeConstraintRow & constraint_row = constraints[my_node].first;
1046 
1047  constraint_row.emplace(their_node, 0.);
1048  }
1049  // To get nodal coordinate constraints right, only
1050  // add non-zero and non-identity values for Lagrange
1051  // basis functions.
1052  else // (1.e-5 <= their_mag <= .999)
1053  {
1054  // since we may be running this method concurrently
1055  // on multiple threads we need to acquire a lock
1056  // before modifying the shared constraint_row object.
1057  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1058 
1059  // A reference to the constraint row.
1060  NodeConstraintRow & constraint_row = constraints[my_node].first;
1061 
1062  constraint_row.emplace(their_node, their_value);
1063  }
1064  }
1065  }
1066  }
1067 }
1068 
1069 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
1070 
1071 #endif // #ifdef LIBMESH_ENABLE_AMR
1072 
1073 
1074 
1075 #ifdef LIBMESH_ENABLE_PERIODIC
1076 
1077 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
1079  const PeriodicBoundaries & boundaries,
1080  const MeshBase & mesh,
1081  const PointLocatorBase * point_locator,
1082  const Elem * elem)
1083 {
1084  // Only bother if we truly have periodic boundaries
1085  if (boundaries.empty())
1086  return;
1087 
1088  libmesh_assert(elem);
1089 
1090  // Only constrain active elements with this method
1091  if (!elem->active())
1092  return;
1093 
1094  const unsigned int Dim = elem->dim();
1095 
1096  const FEFamily mapping_family = FEMap::map_fe_type(*elem);
1097  const FEType fe_type(elem->default_side_order(), mapping_family);
1098 
1099  // Pull objects out of the loop to reduce heap operations
1100  std::vector<const Node *> my_nodes, neigh_nodes;
1101  std::unique_ptr<const Elem> my_side, neigh_side;
1102 
1103  // Look at the element faces. Check to see if we need to
1104  // build constraints.
1105  std::vector<boundary_id_type> bc_ids;
1106  for (auto s : elem->side_index_range())
1107  {
1108  if (elem->neighbor_ptr(s))
1109  continue;
1110 
1111  mesh.get_boundary_info().boundary_ids (elem, s, bc_ids);
1112  for (const auto & boundary_id : bc_ids)
1113  {
1114  const PeriodicBoundaryBase * periodic = boundaries.boundary(boundary_id);
1115  if (periodic)
1116  {
1117  libmesh_assert(point_locator);
1118 
1119  // Get pointers to the element's neighbor.
1120  unsigned int s_neigh;
1121  const Elem * neigh = boundaries.neighbor(boundary_id, *point_locator, elem, s, &s_neigh);
1122 
1123  libmesh_error_msg_if
1124  (!neigh, "PeriodicBoundaries can't find a periodic neighbor for element " <<
1125  elem->id() << " side " << s);
1126 
1127  // h refinement constraints:
1128  // constrain dofs shared between
1129  // this element and ones as coarse
1130  // as or coarser than this element.
1131  if (neigh->level() <= elem->level())
1132  {
1133 #ifdef LIBMESH_ENABLE_AMR
1134  libmesh_assert(neigh->active());
1135 #endif // #ifdef LIBMESH_ENABLE_AMR
1136 
1137  elem->build_side_ptr(my_side, s);
1138  neigh->build_side_ptr(neigh_side, s_neigh);
1139 
1140  const unsigned int n_side_nodes = my_side->n_nodes();
1141 
1142  my_nodes.clear();
1143  my_nodes.reserve (n_side_nodes);
1144  neigh_nodes.clear();
1145  neigh_nodes.reserve (n_side_nodes);
1146 
1147  for (unsigned int n=0; n != n_side_nodes; ++n)
1148  my_nodes.push_back(my_side->node_ptr(n));
1149 
1150  for (unsigned int n=0; n != n_side_nodes; ++n)
1151  neigh_nodes.push_back(neigh_side->node_ptr(n));
1152 
1153  // Make sure we're not adding recursive constraints
1154  // due to the redundancy in the way we add periodic
1155  // boundary constraints, or adding constraints to
1156  // nodes that already have AMR constraints
1157  std::vector<bool> skip_constraint(n_side_nodes, false);
1158 
1159  for (unsigned int my_side_n=0;
1160  my_side_n < n_side_nodes;
1161  my_side_n++)
1162  {
1163  // Do not use the p_level(), if any, that is inherited by the side.
1164  libmesh_assert_less (my_side_n, FEInterface::n_dofs(fe_type, /*extra_order=*/0, my_side.get()));
1165 
1166  const Node * my_node = my_nodes[my_side_n];
1167 
1168  // If we've already got a constraint on this
1169  // node, then the periodic constraint is
1170  // redundant
1171  {
1172  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1173 
1174  if (constraints.count(my_node))
1175  {
1176  skip_constraint[my_side_n] = true;
1177  continue;
1178  }
1179  }
1180 
1181  // Compute the neighbors's side shape function values.
1182  for (unsigned int their_side_n=0;
1183  their_side_n < n_side_nodes;
1184  their_side_n++)
1185  {
1186  // Do not use the p_level(), if any, that is inherited by the side.
1187  libmesh_assert_less (their_side_n, FEInterface::n_dofs(fe_type, /*extra_order=*/0, neigh_side.get()));
1188 
1189  const Node * their_node = neigh_nodes[their_side_n];
1190 
1191  // If there's a constraint on an opposing node,
1192  // we need to see if it's constrained by
1193  // *our side* making any periodic constraint
1194  // on us recursive
1195  {
1196  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1197 
1198  if (!constraints.count(their_node))
1199  continue;
1200 
1201  const NodeConstraintRow & their_constraint_row =
1202  constraints[their_node].first;
1203 
1204  for (unsigned int orig_side_n=0;
1205  orig_side_n < n_side_nodes;
1206  orig_side_n++)
1207  {
1208  // Do not use the p_level(), if any, that is inherited by the side.
1209  libmesh_assert_less (orig_side_n, FEInterface::n_dofs(fe_type, /*extra_order=*/0, my_side.get()));
1210 
1211  const Node * orig_node = my_nodes[orig_side_n];
1212 
1213  if (their_constraint_row.count(orig_node))
1214  skip_constraint[orig_side_n] = true;
1215  }
1216  }
1217  }
1218  }
1219  for (unsigned int my_side_n=0;
1220  my_side_n < n_side_nodes;
1221  my_side_n++)
1222  {
1223  // Do not use the p_level(), if any, that is inherited by the side.
1224  libmesh_assert_less (my_side_n, FEInterface::n_dofs(fe_type, /*extra_order=*/0, my_side.get()));
1225 
1226  if (skip_constraint[my_side_n])
1227  continue;
1228 
1229  const Node * my_node = my_nodes[my_side_n];
1230 
1231  // Figure out where my node lies on their reference element.
1232  const Point neigh_point = periodic->get_corresponding_pos(*my_node);
1233 
1234  // Figure out where my node lies on their reference element.
1235  const Point mapped_point =
1236  FEMap::inverse_map(Dim-1, neigh_side.get(),
1237  neigh_point);
1238 
1239  for (unsigned int their_side_n=0;
1240  their_side_n < n_side_nodes;
1241  their_side_n++)
1242  {
1243  // Do not use the p_level(), if any, that is inherited by the side.
1244  libmesh_assert_less (their_side_n, FEInterface::n_dofs(fe_type, /*extra_order=*/0, neigh_side.get()));
1245 
1246  const Node * their_node = neigh_nodes[their_side_n];
1247  libmesh_assert(their_node);
1248 
1249  // Do not use the p_level(), if any, that is inherited by the side.
1250  const Real their_value = FEInterface::shape(fe_type,
1251  /*extra_order=*/0,
1252  neigh_side.get(),
1253  their_side_n,
1254  mapped_point);
1255 
1256  // since we may be running this method concurrently
1257  // on multiple threads we need to acquire a lock
1258  // before modifying the shared constraint_row object.
1259  {
1260  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
1261 
1262  NodeConstraintRow & constraint_row =
1263  constraints[my_node].first;
1264 
1265  constraint_row.emplace(their_node, their_value);
1266  }
1267  }
1268  }
1269  }
1270  }
1271  }
1272  }
1273 }
1274 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
1275 
1276 #endif // LIBMESH_ENABLE_PERIODIC
1277 
1278 
1279 unsigned int FEAbstract::n_quadrature_points () const
1280 {
1281  if (this->shapes_on_quadrature)
1282  {
1283  libmesh_assert(this->qrule);
1284  libmesh_assert_equal_to(this->qrule->n_points(),
1285  this->_n_total_qp);
1286  }
1287  return this->_n_total_qp;
1288 }
1289 
1290 } // namespace libMesh
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:196
FEFamily family
The type of finite element.
Definition: fe_type.h:221
ElemType
Defines an enum for geometric element types.
virtual void print_phi(std::ostream &os) const =0
Prints the value of each shape function at each quadrature point.
unsigned int _n_total_qp
The total number of quadrature points for the current configuration.
Definition: fe_abstract.h:774
const Elem * parent() const
Definition: elem.h:3030
A Node is like a Point, but with more information.
Definition: node.h:52
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:355
virtual Order default_side_order() const
Definition: elem.h:1012
const unsigned int invalid_uint
A number which is used quite often to represent an invalid or uninitialized value for an unsigned int...
Definition: libmesh.h:310
bool shapes_on_quadrature
A flag indicating if current data structures correspond to quadrature rule points.
Definition: fe_abstract.h:768
static void compute_node_constraints(NodeConstraints &constraints, const Elem *elem)
Computes the nodal constraint contributions (for non-conforming adapted meshes), using Lagrange geome...
Definition: fe_abstract.C:886
IntRange< unsigned short > side_index_range() const
Definition: elem.h:2710
virtual std::unique_ptr< Elem > build_side_ptr(const unsigned int i)=0
unsigned int dim
static Point inverse_map(const unsigned int dim, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true, const bool extra_checks=true)
Definition: fe_map.C:1628
We&#39;re using a class instead of a typedef to allow forward declarations and future flexibility...
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
MeshBase & mesh
PeriodicBoundaryBase * boundary(boundary_id_type id)
FEAbstract(const unsigned int dim, const FEType &fet)
Constructor.
Definition: fe_abstract.C:47
std::ostream & operator<<(std::ostream &os, const OrderWrapper &order)
Overload stream operators.
Definition: fe_type.h:182
The Node constraint storage format.
Definition: dof_map.h:156
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:215
The libMesh namespace provides an interface to certain functionality in the library.
static const unsigned int type_to_n_nodes_map[INVALID_ELEM]
This array maps the integer representation of the ElemType enum to the number of nodes in the element...
Definition: elem.h:650
static unsigned int max_order(const FEType &fe_t, const ElemType &el_t)
static bool on_reference_element(const Point &p, const ElemType t, const Real eps=TOLERANCE)
Definition: fe_abstract.C:637
This is the MeshBase class.
Definition: mesh_base.h:75
virtual Point get_corresponding_pos(const Point &pt) const =0
This function should be overridden by derived classes to define how one finds corresponding nodes on ...
virtual void print_dphi(std::ostream &os) const =0
Prints the value of each shape function&#39;s derivative at each quadrature point.
const dof_id_type n_nodes
Definition: tecplot_io.C:67
const unsigned int dim
The dimensionality of the object.
Definition: fe_abstract.h:660
dof_id_type id() const
Definition: dof_object.h:828
void print_xyz(std::ostream &os) const
Prints the spatial location of each quadrature point (on the physical element).
Definition: fe_abstract.C:853
static Real shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
Definition: fe_interface.C:760
virtual unsigned int n_quadrature_points() const
Definition: fe_abstract.C:1279
static void get_refspace_nodes(const ElemType t, std::vector< Point > &nodes)
Definition: fe_abstract.C:400
libmesh_assert(ctx)
unsigned int n_points() const
Definition: quadrature.h:131
QBase * qrule
A pointer to the quadrature rule employed.
Definition: fe_abstract.h:762
This is the base class for point locators.
void print_JxW(std::ostream &os) const
Prints the Jacobian times the weight for each quadrature point.
Definition: fe_abstract.C:846
std::string enum_to_string(const T e)
void print_info(std::ostream &os) const
Prints all the relevant information about the current element.
Definition: fe_abstract.C:859
const Elem * neighbor_ptr(unsigned int i) const
Definition: elem.h:2598
unsigned int level() const
Definition: elem.h:3074
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual unsigned short dim() const =0
bool subactive() const
Definition: elem.h:2959
static std::unique_ptr< FEAbstract > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
Definition: fe_abstract.C:78
const Elem * neighbor(boundary_id_type boundary_id, const PointLocatorBase &point_locator, const Elem *e, unsigned int side, unsigned int *neigh_side=nullptr) const
virtual ~FEAbstract()
Destructor.
The base class for defining periodic boundaries.
std::map< const Node *, Real, std::less< const Node * >, Threads::scalable_allocator< std::pair< const Node *const, Real > > > NodeConstraintRow
A row of the Node constraint mapping.
Definition: dof_map.h:148
virtual bool infinite() const =0
FEFamily
defines an enum for finite element families.
This class forms the foundation from which generic finite elements may be derived.
Definition: fe_abstract.h:99
FEType fe_type
The finite element type for this object.
Definition: fe_abstract.h:730
std::unique_ptr< FEMap > _fe_map
Definition: fe_abstract.h:654
virtual Order default_order() const =0
bool active() const
Definition: elem.h:2941
Class contained in FE that encapsulates mapping (i.e.
Definition: fe_map.h:47
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
static void compute_periodic_node_constraints(NodeConstraints &constraints, const PeriodicBoundaries &boundaries, const MeshBase &mesh, const PointLocatorBase *point_locator, const Elem *elem)
Computes the node position constraint equation contributions (for meshes with periodic boundary condi...
Definition: fe_abstract.C:1078
static FEFamily map_fe_type(const Elem &elem)
Definition: fe_map.C:46
spin_mutex spin_mtx
A convenient spin mutex object which can be used for obtaining locks.
Definition: threads.C:30
const RemoteElem * remote_elem
Definition: remote_elem.C:57