MFEM  v3.2
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
fe_coll.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010, Lawrence Livermore National Security, LLC. Produced at
2 // the Lawrence Livermore National Laboratory. LLNL-CODE-443211. All Rights
3 // reserved. See file COPYRIGHT for details.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability see http://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the GNU Lesser General Public License (as published by the Free
10 // Software Foundation) version 2.1 dated February 1999.
11 
12 #ifndef MFEM_FE_COLLECTION
13 #define MFEM_FE_COLLECTION
14 
15 #include "../config/config.hpp"
16 #include "geom.hpp"
17 #include "fe.hpp"
18 
19 namespace mfem
20 {
21 
27 {
28 protected:
29  template <Geometry::Type geom>
30  static inline void GetNVE(int &nv, int &ne);
31 
32  template <Geometry::Type geom, typename v_t>
33  static inline void GetEdge(int &nv, v_t &v, int &ne, int &e, int &eo,
34  const int edge_info);
35 
36  template <Geometry::Type geom, Geometry::Type f_geom,
37  typename v_t, typename e_t, typename eo_t>
38  static inline void GetFace(int &nv, v_t &v, int &ne, e_t &e, eo_t &eo,
39  int &nf, int &f, int &fg, int &fo,
40  const int face_info);
41 
42 public:
43  virtual const FiniteElement *
44  FiniteElementForGeometry(int GeomType) const = 0;
45 
46  virtual int DofForGeometry(int GeomType) const = 0;
47 
48  virtual int * DofOrderForOrientation(int GeomType, int Or) const = 0;
49 
50  virtual const char * Name() const { return "Undefined"; }
51 
52  int HasFaceDofs(int GeomType) const;
53 
55  int GeomType) const
56  {
57  return FiniteElementForGeometry(GeomType);
58  }
59 
61 
63 
64  static FiniteElementCollection *New(const char *name);
65 
73  void SubDofOrder(int Geom, int SDim, int Info, Array<int> &dofs) const;
74 };
75 
78 {
79 public:
80  enum BasisType
81  {
82  GaussLobatto = 0, // Nodal basis, with nodes at the Gauss-Lobatto points
83  Positive = 1 // Positive basis, Bernstein polynomials
84  };
85 
86 protected:
88  char h1_name[32];
91  int *SegDofOrd[2], *TriDofOrd[6], *QuadDofOrd[8];
92 
93 public:
94  explicit H1_FECollection(const int p, const int dim = 3,
95  const int type = GaussLobatto);
96 
97  virtual const FiniteElement *FiniteElementForGeometry(int GeomType) const
98  { return H1_Elements[GeomType]; }
99  virtual int DofForGeometry(int GeomType) const
100  { return H1_dof[GeomType]; }
101  virtual int *DofOrderForOrientation(int GeomType, int Or) const;
102  virtual const char *Name() const { return h1_name; }
104 
105  BasisType GetBasisType() const { return m_type; }
106 
107  virtual ~H1_FECollection();
108 };
109 
113 {
114 public:
115  explicit H1Pos_FECollection(const int p, const int dim = 3)
116  : H1_FECollection(p, dim, Positive) { }
117 };
118 
123 {
124 public:
125  H1_Trace_FECollection(const int p, const int dim,
126  const int type = GaussLobatto);
127 };
128 
131 {
132 public:
134  {
135  GaussLegendre = 0, // Nodal basis, with nodes at the Gauss-Legendre points
136  GaussLobatto = 1, // Nodal basis, with nodes at the Gauss-Lobatto points
137  Positive = 2 // Positive basis, Bernstein polynomials
138  };
139 
140 private:
141  BasisType m_type;
142  char d_name[32];
143  FiniteElement *L2_Elements[Geometry::NumGeom];
144  FiniteElement *Tr_Elements[Geometry::NumGeom];
145  int *SegDofOrd[2]; // for rotating segment dofs in 1D
146  int *TriDofOrd[6]; // for rotating triangle dofs in 2D
147 
148 public:
149  L2_FECollection(const int p, const int dim, const int type = GaussLegendre);
150 
151  virtual const FiniteElement *FiniteElementForGeometry(int GeomType) const
152  { return L2_Elements[GeomType]; }
153  virtual int DofForGeometry(int GeomType) const
154  {
155  if (L2_Elements[GeomType])
156  {
157  return L2_Elements[GeomType]->GetDof();
158  }
159  return 0;
160  }
161  virtual int *DofOrderForOrientation(int GeomType, int Or) const;
162  virtual const char *Name() const { return d_name; }
163 
165  int GeomType) const
166  {
167  return Tr_Elements[GeomType];
168  }
169 
170  BasisType GetBasisType() const { return m_type; }
171 
172  virtual ~L2_FECollection();
173 };
174 
175 // Declare an alternative name for L2_FECollection = DG_FECollection
177 
180 {
181 protected:
182  char rt_name[32];
185  int *SegDofOrd[2], *TriDofOrd[6], *QuadDofOrd[8];
186 
187  // Initialize only the face elements
188  void InitFaces(const int p, const int dim, const int map_type,
189  const bool signs);
190 
191  // Constructor used by the constructor of RT_Trace_FECollection
192  RT_FECollection(const int p, const int dim, const int map_type,
193  const bool signs)
194  { InitFaces(p, dim, map_type, signs); }
195 
196 public:
197  RT_FECollection(const int p, const int dim);
198 
199  virtual const FiniteElement *FiniteElementForGeometry(int GeomType) const
200  { return RT_Elements[GeomType]; }
201  virtual int DofForGeometry(int GeomType) const
202  { return RT_dof[GeomType]; }
203  virtual int *DofOrderForOrientation(int GeomType, int Or) const;
204  virtual const char *Name() const { return rt_name; }
206 
207  virtual ~RT_FECollection();
208 };
209 
214 {
215 public:
216  RT_Trace_FECollection(const int p, const int dim,
217  const int map_type = FiniteElement::INTEGRAL);
218 };
219 
224 {
225 public:
226  DG_Interface_FECollection(const int p, const int dim,
227  const int map_type = FiniteElement::VALUE);
228 };
229 
232 {
233 protected:
234  char nd_name[32];
237  int *SegDofOrd[2], *TriDofOrd[6], *QuadDofOrd[8];
238 
239 public:
240  ND_FECollection(const int p, const int dim);
241 
242  virtual const FiniteElement *FiniteElementForGeometry(int GeomType) const
243  { return ND_Elements[GeomType]; }
244  virtual int DofForGeometry(int GeomType) const
245  { return ND_dof[GeomType]; }
246  virtual int *DofOrderForOrientation(int GeomType, int Or) const;
247  virtual const char *Name() const { return nd_name; }
249 
250  virtual ~ND_FECollection();
251 };
252 
257 {
258 public:
259  ND_Trace_FECollection(const int p, const int dim);
260 };
261 
264 {
265 private:
266  NURBS1DFiniteElement *SegmentFE;
267  NURBS2DFiniteElement *QuadrilateralFE;
268  NURBS3DFiniteElement *ParallelepipedFE;
269 
270  char name[16];
271 
272  void Allocate(int Order);
273  void Deallocate();
274 
275 public:
276  explicit NURBSFECollection(int Order) { Allocate(Order); }
277 
278  int GetOrder() const { return SegmentFE->GetOrder(); }
279 
281  void UpdateOrder(int Order) { Deallocate(); Allocate(Order); }
282 
283  void Reset() const
284  {
285  SegmentFE->Reset();
286  QuadrilateralFE->Reset();
287  ParallelepipedFE->Reset();
288  }
289 
290  virtual const FiniteElement *
291  FiniteElementForGeometry(int GeomType) const;
292 
293  virtual int DofForGeometry(int GeomType) const;
294 
295  virtual int *DofOrderForOrientation(int GeomType, int Or) const;
296 
297  virtual const char *Name() const { return name; }
298 
300 
301  virtual ~NURBSFECollection() { Deallocate(); }
302 };
303 
304 
307 {
308 private:
309  const PointFiniteElement PointFE;
310  const Linear1DFiniteElement SegmentFE;
311  const Linear2DFiniteElement TriangleFE;
312  const BiLinear2DFiniteElement QuadrilateralFE;
313  const Linear3DFiniteElement TetrahedronFE;
314  const TriLinear3DFiniteElement ParallelepipedFE;
315 public:
317 
318  virtual const FiniteElement *
319  FiniteElementForGeometry(int GeomType) const;
320 
321  virtual int DofForGeometry(int GeomType) const;
322 
323  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
324 
325  virtual const char * Name() const { return "Linear"; }
326 };
327 
330 {
331 private:
332  const PointFiniteElement PointFE;
333  const Quad1DFiniteElement SegmentFE;
334  const Quad2DFiniteElement TriangleFE;
335  const BiQuad2DFiniteElement QuadrilateralFE;
336  const Quadratic3DFiniteElement TetrahedronFE;
337  const LagrangeHexFiniteElement ParallelepipedFE;
338 
339 public:
340  QuadraticFECollection() : ParallelepipedFE(2) { }
341 
342  virtual const FiniteElement *
343  FiniteElementForGeometry(int GeomType) const;
344 
345  virtual int DofForGeometry(int GeomType) const;
346 
347  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
348 
349  virtual const char * Name() const { return "Quadratic"; }
350 };
351 
354 {
355 private:
356  const QuadPos1DFiniteElement SegmentFE;
357  const BiQuadPos2DFiniteElement QuadrilateralFE;
358 
359 public:
361 
362  virtual const FiniteElement *
363  FiniteElementForGeometry(int GeomType) const;
364 
365  virtual int DofForGeometry(int GeomType) const;
366 
367  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
368 
369  virtual const char * Name() const { return "QuadraticPos"; }
370 };
371 
374 {
375 private:
376  const PointFiniteElement PointFE;
377  const Cubic1DFiniteElement SegmentFE;
378  const Cubic2DFiniteElement TriangleFE;
379  const BiCubic2DFiniteElement QuadrilateralFE;
380  const Cubic3DFiniteElement TetrahedronFE;
381  const LagrangeHexFiniteElement ParallelepipedFE;
382 
383 public:
384  CubicFECollection() : ParallelepipedFE(3) { }
385 
386  virtual const FiniteElement *
387  FiniteElementForGeometry(int GeomType) const;
388 
389  virtual int DofForGeometry(int GeomType) const;
390 
391  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
392 
393  virtual const char * Name() const { return "Cubic"; }
394 };
395 
398 {
399 private:
400  const P0SegmentFiniteElement SegmentFE;
401  const CrouzeixRaviartFiniteElement TriangleFE;
402  const CrouzeixRaviartQuadFiniteElement QuadrilateralFE;
403 public:
404  CrouzeixRaviartFECollection() : SegmentFE(1) { }
405 
406  virtual const FiniteElement *
407  FiniteElementForGeometry(int GeomType) const;
408 
409  virtual int DofForGeometry(int GeomType) const;
410 
411  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
412 
413  virtual const char * Name() const { return "CrouzeixRaviart"; }
414 };
415 
418 {
419 private:
420  const P0TriangleFiniteElement TriangleFE;
421  const P1TetNonConfFiniteElement TetrahedronFE;
422  const P0QuadFiniteElement QuadrilateralFE;
423  const RotTriLinearHexFiniteElement ParallelepipedFE;
424 
425 public:
427 
428  virtual const FiniteElement *
429  FiniteElementForGeometry(int GeomType) const;
430 
431  virtual int DofForGeometry(int GeomType) const;
432 
433  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
434 
435  virtual const char * Name() const { return "LinearNonConf3D"; }
436 };
437 
438 
442 {
443 private:
444  const P0SegmentFiniteElement SegmentFE; // normal component on edge
445  const RT0TriangleFiniteElement TriangleFE;
446  const RT0QuadFiniteElement QuadrilateralFE;
447 public:
448  RT0_2DFECollection() : SegmentFE(0) { }
449 
450  virtual const FiniteElement *
451  FiniteElementForGeometry(int GeomType) const;
452 
453  virtual int DofForGeometry(int GeomType) const;
454 
455  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
456 
457  virtual const char * Name() const { return "RT0_2D"; }
458 };
459 
463 {
464 private:
465  const P1SegmentFiniteElement SegmentFE; // normal component on edge
466  const RT1TriangleFiniteElement TriangleFE;
467  const RT1QuadFiniteElement QuadrilateralFE;
468 public:
470 
471  virtual const FiniteElement *
472  FiniteElementForGeometry(int GeomType) const;
473 
474  virtual int DofForGeometry(int GeomType) const;
475 
476  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
477 
478  virtual const char * Name() const { return "RT1_2D"; }
479 };
480 
484 {
485 private:
486  const P2SegmentFiniteElement SegmentFE; // normal component on edge
487  const RT2TriangleFiniteElement TriangleFE;
488  const RT2QuadFiniteElement QuadrilateralFE;
489 public:
491 
492  virtual const FiniteElement *
493  FiniteElementForGeometry(int GeomType) const;
494 
495  virtual int DofForGeometry(int GeomType) const;
496 
497  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
498 
499  virtual const char * Name() const { return "RT2_2D"; }
500 };
501 
505 {
506 private:
507  const P0TriangleFiniteElement TriangleFE;
508  const P0QuadFiniteElement QuadrilateralFE;
509 public:
511 
512  virtual const FiniteElement *
513  FiniteElementForGeometry(int GeomType) const;
514 
515  virtual int DofForGeometry(int GeomType) const;
516 
517  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
518 
519  virtual const char * Name() const { return "Const2D"; }
520 };
521 
525 {
526 private:
527  const Linear2DFiniteElement TriangleFE;
528  const BiLinear2DFiniteElement QuadrilateralFE;
529 
530 public:
532 
533  virtual const FiniteElement *
534  FiniteElementForGeometry(int GeomType) const;
535 
536  virtual int DofForGeometry(int GeomType) const;
537 
538  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
539 
540  virtual const char * Name() const { return "LinearDiscont2D"; }
541 };
542 
545 {
546 private:
547  // const CrouzeixRaviartFiniteElement TriangleFE;
548  const GaussLinear2DFiniteElement TriangleFE;
549  const GaussBiLinear2DFiniteElement QuadrilateralFE;
550 
551 public:
553 
554  virtual const FiniteElement *
555  FiniteElementForGeometry(int GeomType) const;
556 
557  virtual int DofForGeometry(int GeomType) const;
558 
559  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
560 
561  virtual const char * Name() const { return "GaussLinearDiscont2D"; }
562 };
563 
566 {
567 private:
568  const P1OnQuadFiniteElement QuadrilateralFE;
569 public:
571  virtual const FiniteElement *
572  FiniteElementForGeometry(int GeomType) const;
573  virtual int DofForGeometry(int GeomType) const;
574  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
575  virtual const char * Name() const { return "P1OnQuad"; }
576 };
577 
581 {
582 private:
583  const Quad2DFiniteElement TriangleFE;
584  const BiQuad2DFiniteElement QuadrilateralFE;
585 
586 public:
588 
589  virtual const FiniteElement *
590  FiniteElementForGeometry(int GeomType) const;
591 
592  virtual int DofForGeometry(int GeomType) const;
593 
594  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
595 
596  virtual const char * Name() const { return "QuadraticDiscont2D"; }
597 };
598 
601 {
602 private:
603  const BiQuadPos2DFiniteElement QuadrilateralFE;
604 
605 public:
607  virtual const FiniteElement *
608  FiniteElementForGeometry(int GeomType) const;
609  virtual int DofForGeometry(int GeomType) const;
610  virtual int * DofOrderForOrientation(int GeomType, int Or) const
611  { return NULL; }
612  virtual const char * Name() const { return "QuadraticPosDiscont2D"; }
613 };
614 
617 {
618 private:
619  // const Quad2DFiniteElement TriangleFE;
620  const GaussQuad2DFiniteElement TriangleFE;
621  const GaussBiQuad2DFiniteElement QuadrilateralFE;
622 
623 public:
625 
626  virtual const FiniteElement *
627  FiniteElementForGeometry(int GeomType) const;
628 
629  virtual int DofForGeometry(int GeomType) const;
630 
631  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
632 
633  virtual const char * Name() const { return "GaussQuadraticDiscont2D"; }
634 };
635 
639 {
640 private:
641  const Cubic2DFiniteElement TriangleFE;
642  const BiCubic2DFiniteElement QuadrilateralFE;
643 
644 public:
646 
647  virtual const FiniteElement *
648  FiniteElementForGeometry(int GeomType) const;
649 
650  virtual int DofForGeometry(int GeomType) const;
651 
652  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
653 
654  virtual const char * Name() const { return "CubicDiscont2D"; }
655 };
656 
660 {
661 private:
662  const P0TetFiniteElement TetrahedronFE;
663  const P0HexFiniteElement ParallelepipedFE;
664 
665 public:
667 
668  virtual const FiniteElement *
669  FiniteElementForGeometry(int GeomType) const;
670 
671  virtual int DofForGeometry(int GeomType) const;
672 
673  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
674 
675  virtual const char * Name() const { return "Const3D"; }
676 };
677 
681 {
682 private:
683  const Linear3DFiniteElement TetrahedronFE;
684  const TriLinear3DFiniteElement ParallelepipedFE;
685 
686 public:
688 
689  virtual const FiniteElement *
690  FiniteElementForGeometry(int GeomType) const;
691 
692  virtual int DofForGeometry(int GeomType) const;
693 
694  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
695 
696  virtual const char * Name() const { return "LinearDiscont3D"; }
697 };
698 
702 {
703 private:
704  const Quadratic3DFiniteElement TetrahedronFE;
705  const LagrangeHexFiniteElement ParallelepipedFE;
706 
707 public:
708  QuadraticDiscont3DFECollection () : ParallelepipedFE(2) { }
709 
710  virtual const FiniteElement *
711  FiniteElementForGeometry(int GeomType) const;
712 
713  virtual int DofForGeometry(int GeomType) const;
714 
715  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
716 
717  virtual const char * Name() const { return "QuadraticDiscont3D"; }
718 };
719 
722 {
723 private:
724  const PointFiniteElement PointFE;
725  const RefinedLinear1DFiniteElement SegmentFE;
726  const RefinedLinear2DFiniteElement TriangleFE;
727  const RefinedBiLinear2DFiniteElement QuadrilateralFE;
728  const RefinedLinear3DFiniteElement TetrahedronFE;
729  const RefinedTriLinear3DFiniteElement ParallelepipedFE;
730 
731 public:
733 
734  virtual const FiniteElement *
735  FiniteElementForGeometry(int GeomType) const;
736 
737  virtual int DofForGeometry(int GeomType) const;
738 
739  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
740 
741  virtual const char * Name() const { return "RefinedLinear"; }
742 };
743 
747 {
748 private:
749  const Nedelec1HexFiniteElement HexahedronFE;
750  const Nedelec1TetFiniteElement TetrahedronFE;
751 
752 public:
754 
755  virtual const FiniteElement *
756  FiniteElementForGeometry(int GeomType) const;
757 
758  virtual int DofForGeometry(int GeomType) const;
759 
760  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
761 
762  virtual const char * Name() const { return "ND1_3D"; }
763 };
764 
768 {
769 private:
770  const P0TriangleFiniteElement TriangleFE;
771  const P0QuadFiniteElement QuadrilateralFE;
772  const RT0HexFiniteElement HexahedronFE;
773  const RT0TetFiniteElement TetrahedronFE;
774 public:
776 
777  virtual const FiniteElement *
778  FiniteElementForGeometry(int GeomType) const;
779 
780  virtual int DofForGeometry(int GeomType) const;
781 
782  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
783 
784  virtual const char * Name() const { return "RT0_3D"; }
785 };
786 
790 {
791 private:
792  const Linear2DFiniteElement TriangleFE;
793  const BiLinear2DFiniteElement QuadrilateralFE;
794  const RT1HexFiniteElement HexahedronFE;
795 public:
797 
798  virtual const FiniteElement *
799  FiniteElementForGeometry(int GeomType) const;
800 
801  virtual int DofForGeometry(int GeomType) const;
802 
803  virtual int * DofOrderForOrientation(int GeomType, int Or) const;
804 
805  virtual const char * Name() const { return "RT1_3D"; }
806 };
807 
810 {
811 private:
812  char d_name[32];
813  int GeomType;
814  FiniteElement *Local_Element;
815 
816 public:
817  Local_FECollection(const char *fe_name);
818 
819  virtual const FiniteElement *FiniteElementForGeometry(int _GeomType) const
820  { return (GeomType == _GeomType) ? Local_Element : NULL; }
821  virtual int DofForGeometry(int _GeomType) const
822  { return (GeomType == _GeomType) ? Local_Element->GetDof() : 0; }
823  virtual int *DofOrderForOrientation(int GeomType, int Or) const
824  { return NULL; }
825  virtual const char *Name() const { return d_name; }
826 
827  virtual ~Local_FECollection() { delete Local_Element; }
828 };
829 
830 }
831 
832 #endif
Abstract class for Finite Elements.
Definition: fe.hpp:44
int GetOrder() const
Definition: fe_coll.hpp:278
Arbitrary order non-uniform rational B-splines (NURBS) finite elements.
Definition: fe_coll.hpp:263
virtual const char * Name() const
Definition: fe_coll.hpp:393
virtual const char * Name() const
Definition: fe_coll.hpp:717
virtual const char * Name() const
Definition: fe_coll.hpp:633
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:1174
Version of QuadraticDiscont2DFECollection with dofs in the Gaussian points.
Definition: fe_coll.hpp:616
int RT_dof[Geometry::NumGeom]
Definition: fe_coll.hpp:184
FiniteElement * ND_Elements[Geometry::NumGeom]
Definition: fe_coll.hpp:235
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1279
void UpdateOrder(int Order)
Change the order of the collection.
Definition: fe_coll.hpp:281
Linear 1D element with nodes 1/3 and 2/3 (trace of RT1)
Definition: fe.hpp:901
Quadratic 1D element with nodes the Gaussian points in [0,1] (trace of RT2)
Definition: fe.hpp:911
virtual const char * Name() const
Definition: fe_coll.hpp:741
virtual FiniteElementCollection * GetTraceCollection() const
Definition: fe_coll.cpp:38
virtual const char * Name() const
Definition: fe_coll.hpp:675
H1_Trace_FECollection(const int p, const int dim, const int type=GaussLobatto)
Definition: fe_coll.cpp:1569
void SubDofOrder(int Geom, int SDim, int Info, Array< int > &dofs) const
Get the local dofs for a given sub-manifold.
Definition: fe_coll.cpp:284
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:546
virtual const char * Name() const
Definition: fe_coll.hpp:540
FiniteElementCollection * GetTraceCollection() const
Definition: fe_coll.cpp:2185
static const int NumGeom
Definition: geom.hpp:34
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:2179
virtual const char * Name() const
Definition: fe_coll.hpp:349
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:1342
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:1102
const Geometry::Type geom
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1020
virtual ~L2_FECollection()
Definition: fe_coll.cpp:1707
virtual const char * Name() const
Definition: fe_coll.hpp:575
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:635
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1263
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:865
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:923
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:901
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:831
int GetOrder() const
Returns the order of the finite element.
Definition: fe.hpp:91
Piecewise-(bi)linear continuous finite elements.
Definition: fe_coll.hpp:306
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1370
void Reset() const
Definition: fe.hpp:2073
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1202
ND_FECollection(const int p, const int dim)
Definition: fe_coll.cpp:1926
Class for quadratic FE on triangle.
Definition: fe.hpp:485
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:741
Class for quadratic FE on interval.
Definition: fe.hpp:456
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:1209
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1308
virtual const char * Name() const
Definition: fe_coll.hpp:204
FiniteElementCollection * GetTraceCollection() const
Definition: fe_coll.cpp:2081
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.hpp:201
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1054
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1040
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:823
BasisType GetBasisType() const
Definition: fe_coll.hpp:170
void InitFaces(const int p, const int dim, const int map_type, const bool signs)
Definition: fe_coll.cpp:1749
FiniteElementCollection * GetTraceCollection() const
Definition: fe_coll.cpp:1878
virtual const char * Name() const
Definition: fe_coll.hpp:297
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1241
virtual const char * Name() const
Definition: fe_coll.hpp:457
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:1138
virtual const char * Name() const
Definition: fe_coll.hpp:654
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:768
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:991
Finite element collection on a macro-element.
Definition: fe_coll.hpp:721
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:714
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:470
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:493
Class for refined bi-linear FE on quadrilateral.
Definition: fe.hpp:1040
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.hpp:242
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:755
virtual const char * Name() const
Definition: fe_coll.hpp:247
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:893
virtual ~RT_FECollection()
Definition: fe_coll.cpp:1883
virtual int DofForGeometry(int _GeomType) const
Definition: fe_coll.hpp:821
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:809
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1857
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1324
Piecewise-(bi)cubic continuous finite elements.
Definition: fe_coll.hpp:373
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:857
virtual const char * Name() const
Definition: fe_coll.hpp:478
int dim
Definition: ex3.cpp:47
Class for bilinear FE on quad with nodes at the 4 Gaussian points.
Definition: fe.hpp:431
Class for refined linear FE on interval.
Definition: fe.hpp:987
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1355
int H1_dof[Geometry::NumGeom]
Definition: fe_coll.hpp:90
Class for refined linear FE on triangle.
Definition: fe.hpp:1007
Class for refined linear FE on tetrahedron.
Definition: fe.hpp:1027
FiniteElement * RT_Elements[Geometry::NumGeom]
Definition: fe_coll.hpp:183
Class for constant FE on triangle.
Definition: fe.hpp:624
static void GetNVE(int &nv, int &ne)
Definition: fe_coll.cpp:209
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:559
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:843
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:486
L2_FECollection(const int p, const int dim, const int type=GaussLegendre)
Definition: fe_coll.cpp:1584
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:910
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1114
ND_Trace_FECollection(const int p, const int dim)
Definition: fe_coll.cpp:2099
Class for refined trilinear FE on a hexahedron.
Definition: fe.hpp:1060
Discontinuous collection defined locally by a given finite element.
Definition: fe_coll.hpp:809
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:534
Version of QuadraticDiscont2DFECollection with positive basis functions.
Definition: fe_coll.hpp:600
virtual ~ND_FECollection()
Definition: fe_coll.cpp:2087
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1130
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:672
Class for quadratic FE on tetrahedron.
Definition: fe.hpp:679
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.hpp:823
RT_Trace_FECollection(const int p, const int dim, const int map_type=FiniteElement::INTEGRAL)
Definition: fe_coll.cpp:1894
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:943
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:797
Version of QuadraticFECollection with positive basis functions.
Definition: fe_coll.hpp:353
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:1062
virtual ~Local_FECollection()
Definition: fe_coll.hpp:827
virtual const char * Name() const
Definition: fe_coll.hpp:561
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1225
virtual ~H1_FECollection()
Definition: fe_coll.cpp:1557
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:454
Class for linear FE on tetrahedron.
Definition: fe.hpp:654
Class for linear FE on interval.
Definition: fe.hpp:354
Crouzeix-Raviart finite element on quadrilateral.
Definition: fe.hpp:727
Class for linear FE on triangle.
Definition: fe.hpp:374
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:1294
Class for quadratic FE on triangle with nodes at the &quot;Gaussian&quot; points.
Definition: fe.hpp:509
virtual const char * Name() const
Definition: fe_coll.hpp:102
virtual const char * Name() const
Definition: fe_coll.hpp:696
FiniteElementCollection * GetTraceCollection() const
Definition: fe_coll.cpp:1543
FiniteElement * H1_Elements[Geometry::NumGeom]
Definition: fe_coll.hpp:89
Piecewise-linear nonconforming finite elements in 3D.
Definition: fe_coll.hpp:417
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:601
DG_Interface_FECollection(const int p, const int dim, const int map_type=FiniteElement::VALUE)
Definition: fe_coll.cpp:1910
Crouzeix-Raviart nonconforming elements in 2D.
Definition: fe_coll.hpp:397
NURBSFECollection(int Order)
Definition: fe_coll.hpp:276
Arbitrary order H(div)-conforming Raviart-Thomas finite elements.
Definition: fe_coll.hpp:179
L2_FECollection DG_FECollection
Definition: fe_coll.hpp:176
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:662
H1Pos_FECollection(const int p, const int dim=3)
Definition: fe_coll.hpp:115
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:648
static void GetEdge(int &nv, v_t &v, int &ne, int &e, int &eo, const int edge_info)
Definition: fe_coll.cpp:219
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:727
virtual const char * Name() const
Definition: fe_coll.hpp:413
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.hpp:153
virtual const char * Name() const
Definition: fe_coll.hpp:325
virtual const char * Name() const
Definition: fe_coll.hpp:499
virtual const FiniteElement * TraceFiniteElementForGeometry(int GeomType) const
Definition: fe_coll.hpp:164
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1150
virtual const char * Name() const
Definition: fe_coll.hpp:435
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:957
virtual const char * Name() const
Definition: fe_coll.hpp:762
int GetDof() const
Returns the degrees of freedom in the FE space.
Definition: fe.hpp:88
Class for bi-quadratic FE on quadrilateral.
Definition: fe.hpp:525
Class for tri-linear FE on cube.
Definition: fe.hpp:692
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:965
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:525
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:878
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.hpp:199
virtual const char * Name() const
Definition: fe_coll.hpp:596
Class for bilinear FE on quadrilateral.
Definition: fe.hpp:396
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1186
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:1028
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1005
Version of LinearDiscont2DFECollection with dofs in the Gaussian points.
Definition: fe_coll.hpp:544
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:2160
virtual int * DofOrderForOrientation(int GeomType, int Or) const =0
virtual const char * Name() const
Definition: fe_coll.hpp:162
static FiniteElementCollection * New(const char *name)
Definition: fe_coll.cpp:44
virtual const FiniteElement * TraceFiniteElementForGeometry(int GeomType) const
Definition: fe_coll.hpp:54
Class for linear FE on triangle with nodes at the 3 &quot;Gaussian&quot; points.
Definition: fe.hpp:420
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1166
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:1251
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.hpp:244
virtual const char * Name() const
Definition: fe_coll.hpp:50
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:509
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.hpp:97
Bi-quadratic element on quad with nodes at the 9 Gaussian points.
Definition: fe.hpp:564
virtual const char * Name() const
Definition: fe_coll.hpp:519
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:782
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:2054
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1522
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:699
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.hpp:99
BasisType GetBasisType() const
Definition: fe_coll.hpp:105
RT_FECollection(const int p, const int dim, const int map_type, const bool signs)
Definition: fe_coll.hpp:192
virtual const char * Name() const
Definition: fe_coll.hpp:784
H1_FECollection(const int p, const int dim=3, const int type=GaussLobatto)
Definition: fe_coll.cpp:1392
Tensor products of 1D FEs (only degree 2 is functional)
Definition: fe.hpp:966
Piecewise-(bi)quadratic continuous finite elements.
Definition: fe_coll.hpp:329
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.hpp:151
int HasFaceDofs(int GeomType) const
Definition: fe_coll.cpp:25
Arbitrary order H(curl)-conforming Nedelec finite elements.
Definition: fe_coll.hpp:231
void Reset() const
Definition: fe_coll.hpp:283
Class for cubic FE on tetrahedron.
Definition: fe.hpp:611
static void GetFace(int &nv, v_t &v, int &ne, e_t &e, eo_t &eo, int &nf, int &f, int &fg, int &fo, const int face_info)
Definition: fe_coll.cpp:239
virtual const char * Name() const
Definition: fe_coll.hpp:825
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1092
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.cpp:1690
Arbitrary order H1-conforming (continuous) finite elements.
Definition: fe_coll.hpp:77
virtual const char * Name() const
Definition: fe_coll.hpp:805
virtual const char * Name() const
Definition: fe_coll.hpp:612
virtual const FiniteElement * FiniteElementForGeometry(int _GeomType) const
Definition: fe_coll.hpp:819
virtual int DofForGeometry(int GeomType) const =0
virtual int * DofOrderForOrientation(int GeomType, int Or) const
Definition: fe_coll.hpp:610
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:2173
Local_FECollection(const char *fe_name)
Definition: fe_coll.cpp:2106
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:976
virtual const char * Name() const
Definition: fe_coll.hpp:369
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:1076
Linear (P1) finite elements on quadrilaterals.
Definition: fe_coll.hpp:565
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:585
virtual ~FiniteElementCollection()
Definition: fe_coll.hpp:62
virtual ~NURBSFECollection()
Definition: fe_coll.hpp:301
int ND_dof[Geometry::NumGeom]
Definition: fe_coll.hpp:236
virtual int DofForGeometry(int GeomType) const
Definition: fe_coll.cpp:685
Crouzeix-Raviart finite element on triangle.
Definition: fe.hpp:715
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:569
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const =0
Arbitrary order &quot;L2-conforming&quot; discontinuous finite elements.
Definition: fe_coll.hpp:130
virtual const FiniteElement * FiniteElementForGeometry(int GeomType) const
Definition: fe_coll.cpp:931