libMesh
fe_abstract.h
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 #ifndef LIBMESH_FE_ABSTRACT_H
21 #define LIBMESH_FE_ABSTRACT_H
22 
23 // Local includes
24 #include "libmesh/reference_counted_object.h"
25 #include "libmesh/point.h"
26 #include "libmesh/vector_value.h"
27 #include "libmesh/enum_elem_type.h"
28 #include "libmesh/fe_type.h"
29 #include "libmesh/auto_ptr.h"
30 #include "libmesh/fe_map.h"
31 
32 // C++ includes
33 #include <cstddef>
34 #include <vector>
35 
36 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
37 #include "libmesh/tensor_value.h"
38 #endif
39 
40 namespace libMesh
41 {
42 
43 
44 // forward declarations
45 template <typename T> class DenseMatrix;
46 template <typename T> class DenseVector;
47 class BoundaryInfo;
48 class DofConstraints;
49 class DofMap;
50 class Elem;
51 class MeshBase;
52 template <typename T> class NumericVector;
53 class QBase;
54 
55 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
56 class NodeConstraints;
57 #endif
58 
59 #ifdef LIBMESH_ENABLE_PERIODIC
60 class PeriodicBoundaries;
61 class PointLocatorBase;
62 #endif
63 
64 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
65 template <unsigned int Dim, FEFamily T_radial, InfMapType T_map>
66 class InfFE;
67 #endif
68 
69 
70 
93 class FEAbstract : public ReferenceCountedObject<FEAbstract>
94 {
95 protected:
96 
102  FEAbstract (const unsigned int dim,
103  const FEType & fet);
104 
105 public:
106 
110  virtual ~FEAbstract();
111 
118  static UniquePtr<FEAbstract> build (const unsigned int dim,
119  const FEType & type);
120 
135  virtual void reinit (const Elem * elem,
136  const std::vector<Point> * const pts = libmesh_nullptr,
137  const std::vector<Real> * const weights = libmesh_nullptr) = 0;
138 
147  virtual void reinit (const Elem * elem,
148  const unsigned int side,
149  const Real tolerance = TOLERANCE,
150  const std::vector<Point> * const pts = libmesh_nullptr,
151  const std::vector<Real> * const weights = libmesh_nullptr) = 0;
152 
161  virtual void edge_reinit (const Elem * elem,
162  const unsigned int edge,
163  const Real tolerance = TOLERANCE,
164  const std::vector<Point> * pts = libmesh_nullptr,
165  const std::vector<Real> * weights = libmesh_nullptr) = 0;
166 
171  virtual void side_map (const Elem * elem,
172  const Elem * side,
173  const unsigned int s,
174  const std::vector<Point> & reference_side_points,
175  std::vector<Point> & reference_points) = 0;
176 
184  static bool on_reference_element(const Point & p,
185  const ElemType t,
186  const Real eps = TOLERANCE);
191  static void get_refspace_nodes(const ElemType t,
192  std::vector<Point> & nodes);
193 
194 
195 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
196 
200  static void compute_node_constraints (NodeConstraints & constraints,
201  const Elem * elem);
202 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
203 
204 #ifdef LIBMESH_ENABLE_PERIODIC
205 
206 #ifdef LIBMESH_ENABLE_NODE_CONSTRAINTS
207 
211  static void compute_periodic_node_constraints (NodeConstraints & constraints,
212  const PeriodicBoundaries & boundaries,
213  const MeshBase & mesh,
214  const PointLocatorBase * point_locator,
215  const Elem * elem);
216 #endif // LIBMESH_ENABLE_NODE_CONSTRAINTS
217 
218 #endif // LIBMESH_ENABLE_PERIODIC
219 
223  unsigned int get_dim() const
224  { return dim; }
225 
230  const std::vector<Point> & get_xyz() const
231  { return this->_fe_map->get_xyz(); }
232 
237  const std::vector<Real> & get_JxW() const
238  { return this->_fe_map->get_JxW(); }
239 
244  const std::vector<RealGradient> & get_dxyzdxi() const
245  { return this->_fe_map->get_dxyzdxi(); }
246 
251  const std::vector<RealGradient> & get_dxyzdeta() const
252  { return this->_fe_map->get_dxyzdeta(); }
253 
258  const std::vector<RealGradient> & get_dxyzdzeta() const
259  { return _fe_map->get_dxyzdzeta(); }
260 
264  const std::vector<RealGradient> & get_d2xyzdxi2() const
265  { return this->_fe_map->get_d2xyzdxi2(); }
266 
270  const std::vector<RealGradient> & get_d2xyzdeta2() const
271  { return this->_fe_map->get_d2xyzdeta2(); }
272 
273 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
274 
278  const std::vector<RealGradient> & get_d2xyzdzeta2() const
279  { return this->_fe_map->get_d2xyzdzeta2(); }
280 
281 #endif
282 
286  const std::vector<RealGradient> & get_d2xyzdxideta() const
287  { return this->_fe_map->get_d2xyzdxideta(); }
288 
289 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
290 
294  const std::vector<RealGradient> & get_d2xyzdxidzeta() const
295  { return this->_fe_map->get_d2xyzdxidzeta(); }
296 
300  const std::vector<RealGradient> & get_d2xyzdetadzeta() const
301  { return this->_fe_map->get_d2xyzdetadzeta(); }
302 
303 #endif
304 
309  const std::vector<Real> & get_dxidx() const
310  { return this->_fe_map->get_dxidx(); }
311 
316  const std::vector<Real> & get_dxidy() const
317  { return this->_fe_map->get_dxidy(); }
318 
323  const std::vector<Real> & get_dxidz() const
324  { return this->_fe_map->get_dxidz(); }
325 
330  const std::vector<Real> & get_detadx() const
331  { return this->_fe_map->get_detadx(); }
332 
337  const std::vector<Real> & get_detady() const
338  { return this->_fe_map->get_detady(); }
339 
344  const std::vector<Real> & get_detadz() const
345  { return this->_fe_map->get_detadz(); }
346 
351  const std::vector<Real> & get_dzetadx() const
352  { return this->_fe_map->get_dzetadx(); }
353 
358  const std::vector<Real> & get_dzetady() const
359  { return this->_fe_map->get_dzetady(); }
360 
365  const std::vector<Real> & get_dzetadz() const
366  { return this->_fe_map->get_dzetadz(); }
367 
371  const std::vector<std::vector<Point>> & get_tangents() const
372  { return this->_fe_map->get_tangents(); }
373 
377  const std::vector<Point> & get_normals() const
378  { return this->_fe_map->get_normals(); }
379 
383  const std::vector<Real> & get_curvatures() const
384  { return this->_fe_map->get_curvatures();}
385 
390  virtual void attach_quadrature_rule (QBase * q) = 0;
391 
397  virtual unsigned int n_shape_functions () const = 0;
398 
403  virtual unsigned int n_quadrature_points () const = 0;
404 
410  ElemType get_type() const { return elem_type; }
411 
416  unsigned int get_p_level() const { return _p_level; }
417 
421  FEType get_fe_type() const { return fe_type; }
422 
426  Order get_order() const { return static_cast<Order>(fe_type.order + _p_level); }
427 
431  void set_fe_order(int new_order) { fe_type.order = new_order; }
432 
436  virtual FEContinuity get_continuity() const = 0;
437 
442  virtual bool is_hierarchic() const = 0;
443 
447  FEFamily get_family() const { return fe_type.family; }
448 
452  const FEMap & get_fe_map() const { return *_fe_map.get(); }
453 
457  void print_JxW(std::ostream & os) const;
458 
464  virtual void print_phi(std::ostream & os) const =0;
465 
471  virtual void print_dphi(std::ostream & os) const =0;
472 
473 #ifdef LIBMESH_ENABLE_SECOND_DERIVATIVES
474 
480  virtual void print_d2phi(std::ostream & os) const =0;
481 
482 #endif
483 
488  void print_xyz(std::ostream & os) const;
489 
493  void print_info(std::ostream & os) const;
494 
498  friend std::ostream & operator << (std::ostream & os, const FEAbstract & fe);
499 
500 
501 protected:
502 
515  virtual void compute_shape_functions(const Elem *, const std::vector<Point> & ) =0;
516 
518 
519 
523  const unsigned int dim;
524 
529  mutable bool calculations_started;
530 
534  mutable bool calculate_phi;
535 
539  mutable bool calculate_dphi;
540 
544  mutable bool calculate_d2phi;
545 
549  mutable bool calculate_curl_phi;
550 
554  mutable bool calculate_div_phi;
555 
559  mutable bool calculate_dphiref;
560 
561 
568 
574 
579  unsigned int _p_level;
580 
585 
591 
598  virtual bool shapes_need_reinit() const = 0;
599 
600 };
601 
602 
603 
604 
605 // ------------------------------------------------------------
606 // FEAbstract class inline members
607 inline
608 FEAbstract::FEAbstract(const unsigned int d,
609  const FEType & fet) :
610  _fe_map( FEMap::build(fet) ),
611  dim(d),
612  calculations_started(false),
613  calculate_phi(false),
614  calculate_dphi(false),
615  calculate_d2phi(false),
616  calculate_curl_phi(false),
617  calculate_div_phi(false),
618  calculate_dphiref(false),
619  fe_type(fet),
621  _p_level(0),
623  shapes_on_quadrature(false)
624 {
625 }
626 
627 
628 inline
630 {
631 }
632 
633 } // namespace libMesh
634 
635 #endif // LIBMESH_FE_ABSTRACT_H
bool calculate_d2phi
Should we calculate shape function hessians?
Definition: fe_abstract.h:544
bool calculate_curl_phi
Should we calculate shape function curls?
Definition: fe_abstract.h:549
class FEType hides (possibly multiple) FEFamily and approximation orders, thereby enabling specialize...
Definition: fe_type.h:178
FEFamily family
The type of finite element.
Definition: fe_type.h:203
const std::vector< Real > & get_dxidx() const
Definition: fe_abstract.h:309
bool calculations_started
Have calculations with this object already been started? Then all get_* functions should already have...
Definition: fe_abstract.h:529
bool calculate_phi
Should we calculate shape functions?
Definition: fe_abstract.h:534
unsigned int _p_level
The p refinement level the current data structures are set up for.
Definition: fe_abstract.h:579
const std::vector< RealGradient > & get_dxyzdxi() const
Definition: fe_abstract.h:244
bool shapes_on_quadrature
A flag indicating if current data structures correspond to quadrature rule points.
Definition: fe_abstract.h:590
const std::vector< RealGradient > & get_dxyzdzeta() const
Definition: fe_abstract.h:258
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:796
const std::vector< Real > & get_dzetadx() const
Definition: fe_abstract.h:351
virtual void print_d2phi(std::ostream &os) const =0
Prints the value of each shape function&#39;s second derivatives at each quadrature point.
unsigned int get_dim() const
Definition: fe_abstract.h:223
void print_JxW(std::ostream &os) const
Prints the Jacobian times the weight for each quadrature point.
Definition: fe_abstract.C:756
virtual unsigned int n_quadrature_points() const =0
We&#39;re using a class instead of a typedef to allow forward declarations and future flexibility...
const std::vector< RealGradient > & get_d2xyzdzeta2() const
Definition: fe_abstract.h:278
virtual void edge_reinit(const Elem *elem, const unsigned int edge, const Real tolerance=TOLERANCE, const std::vector< Point > *pts=libmesh_nullptr, const std::vector< Real > *weights=libmesh_nullptr)=0
Reinitializes all the physical element-dependent data based on the edge of the element elem...
ElemType
Defines an enum for geometric element types.
void print_xyz(std::ostream &os) const
Prints the spatial location of each quadrature point (on the physical element).
Definition: fe_abstract.C:763
const FEMap & get_fe_map() const
Definition: fe_abstract.h:452
unsigned short int side
Definition: xdr_io.C:49
const std::vector< Real > & get_dzetadz() const
Definition: fe_abstract.h:365
This is the base class from which all geometric element types are derived.
Definition: elem.h:89
MeshBase & mesh
FEAbstract(const unsigned int dim, const FEType &fet)
Constructor.
Definition: fe_abstract.h:608
const std::vector< RealGradient > & get_d2xyzdeta2() const
Definition: fe_abstract.h:270
const class libmesh_nullptr_t libmesh_nullptr
Order get_order() const
Definition: fe_abstract.h:426
The Node constraint storage format.
Definition: dof_map.h:144
OrderWrapper order
The approximation order of the element.
Definition: fe_type.h:197
static const Real TOLERANCE
FEFamily get_family() const
Definition: fe_abstract.h:447
The libMesh namespace provides an interface to certain functionality in the library.
FEType get_fe_type() const
Definition: fe_abstract.h:421
bool calculate_div_phi
Should we calculate shape function divergences?
Definition: fe_abstract.h:554
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
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
ElemType get_type() const
Definition: fe_abstract.h:410
const std::vector< Point > & get_normals() const
Definition: fe_abstract.h:377
const std::vector< RealGradient > & get_d2xyzdxidzeta() const
Definition: fe_abstract.h:294
FEFamily
defines an enum for finite element families.
virtual unsigned int n_shape_functions() const =0
const std::vector< Real > & get_JxW() const
Definition: fe_abstract.h:237
const unsigned int dim
The dimensionality of the object.
Definition: fe_abstract.h:523
const std::vector< RealGradient > & get_d2xyzdxi2() const
Definition: fe_abstract.h:264
static void get_refspace_nodes(const ElemType t, std::vector< Point > &nodes)
Definition: fe_abstract.C:259
const std::vector< Real > & get_detadz() const
Definition: fe_abstract.h:344
const std::vector< Real > & get_detady() const
Definition: fe_abstract.h:337
QBase * qrule
A pointer to the quadrature rule employed.
Definition: fe_abstract.h:584
This is the base class for point locators.
This class implements reference counting.
bool calculate_dphiref
Should we calculate reference shape function gradients?
Definition: fe_abstract.h:559
FEContinuity
defines an enum for finite element types to libmesh_assert a certain level (or type? Hcurl?) of continuity.
UniquePtr< FEMap > _fe_map
Definition: fe_abstract.h:517
virtual void compute_shape_functions(const Elem *, const std::vector< Point > &)=0
After having updated the jacobian and the transformation from local to global coordinates in FEMap::c...
const std::vector< RealGradient > & get_d2xyzdetadzeta() const
Definition: fe_abstract.h:300
const std::vector< std::vector< Point > > & get_tangents() const
Definition: fe_abstract.h:371
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const std::vector< Real > & get_curvatures() const
Definition: fe_abstract.h:383
bool calculate_dphi
Should we calculate shape function gradients?
Definition: fe_abstract.h:539
virtual void print_phi(std::ostream &os) const =0
Prints the value of each shape function at each quadrature point.
virtual bool shapes_need_reinit() const =0
void print_info(std::ostream &os) const
Prints all the relevant information about the current element.
Definition: fe_abstract.C:769
const std::vector< Real > & get_dzetady() const
Definition: fe_abstract.h:358
Order
defines an enum for polynomial orders.
Definition: enum_order.h:32
const std::vector< RealGradient > & get_dxyzdeta() const
Definition: fe_abstract.h:251
const std::vector< Real > & get_detadx() const
Definition: fe_abstract.h:330
friend std::ostream & operator<<(std::ostream &os, const FEAbstract &fe)
Same as above, but allows you to print to a stream.
Definition: fe_abstract.C:785
virtual ~FEAbstract()
Destructor.
Definition: fe_abstract.h:629
const std::vector< Real > & get_dxidy() const
Definition: fe_abstract.h:316
virtual FEContinuity get_continuity() const =0
void set_fe_order(int new_order)
Sets the base FE order of the finite element.
Definition: fe_abstract.h:431
const std::vector< RealGradient > & get_d2xyzdxideta() const
Definition: fe_abstract.h:286
static UniquePtr< FEAbstract > build(const unsigned int dim, const FEType &type)
Builds a specific finite element type.
Definition: fe_abstract.C:44
This class forms the foundation from which generic finite elements may be derived.
Definition: fe_abstract.h:93
FEType fe_type
The finite element type for this object.
Definition: fe_abstract.h:567
virtual bool is_hierarchic() const =0
Class contained in FE that encapsulates mapping (i.e.
Definition: fe_map.h:45
ElemType elem_type
The element type the current data structures are set up for.
Definition: fe_abstract.h:573
A Point defines a location in LIBMESH_DIM dimensional Real space.
Definition: point.h:38
virtual void attach_quadrature_rule(QBase *q)=0
Provides the class with the quadrature rule.
const std::vector< Real > & get_dxidz() const
Definition: fe_abstract.h:323
The QBase class provides the basic functionality from which various quadrature rules can be derived...
Definition: quadrature.h:53
const std::vector< Point > & get_xyz() const
Definition: fe_abstract.h:230
virtual void print_dphi(std::ostream &os) const =0
Prints the value of each shape function&#39;s derivative at each quadrature point.
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:939
virtual void side_map(const Elem *elem, const Elem *side, const unsigned int s, const std::vector< Point > &reference_side_points, std::vector< Point > &reference_points)=0
Computes the reference space quadrature points on the side of an element based on the side quadrature...
virtual void reinit(const Elem *elem, const std::vector< Point > *const pts=libmesh_nullptr, const std::vector< Real > *const weights=libmesh_nullptr)=0
This is at the core of this class.
unsigned int get_p_level() const
Definition: fe_abstract.h:416