libMesh
fe_interface.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 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 
19 
20 // Local includes
21 #include "libmesh/fe_interface.h"
22 #include "libmesh/elem.h"
23 #include "libmesh/fe.h"
24 #include "libmesh/fe_compute_data.h"
25 #include "libmesh/dof_map.h"
26 
27 namespace libMesh
28 {
29 
30 //------------------------------------------------------------
31 //FEInterface class members
33 {
34  libmesh_error_msg("ERROR: Do not define an object of this type.");
35 }
36 
37 
38 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
39 #define fe_family_switch(dim, func_and_args, prefix, suffix) \
40  do { \
41  switch (fe_t.family) \
42  { \
43  case CLOUGH: \
44  prefix FE<dim,CLOUGH>::func_and_args suffix \
45  case HERMITE: \
46  prefix FE<dim,HERMITE>::func_and_args suffix \
47  case HIERARCHIC: \
48  prefix FE<dim,HIERARCHIC>::func_and_args suffix \
49  case L2_HIERARCHIC: \
50  prefix FE<dim,L2_HIERARCHIC>::func_and_args suffix \
51  case LAGRANGE: \
52  prefix FE<dim,LAGRANGE>::func_and_args suffix \
53  case L2_LAGRANGE: \
54  prefix FE<dim,L2_LAGRANGE>::func_and_args suffix \
55  case MONOMIAL: \
56  prefix FE<dim,MONOMIAL>::func_and_args suffix \
57  case SCALAR: \
58  prefix FE<dim,SCALAR>::func_and_args suffix \
59  case BERNSTEIN: \
60  prefix FE<dim,BERNSTEIN>::func_and_args suffix \
61  case SZABAB: \
62  prefix FE<dim,SZABAB>::func_and_args suffix \
63  case XYZ: \
64  prefix FEXYZ<dim>::func_and_args suffix \
65  case SUBDIVISION: \
66  libmesh_assert_equal_to (dim, 2); \
67  prefix FE<2,SUBDIVISION>::func_and_args suffix \
68  default: \
69  libmesh_error_msg("Unsupported family = " << fe_t.family); \
70  } \
71  } while (0)
72 
73 #define fe_family_with_vec_switch(dim, func_and_args, prefix, suffix) \
74  do { \
75  switch (fe_t.family) \
76  { \
77  case CLOUGH: \
78  prefix FE<dim,CLOUGH>::func_and_args suffix \
79  case HERMITE: \
80  prefix FE<dim,HERMITE>::func_and_args suffix \
81  case HIERARCHIC: \
82  prefix FE<dim,HIERARCHIC>::func_and_args suffix \
83  case L2_HIERARCHIC: \
84  prefix FE<dim,L2_HIERARCHIC>::func_and_args suffix \
85  case LAGRANGE: \
86  prefix FE<dim,LAGRANGE>::func_and_args suffix \
87  case LAGRANGE_VEC: \
88  prefix FELagrangeVec<dim>::func_and_args suffix \
89  case L2_LAGRANGE: \
90  prefix FE<dim,L2_LAGRANGE>::func_and_args suffix \
91  case MONOMIAL: \
92  prefix FE<dim,MONOMIAL>::func_and_args suffix \
93  case SCALAR: \
94  prefix FE<dim,SCALAR>::func_and_args suffix \
95  case BERNSTEIN: \
96  prefix FE<dim,BERNSTEIN>::func_and_args suffix \
97  case SZABAB: \
98  prefix FE<dim,SZABAB>::func_and_args suffix \
99  case XYZ: \
100  prefix FEXYZ<dim>::func_and_args suffix \
101  case SUBDIVISION: \
102  libmesh_assert_equal_to (dim, 2); \
103  prefix FE<2,SUBDIVISION>::func_and_args suffix \
104  case NEDELEC_ONE: \
105  prefix FENedelecOne<dim>::func_and_args suffix \
106  default: \
107  libmesh_error_msg("Unsupported family = " << fe_t.family); \
108  } \
109  } while (0)
110 
111 #define fe_scalar_vec_error_switch(dim, func_and_args, prefix, suffix) \
112  do { \
113  switch (fe_t.family) \
114  { \
115  case CLOUGH: \
116  prefix FE<dim,CLOUGH>::func_and_args suffix \
117  case HERMITE: \
118  prefix FE<dim,HERMITE>::func_and_args suffix \
119  case HIERARCHIC: \
120  prefix FE<dim,HIERARCHIC>::func_and_args suffix \
121  case L2_HIERARCHIC: \
122  prefix FE<dim,L2_HIERARCHIC>::func_and_args suffix \
123  case LAGRANGE: \
124  prefix FE<dim,LAGRANGE>::func_and_args suffix \
125  case L2_LAGRANGE: \
126  prefix FE<dim,L2_LAGRANGE>::func_and_args suffix \
127  case MONOMIAL: \
128  prefix FE<dim,MONOMIAL>::func_and_args suffix \
129  case SCALAR: \
130  prefix FE<dim,SCALAR>::func_and_args suffix \
131  case BERNSTEIN: \
132  prefix FE<dim,BERNSTEIN>::func_and_args suffix \
133  case SZABAB: \
134  prefix FE<dim,SZABAB>::func_and_args suffix \
135  case XYZ: \
136  prefix FEXYZ<dim>::func_and_args suffix \
137  case SUBDIVISION: \
138  libmesh_assert_equal_to (dim, 2); \
139  prefix FE<2,SUBDIVISION>::func_and_args suffix \
140  case LAGRANGE_VEC: \
141  case NEDELEC_ONE: \
142  libmesh_error_msg("Error: Can only request scalar valued elements for Real FEInterface::func_and_args"); \
143  default: \
144  libmesh_error_msg("Unsupported family = " << fe_t.family); \
145  } \
146  } while (0)
147 
148 
149 #define fe_vector_scalar_error_switch(dim, func_and_args, prefix, suffix) \
150  do { \
151  switch (fe_t.family) \
152  { \
153  case LAGRANGE_VEC: \
154  prefix FELagrangeVec<dim>::func_and_args suffix \
155  case NEDELEC_ONE: \
156  prefix FENedelecOne<dim>::func_and_args suffix \
157  case HERMITE: \
158  case HIERARCHIC: \
159  case L2_HIERARCHIC: \
160  case LAGRANGE: \
161  case L2_LAGRANGE: \
162  case MONOMIAL: \
163  case SCALAR: \
164  case BERNSTEIN: \
165  case SZABAB: \
166  case XYZ: \
167  case SUBDIVISION: \
168  libmesh_error_msg("Error: Can only request vector valued elements for RealGradient FEInterface::shape"); \
169  default: \
170  libmesh_error_msg("Unsupported family = " << fe_t.family); \
171  } \
172  } while (0)
173 
174 #else
175 #define fe_family_switch(dim, func_and_args, prefix, suffix) \
176  do { \
177  switch (fe_t.family) \
178  { \
179  case CLOUGH: \
180  prefix FE<dim,CLOUGH>::func_and_args suffix \
181  case HERMITE: \
182  prefix FE<dim,HERMITE>::func_and_args suffix \
183  case HIERARCHIC: \
184  prefix FE<dim,HIERARCHIC>::func_and_args suffix \
185  case L2_HIERARCHIC: \
186  prefix FE<dim,L2_HIERARCHIC>::func_and_args suffix \
187  case LAGRANGE: \
188  prefix FE<dim,LAGRANGE>::func_and_args suffix \
189  case L2_LAGRANGE: \
190  prefix FE<dim,L2_LAGRANGE>::func_and_args suffix \
191  case MONOMIAL: \
192  prefix FE<dim,MONOMIAL>::func_and_args suffix \
193  case SCALAR: \
194  prefix FE<dim,SCALAR>::func_and_args suffix \
195  case XYZ: \
196  prefix FEXYZ<dim>::func_and_args suffix \
197  case SUBDIVISION: \
198  libmesh_assert_equal_to (dim, 2); \
199  prefix FE<2,SUBDIVISION>::func_and_args suffix \
200  default: \
201  libmesh_error_msg("Unsupported family = " << fe_t.family); \
202  } \
203  } while (0)
204 
205 #define fe_family_with_vec_switch(dim, func_and_args, prefix, suffix) \
206  do { \
207  switch (fe_t.family) \
208  { \
209  case CLOUGH: \
210  prefix FE<dim,CLOUGH>::func_and_args suffix \
211  case HERMITE: \
212  prefix FE<dim,HERMITE>::func_and_args suffix \
213  case HIERARCHIC: \
214  prefix FE<dim,HIERARCHIC>::func_and_args suffix \
215  case L2_HIERARCHIC: \
216  prefix FE<dim,L2_HIERARCHIC>::func_and_args suffix \
217  case LAGRANGE: \
218  prefix FE<dim,LAGRANGE>::func_and_args suffix \
219  case LAGRANGE_VEC: \
220  prefix FELagrangeVec<dim>::func_and_args suffix \
221  case L2_LAGRANGE: \
222  prefix FE<dim,L2_LAGRANGE>::func_and_args suffix \
223  case MONOMIAL: \
224  prefix FE<dim,MONOMIAL>::func_and_args suffix \
225  case SCALAR: \
226  prefix FE<dim,SCALAR>::func_and_args suffix \
227  case XYZ: \
228  prefix FEXYZ<dim>::func_and_args suffix \
229  case SUBDIVISION: \
230  libmesh_assert_equal_to (dim, 2); \
231  prefix FE<2,SUBDIVISION>::func_and_args suffix \
232  case NEDELEC_ONE: \
233  prefix FENedelecOne<dim>::func_and_args suffix \
234  default: \
235  libmesh_error_msg("Unsupported family = " << fe_t.family); \
236  } \
237  } while (0)
238 
239 #define fe_scalar_vec_error_switch(dim, func_and_args, prefix, suffix) \
240  do { \
241  switch (fe_t.family) \
242  { \
243  case CLOUGH: \
244  prefix FE<dim,CLOUGH>::func_and_args suffix \
245  case HERMITE: \
246  prefix FE<dim,HERMITE>::func_and_args suffix \
247  case HIERARCHIC: \
248  prefix FE<dim,HIERARCHIC>::func_and_args suffix \
249  case L2_HIERARCHIC: \
250  prefix FE<dim,L2_HIERARCHIC>::func_and_args suffix \
251  case LAGRANGE: \
252  prefix FE<dim,LAGRANGE>::func_and_args suffix \
253  case L2_LAGRANGE: \
254  prefix FE<dim,L2_LAGRANGE>::func_and_args suffix \
255  case MONOMIAL: \
256  prefix FE<dim,MONOMIAL>::func_and_args suffix \
257  case SCALAR: \
258  prefix FE<dim,SCALAR>::func_and_args suffix \
259  case XYZ: \
260  prefix FEXYZ<dim>::func_and_args suffix \
261  case SUBDIVISION: \
262  libmesh_assert_equal_to (dim, 2); \
263  prefix FE<2,SUBDIVISION>::func_and_args suffix \
264  case LAGRANGE_VEC: \
265  case NEDELEC_ONE: \
266  libmesh_error_msg("Error: Can only request scalar valued elements for Real FEInterface::func_and_args"); \
267  default: \
268  libmesh_error_msg("Unsupported family = " << fe_t.family); \
269  } \
270  } while (0)
271 
272 
273 #define fe_vector_scalar_error_switch(dim, func_and_args, prefix, suffix) \
274  do { \
275  switch (fe_t.family) \
276  { \
277  case LAGRANGE_VEC: \
278  prefix FELagrangeVec<dim>::func_and_args suffix \
279  case NEDELEC_ONE: \
280  prefix FENedelecOne<dim>::func_and_args suffix \
281  case HERMITE: \
282  case HIERARCHIC: \
283  case L2_HIERARCHIC: \
284  case LAGRANGE: \
285  case L2_LAGRANGE: \
286  case MONOMIAL: \
287  case SCALAR: \
288  case XYZ: \
289  case SUBDIVISION: \
290  libmesh_error_msg("Error: Can only request vector valued elements for RealGradient FEInterface::func_and_args"); \
291  default: \
292  libmesh_error_msg("Unsupported family = " << fe_t.family); \
293  } \
294  } while (0)
295 #endif
296 
297 
298 #define fe_switch(func_and_args) \
299  do { \
300  switch (dim) \
301  { \
302  /* 0D */ \
303  case 0: \
304  fe_family_switch (0, func_and_args, return, ;); \
305  /* 1D */ \
306  case 1: \
307  fe_family_switch (1, func_and_args, return, ;); \
308  /* 2D */ \
309  case 2: \
310  fe_family_switch (2, func_and_args, return, ;); \
311  /* 3D */ \
312  case 3: \
313  fe_family_switch (3, func_and_args, return, ;); \
314  default: \
315  libmesh_error_msg("Invalid dim = " << dim); \
316  } \
317  } while (0)
318 
319 #define fe_with_vec_switch(func_and_args) \
320  do { \
321  switch (dim) \
322  { \
323  /* 0D */ \
324  case 0: \
325  fe_family_with_vec_switch (0, func_and_args, return, ;); \
326  /* 1D */ \
327  case 1: \
328  fe_family_with_vec_switch (1, func_and_args, return, ;); \
329  /* 2D */ \
330  case 2: \
331  fe_family_with_vec_switch (2, func_and_args, return, ;); \
332  /* 3D */ \
333  case 3: \
334  fe_family_with_vec_switch (3, func_and_args, return, ;); \
335  default: \
336  libmesh_error_msg("Invalid dim = " << dim); \
337  } \
338  } while (0)
339 
340 
341 #define void_fe_switch(func_and_args) \
342  do { \
343  switch (dim) \
344  { \
345  /* 0D */ \
346  case 0: \
347  fe_family_switch (0, func_and_args, ;, ; return;); \
348  /* 1D */ \
349  case 1: \
350  fe_family_switch (1, func_and_args, ;, ; return;); \
351  /* 2D */ \
352  case 2: \
353  fe_family_switch (2, func_and_args, ;, ; return;); \
354  /* 3D */ \
355  case 3: \
356  fe_family_switch (3, func_and_args, ;, ; return;); \
357  default: \
358  libmesh_error_msg("Invalid dim = " << dim); \
359  } \
360  } while (0)
361 
362 #define void_fe_with_vec_switch(func_and_args) \
363  do { \
364  switch (dim) \
365  { \
366  /* 0D */ \
367  case 0: \
368  fe_family_with_vec_switch (0, func_and_args, ;, ; return;); \
369  /* 1D */ \
370  case 1: \
371  fe_family_with_vec_switch (1, func_and_args, ;, ; return;); \
372  /* 2D */ \
373  case 2: \
374  fe_family_with_vec_switch (2, func_and_args, ;, ; return;); \
375  /* 3D */ \
376  case 3: \
377  fe_family_with_vec_switch (3, func_and_args, ;, ; return;); \
378  default: \
379  libmesh_error_msg("Invalid dim = " << dim); \
380  } \
381  } while (0)
382 
383 
384 
385 unsigned int FEInterface::n_shape_functions(const unsigned int dim,
386  const FEType & fe_t,
387  const ElemType t)
388 {
389 
390 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
391  /*
392  * Since the FEType, stored in DofMap/(some System child), has to
393  * be the _same_ for InfFE and FE, we have to catch calls
394  * to infinite elements through the element type.
395  */
396 
397  if (is_InfFE_elem(t))
398  return ifem_n_shape_functions(dim, fe_t, t);
399 
400 #endif
401 
402  const Order o = fe_t.order;
403 
404  fe_with_vec_switch(n_shape_functions(t, o));
405 
406  libmesh_error_msg("We'll never get here!");
407  return 0;
408 }
409 
410 
411 
412 
413 
414 unsigned int FEInterface::n_dofs(const unsigned int dim,
415  const FEType & fe_t,
416  const ElemType t)
417 {
418 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
419 
420  if (is_InfFE_elem(t))
421  return ifem_n_dofs(dim, fe_t, t);
422 
423 #endif
424 
425  const Order o = fe_t.order;
426 
427  fe_with_vec_switch(n_dofs(t, o));
428 
429  libmesh_error_msg("We'll never get here!");
430  return 0;
431 }
432 
433 
434 
435 
436 unsigned int FEInterface::n_dofs_at_node(const unsigned int dim,
437  const FEType & fe_t,
438  const ElemType t,
439  const unsigned int n)
440 {
441 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
442 
443  if (is_InfFE_elem(t))
444  return ifem_n_dofs_at_node(dim, fe_t, t, n);
445 
446 #endif
447 
448  const Order o = fe_t.order;
449 
450  fe_with_vec_switch(n_dofs_at_node(t, o, n));
451 
452  libmesh_error_msg("We'll never get here!");
453  return 0;
454 }
455 
456 
457 
460  const FEType & fe_t)
461 {
462  fe_with_vec_switch(n_dofs_at_node);
463 
464  libmesh_error_msg("We'll never get here!");
465  return libmesh_nullptr;
466 }
467 
468 
469 
470 
471 
472 
473 unsigned int FEInterface::n_dofs_per_elem(const unsigned int dim,
474  const FEType & fe_t,
475  const ElemType t)
476 {
477 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
478 
479  if (is_InfFE_elem(t))
480  return ifem_n_dofs_per_elem(dim, fe_t, t);
481 
482 #endif
483 
484  const Order o = fe_t.order;
485 
486  fe_with_vec_switch(n_dofs_per_elem(t, o));
487 
488  libmesh_error_msg("We'll never get here!");
489  return 0;
490 }
491 
492 
493 
494 
495 void FEInterface::dofs_on_side(const Elem * const elem,
496  const unsigned int dim,
497  const FEType & fe_t,
498  unsigned int s,
499  std::vector<unsigned int> & di)
500 {
501  const Order o = fe_t.order;
502 
503  void_fe_with_vec_switch(dofs_on_side(elem, o, s, di));
504 
505  libmesh_error_msg("We'll never get here!");
506 }
507 
508 
509 
510 void FEInterface::dofs_on_edge(const Elem * const elem,
511  const unsigned int dim,
512  const FEType & fe_t,
513  unsigned int e,
514  std::vector<unsigned int> & di)
515 {
516  const Order o = fe_t.order;
517 
518  void_fe_with_vec_switch(dofs_on_edge(elem, o, e, di));
519 
520  libmesh_error_msg("We'll never get here!");
521 }
522 
523 
524 
525 
526 void FEInterface::nodal_soln(const unsigned int dim,
527  const FEType & fe_t,
528  const Elem * elem,
529  const std::vector<Number> & elem_soln,
530  std::vector<Number> & nodal_soln)
531 {
532 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
533 
534  if (is_InfFE_elem(elem->type()))
535  {
536  ifem_nodal_soln(dim, fe_t, elem, elem_soln, nodal_soln);
537  return;
538  }
539 
540 #endif
541 
542  const Order order = fe_t.order;
543 
544  void_fe_with_vec_switch(nodal_soln(elem, order, elem_soln, nodal_soln));
545 }
546 
547 
548 
549 
551  const FEType & fe_t,
552  const Elem * elem,
553  const Point & p)
554 {
555 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
556  if (is_InfFE_elem(elem->type()))
557  return ifem_map(dim, fe_t, elem, p);
558 #endif
559  fe_with_vec_switch(map(elem, p));
560 
561  libmesh_error_msg("We'll never get here!");
562  return Point();
563 }
564 
565 
566 
567 
568 
569 Point FEInterface::inverse_map (const unsigned int dim,
570  const FEType & fe_t,
571  const Elem * elem,
572  const Point & p,
573  const Real tolerance,
574  const bool secure)
575 {
576 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
577 
578  if (is_InfFE_elem(elem->type()))
579  return ifem_inverse_map(dim, fe_t, elem, p,tolerance, secure);
580 
581 #endif
582 
583  fe_with_vec_switch(inverse_map(elem, p, tolerance, secure));
584 
585  libmesh_error_msg("We'll never get here!");
586  return Point();
587 }
588 
589 
590 
591 
592 void FEInterface::inverse_map (const unsigned int dim,
593  const FEType & fe_t,
594  const Elem * elem,
595  const std::vector<Point> & physical_points,
596  std::vector<Point> & reference_points,
597  const Real tolerance,
598  const bool secure)
599 {
600  const std::size_t n_pts = physical_points.size();
601 
602  // Resize the vector
603  reference_points.resize(n_pts);
604 
605  if (n_pts == 0)
606  {
607  libMesh::err << "WARNING: empty vector physical_points!"
608  << std::endl;
609  libmesh_here();
610  return;
611  }
612 
613 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
614 
615  if (is_InfFE_elem(elem->type()))
616  {
617  ifem_inverse_map(dim, fe_t, elem, physical_points, reference_points, tolerance, secure);
618  return;
619  // libmesh_not_implemented();
620  }
621 
622 #endif
623 
624  void_fe_with_vec_switch(inverse_map(elem, physical_points, reference_points, tolerance, secure));
625 
626  libmesh_error_msg("We'll never get here!");
627 }
628 
629 
630 
632  const ElemType t,
633  const Real eps)
634 {
635  return FEBase::on_reference_element(p,t,eps);
636 }
637 
638 
639 
640 
641 Real FEInterface::shape(const unsigned int dim,
642  const FEType & fe_t,
643  const ElemType t,
644  const unsigned int i,
645  const Point & p)
646 {
647 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
648 
649  if (is_InfFE_elem(t))
650  return ifem_shape(dim, fe_t, t, i, p);
651 
652 #endif
653 
654  const Order o = fe_t.order;
655 
656  fe_switch(shape(t,o,i,p));
657 
658  libmesh_error_msg("We'll never get here!");
659  return 0.;
660 }
661 
662 Real FEInterface::shape(const unsigned int dim,
663  const FEType & fe_t,
664  const Elem * elem,
665  const unsigned int i,
666  const Point & p)
667 {
668 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
669 
670  if (elem && is_InfFE_elem(elem->type()))
671  return ifem_shape(dim, fe_t, elem, i, p);
672 
673 #endif
674 
675  const Order o = fe_t.order;
676 
677  fe_switch(shape(elem,o,i,p));
678 
679  libmesh_error_msg("We'll never get here!");
680  return 0.;
681 }
682 
683 template<>
684 void FEInterface::shape<Real>(const unsigned int dim,
685  const FEType & fe_t,
686  const ElemType t,
687  const unsigned int i,
688  const Point & p,
689  Real & phi)
690 {
691 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
692 
693  if (is_InfFE_elem(t))
694  {
695  phi = ifem_shape(dim, fe_t, t, i, p);
696  return;
697  }
698 
699 #endif
700 
701  const Order o = fe_t.order;
702 
703  switch(dim)
704  {
705  case 0:
706  fe_scalar_vec_error_switch(0, shape(t,o,i,p), phi = , ; break;);
707  break;
708  case 1:
709  fe_scalar_vec_error_switch(1, shape(t,o,i,p), phi = , ; break;);
710  break;
711  case 2:
712  fe_scalar_vec_error_switch(2, shape(t,o,i,p), phi = , ; break;);
713  break;
714  case 3:
715  fe_scalar_vec_error_switch(3, shape(t,o,i,p), phi = , ; break;);
716  break;
717  default:
718  libmesh_error_msg("Invalid dimension = " << dim);
719  }
720 
721  return;
722 }
723 
724 template<>
725 void FEInterface::shape<Real>(const unsigned int dim,
726  const FEType & fe_t,
727  const Elem * elem,
728  const unsigned int i,
729  const Point & p,
730  Real & phi)
731 {
732 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
733 
734  if (elem && is_InfFE_elem(elem->type()))
735  {
736  phi = ifem_shape(dim, fe_t, elem, i, p);
737  return;
738  }
739 #endif
740 
741  const Order o = fe_t.order;
742 
743  switch(dim)
744  {
745  case 0:
746  fe_scalar_vec_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
747  break;
748  case 1:
749  fe_scalar_vec_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
750  break;
751  case 2:
752  fe_scalar_vec_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
753  break;
754  case 3:
755  fe_scalar_vec_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
756  break;
757  default:
758  libmesh_error_msg("Invalid dimension = " << dim);
759  }
760 
761  return;
762 }
763 
764 template<>
765 void FEInterface::shape<RealGradient>(const unsigned int dim,
766  const FEType & fe_t,
767  const ElemType t,
768  const unsigned int i,
769  const Point & p,
770  RealGradient & phi)
771 {
772  // This even does not handle infinite elements at all!?
773  const Order o = fe_t.order;
774 
775  switch(dim)
776  {
777  case 0:
778  fe_vector_scalar_error_switch(0, shape(t,o,i,p), phi = , ; break;);
779  break;
780  case 1:
781  fe_vector_scalar_error_switch(1, shape(t,o,i,p), phi = , ; break;);
782  break;
783  case 2:
784  fe_vector_scalar_error_switch(2, shape(t,o,i,p), phi = , ; break;);
785  break;
786  case 3:
787  fe_vector_scalar_error_switch(3, shape(t,o,i,p), phi = , ; break;);
788  break;
789  default:
790  libmesh_error_msg("Invalid dimension = " << dim);
791  }
792 
793  return;
794 }
795 
796 template<>
797 void FEInterface::shape<RealGradient>(const unsigned int dim,
798  const FEType & fe_t,
799  const Elem * elem,
800  const unsigned int i,
801  const Point & p,
802  RealGradient & phi)
803 {
804  const Order o = fe_t.order;
805 
806  switch(dim)
807  {
808  case 0:
809  fe_vector_scalar_error_switch(0, shape(elem,o,i,p), phi = , ; break;);
810  break;
811  case 1:
812  fe_vector_scalar_error_switch(1, shape(elem,o,i,p), phi = , ; break;);
813  break;
814  case 2:
815  fe_vector_scalar_error_switch(2, shape(elem,o,i,p), phi = , ; break;);
816  break;
817  case 3:
818  fe_vector_scalar_error_switch(3, shape(elem,o,i,p), phi = , ; break;);
819  break;
820  default:
821  libmesh_error_msg("Invalid dimension = " << dim);
822  }
823 
824  return;
825 }
826 
827 void FEInterface::compute_data(const unsigned int dim,
828  const FEType & fe_t,
829  const Elem * elem,
831 {
832 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
833 
834  if (elem && is_InfFE_elem(elem->type()))
835  {
836  data.init();
837  ifem_compute_data(dim, fe_t, elem, data);
838  return;
839  }
840 
841 #endif
842 
843  FEType p_refined = fe_t;
844  p_refined.order = static_cast<Order>(p_refined.order + elem->p_level());
845 
846  const unsigned int n_dof = n_dofs (dim, p_refined, elem->type());
847  const Point & p = data.p;
848  data.shape.resize(n_dof);
849 
850  // set default values for all the output fields
851  data.init();
852 
853  for (unsigned int n=0; n<n_dof; n++)
854  data.shape[n] = shape(dim, p_refined, elem, n, p);
855 
856  return;
857 }
858 
859 
860 
861 #ifdef LIBMESH_ENABLE_AMR
862 
864  DofMap & dof_map,
865  const unsigned int variable_number,
866  const Elem * elem)
867 {
868  libmesh_assert(elem);
869 
870  const FEType & fe_t = dof_map.variable_type(variable_number);
871 
872  switch (elem->dim())
873  {
874  case 0:
875  case 1:
876  {
877  // No constraints in 0D/1D.
878  return;
879  }
880 
881 
882  case 2:
883  {
884  switch (fe_t.family)
885  {
886  case CLOUGH:
888  dof_map,
889  variable_number,
890  elem); return;
891 
892  case HERMITE:
894  dof_map,
895  variable_number,
896  elem); return;
897 
898  case LAGRANGE:
900  dof_map,
901  variable_number,
902  elem); return;
903 
904  case HIERARCHIC:
906  dof_map,
907  variable_number,
908  elem); return;
909 
910  case L2_HIERARCHIC:
912  dof_map,
913  variable_number,
914  elem); return;
915 
916  case LAGRANGE_VEC:
918  dof_map,
919  variable_number,
920  elem); return;
921 
922 
923 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
924  case SZABAB:
926  dof_map,
927  variable_number,
928  elem); return;
929 
930  case BERNSTEIN:
932  dof_map,
933  variable_number,
934  elem); return;
935 
936 #endif
937  default:
938  return;
939  }
940  }
941 
942 
943  case 3:
944  {
945  switch (fe_t.family)
946  {
947  case HERMITE:
949  dof_map,
950  variable_number,
951  elem); return;
952 
953  case LAGRANGE:
955  dof_map,
956  variable_number,
957  elem); return;
958 
959  case HIERARCHIC:
961  dof_map,
962  variable_number,
963  elem); return;
964 
965  case L2_HIERARCHIC:
967  dof_map,
968  variable_number,
969  elem); return;
970 
971  case LAGRANGE_VEC:
973  dof_map,
974  variable_number,
975  elem); return;
976 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
977  case SZABAB:
979  dof_map,
980  variable_number,
981  elem); return;
982 
983  case BERNSTEIN:
985  dof_map,
986  variable_number,
987  elem); return;
988 
989 #endif
990  default:
991  return;
992  }
993  }
994 
995 
996  default:
997  libmesh_error_msg("Invalid dimension = " << elem->dim());
998  }
999 }
1000 
1001 #endif // #ifdef LIBMESH_ENABLE_AMR
1002 
1003 
1004 
1005 #ifdef LIBMESH_ENABLE_PERIODIC
1006 
1008  DofMap & dof_map,
1009  const PeriodicBoundaries & boundaries,
1010  const MeshBase & mesh,
1011  const PointLocatorBase * point_locator,
1012  const unsigned int variable_number,
1013  const Elem * elem)
1014 {
1015  // No element-specific optimizations currently exist
1017  dof_map,
1018  boundaries,
1019  mesh,
1020  point_locator,
1021  variable_number,
1022  elem);
1023 }
1024 
1025 #endif // #ifdef LIBMESH_ENABLE_PERIODIC
1026 
1027 
1028 
1029 unsigned int FEInterface::max_order(const FEType & fe_t,
1030  const ElemType & el_t)
1031 {
1032  // Yeah, I know, infinity is much larger than 11, but our
1033  // solvers don't seem to like high degree polynomials, and our
1034  // quadrature rules and number_lookups tables
1035  // need to go up higher.
1036  const unsigned int unlimited = 11;
1037 
1038  // If we used 0 as a default, then elements missing from this
1039  // table (e.g. infinite elements) would be considered broken.
1040  const unsigned int unknown = unlimited;
1041 
1042  switch (fe_t.family)
1043  {
1044  case LAGRANGE:
1045  case L2_LAGRANGE: // TODO: L2_LAGRANGE can have higher "max_order" than LAGRANGE
1046  case LAGRANGE_VEC:
1047  switch (el_t)
1048  {
1049  case EDGE2:
1050  case EDGE3:
1051  case EDGE4:
1052  return 3;
1053  case TRI3:
1054  case TRISHELL3:
1055  return 1;
1056  case TRI6:
1057  return 2;
1058  case QUAD4:
1059  case QUADSHELL4:
1060  return 1;
1061  case QUAD8:
1062  case QUADSHELL8:
1063  case QUAD9:
1064  return 2;
1065  case TET4:
1066  return 1;
1067  case TET10:
1068  return 2;
1069  case HEX8:
1070  return 1;
1071  case HEX20:
1072  case HEX27:
1073  return 2;
1074  case PRISM6:
1075  return 1;
1076  case PRISM15:
1077  case PRISM18:
1078  return 2;
1079  case PYRAMID5:
1080  return 1;
1081  case PYRAMID13:
1082  case PYRAMID14:
1083  return 2;
1084  default:
1085  return unknown;
1086  }
1087  break;
1088  case MONOMIAL:
1089  switch (el_t)
1090  {
1091  case EDGE2:
1092  case EDGE3:
1093  case EDGE4:
1094  case TRI3:
1095  case TRISHELL3:
1096  case TRI6:
1097  case QUAD4:
1098  case QUADSHELL4:
1099  case QUAD8:
1100  case QUADSHELL8:
1101  case QUAD9:
1102  case TET4:
1103  case TET10:
1104  case HEX8:
1105  case HEX20:
1106  case HEX27:
1107  case PRISM6:
1108  case PRISM15:
1109  case PRISM18:
1110  case PYRAMID5:
1111  case PYRAMID13:
1112  case PYRAMID14:
1113  return unlimited;
1114  default:
1115  return unknown;
1116  }
1117  break;
1118 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
1119  case BERNSTEIN:
1120  switch (el_t)
1121  {
1122  case EDGE2:
1123  case EDGE3:
1124  case EDGE4:
1125  return unlimited;
1126  case TRI3:
1127  case TRISHELL3:
1128  return 1;
1129  case TRI6:
1130  return 6;
1131  case QUAD4:
1132  case QUADSHELL4:
1133  return 1;
1134  case QUAD8:
1135  case QUADSHELL8:
1136  case QUAD9:
1137  return unlimited;
1138  case TET4:
1139  return 1;
1140  case TET10:
1141  return 2;
1142  case HEX8:
1143  return 1;
1144  case HEX20:
1145  return 2;
1146  case HEX27:
1147  return 4;
1148  case PRISM6:
1149  case PRISM15:
1150  case PRISM18:
1151  case PYRAMID5:
1152  case PYRAMID13:
1153  case PYRAMID14:
1154  return 0;
1155  default:
1156  return unknown;
1157  }
1158  break;
1159  case SZABAB:
1160  switch (el_t)
1161  {
1162  case EDGE2:
1163  return 1;
1164  case EDGE3:
1165  case EDGE4:
1166  return 7;
1167  case TRI3:
1168  case TRISHELL3:
1169  return 1;
1170  case TRI6:
1171  return 7;
1172  case QUAD4:
1173  case QUADSHELL4:
1174  return 1;
1175  case QUAD8:
1176  case QUADSHELL8:
1177  case QUAD9:
1178  return 7;
1179  case TET4:
1180  case TET10:
1181  case HEX8:
1182  case HEX20:
1183  case HEX27:
1184  case PRISM6:
1185  case PRISM15:
1186  case PRISM18:
1187  case PYRAMID5:
1188  case PYRAMID13:
1189  case PYRAMID14:
1190  return 0;
1191  default:
1192  return unknown;
1193  }
1194  break;
1195 #endif
1196  case XYZ:
1197  switch (el_t)
1198  {
1199  case EDGE2:
1200  case EDGE3:
1201  case EDGE4:
1202  case TRI3:
1203  case TRISHELL3:
1204  case TRI6:
1205  case QUAD4:
1206  case QUADSHELL4:
1207  case QUAD8:
1208  case QUADSHELL8:
1209  case QUAD9:
1210  case TET4:
1211  case TET10:
1212  case HEX8:
1213  case HEX20:
1214  case HEX27:
1215  case PRISM6:
1216  case PRISM15:
1217  case PRISM18:
1218  case PYRAMID5:
1219  case PYRAMID13:
1220  case PYRAMID14:
1221  return unlimited;
1222  default:
1223  return unknown;
1224  }
1225  break;
1226  case CLOUGH:
1227  switch (el_t)
1228  {
1229  case EDGE2:
1230  case EDGE3:
1231  return 3;
1232  case EDGE4:
1233  case TRI3:
1234  case TRISHELL3:
1235  return 0;
1236  case TRI6:
1237  return 3;
1238  case QUAD4:
1239  case QUADSHELL4:
1240  case QUAD8:
1241  case QUADSHELL8:
1242  case QUAD9:
1243  case TET4:
1244  case TET10:
1245  case HEX8:
1246  case HEX20:
1247  case HEX27:
1248  case PRISM6:
1249  case PRISM15:
1250  case PRISM18:
1251  case PYRAMID5:
1252  case PYRAMID13:
1253  case PYRAMID14:
1254  return 0;
1255  default:
1256  return unknown;
1257  }
1258  break;
1259  case HERMITE:
1260  switch (el_t)
1261  {
1262  case EDGE2:
1263  case EDGE3:
1264  return unlimited;
1265  case EDGE4:
1266  case TRI3:
1267  case TRISHELL3:
1268  case TRI6:
1269  return 0;
1270  case QUAD4:
1271  case QUADSHELL4:
1272  return 3;
1273  case QUAD8:
1274  case QUADSHELL8:
1275  case QUAD9:
1276  return unlimited;
1277  case TET4:
1278  case TET10:
1279  return 0;
1280  case HEX8:
1281  return 3;
1282  case HEX20:
1283  case HEX27:
1284  return unlimited;
1285  case PRISM6:
1286  case PRISM15:
1287  case PRISM18:
1288  case PYRAMID5:
1289  case PYRAMID13:
1290  case PYRAMID14:
1291  return 0;
1292  default:
1293  return unknown;
1294  }
1295  break;
1296  case HIERARCHIC:
1297  switch (el_t)
1298  {
1299  case EDGE2:
1300  case EDGE3:
1301  case EDGE4:
1302  return unlimited;
1303  case TRI3:
1304  case TRISHELL3:
1305  return 1;
1306  case TRI6:
1307  return unlimited;
1308  case QUAD4:
1309  case QUADSHELL4:
1310  return 1;
1311  case QUAD8:
1312  case QUADSHELL8:
1313  case QUAD9:
1314  return unlimited;
1315  case TET4:
1316  case TET10:
1317  return 0;
1318  case HEX8:
1319  case HEX20:
1320  return 1;
1321  case HEX27:
1322  return unlimited;
1323  case PRISM6:
1324  case PRISM15:
1325  case PRISM18:
1326  case PYRAMID5:
1327  case PYRAMID13:
1328  case PYRAMID14:
1329  return 0;
1330  default:
1331  return unknown;
1332  }
1333  break;
1334  case L2_HIERARCHIC:
1335  switch (el_t)
1336  {
1337  case EDGE2:
1338  case EDGE3:
1339  case EDGE4:
1340  return unlimited;
1341  case TRI3:
1342  case TRISHELL3:
1343  return 1;
1344  case TRI6:
1345  return unlimited;
1346  case QUAD4:
1347  case QUADSHELL4:
1348  return 1;
1349  case QUAD8:
1350  case QUADSHELL8:
1351  case QUAD9:
1352  return unlimited;
1353  case TET4:
1354  case TET10:
1355  return 0;
1356  case HEX8:
1357  case HEX20:
1358  return 1;
1359  case HEX27:
1360  return unlimited;
1361  case PRISM6:
1362  case PRISM15:
1363  case PRISM18:
1364  case PYRAMID5:
1365  case PYRAMID13:
1366  case PYRAMID14:
1367  return 0;
1368  default:
1369  return unknown;
1370  }
1371  break;
1372  case SUBDIVISION:
1373  switch (el_t)
1374  {
1375  case TRI3SUBDIVISION:
1376  return unlimited;
1377  default:
1378  return unknown;
1379  }
1380  break;
1381  case NEDELEC_ONE:
1382  switch (el_t)
1383  {
1384  case TRI6:
1385  case QUAD8:
1386  case QUAD9:
1387  case HEX20:
1388  case HEX27:
1389  return 1;
1390  default:
1391  return 0;
1392  }
1393  break;
1394  default:
1395  return 0;
1396  break;
1397  }
1398 }
1399 
1400 
1401 
1403 {
1404  switch (fe_t.family)
1405  {
1406  case LAGRANGE:
1407  case L2_LAGRANGE:
1408  case MONOMIAL:
1409  case L2_HIERARCHIC:
1410  case XYZ:
1411  case SUBDIVISION:
1412  case LAGRANGE_VEC:
1413  case NEDELEC_ONE:
1414  return false;
1415  case CLOUGH:
1416  case HERMITE:
1417  case HIERARCHIC:
1418 #ifdef LIBMESH_ENABLE_HIGHER_ORDER_SHAPES
1419  case BERNSTEIN:
1420  case SZABAB:
1421 #endif
1422  default:
1423  return true;
1424  }
1425 }
1426 
1428 {
1429  return FEInterface::field_type(fe_type.family);
1430 }
1431 
1433 {
1434  switch (fe_family)
1435  {
1436  case LAGRANGE_VEC:
1437  case NEDELEC_ONE:
1438  return TYPE_VECTOR;
1439  default:
1440  return TYPE_SCALAR;
1441  }
1442 }
1443 
1444 unsigned int FEInterface::n_vec_dim (const MeshBase & mesh,
1445  const FEType & fe_type)
1446 {
1447  switch (fe_type.family)
1448  {
1449  //FIXME: We currently assume that the number of vector components is tied
1450  // to the mesh dimension. This will break for mixed-dimension meshes.
1451  case LAGRANGE_VEC:
1452  case NEDELEC_ONE:
1453  return mesh.mesh_dimension();
1454  default:
1455  return 1;
1456  }
1457 }
1458 
1459 } // namespace libMesh
static unsigned int n_dofs_per_elem(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:473
static void dofs_on_edge(const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int e, std::vector< unsigned int > &di)
Fills the vector di with the local degree of freedom indices associated with edge e of element elem A...
Definition: fe_interface.C:510
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
OStreamProxy err
FEFamily family
The type of finite element.
Definition: fe_type.h:203
static Point map(unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
Definition: fe_interface.C:550
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:414
virtual ElemType type() const =0
static unsigned int ifem_n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
unsigned int p_level() const
Definition: elem.h:2422
const FEType & variable_type(const unsigned int c) const
Definition: dof_map.h:1697
static Real ifem_shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
static void dofs_on_side(const Elem *const elem, const unsigned int dim, const FEType &fe_t, unsigned int s, std::vector< unsigned int > &di)
Fills the vector di with the local degree of freedom indices associated with side s of element elem A...
Definition: fe_interface.C:495
unsigned int dim
class FEComputeData hides arbitrary data to be passed to and from children of FEBase through the FEIn...
We&#39;re using a class instead of a typedef to allow forward declarations and future flexibility...
const Point & p
Holds the point where the data are to be computed.
ElemType
Defines an enum for geometric element types.
FEFieldType
defines an enum for finite element field types - i.e.
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
MeshBase & mesh
static FEFieldType field_type(const FEType &fe_type)
const class libmesh_nullptr_t libmesh_nullptr
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:197
The libMesh namespace provides an interface to certain functionality in the library.
static Point ifem_inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
void init()
Inits the output data to default values, provided the fields are correctly resized.
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:557
This is the MeshBase class.
Definition: mesh_base.h:68
libmesh_assert(j)
static void compute_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to var...
Definition: fe_interface.C:863
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:167
FEFamily
defines an enum for finite element families.
FEInterface()
Empty constructor.
Definition: fe_interface.C:32
static void ifem_compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
static bool extra_hanging_dofs(const FEType &fe_t)
static unsigned int n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:385
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:641
static void compute_data(const unsigned int dim, const FEType &fe_t, const Elem *elem, FEComputeData &data)
Lets the appropriate child of FEBase compute the requested data for the input specified in data...
Definition: fe_interface.C:827
static Point inverse_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p, const Real tolerance=TOLERANCE, const bool secure=true)
Definition: fe_interface.C:569
static Point ifem_map(const unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
This is the base class for point locators.
static void ifem_nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
static bool on_reference_element(const Point &p, const ElemType t, const Real eps=TOLERANCE)
Definition: fe_interface.C:631
static unsigned int ifem_n_dofs_per_elem(const unsigned int dim, const FEType &fe_t, const ElemType t)
static unsigned int n_vec_dim(const MeshBase &mesh, const FEType &fe_type)
std::vector< Number > shape
Storage for the computed shape function values.
static n_dofs_at_node_ptr n_dofs_at_node_function(const unsigned int dim, const FEType &fe_t)
Definition: fe_interface.C:459
static void compute_periodic_constraints(DofConstraints &constraints, DofMap &dof_map, const PeriodicBoundaries &boundaries, const MeshBase &mesh, const PointLocatorBase *point_locator, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for meshes with periodic boundary conditions) correspon...
Definition: fe_base.C:1655
static void compute_constraints(DofConstraints &constraints, DofMap &dof_map, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for non-conforming adapted meshes) corresponding to var...
static unsigned int n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
Definition: fe_interface.C:436
static unsigned int ifem_n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
unsigned int(* n_dofs_at_node_ptr)(const ElemType, const Order, const unsigned int)
Definition: fe_interface.h:113
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
unsigned int mesh_dimension() const
Definition: mesh_base.C:148
virtual unsigned int dim() const =0
Order
defines an enum for polynomial orders.
Definition: enum_order.h:32
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.h:468
static void nodal_soln(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln)
Build the nodal soln from the element soln.
Definition: fe_interface.C:526
Real size() const
Definition: type_vector.h:898
IterBase * data
Ideally this private member data should have protected access.
static void compute_periodic_constraints(DofConstraints &constraints, DofMap &dof_map, const PeriodicBoundaries &boundaries, const MeshBase &mesh, const PointLocatorBase *point_locator, const unsigned int variable_number, const Elem *elem)
Computes the constraint matrix contributions (for periodic boundary conditions) corresponding to vari...
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
The constraint matrix storage format.
Definition: dof_map.h:96
static unsigned int ifem_n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)