fe.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2014 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_H
21 #define LIBMESH_FE_H
22 
23 // Local includes
24 #include "libmesh/fe_base.h"
25 #include "libmesh/libmesh.h"
26 
27 // C++ includes
28 #include <cstddef>
29 
30 namespace libMesh
31 {
32 
33 // forward declarations
34 class DofConstraints;
35 class DofMap;
36 
37 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
38 
39 template <unsigned int friend_Dim, FEFamily friend_T_radial, InfMapType friend_T_map>
40 class InfFE;
41 
42 #endif
43 
44 
48 template <FEFamily T>
50 {
51  typedef Real type;
52 };
53 
54 
58 template<>
60 {
61  typedef RealGradient type;
62 };
63 
64 template<>
66 {
67  typedef RealGradient type;
68 };
69 
70 
88 //-------------------------------------------------------------
89 // FE class definition
90 template <unsigned int Dim, FEFamily T>
91 class FE : public FEGenericBase<typename FEOutputType<T>::type>
92 {
93 public:
94 
98  explicit
99  FE(const FEType& fet);
100 
101  typedef typename
104 
113  static OutputShape shape(const ElemType t,
114  const Order o,
115  const unsigned int i,
116  const Point& p);
117 
126  static OutputShape shape(const Elem* elem,
127  const Order o,
128  const unsigned int i,
129  const Point& p);
130 
138  static OutputShape shape_deriv(const ElemType t,
139  const Order o,
140  const unsigned int i,
141  const unsigned int j,
142  const Point& p);
143 
150  static OutputShape shape_deriv(const Elem* elem,
151  const Order o,
152  const unsigned int i,
153  const unsigned int j,
154  const Point& p);
155 
174  static OutputShape shape_second_deriv(const ElemType t,
175  const Order o,
176  const unsigned int i,
177  const unsigned int j,
178  const Point& p);
179 
198  static OutputShape shape_second_deriv(const Elem* elem,
199  const Order o,
200  const unsigned int i,
201  const unsigned int j,
202  const Point& p);
203 
210  static void nodal_soln(const Elem* elem, const Order o,
211  const std::vector<Number>& elem_soln,
212  std::vector<Number>& nodal_soln);
213 
218  virtual unsigned int n_shape_functions () const;
219 
226  static unsigned int n_shape_functions (const ElemType t,
227  const Order o)
228  { return FE<Dim,T>::n_dofs (t,o); }
229 
236  static unsigned int n_dofs(const ElemType t,
237  const Order o);
238 
245  static unsigned int n_dofs_at_node(const ElemType t,
246  const Order o,
247  const unsigned int n);
248 
255  static unsigned int n_dofs_per_elem(const ElemType t,
256  const Order o);
257 
261  virtual FEContinuity get_continuity() const;
262 
267  virtual bool is_hierarchic() const;
268 
275  static void dofs_on_side(const Elem* const elem,
276  const Order o,
277  unsigned int s,
278  std::vector<unsigned int>& di);
285  static void dofs_on_edge(const Elem* const elem,
286  const Order o,
287  unsigned int e,
288  std::vector<unsigned int>& di);
289 
299  static Point inverse_map (const Elem* elem,
300  const Point& p,
301  const Real tolerance = TOLERANCE,
302  const bool secure = true);
303 
314  static void inverse_map (const Elem* elem,
315  const std::vector<Point>& physical_points,
316  std::vector<Point>& reference_points,
317  const Real tolerance = TOLERANCE,
318  const bool secure = true);
319 
330  virtual void reinit (const Elem* elem,
331  const std::vector<Point>* const pts = NULL,
332  const std::vector<Real>* const weights = NULL);
333 
343  virtual void reinit (const Elem* elem,
344  const unsigned int side,
345  const Real tolerance = TOLERANCE,
346  const std::vector<Point>* const pts = NULL,
347  const std::vector<Real>* const weights = NULL);
348 
358  virtual void edge_reinit (const Elem* elem,
359  const unsigned int edge,
360  const Real tolerance = TOLERANCE,
361  const std::vector<Point>* const pts = NULL,
362  const std::vector<Real>* const weights = NULL);
363 
368  virtual void side_map (const Elem* elem,
369  const Elem* side,
370  const unsigned int s,
371  const std::vector<Point>& reference_side_points,
372  std::vector<Point>& reference_points);
373 
379  virtual void attach_quadrature_rule (QBase* q);
380 
386  virtual unsigned int n_quadrature_points () const;
387 
388 #ifdef LIBMESH_ENABLE_AMR
389 
395  static void compute_constraints (DofConstraints &constraints,
396  DofMap &dof_map,
397  const unsigned int variable_number,
398  const Elem* elem);
399 #endif // #ifdef LIBMESH_ENABLE_AMR
400 
407  virtual bool shapes_need_reinit() const;
408 
413  static Point map (const Elem* elem,
414  const Point& reference_point);
415 
420  static Point map_xi (const Elem* elem,
421  const Point& reference_point);
422 
427  static Point map_eta (const Elem* elem,
428  const Point& reference_point);
429 
434  static Point map_zeta (const Elem* elem,
435  const Point& reference_point);
436 
437 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
438 
442  template <unsigned int friend_Dim, FEFamily friend_T_radial, InfMapType friend_T_map>
443  friend class InfFE;
444 #endif
445 
446 protected:
447 
455  virtual void init_shape_functions(const std::vector<Point>& qp,
456  const Elem* e);
457 
458 #ifdef LIBMESH_ENABLE_INFINITE_ELEMENTS
459 
464  virtual void init_base_shape_functions(const std::vector<Point>& qp,
465  const Elem* e);
466 
467 #endif
468 
473  std::vector<Point> cached_nodes;
474 
479 
480  unsigned int last_edge;
481 };
482 
483 
484 
493 //-------------------------------------------------------------
494 // FEHierarchic class definition
495 template <unsigned int Dim>
496 class FEClough : public FE<Dim,CLOUGH>
497 {
498 public:
499 
504  explicit
505  FEClough(const FEType& fet);
506 };
507 
508 
509 
518 //-------------------------------------------------------------
519 // FEHierarchic class definition
520 template <unsigned int Dim>
521 class FEHermite : public FE<Dim,HERMITE>
522 {
523 public:
524 
529  explicit
530  FEHermite(const FEType& fet);
531 
535  static Real hermite_raw_shape_second_deriv(const unsigned int basis_num,
536  const Real xi);
537  static Real hermite_raw_shape_deriv(const unsigned int basis_num,
538  const Real xi);
539  static Real hermite_raw_shape(const unsigned int basis_num,
540  const Real xi);
541 };
542 
543 
544 
553 //-------------------------------------------------------------
554 // FEHierarchic class definition
555 template <unsigned int Dim>
556 class FEHierarchic : public FE<Dim,HIERARCHIC>
557 {
558 public:
559 
564  explicit
565  FEHierarchic(const FEType& fet);
566 };
567 
568 
569 
578 //-------------------------------------------------------------
579 // FEL2Hierarchic class definition
580 template <unsigned int Dim>
581 class FEL2Hierarchic : public FE<Dim,L2_HIERARCHIC>
582 {
583 public:
584 
589  explicit
590  FEL2Hierarchic(const FEType& fet);
591 };
592 
593 
594 
603 //-------------------------------------------------------------
604 // FELagrange class definition
605 template <unsigned int Dim>
606 class FELagrange : public FE<Dim,LAGRANGE>
607 {
608 public:
609 
614  explicit
615  FELagrange(const FEType& fet);
616 };
617 
618 
622 //-------------------------------------------------------------
623 // FEL2Lagrange class definition
624 template <unsigned int Dim>
625 class FEL2Lagrange : public FE<Dim,L2_LAGRANGE>
626 {
627 public:
628 
633  explicit
634  FEL2Lagrange(const FEType& fet);
635 };
636 
637 
646 //-------------------------------------------------------------
647 // FEMonomial class definition
648 template <unsigned int Dim>
649 class FEMonomial : public FE<Dim,MONOMIAL>
650 {
651 public:
652 
657  explicit
658  FEMonomial(const FEType& fet);
659 };
660 
661 
662 //-------------------------------------------------------------
663 // FEScalar class definition
664 template <unsigned int Dim>
665 class FEScalar : public FE<Dim,SCALAR>
666 {
667 public:
668 
675  explicit
676  FEScalar(const FEType& fet);
677 };
678 
679 
689 //-------------------------------------------------------------
690 // FEXYZ class definition
691 template <unsigned int Dim>
692 class FEXYZ : public FE<Dim,XYZ>
693 {
694 public:
695 
700  explicit
701  FEXYZ(const FEType& fet);
702 
707  virtual void reinit (const Elem* elem,
708  const std::vector<Point>* const pts = NULL,
709  const std::vector<Real>* const weights = NULL)
710  { FE<Dim,XYZ>::reinit (elem, pts, weights); }
711 
716  virtual void reinit (const Elem* elem,
717  const unsigned int side,
718  const Real tolerance = TOLERANCE,
719  const std::vector<Point>* const pts = NULL,
720  const std::vector<Real>* const weights = NULL);
721 
722 
723 protected:
724 
732  virtual void init_shape_functions(const std::vector<Point>& qp,
733  const Elem* e);
734 
745  virtual void compute_shape_functions(const Elem* elem, const std::vector<Point>& qp);
746 
750  void compute_face_values (const Elem* elem,
751  const Elem* side,
752  const std::vector<Real>& weights);
753 };
754 
755 //-------------------------------------------------------------
756 // FELagrangeVec class definition
757 template <unsigned int Dim>
758 class FELagrangeVec : public FE<Dim,LAGRANGE_VEC>
759 {
760 public:
761 
766  explicit
767  FELagrangeVec(const FEType& fet);
768 
769 };
770 
771 //-------------------------------------------------------------
772 // FENedelecOne class definition
773 template <unsigned int Dim>
774 class FENedelecOne : public FE<Dim,NEDELEC_ONE>
775 {
776 public:
777 
782  explicit
783  FENedelecOne(const FEType& fet);
784 
785 };
786 
787 
788 
789 
793 namespace FiniteElements
794 {
800 
806 
812 
818 
819 
825 
831 
837 
838 
844 
850 
856 
857 
863 
869 
875 
876 
882 
888 
894 
895 }
896 
897 
898 
899 
900 // ------------------------------------------------------------
901 // FE class inline members
902 template <unsigned int Dim, FEFamily T>
903 inline
904 FE<Dim,T>::FE (const FEType& fet) :
905  FEGenericBase<typename FEOutputType<T>::type> (Dim,fet),
906  last_side(INVALID_ELEM),
907  last_edge(libMesh::invalid_uint)
908 {
909  // Sanity check. Make sure the
910  // Family specified in the template instantiation
911  // matches the one in the FEType object
912  libmesh_assert_equal_to (T, this->get_family());
913 }
914 
915 
916 
917 // ------------------------------------------------------------
918 // FEClough class inline members
919 template <unsigned int Dim>
920 inline
922  FE<Dim,CLOUGH> (fet)
923 {
924 }
925 
926 
927 
928 // ------------------------------------------------------------
929 // FEHermite class inline members
930 template <unsigned int Dim>
931 inline
933  FE<Dim,HERMITE> (fet)
934 {
935 }
936 
937 
938 
939 // ------------------------------------------------------------
940 // FEHierarchic class inline members
941 template <unsigned int Dim>
942 inline
944  FE<Dim,HIERARCHIC> (fet)
945 {
946 }
947 
948 
949 
950 // ------------------------------------------------------------
951 // FEL2Hierarchic class inline members
952 template <unsigned int Dim>
953 inline
955  FE<Dim,L2_HIERARCHIC> (fet)
956 {
957 }
958 
959 
960 
961 // ------------------------------------------------------------
962 // FELagrange class inline members
963 template <unsigned int Dim>
964 inline
966  FE<Dim,LAGRANGE> (fet)
967 {
968 }
969 
970 // ------------------------------------------------------------
971 // FELagrangeVec class inline members
972 template <unsigned int Dim>
973 inline
975  FE<Dim,LAGRANGE_VEC> (fet)
976 {
977 }
978 
979 // ------------------------------------------------------------
980 // FEL2Lagrange class inline members
981 template <unsigned int Dim>
982 inline
984  FE<Dim,L2_LAGRANGE> (fet)
985 {
986 }
987 
988 
989 
990 // ------------------------------------------------------------
991 // FEMonomial class inline members
992 template <unsigned int Dim>
993 inline
995  FE<Dim,MONOMIAL> (fet)
996 {
997 }
998 
999 
1000 
1001 
1002 // ------------------------------------------------------------
1003 // FEXYZ class inline members
1004 template <unsigned int Dim>
1005 inline
1007  FE<Dim,XYZ> (fet)
1008 {
1009 }
1010 
1011 // ------------------------------------------------------------
1012 // FEScalar class inline members
1013 template <unsigned int Dim>
1014 inline
1016  FE<Dim,SCALAR> (fet)
1017 {
1018 }
1019 
1020 // ------------------------------------------------------------
1021 // FENedelecOne class inline members
1022 template <unsigned int Dim>
1023 inline
1025  FE<Dim,NEDELEC_ONE> (fet)
1026 {
1027 }
1028 
1029 } // namespace libMesh
1030 
1031 #endif // LIBMESH_FE_H

Site Created By: libMesh Developers
Last modified: February 07 2014 16:57:04 UTC

Hosted By:
SourceForge.net Logo