libMesh
fe_interface.h
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 
19 
20 #ifndef LIBMESH_FE_INTERFACE_H
21 #define LIBMESH_FE_INTERFACE_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/vector_value.h"
26 
27 // C++ includes
28 #include <map>
29 #include <vector>
30 
31 namespace libMesh
32 {
33 
34 
35 // forward declarations
36 class BoundaryInfo;
37 class DofConstraints;
38 class DofMap;
39 class Elem;
40 class FEType;
41 class FEComputeData;
42 class Point;
43 class MeshBase;
44 enum FEFamily : int;
45 enum Order : int;
46 enum FEFieldType : int;
47 enum ElemType : int;
48 enum FEContinuity : int;
49 
50 #ifdef LIBMESH_ENABLE_PERIODIC
51 class PeriodicBoundaries;
52 class PointLocatorBase;
53 #endif
54 
66 {
67 private:
68 
72  FEInterface();
73 
74 public:
75 
79  virtual ~FEInterface() = default;
80 
81 #ifdef LIBMESH_ENABLE_DEPRECATED
82 
85  static unsigned int n_shape_functions(const unsigned int dim,
86  const FEType & fe_t,
87  const ElemType t);
88 #endif // LIBMESH_ENABLE_DEPRECATED
89 
97  static unsigned int n_shape_functions(const FEType & fe_t,
98  const Elem * elem,
99  const bool add_p_level = true);
100 
105  static unsigned int n_shape_functions(const FEType & fe_t,
106  const int extra_order,
107  const Elem * elem);
108 
109 #ifdef LIBMESH_ENABLE_DEPRECATED
110 
113  static unsigned int n_dofs(const unsigned int dim,
114  const FEType & fe_t,
115  const ElemType t);
116 
126  static unsigned int n_dofs(const unsigned int dim,
127  const FEType & fe_t,
128  const Elem * elem);
129 #endif // LIBMESH_ENABLE_DEPRECATED
130 
137  static unsigned int n_dofs(const FEType & fe_t,
138  const Elem * elem,
139  const bool add_p_level = true);
140 
147  static unsigned int n_dofs(const FEType & fe_t,
148  int extra_order,
149  const Elem * elem);
150 
151  typedef unsigned int (*n_dofs_at_node_ptr) (const ElemType,
152  const Order,
153  const unsigned int);
154 
155 #ifdef LIBMESH_ENABLE_DEPRECATED
156 
167  static unsigned int n_dofs_at_node(const unsigned int dim,
168  const FEType & fe_t,
169  const ElemType t,
170  const unsigned int n);
171 
177  static n_dofs_at_node_ptr
178  n_dofs_at_node_function(const unsigned int dim,
179  const FEType & fe_t);
180 #endif // LIBMESH_ENABLE_DEPRECATED
181 
186  static n_dofs_at_node_ptr
187  n_dofs_at_node_function(const FEType & fe_t,
188  const Elem * elem);
189 
195  static unsigned int n_dofs_at_node(const FEType & fe_t,
196  const Elem * elem,
197  const unsigned int n,
198  const bool add_p_level = true);
199 
205  static unsigned int n_dofs_at_node(const FEType & fe_t,
206  const int extra_order,
207  const Elem * elem,
208  const unsigned int n);
209 
210 #ifdef LIBMESH_ENABLE_DEPRECATED
211 
214  static unsigned int n_dofs_per_elem(const unsigned int dim,
215  const FEType & fe_t,
216  const ElemType t);
217 #endif // LIBMESH_ENABLE_DEPRECATED
218 
226  static unsigned int n_dofs_per_elem(const FEType & fe_t,
227  const Elem * elem,
228  const bool add_p_level = true);
229 
233  static unsigned int n_dofs_per_elem(const FEType & fe_t,
234  const int extra_order,
235  const Elem * elem);
236 
247  static void dofs_on_side(const Elem * const elem,
248  const unsigned int dim,
249  const FEType & fe_t,
250  unsigned int s,
251  std::vector<unsigned int> & di,
252  const bool add_p_level = true);
253 
264  static void dofs_on_edge(const Elem * const elem,
265  const unsigned int dim,
266  const FEType & fe_t,
267  unsigned int e,
268  std::vector<unsigned int> & di,
269  const bool add_p_level = true);
270 
287  static void nodal_soln(const unsigned int dim,
288  const FEType & fe_t,
289  const Elem * elem,
290  const std::vector<Number> & elem_soln,
291  std::vector<Number> & nodal_soln,
292  const bool add_p_level = true,
293  const unsigned int vdim = 1);
294 
295 
304  static void side_nodal_soln(const FEType & fe_t,
305  const Elem * elem,
306  const unsigned int side,
307  const std::vector<Number> & elem_soln,
308  std::vector<Number> & nodal_soln,
309  const bool add_p_level = true,
310  const unsigned int vdim = 1);
311 
312 #ifdef LIBMESH_ENABLE_DEPRECATED
313 
316  static Point map(unsigned int dim,
317  const FEType & fe_t,
318  const Elem * elem,
319  const Point & p);
320 
324  static Point inverse_map (const unsigned int dim,
325  const FEType & fe_t,
326  const Elem * elem,
327  const Point & p,
328  const Real tolerance = TOLERANCE,
329  const bool secure = true);
330 
334  static void inverse_map (const unsigned int dim,
335  const FEType & fe_t,
336  const Elem * elem,
337  const std::vector<Point> & physical_points,
338  std::vector<Point> & reference_points,
339  const Real tolerance = TOLERANCE,
340  const bool secure = true);
341 
355  static bool on_reference_element(const Point & p,
356  const ElemType t,
357  const Real eps=TOLERANCE);
358 
372  static Real shape(const unsigned int dim,
373  const FEType & fe_t,
374  const ElemType t,
375  const unsigned int i,
376  const Point & p);
377 
391  static Real shape(const unsigned int dim,
392  const FEType & fe_t,
393  const Elem * elem,
394  const unsigned int i,
395  const Point & p);
396 #endif // LIBMESH_ENABLE_DEPRECATED
397 
405  static Real shape(const FEType & fe_t,
406  const Elem * elem,
407  const unsigned int i,
408  const Point & p,
409  const bool add_p_level = true);
410 
418  static Real shape(const FEType & fe_t,
419  int extra_order,
420  const Elem * elem,
421  const unsigned int i,
422  const Point & p);
423 
424 #ifdef LIBMESH_ENABLE_DEPRECATED
425 
438  template<typename OutputType>
439  static void shape(const unsigned int dim,
440  const FEType & fe_t,
441  const ElemType t,
442  const unsigned int i,
443  const Point & p,
444  OutputType & phi);
445 
459  template<typename OutputType>
460  static void shape(const unsigned int dim,
461  const FEType & fe_t,
462  const Elem * elem,
463  const unsigned int i,
464  const Point & p,
465  OutputType & phi);
466 #endif // LIBMESH_ENABLE_DEPRECATED
467 
477  template<typename OutputType>
478  static void shape(const FEType & fe_t,
479  const Elem * elem,
480  const unsigned int i,
481  const Point & p,
482  OutputType & phi);
483 
493  template<typename OutputType>
494  static void shape(const FEType & fe_t,
495  int extra_order,
496  const Elem * elem,
497  const unsigned int i,
498  const Point & p,
499  OutputType & phi);
500 
519  template<typename OutputType>
520  static void shapes(const unsigned int dim,
521  const FEType & fe_t,
522  const Elem * elem,
523  const unsigned int i,
524  const std::vector<Point> & p,
525  std::vector<OutputType> & phi,
526  const bool add_p_level = true);
527 
528  template<typename OutputType>
529  static void all_shapes(const unsigned int dim,
530  const FEType & fe_t,
531  const Elem * elem,
532  const std::vector<Point> & p,
533  std::vector<std::vector<OutputType>> & phi,
534  const bool add_p_level = true);
535 
541  typedef Real (*shape_ptr) (const FEType fe_t,
542  const Elem * elem,
543  const unsigned int i,
544  const Point & p,
545  const bool add_p_level);
546 
547 #ifdef LIBMESH_ENABLE_DEPRECATED
548 
551  static shape_ptr
552  shape_function(const unsigned int dim,
553  const FEType & fe_t,
554  const ElemType t);
555 #endif // LIBMESH_ENABLE_DEPRECATED
556 
561  static shape_ptr
562  shape_function(const FEType & fe_t,
563  const Elem * elem);
564 
565 #ifdef LIBMESH_ENABLE_DEPRECATED
566 
578  static Real shape_deriv(const unsigned int dim,
579  const FEType & fe_t,
580  const ElemType t,
581  const unsigned int i,
582  const unsigned int j,
583  const Point & p);
584 
588  static Real shape_deriv (const unsigned int dim,
589  const FEType & fe_t,
590  const Elem *elem,
591  const unsigned int i,
592  const unsigned int j,
593  const Point & p);
594 #endif // LIBMESH_ENABLE_DEPRECATED
595 
606  static Real shape_deriv(const FEType & fe_t,
607  const Elem * elem,
608  const unsigned int i,
609  const unsigned int j,
610  const Point & p);
611 
615  static Real shape_deriv(const FEType & fe_t,
616  int extra_order,
617  const Elem * elem,
618  const unsigned int i,
619  const unsigned int j,
620  const Point & p);
621 
626  template<typename OutputType>
627  static void shape_derivs(const FEType & fe_t,
628  const Elem * elem,
629  const unsigned int i,
630  const unsigned int j,
631  const std::vector<Point> & p,
632  std::vector<OutputType> & dphi,
633  const bool add_p_level = true);
634 
635  template<typename OutputType>
636  static void all_shape_derivs(const unsigned int dim,
637  const FEType & fe_t,
638  const Elem * elem,
639  const std::vector<Point> & p,
640  std::vector<std::vector<OutputType>> * comps[3],
641  const bool add_p_level = true);
642 
648  typedef Real (*shape_deriv_ptr) (const FEType fet,
649  const Elem * elem,
650  const unsigned int i,
651  const unsigned int j,
652  const Point & p,
653  const bool add_p_level);
654 
661  static shape_deriv_ptr
662  shape_deriv_function(const unsigned int dim,
663  const FEType & fe_t,
664  const ElemType t);
665 
669  static shape_deriv_ptr
670  shape_deriv_function(const FEType & fe_t,
671  const Elem * elem);
672 
673 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
674 
675 #ifdef LIBMESH_ENABLE_DEPRECATED
676 
698  static Real shape_second_deriv(const unsigned int dim,
699  const FEType & fe_t,
700  const ElemType t,
701  const unsigned int i,
702  const unsigned int j,
703  const Point & p);
704 
709  static Real shape_second_deriv (const unsigned int dim,
710  const FEType & fe_t,
711  const Elem *elem,
712  const unsigned int i,
713  const unsigned int j,
714  const Point & p);
715 #endif // LIBMESH_ENABLE_DEPRECATED
716 
732  static Real shape_second_deriv(const FEType & fe_t,
733  const Elem * elem,
734  const unsigned int i,
735  const unsigned int j,
736  const Point & p);
737 
741  static Real shape_second_deriv(const FEType & fe_t,
742  int extra_order,
743  const Elem * elem,
744  const unsigned int i,
745  const unsigned int j,
746  const Point & p);
747 
753  typedef Real (*shape_second_deriv_ptr) (const FEType fet,
754  const Elem * elem,
755  const unsigned int i,
756  const unsigned int j,
757  const Point & p,
758  const bool add_p_level);
759 
760 #ifdef LIBMESH_ENABLE_DEPRECATED
761 
765  shape_second_deriv_function(const unsigned int dim,
766  const FEType & fe_t,
767  const ElemType t);
768 #endif // LIBMESH_ENABLE_DEPRECATED
769 
775  shape_second_deriv_function(const FEType & fe_t,
776  const Elem * elem);
777 
778 #endif
779 
793  static void compute_data(const unsigned int dim,
794  const FEType & fe_t,
795  const Elem * elem,
796  FEComputeData & data);
797 
798 #ifdef LIBMESH_ENABLE_AMR
799 
804  static void compute_constraints (DofConstraints & constraints,
805  DofMap & dof_map,
806  const unsigned int variable_number,
807  const Elem * elem);
808 #endif // #ifdef LIBMESH_ENABLE_AMR
809 
810 #ifdef LIBMESH_ENABLE_PERIODIC
811 
816  static void compute_periodic_constraints (DofConstraints & constraints,
817  DofMap & dof_map,
818  const PeriodicBoundaries & boundaries,
819  const MeshBase & mesh,
820  const PointLocatorBase * point_locator,
821  const unsigned int variable_number,
822  const Elem * elem);
823 #endif // #ifdef LIBMESH_ENABLE_PERIODIC
824 
829  static unsigned int max_order (const FEType & fe_t,
830  const ElemType & el_t);
831 
836  static bool extra_hanging_dofs (const FEType & fe_t);
837 
841  static FEFieldType field_type (const FEType & fe_type);
842 
846  static FEFieldType field_type (const FEFamily & fe_family);
847 
851  static bool orientation_dependent (const FEFamily & fe_family);
852 
857  static unsigned int n_vec_dim (const MeshBase & mesh,
858  const FEType & fe_type);
859 
867  static FEContinuity get_continuity(const FEType & fe_type);
868 
873  static bool is_hierarchic(const FEType & fe_type);
874 
875 private:
876 
877 
882  static bool is_InfFE_elem(const ElemType et);
883 
884 
885 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
886 
887  /*
888  * All these private members do the same as their public
889  * counterparts, except for infinite elements. This disentangles
890  * the calls to \p FE and \p InfFE.
891  */
892 
893 #ifdef LIBMESH_ENABLE_DEPRECATED
894 
898  static unsigned int ifem_n_shape_functions(const unsigned int dim,
899  const FEType & fe_t,
900  const ElemType t);
901 #endif // LIBMESH_ENABLE_DEPRECATED
902 
903  static unsigned int ifem_n_shape_functions(const FEType & fe_t,
904  const Elem * elem);
905 
906 #ifdef LIBMESH_ENABLE_DEPRECATED
907 
911  static unsigned int ifem_n_dofs(const unsigned int dim,
912  const FEType & fe_t,
913  const ElemType t);
914 #endif // LIBMESH_ENABLE_DEPRECATED
915 
916  static unsigned int ifem_n_dofs(const FEType & fe_t,
917  const Elem * elem);
918 
919 #ifdef LIBMESH_ENABLE_DEPRECATED
920 
924  static unsigned int ifem_n_dofs_at_node(const unsigned int dim,
925  const FEType & fe_t,
926  const ElemType t,
927  const unsigned int n);
928 #endif // LIBMESH_ENABLE_DEPRECATED
929 
930  static unsigned int ifem_n_dofs_at_node(const FEType & fe_t,
931  const Elem * elem,
932  const unsigned int n);
933 
934 #ifdef LIBMESH_ENABLE_DEPRECATED
935 
939  static unsigned int ifem_n_dofs_per_elem(const unsigned int dim,
940  const FEType & fe_t,
941  const ElemType t);
942 #endif // LIBMESH_ENABLE_DEPRECATED
943 
944  static unsigned int ifem_n_dofs_per_elem(const FEType & fe_t,
945  const Elem * elem);
946 
947  static void ifem_nodal_soln(const unsigned int dim,
948  const FEType & fe_t,
949  const Elem * elem,
950  const std::vector<Number> & elem_soln,
951  std::vector<Number> & nodal_soln);
952 
953  static Point ifem_map (const unsigned int dim,
954  const FEType & fe_t,
955  const Elem * elem,
956  const Point & p);
957 
958  static Point ifem_inverse_map (const unsigned int dim,
959  const FEType & fe_t,
960  const Elem * elem,
961  const Point & p,
962  const Real tolerance = TOLERANCE,
963  const bool secure = true);
964 
965  static void ifem_inverse_map (const unsigned int dim,
966  const FEType & fe_t,
967  const Elem * elem,
968  const std::vector<Point> & physical_points,
969  std::vector<Point> & reference_points,
970  const Real tolerance = TOLERANCE,
971  const bool secure = true);
972 
973 #ifdef LIBMESH_ENABLE_DEPRECATED
974  /*
975  * \deprecated This method overload may not support all finite
976  * element types. Use \p Elem::on_reference_element() instead.
977  */
978  static bool ifem_on_reference_element(const Point & p,
979  const ElemType t,
980  const Real eps);
981 
986  static Real ifem_shape(const unsigned int dim,
987  const FEType & fe_t,
988  const ElemType t,
989  const unsigned int i,
990  const Point & p);
991 
996  static Real ifem_shape(const unsigned int dim,
997  const FEType & fe_t,
998  const Elem * elem,
999  const unsigned int i,
1000  const Point & p);
1001 #endif // LIBMESH_ENABLE_DEPRECATED
1002 
1003  static Real ifem_shape(const FEType & fe_t,
1004  const Elem * t,
1005  const unsigned int i,
1006  const Point & p);
1007 
1008 #ifdef LIBMESH_ENABLE_DEPRECATED
1009 
1013  static Real ifem_shape_deriv(const unsigned int dim,
1014  const FEType & fe_t,
1015  const ElemType t,
1016  const unsigned int i,
1017  const unsigned int j,
1018  const Point & p);
1019 
1024  static Real ifem_shape_deriv(const unsigned int dim,
1025  const FEType & fe_t,
1026  const Elem * elem,
1027  const unsigned int i,
1028  const unsigned int j,
1029  const Point & p);
1030 #endif // LIBMESH_ENABLE_DEPRECATED
1031 
1032  static Real ifem_shape_deriv(const FEType & fe_t,
1033  const Elem * elem,
1034  const unsigned int i,
1035  const unsigned int j,
1036  const Point & p);
1037 
1038  static void ifem_compute_data(const unsigned int dim,
1039  const FEType & fe_t,
1040  const Elem * elem,
1041  FEComputeData & data);
1042 
1043 #endif
1044 
1045 
1046 };
1047 
1048 } // namespace libMesh
1049 
1050 #endif // LIBMESH_FE_INTERFACE_H
Real(* shape_deriv_ptr)(const FEType fet, const Elem *elem, const unsigned int i, const unsigned int j, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function derivative values.
Definition: fe_interface.h:648
static unsigned int n_dofs_per_elem(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:530
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:196
static Point map(unsigned int dim, const FEType &fe_t, const Elem *elem, const Point &p)
This is now deprecated; use FEMap::map instead.
Definition: fe_interface.C:677
ElemType
Defines an enum for geometric element types.
Order
defines an enum for polynomial orders.
Definition: enum_order.h:40
static unsigned int n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Definition: fe_interface.C:355
static shape_ptr shape_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
static unsigned int ifem_n_shape_functions(const unsigned int dim, const FEType &fe_t, const ElemType t)
static bool ifem_on_reference_element(const Point &p, const ElemType t, const Real eps)
static Real shape_deriv(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const unsigned int j, const Point &p)
static Real ifem_shape(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const Point &p)
static constexpr Real TOLERANCE
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...
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, const bool add_p_level=true)
Fills the vector di with the local degree of freedom indices associated with edge e of element elem A...
Definition: fe_interface.C:611
This is the base class from which all geometric element types are derived.
Definition: elem.h:94
MeshBase & mesh
static FEFieldType field_type(const FEType &fe_type)
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)
static unsigned int max_order(const FEType &fe_t, const ElemType &el_t)
static void side_nodal_soln(const FEType &fe_t, const Elem *elem, const unsigned int side, const std::vector< Number > &elem_soln, std::vector< Number > &nodal_soln, const bool add_p_level=true, const unsigned int vdim=1)
Build the nodal soln on one side from the (full) element soln.
Definition: fe_interface.C:650
This is the MeshBase class.
Definition: mesh_base.h:75
static bool orientation_dependent(const FEFamily &fe_family)
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...
virtual ~FEInterface()=default
Destructor.
This class handles the numbering of degrees of freedom on a mesh.
Definition: dof_map.h:179
FEInterface()
Empty constructor.
Definition: fe_interface.C:42
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:272
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
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...
static bool is_hierarchic(const FEType &fe_type)
Returns whether or not the input FEType&#39;s higher-order shape functions are always hierarchic...
static void shapes(const unsigned int dim, const FEType &fe_t, const Elem *elem, const unsigned int i, const std::vector< Point > &p, std::vector< OutputType > &phi, const bool add_p_level=true)
Fills phi with the values of the shape function at point p.
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)
This is now deprecated; use FEMap::inverse_map instead.
Definition: fe_interface.C:692
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:751
This class provides an encapsulated access to all static public member functions of finite element cl...
Definition: fe_interface.h:65
static shape_deriv_ptr shape_deriv_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
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)
static n_dofs_at_node_ptr n_dofs_at_node_function(const unsigned int dim, const FEType &fe_t)
Definition: fe_interface.C:458
static FEContinuity get_continuity(const FEType &fe_type)
Returns the input FEType&#39;s FEContinuity based on the underlying FEFamily and potentially the Order...
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 Real shape_second_deriv(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const unsigned int j, const Point &p)
static unsigned int ifem_n_dofs(const unsigned int dim, const FEType &fe_t, const ElemType t)
Real(* shape_ptr)(const FEType fe_t, const Elem *elem, const unsigned int i, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function values.
Definition: fe_interface.h:541
unsigned int(* n_dofs_at_node_ptr)(const ElemType, const Order, const unsigned int)
Definition: fe_interface.h:151
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static shape_second_deriv_ptr shape_second_deriv_function(const unsigned int dim, const FEType &fe_t, const ElemType t)
static void shape_derivs(const FEType &fe_t, const Elem *elem, const unsigned int i, const unsigned int j, const std::vector< Point > &p, std::vector< OutputType > &dphi, const bool add_p_level=true)
Fills dphi with the derivatives of the shape function at point p in direction j. ...
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
static bool is_InfFE_elem(const ElemType et)
Definition: fe_interface.C:52
Real(* shape_second_deriv_ptr)(const FEType fet, const Elem *elem, const unsigned int i, const unsigned int j, const Point &p, const bool add_p_level)
Typedef for pointer to a function that returns FE shape function second derivative values...
Definition: fe_interface.h:753
FEFamily
defines an enum for finite element families.
static void all_shape_derivs(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Point > &p, std::vector< std::vector< OutputType >> *comps[3], const bool add_p_level=true)
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, const bool add_p_level=true, const unsigned int vdim=1)
Build the nodal soln from the element soln.
Definition: fe_interface.C:625
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...
static Real ifem_shape_deriv(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int i, const unsigned int j, 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, const bool add_p_level=true)
Fills the vector di with the local degree of freedom indices associated with side s of element elem A...
Definition: fe_interface.C:597
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:39
The constraint matrix storage format.
Definition: dof_map.h:108
void ErrorVector unsigned int
Definition: adjoints_ex3.C:360
static void all_shapes(const unsigned int dim, const FEType &fe_t, const Elem *elem, const std::vector< Point > &p, std::vector< std::vector< OutputType >> &phi, const bool add_p_level=true)
static unsigned int ifem_n_dofs_at_node(const unsigned int dim, const FEType &fe_t, const ElemType t, const unsigned int n)
FEFieldType
defines an enum for finite element field types - i.e.