MFEM  v3.2
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
fe.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
13 #define MFEM_FE
14 
15 #include "../config/config.hpp"
16 #include "../general/array.hpp"
17 #include "../linalg/linalg.hpp"
18 #include "intrules.hpp"
19 #include "geom.hpp"
20 
21 namespace mfem
22 {
23 
24 // Base and derived classes for finite elements
25 
28 {
29 public:
30  enum
31  {
32  Pk, // polynomials of order k
33  Qk, // tensor products of polynomials of order k
34  rQk // refined tensor products of polynomials of order k
35  };
36 };
37 
38 class ElementTransformation;
39 class Coefficient;
40 class VectorCoefficient;
41 class KnotVector;
42 
45 {
46 protected:
49 
50 public:
52  enum { SCALAR, VECTOR };
53 
71  enum { VALUE, INTEGRAL, H_DIV, H_CURL };
72 
79  FiniteElement(int D, int G, int Do, int O, int F = FunctionSpace::Pk);
80 
82  int GetDim() const { return Dim; }
83 
85  int GetGeomType() const { return GeomType; }
86 
88  int GetDof() const { return Dof; }
89 
91  int GetOrder() const { return Order; }
92 
94  int Space() const { return FuncSpace; }
95 
96  int GetRangeType() const { return RangeType; }
97 
98  int GetMapType() const { return MapType; }
99 
103  virtual void CalcShape(const IntegrationPoint &ip,
104  Vector &shape) const = 0;
105 
110  virtual void CalcDShape(const IntegrationPoint &ip,
111  DenseMatrix &dshape) const = 0;
112 
113  const IntegrationRule & GetNodes() const { return Nodes; }
114 
115  // virtual functions for finite elements on vector spaces
116 
121  virtual void CalcVShape(const IntegrationPoint &ip,
122  DenseMatrix &shape) const;
123 
124  virtual void CalcVShape(ElementTransformation &Trans,
125  DenseMatrix &shape) const;
126 
130  virtual void CalcDivShape(const IntegrationPoint &ip,
131  Vector &divshape) const;
132 
137  virtual void CalcCurlShape(const IntegrationPoint &ip,
138  DenseMatrix &curl_shape) const;
139 
140  virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const;
141 
144  virtual void CalcHessian (const IntegrationPoint &ip,
145  DenseMatrix &h) const;
146 
151  DenseMatrix &I) const;
152 
156  virtual void Project (Coefficient &coeff,
157  ElementTransformation &Trans, Vector &dofs) const;
158 
162  virtual void Project (VectorCoefficient &vc,
163  ElementTransformation &Trans, Vector &dofs) const;
164 
167  virtual void ProjectDelta(int vertex, Vector &dofs) const;
168 
172  virtual void Project(const FiniteElement &fe, ElementTransformation &Trans,
173  DenseMatrix &I) const;
174 
178  virtual void ProjectGrad(const FiniteElement &fe,
180  DenseMatrix &grad) const;
181 
185  virtual void ProjectCurl(const FiniteElement &fe,
187  DenseMatrix &curl) const;
188 
192  virtual void ProjectDiv(const FiniteElement &fe,
194  DenseMatrix &div) const;
195 
196  virtual ~FiniteElement () { }
197 };
198 
200 {
201 protected:
203  DenseMatrix &I,
204  const NodalFiniteElement &fine_fe) const;
205 
206  void ProjectCurl_2D(const FiniteElement &fe,
207  ElementTransformation &Trans,
208  DenseMatrix &curl) const;
209 
210 #ifndef MFEM_THREAD_SAFE
211  mutable Vector c_shape;
212 #endif
213 
214 public:
215  NodalFiniteElement(int D, int G, int Do, int O,
216  int F = FunctionSpace::Pk) :
217 #ifdef MFEM_THREAD_SAFE
218  FiniteElement(D, G, Do, O, F)
219 #else
220  FiniteElement(D, G, Do, O, F), c_shape(Do)
221 #endif
222  { }
223 
224  void SetMapType(int M)
225  {
226  MFEM_VERIFY(M == VALUE || M == INTEGRAL, "unknown MapType");
227  MapType = M;
228  }
229 
231  DenseMatrix &I) const
232  { NodalLocalInterpolation (Trans, I, *this); }
233 
234  virtual void Project (Coefficient &coeff,
235  ElementTransformation &Trans, Vector &dofs) const;
236 
237  virtual void Project (VectorCoefficient &vc,
238  ElementTransformation &Trans, Vector &dofs) const;
239 
240  virtual void Project(const FiniteElement &fe, ElementTransformation &Trans,
241  DenseMatrix &I) const;
242 
243  virtual void ProjectGrad(const FiniteElement &fe,
244  ElementTransformation &Trans,
245  DenseMatrix &grad) const;
246 
247  virtual void ProjectDiv(const FiniteElement &fe,
248  ElementTransformation &Trans,
249  DenseMatrix &div) const;
250 };
251 
252 
254 {
255 public:
256  PositiveFiniteElement(int D, int G, int Do, int O, int F = FunctionSpace::Pk)
257  : FiniteElement(D, G, Do, O, F) { }
259  virtual void Project(Coefficient &coeff,
260  ElementTransformation &Trans, Vector &dofs) const;
261  virtual void Project(const FiniteElement &fe, ElementTransformation &Trans,
262  DenseMatrix &I) const;
263 };
264 
266 {
267  // Hide the scalar functions CalcShape and CalcDShape.
268 private:
270  virtual void CalcShape(const IntegrationPoint &ip,
271  Vector &shape) const;
272 
274  virtual void CalcDShape(const IntegrationPoint &ip,
275  DenseMatrix &dshape) const;
276 
277 protected:
278 #ifndef MFEM_THREAD_SAFE
279  mutable DenseMatrix J, Jinv;
281 #endif
282 
284  DenseMatrix &shape) const;
285 
287  DenseMatrix &shape) const;
288 
289  void Project_RT(const double *nk, const Array<int> &d2n,
291  Vector &dofs) const;
292 
293  void Project_RT(const double *nk, const Array<int> &d2n,
295  DenseMatrix &I) const;
296 
297  // rotated gradient in 2D
298  void ProjectGrad_RT(const double *nk, const Array<int> &d2n,
300  DenseMatrix &grad) const;
301 
302  void ProjectCurl_ND(const double *tk, const Array<int> &d2t,
304  DenseMatrix &curl) const;
305 
306  void ProjectCurl_RT(const double *nk, const Array<int> &d2n,
308  DenseMatrix &curl) const;
309 
310  void Project_ND(const double *tk, const Array<int> &d2t,
312  Vector &dofs) const;
313 
314  void Project_ND(const double *tk, const Array<int> &d2t,
316  DenseMatrix &I) const;
317 
318  void ProjectGrad_ND(const double *tk, const Array<int> &d2t,
320  DenseMatrix &grad) const;
321 
322  void LocalInterpolation_RT(const double *nk, const Array<int> &d2n,
324  DenseMatrix &I) const;
325 
326  void LocalInterpolation_ND(const double *tk, const Array<int> &d2t,
328  DenseMatrix &I) const;
329 
330 public:
331  VectorFiniteElement (int D, int G, int Do, int O, int M,
332  int F = FunctionSpace::Pk) :
333 #ifdef MFEM_THREAD_SAFE
334  FiniteElement(D, G, Do, O, F)
335  { RangeType = VECTOR; MapType = M; }
336 #else
337  FiniteElement(D, G, Do, O, F), Jinv(D), vshape(Do, D)
338  { RangeType = VECTOR; MapType = M; }
339 #endif
340 };
341 
343 {
344 public:
346 
347  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
348 
349  virtual void CalcDShape(const IntegrationPoint &ip,
350  DenseMatrix &dshape) const;
351 };
352 
355 {
356 public:
359 
363  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
364 
369  virtual void CalcDShape(const IntegrationPoint &ip,
370  DenseMatrix &dshape) const;
371 };
372 
375 {
376 public:
379 
383  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
384 
389  virtual void CalcDShape(const IntegrationPoint &ip,
390  DenseMatrix &dshape) const;
391  virtual void ProjectDelta(int vertex, Vector &dofs) const
392  { dofs = 0.0; dofs(vertex) = 1.0; }
393 };
394 
397 {
398 public:
401 
405  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
406 
411  virtual void CalcDShape(const IntegrationPoint &ip,
412  DenseMatrix &dshape) const;
413  virtual void CalcHessian (const IntegrationPoint &ip,
414  DenseMatrix &h) const;
415  virtual void ProjectDelta(int vertex, Vector &dofs) const
416  { dofs = 0.0; dofs(vertex) = 1.0; } // { dofs = 1.0; }
417 };
418 
421 {
422 public:
424  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
425  virtual void CalcDShape(const IntegrationPoint &ip,
426  DenseMatrix &dshape) const;
427  virtual void ProjectDelta(int vertex, Vector &dofs) const;
428 };
429 
432 {
433 private:
434  static const double p[2];
435 
436 public:
438  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
439  virtual void CalcDShape(const IntegrationPoint &ip,
440  DenseMatrix &dshape) const;
441  virtual void ProjectDelta(int vertex, Vector &dofs) const;
442 };
443 
445 {
446 public:
448  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
449  virtual void CalcDShape(const IntegrationPoint &ip,
450  DenseMatrix &dshape) const;
451  virtual void ProjectDelta(int vertex, Vector &dofs) const
452  { dofs = 1.0; }
453 };
454 
457 {
458 public:
461 
465  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
466 
471  virtual void CalcDShape(const IntegrationPoint &ip,
472  DenseMatrix &dshape) const;
473 };
474 
476 {
477 public:
479  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
480  virtual void CalcDShape(const IntegrationPoint &ip,
481  DenseMatrix &dshape) const;
482 };
483 
486 {
487 public:
490 
494  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
495 
500  virtual void CalcDShape(const IntegrationPoint &ip,
501  DenseMatrix &dshape) const;
502 
503  virtual void CalcHessian (const IntegrationPoint &ip,
504  DenseMatrix &h) const;
505  virtual void ProjectDelta(int vertex, Vector &dofs) const;
506 };
507 
510 {
511 private:
512  static const double p[2];
513  DenseMatrix A;
514  mutable DenseMatrix D;
515  mutable Vector pol;
516 public:
518  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
519  virtual void CalcDShape(const IntegrationPoint &ip,
520  DenseMatrix &dshape) const;
521  // virtual void ProjectDelta(int vertex, Vector &dofs) const;
522 };
523 
526 {
527 public:
530 
534  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
535 
540  virtual void CalcDShape(const IntegrationPoint &ip,
541  DenseMatrix &dshape) const;
542  virtual void ProjectDelta(int vertex, Vector &dofs) const;
543 };
544 
546 {
547 public:
549  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
550  virtual void CalcDShape(const IntegrationPoint &ip,
551  DenseMatrix &dshape) const;
553  DenseMatrix &I) const;
555  virtual void Project(Coefficient &coeff, ElementTransformation &Trans,
556  Vector &dofs) const;
557  virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans,
558  Vector &dofs) const;
559  virtual void ProjectDelta(int vertex, Vector &dofs) const
560  { dofs = 0.; dofs(vertex) = 1.; }
561 };
562 
565 {
566 public:
568  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
569  virtual void CalcDShape(const IntegrationPoint &ip,
570  DenseMatrix &dshape) const;
571  // virtual void ProjectDelta(int vertex, Vector &dofs) const { dofs = 1.; }
572 };
573 
575 {
576 public:
578  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
579  virtual void CalcDShape(const IntegrationPoint &ip,
580  DenseMatrix &dshape) const;
581  virtual void CalcHessian (const IntegrationPoint &ip,
582  DenseMatrix &h) const;
583 };
584 
586 {
587 public:
589 
590  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
591 
592  virtual void CalcDShape(const IntegrationPoint &ip,
593  DenseMatrix &dshape) const;
594 };
595 
597 {
598 public:
600 
601  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
602 
603  virtual void CalcDShape(const IntegrationPoint &ip,
604  DenseMatrix &dshape) const;
605 
606  virtual void CalcHessian (const IntegrationPoint &ip,
607  DenseMatrix &h) const;
608 };
609 
612 {
613 public:
616 
617  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
618 
619  virtual void CalcDShape(const IntegrationPoint &ip,
620  DenseMatrix &dshape) const;
621 };
622 
625 {
626 public:
629 
631  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
632 
634  virtual void CalcDShape(const IntegrationPoint &ip,
635  DenseMatrix &dshape) const;
636  virtual void ProjectDelta(int vertex, Vector &dofs) const
637  { dofs(0) = 1.0; }
638 };
639 
640 
642 {
643 public:
645  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
646  virtual void CalcDShape(const IntegrationPoint &ip,
647  DenseMatrix &dshape) const;
648  virtual void ProjectDelta(int vertex, Vector &dofs) const
649  { dofs(0) = 1.0; }
650 };
651 
652 
655 {
656 public:
659 
663  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
664 
669  virtual void CalcDShape(const IntegrationPoint &ip,
670  DenseMatrix &dshape) const;
671 
672  virtual void ProjectDelta(int vertex, Vector &dofs) const
673  { dofs = 0.0; dofs(vertex) = 1.0; }
674 
675  virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const;
676 };
677 
680 {
681 public:
684 
685  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
686 
687  virtual void CalcDShape(const IntegrationPoint &ip,
688  DenseMatrix &dshape) const;
689 };
690 
693 {
694 public:
697 
701  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
702 
707  virtual void CalcDShape(const IntegrationPoint &ip,
708  DenseMatrix &dshape) const;
709 
710  virtual void ProjectDelta(int vertex, Vector &dofs) const
711  { dofs = 0.0; dofs(vertex) = 1.0; }
712 };
713 
716 {
717 public:
719  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
720  virtual void CalcDShape(const IntegrationPoint &ip,
721  DenseMatrix &dshape) const;
722  virtual void ProjectDelta(int vertex, Vector &dofs) const
723  { dofs = 1.0; }
724 };
725 
728 {
729 public:
731  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
732  virtual void CalcDShape(const IntegrationPoint &ip,
733  DenseMatrix &dshape) const;
734 };
735 
737 {
738 public:
739  P0SegmentFiniteElement(int Ord = 0);
740  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
741  virtual void CalcDShape(const IntegrationPoint &ip,
742  DenseMatrix &dshape) const;
743 };
744 
746 {
747 private:
748  static const double nk[3][2];
749 
750 public:
752 
753  virtual void CalcVShape(const IntegrationPoint &ip,
754  DenseMatrix &shape) const;
755 
757  DenseMatrix &shape) const
758  { CalcVShape_RT(Trans, shape); }
759 
760  virtual void CalcDivShape(const IntegrationPoint &ip,
761  Vector &divshape) const;
762 
764  DenseMatrix &I) const;
765 
767 
768  virtual void Project (VectorCoefficient &vc,
769  ElementTransformation &Trans, Vector &dofs) const;
770 };
771 
773 {
774 private:
775  static const double nk[4][2];
776 
777 public:
779 
780  virtual void CalcVShape(const IntegrationPoint &ip,
781  DenseMatrix &shape) const;
782 
784  DenseMatrix &shape) const
785  { CalcVShape_RT(Trans, shape); }
786 
787  virtual void CalcDivShape(const IntegrationPoint &ip,
788  Vector &divshape) const;
789 
791  DenseMatrix &I) const;
792 
794 
795  virtual void Project (VectorCoefficient &vc,
796  ElementTransformation &Trans, Vector &dofs) const;
797 };
798 
800 {
801 private:
802  static const double nk[8][2];
803 
804 public:
806 
807  virtual void CalcVShape(const IntegrationPoint &ip,
808  DenseMatrix &shape) const;
809 
811  DenseMatrix &shape) const
812  { CalcVShape_RT(Trans, shape); }
813 
814  virtual void CalcDivShape(const IntegrationPoint &ip,
815  Vector &divshape) const;
816 
818  DenseMatrix &I) const;
819 
821 
822  virtual void Project (VectorCoefficient &vc,
823  ElementTransformation &Trans, Vector &dofs) const;
824 };
825 
827 {
828 private:
829  static const double nk[12][2];
830 
831 public:
833 
834  virtual void CalcVShape(const IntegrationPoint &ip,
835  DenseMatrix &shape) const;
836 
838  DenseMatrix &shape) const
839  { CalcVShape_RT(Trans, shape); }
840 
841  virtual void CalcDivShape(const IntegrationPoint &ip,
842  Vector &divshape) const;
843 
845  DenseMatrix &I) const;
846 
848 
849  virtual void Project (VectorCoefficient &vc,
850  ElementTransformation &Trans, Vector &dofs) const;
851 };
852 
854 {
855 private:
856  static const double M[15][15];
857 public:
859 
860  virtual void CalcVShape(const IntegrationPoint &ip,
861  DenseMatrix &shape) const;
862 
864  DenseMatrix &shape) const
865  { CalcVShape_RT(Trans, shape); }
866 
867  virtual void CalcDivShape(const IntegrationPoint &ip,
868  Vector &divshape) const;
869 };
870 
872 {
873 private:
874  static const double nk[24][2];
875  static const double pt[4];
876  static const double dpt[3];
877 
878 public:
880 
881  virtual void CalcVShape(const IntegrationPoint &ip,
882  DenseMatrix &shape) const;
883 
885  DenseMatrix &shape) const
886  { CalcVShape_RT(Trans, shape); }
887 
888  virtual void CalcDivShape(const IntegrationPoint &ip,
889  Vector &divshape) const;
890 
892  DenseMatrix &I) const;
893 
895 
896  virtual void Project (VectorCoefficient &vc,
897  ElementTransformation &Trans, Vector &dofs) const;
898 };
899 
902 {
903 public:
905  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
906  virtual void CalcDShape(const IntegrationPoint &ip,
907  DenseMatrix &dshape) const;
908 };
909 
912 {
913 public:
915  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
916  virtual void CalcDShape(const IntegrationPoint &ip,
917  DenseMatrix &dshape) const;
918 };
919 
921 {
922 private:
923  Vector rwk;
924 #ifndef MFEM_THREAD_SAFE
925  mutable Vector rxxk;
926 #endif
927 public:
928  Lagrange1DFiniteElement (int degree);
929  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
930  virtual void CalcDShape(const IntegrationPoint &ip,
931  DenseMatrix &dshape) const;
932 };
933 
935 {
936 public:
938  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
939  virtual void CalcDShape(const IntegrationPoint &ip,
940  DenseMatrix &dshape) const;
941 };
942 
944 {
945 public:
947  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
948  virtual void CalcDShape(const IntegrationPoint &ip,
949  DenseMatrix &dshape) const;
950  virtual void ProjectDelta(int vertex, Vector &dofs) const
951  { dofs(0) = 1.0; }
952 };
953 
955 {
956 public:
958  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
959  virtual void CalcDShape(const IntegrationPoint &ip,
960  DenseMatrix &dshape) const;
961  virtual void ProjectDelta(int vertex, Vector &dofs) const
962  { dofs(0) = 1.0; }
963 };
964 
967 {
968 private:
970  int dof1d;
971  int *I, *J, *K;
972 #ifndef MFEM_THREAD_SAFE
973  mutable Vector shape1dx, shape1dy, shape1dz;
974  mutable DenseMatrix dshape1dx, dshape1dy, dshape1dz;
975 #endif
976 
977 public:
978  LagrangeHexFiniteElement (int degree);
979  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
980  virtual void CalcDShape(const IntegrationPoint &ip,
981  DenseMatrix &dshape) const;
983 };
984 
985 
988 {
989 public:
992 
996  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
997 
1002  virtual void CalcDShape(const IntegrationPoint &ip,
1003  DenseMatrix &dshape) const;
1004 };
1005 
1008 {
1009 public:
1012 
1016  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1017 
1022  virtual void CalcDShape(const IntegrationPoint &ip,
1023  DenseMatrix &dshape) const;
1024 };
1025 
1028 {
1029 public:
1032 
1033  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1034 
1035  virtual void CalcDShape(const IntegrationPoint &ip,
1036  DenseMatrix &dshape) const;
1037 };
1038 
1041 {
1042 public:
1045 
1049  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1050 
1055  virtual void CalcDShape(const IntegrationPoint &ip,
1056  DenseMatrix &dshape) const;
1057 };
1058 
1061 {
1062 public:
1065 
1069  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1070 
1075  virtual void CalcDShape(const IntegrationPoint &ip,
1076  DenseMatrix &dshape) const;
1077 };
1078 
1079 
1081 {
1082 private:
1083  static const double tk[12][3];
1084 
1085 public:
1087  virtual void CalcVShape(const IntegrationPoint &ip,
1088  DenseMatrix &shape) const;
1090  DenseMatrix &shape) const
1091  { CalcVShape_ND(Trans, shape); }
1092  virtual void CalcCurlShape(const IntegrationPoint &ip,
1093  DenseMatrix &curl_shape) const;
1095  DenseMatrix &I) const;
1096  using FiniteElement::Project;
1097  virtual void Project (VectorCoefficient &vc,
1098  ElementTransformation &Trans, Vector &dofs) const;
1099 };
1100 
1101 
1103 {
1104 private:
1105  static const double tk[6][3];
1106 
1107 public:
1109  virtual void CalcVShape(const IntegrationPoint &ip,
1110  DenseMatrix &shape) const;
1112  DenseMatrix &shape) const
1113  { CalcVShape_ND(Trans, shape); }
1114  virtual void CalcCurlShape(const IntegrationPoint &ip,
1115  DenseMatrix &curl_shape) const;
1117  DenseMatrix &I) const;
1118  using FiniteElement::Project;
1119  virtual void Project (VectorCoefficient &vc,
1120  ElementTransformation &Trans, Vector &dofs) const;
1121 };
1122 
1123 
1125 {
1126 private:
1127  static const double nk[6][3];
1128 
1129 public:
1131 
1132  virtual void CalcVShape(const IntegrationPoint &ip,
1133  DenseMatrix &shape) const;
1134 
1136  DenseMatrix &shape) const
1137  { CalcVShape_RT(Trans, shape); }
1138 
1139  virtual void CalcDivShape(const IntegrationPoint &ip,
1140  Vector &divshape) const;
1141 
1143  DenseMatrix &I) const;
1144 
1145  using FiniteElement::Project;
1146 
1147  virtual void Project (VectorCoefficient &vc,
1148  ElementTransformation &Trans, Vector &dofs) const;
1149 };
1150 
1151 
1153 {
1154 private:
1155  static const double nk[36][3];
1156 
1157 public:
1159 
1160  virtual void CalcVShape(const IntegrationPoint &ip,
1161  DenseMatrix &shape) const;
1162 
1164  DenseMatrix &shape) const
1165  { CalcVShape_RT(Trans, shape); }
1166 
1167  virtual void CalcDivShape(const IntegrationPoint &ip,
1168  Vector &divshape) const;
1169 
1171  DenseMatrix &I) const;
1172 
1173  using FiniteElement::Project;
1174 
1175  virtual void Project (VectorCoefficient &vc,
1176  ElementTransformation &Trans, Vector &dofs) const;
1177 };
1178 
1179 
1181 {
1182 private:
1183  static const double nk[4][3];
1184 
1185 public:
1187 
1188  virtual void CalcVShape(const IntegrationPoint &ip,
1189  DenseMatrix &shape) const;
1190 
1192  DenseMatrix &shape) const
1193  { CalcVShape_RT(Trans, shape); }
1194 
1195  virtual void CalcDivShape(const IntegrationPoint &ip,
1196  Vector &divshape) const;
1197 
1199  DenseMatrix &I) const;
1200 
1201  using FiniteElement::Project;
1202 
1203  virtual void Project (VectorCoefficient &vc,
1204  ElementTransformation &Trans, Vector &dofs) const;
1205 };
1206 
1207 
1209 {
1210 public:
1212  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1213  virtual void CalcDShape(const IntegrationPoint &ip,
1214  DenseMatrix &dshape) const;
1215 };
1216 
1217 
1218 class Poly_1D
1219 {
1220 public:
1221  class Basis
1222  {
1223  private:
1224  int mode; // 0 - use change of basis, O(p^2) Evals
1225  // 1 - use barycentric Lagrangian interpolation, O(p) Evals
1226  DenseMatrixInverse Ai;
1227  mutable Vector x, w;
1228 
1229  public:
1230  Basis(const int p, const double *nodes, const int _mode = 1);
1231  void Eval(const double x, Vector &u) const;
1232  void Eval(const double x, Vector &u, Vector &d) const;
1233  };
1234 
1235 private:
1236  Array<double *> open_pts, closed_pts;
1237  Array<Basis *> open_basis, closed_basis;
1238 
1239  static Array2D<int> binom;
1240 
1241  static void CalcMono(const int p, const double x, double *u);
1242  static void CalcMono(const int p, const double x, double *u, double *d);
1243 
1244  static void CalcLegendre(const int p, const double x, double *u);
1245  static void CalcLegendre(const int p, const double x, double *u, double *d);
1246 
1247  static void CalcChebyshev(const int p, const double x, double *u);
1248  static void CalcChebyshev(const int p, const double x, double *u, double *d);
1249 
1250 public:
1251  Poly_1D() { }
1252 
1253  static const int *Binom(const int p);
1254 
1255  const double *OpenPoints(const int p);
1256  const double *ClosedPoints(const int p);
1257 
1258  Basis &OpenBasis(const int p);
1259  Basis &ClosedBasis(const int p);
1260 
1261  // Evaluate the values of a hierarchical 1D basis at point x
1262  // hierarchical = k-th basis function is degree k polynomial
1263  static void CalcBasis(const int p, const double x, double *u)
1264  // { CalcMono(p, x, u); }
1265  // Bernstein basis is not hierarchical --> does not work for triangles
1266  // and tetrahedra
1267  // { CalcBernstein(p, x, u); }
1268  // { CalcLegendre(p, x, u); }
1269  { CalcChebyshev(p, x, u); }
1270 
1271  // Evaluate the values and derivatives of a hierarchical 1D basis at point x
1272  static void CalcBasis(const int p, const double x, double *u, double *d)
1273  // { CalcMono(p, x, u, d); }
1274  // { CalcBernstein(p, x, u, d); }
1275  // { CalcLegendre(p, x, u, d); }
1276  { CalcChebyshev(p, x, u, d); }
1277 
1278  // Evaluate a representation of a Delta function at point x
1279  static double CalcDelta(const int p, const double x)
1280  { return pow(x, (double) p); }
1281 
1282  static void UniformPoints(const int p, double *x);
1283  static void GaussPoints(const int p, double *x);
1284  static void GaussLobattoPoints(const int p, double *x);
1285  static void ChebyshevPoints(const int p, double *x);
1286 
1288  static void CalcBinomTerms(const int p, const double x, const double y,
1289  double *u);
1292  static void CalcBinomTerms(const int p, const double x, const double y,
1293  double *u, double *d);
1296  static void CalcDBinomTerms(const int p, const double x, const double y,
1297  double *d);
1298  static void CalcBernstein(const int p, const double x, double *u)
1299  { CalcBinomTerms(p, x, 1. - x, u); }
1300  static void CalcBernstein(const int p, const double x, double *u, double *d)
1301  { CalcBinomTerms(p, x, 1. - x, u, d); }
1302 
1303  ~Poly_1D();
1304 };
1305 
1306 extern Poly_1D poly1d;
1307 
1308 
1310 {
1311 private:
1312  Poly_1D::Basis &basis1d;
1313 #ifndef MFEM_THREAD_SAFE
1314  mutable Vector shape_x, dshape_x;
1315 #endif
1316  Array<int> dof_map;
1317 
1318 public:
1319  H1_SegmentElement(const int p);
1320  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1321  virtual void CalcDShape(const IntegrationPoint &ip,
1322  DenseMatrix &dshape) const;
1323  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1324  const Array<int> &GetDofMap() const { return dof_map; }
1325 };
1326 
1327 
1329 {
1330 private:
1331  Poly_1D::Basis &basis1d;
1332 #ifndef MFEM_THREAD_SAFE
1333  mutable Vector shape_x, shape_y, dshape_x, dshape_y;
1334 #endif
1335  Array<int> dof_map;
1336 
1337 public:
1338  H1_QuadrilateralElement(const int p);
1339  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1340  virtual void CalcDShape(const IntegrationPoint &ip,
1341  DenseMatrix &dshape) const;
1342  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1343  const Array<int> &GetDofMap() const { return dof_map; }
1344 };
1345 
1346 
1348 {
1349 private:
1350  Poly_1D::Basis &basis1d;
1351 #ifndef MFEM_THREAD_SAFE
1352  mutable Vector shape_x, shape_y, shape_z, dshape_x, dshape_y, dshape_z;
1353 #endif
1354  Array<int> dof_map;
1355 
1356 public:
1357  H1_HexahedronElement(const int p);
1358  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1359  virtual void CalcDShape(const IntegrationPoint &ip,
1360  DenseMatrix &dshape) const;
1361  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1362  const Array<int> &GetDofMap() const { return dof_map; }
1363 };
1364 
1366 {
1367 private:
1368 #ifndef MFEM_THREAD_SAFE
1369  // This is to share scratch space between invocations, which helps
1370  // speed things up, but with OpenMP, we need one copy per thread.
1371  // Right now, we solve this by allocating this space within each function
1372  // call every time we call it. Alternatively, we should do some sort
1373  // thread private thing. Brunner, Jan 2014
1374  mutable Vector shape_x, dshape_x;
1375 #endif
1376  Array<int> dof_map;
1377 
1378 public:
1379  H1Pos_SegmentElement(const int p);
1380  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1381  virtual void CalcDShape(const IntegrationPoint &ip,
1382  DenseMatrix &dshape) const;
1383  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1384  const Array<int> &GetDofMap() const { return dof_map; }
1385 };
1386 
1387 
1389 {
1390 private:
1391 #ifndef MFEM_THREAD_SAFE
1392  // See comment in H1Pos_SegmentElement
1393  mutable Vector shape_x, shape_y, dshape_x, dshape_y;
1394 #endif
1395  Array<int> dof_map;
1396 
1397 public:
1398  H1Pos_QuadrilateralElement(const int p);
1399  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1400  virtual void CalcDShape(const IntegrationPoint &ip,
1401  DenseMatrix &dshape) const;
1402  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1403  const Array<int> &GetDofMap() const { return dof_map; }
1404 };
1405 
1406 
1408 {
1409 private:
1410 #ifndef MFEM_THREAD_SAFE
1411  // See comment in H1Pos_SegementElement.
1412  mutable Vector shape_x, shape_y, shape_z, dshape_x, dshape_y, dshape_z;
1413 #endif
1414  Array<int> dof_map;
1415 
1416 public:
1417  H1Pos_HexahedronElement(const int p);
1418  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1419  virtual void CalcDShape(const IntegrationPoint &ip,
1420  DenseMatrix &dshape) const;
1421  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1422  const Array<int> &GetDofMap() const { return dof_map; }
1423 };
1424 
1425 
1427 {
1428 private:
1429 #ifndef MFEM_THREAD_SAFE
1430  mutable Vector shape_x, shape_y, shape_l, dshape_x, dshape_y, dshape_l, u;
1431  mutable DenseMatrix du;
1432 #endif
1433  DenseMatrixInverse Ti;
1434 
1435 public:
1436  H1_TriangleElement(const int p);
1437  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1438  virtual void CalcDShape(const IntegrationPoint &ip,
1439  DenseMatrix &dshape) const;
1440 };
1441 
1442 
1444 {
1445 private:
1446 #ifndef MFEM_THREAD_SAFE
1447  mutable Vector shape_x, shape_y, shape_z, shape_l;
1448  mutable Vector dshape_x, dshape_y, dshape_z, dshape_l, u;
1449  mutable DenseMatrix du;
1450 #endif
1451  DenseMatrixInverse Ti;
1452 
1453 public:
1454  H1_TetrahedronElement(const int p);
1455  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1456  virtual void CalcDShape(const IntegrationPoint &ip,
1457  DenseMatrix &dshape) const;
1458 };
1459 
1460 
1462 {
1463 protected:
1464 #ifndef MFEM_THREAD_SAFE
1467 #endif
1469 
1470 public:
1471  H1Pos_TriangleElement(const int p);
1472 
1473  // The size of shape is (p+1)(p+2)/2 (dof).
1474  static void CalcShape(const int p, const double x, const double y,
1475  double *shape);
1476 
1477  // The size of dshape_1d is p+1; the size of dshape is (dof x dim).
1478  static void CalcDShape(const int p, const double x, const double y,
1479  double *dshape_1d, double *dshape);
1480 
1481  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1482  virtual void CalcDShape(const IntegrationPoint &ip,
1483  DenseMatrix &dshape) const;
1484 };
1485 
1486 
1488 {
1489 protected:
1490 #ifndef MFEM_THREAD_SAFE
1493 #endif
1495 
1496 public:
1497  H1Pos_TetrahedronElement(const int p);
1498 
1499  // The size of shape is (p+1)(p+2)(p+3)/6 (dof).
1500  static void CalcShape(const int p, const double x, const double y,
1501  const double z, double *shape);
1502 
1503  // The size of dshape_1d is p+1; the size of dshape is (dof x dim).
1504  static void CalcDShape(const int p, const double x, const double y,
1505  const double z, double *dshape_1d, double *dshape);
1506 
1507  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1508  virtual void CalcDShape(const IntegrationPoint &ip,
1509  DenseMatrix &dshape) const;
1510 };
1511 
1512 
1514 {
1515 private:
1516  int type;
1517  Poly_1D::Basis *basis1d;
1518 #ifndef MFEM_THREAD_SAFE
1519  mutable Vector shape_x, dshape_x;
1520 #endif
1521 
1522 public:
1523  L2_SegmentElement(const int p, const int _type = 0);
1524  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1525  virtual void CalcDShape(const IntegrationPoint &ip,
1526  DenseMatrix &dshape) const;
1527  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1528 };
1529 
1530 
1532 {
1533 private:
1534 #ifndef MFEM_THREAD_SAFE
1535  mutable Vector shape_x, dshape_x;
1536 #endif
1537 
1538 public:
1539  L2Pos_SegmentElement(const int p);
1540  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1541  virtual void CalcDShape(const IntegrationPoint &ip,
1542  DenseMatrix &dshape) const;
1543  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1544 };
1545 
1546 
1548 {
1549 private:
1550  int type;
1551  Poly_1D::Basis *basis1d;
1552 #ifndef MFEM_THREAD_SAFE
1553  mutable Vector shape_x, shape_y, dshape_x, dshape_y;
1554 #endif
1555 
1556 public:
1557  L2_QuadrilateralElement(const int p, const int _type = 0);
1558  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1559  virtual void CalcDShape(const IntegrationPoint &ip,
1560  DenseMatrix &dshape) const;
1561  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1562  virtual void ProjectCurl(const FiniteElement &fe,
1564  DenseMatrix &curl) const
1565  { ProjectCurl_2D(fe, Trans, curl); }
1566 };
1567 
1568 
1570 {
1571 private:
1572 #ifndef MFEM_THREAD_SAFE
1573  mutable Vector shape_x, shape_y, dshape_x, dshape_y;
1574 #endif
1575 
1576 public:
1577  L2Pos_QuadrilateralElement(const int p);
1578  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1579  virtual void CalcDShape(const IntegrationPoint &ip,
1580  DenseMatrix &dshape) const;
1581  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1582 };
1583 
1584 
1586 {
1587 private:
1588  int type;
1589  Poly_1D::Basis *basis1d;
1590 #ifndef MFEM_THREAD_SAFE
1591  mutable Vector shape_x, shape_y, shape_z, dshape_x, dshape_y, dshape_z;
1592 #endif
1593 
1594 public:
1595  L2_HexahedronElement(const int p, const int _type = 0);
1596  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1597  virtual void CalcDShape(const IntegrationPoint &ip,
1598  DenseMatrix &dshape) const;
1599  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1600 };
1601 
1602 
1604 {
1605 private:
1606 #ifndef MFEM_THREAD_SAFE
1607  mutable Vector shape_x, shape_y, shape_z, dshape_x, dshape_y, dshape_z;
1608 #endif
1609 
1610 public:
1611  L2Pos_HexahedronElement(const int p);
1612  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1613  virtual void CalcDShape(const IntegrationPoint &ip,
1614  DenseMatrix &dshape) const;
1615  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1616 };
1617 
1618 
1620 {
1621 private:
1622  int type;
1623 #ifndef MFEM_THREAD_SAFE
1624  mutable Vector shape_x, shape_y, shape_l, dshape_x, dshape_y, dshape_l, u;
1625  mutable DenseMatrix du;
1626 #endif
1627  DenseMatrixInverse Ti;
1628 
1629 public:
1630  L2_TriangleElement(const int p, const int _type = 0);
1631  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1632  virtual void CalcDShape(const IntegrationPoint &ip,
1633  DenseMatrix &dshape) const;
1634  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1635  virtual void ProjectCurl(const FiniteElement &fe,
1637  DenseMatrix &curl) const
1638  { ProjectCurl_2D(fe, Trans, curl); }
1639 };
1640 
1641 
1643 {
1644 private:
1645 #ifndef MFEM_THREAD_SAFE
1646  mutable Vector dshape_1d;
1647 #endif
1648 
1649 public:
1650  L2Pos_TriangleElement(const int p);
1651  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1652  virtual void CalcDShape(const IntegrationPoint &ip,
1653  DenseMatrix &dshape) const;
1654  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1655 };
1656 
1657 
1659 {
1660 private:
1661  int type;
1662 #ifndef MFEM_THREAD_SAFE
1663  mutable Vector shape_x, shape_y, shape_z, shape_l;
1664  mutable Vector dshape_x, dshape_y, dshape_z, dshape_l, u;
1665  mutable DenseMatrix du;
1666 #endif
1667  DenseMatrixInverse Ti;
1668 
1669 public:
1670  L2_TetrahedronElement(const int p, const int _type = 0);
1671  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1672  virtual void CalcDShape(const IntegrationPoint &ip,
1673  DenseMatrix &dshape) const;
1674  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1675 };
1676 
1677 
1679 {
1680 private:
1681 #ifndef MFEM_THREAD_SAFE
1682  mutable Vector dshape_1d;
1683 #endif
1684 
1685 public:
1686  L2Pos_TetrahedronElement(const int p);
1687  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1688  virtual void CalcDShape(const IntegrationPoint &ip,
1689  DenseMatrix &dshape) const;
1690  virtual void ProjectDelta(int vertex, Vector &dofs) const;
1691 };
1692 
1693 
1695 {
1696 private:
1697  static const double nk[8];
1698 
1699  Poly_1D::Basis &cbasis1d, &obasis1d;
1700 #ifndef MFEM_THREAD_SAFE
1701  mutable Vector shape_cx, shape_ox, shape_cy, shape_oy;
1702  mutable Vector dshape_cx, dshape_cy;
1703 #endif
1704  Array<int> dof_map, dof2nk;
1705 
1706 public:
1707  RT_QuadrilateralElement(const int p);
1708  virtual void CalcVShape(const IntegrationPoint &ip,
1709  DenseMatrix &shape) const;
1711  DenseMatrix &shape) const
1712  { CalcVShape_RT(Trans, shape); }
1713  virtual void CalcDivShape(const IntegrationPoint &ip,
1714  Vector &divshape) const;
1716  DenseMatrix &I) const
1717  { LocalInterpolation_RT(nk, dof2nk, Trans, I); }
1718  using FiniteElement::Project;
1719  virtual void Project(VectorCoefficient &vc,
1720  ElementTransformation &Trans, Vector &dofs) const
1721  { Project_RT(nk, dof2nk, vc, Trans, dofs); }
1723  DenseMatrix &I) const
1724  { Project_RT(nk, dof2nk, fe, Trans, I); }
1725  virtual void ProjectGrad(const FiniteElement &fe,
1727  DenseMatrix &grad) const
1728  { ProjectGrad_RT(nk, dof2nk, fe, Trans, grad); }
1729 };
1730 
1731 
1733 {
1734  static const double nk[18];
1735 
1736  Poly_1D::Basis &cbasis1d, &obasis1d;
1737 #ifndef MFEM_THREAD_SAFE
1738  mutable Vector shape_cx, shape_ox, shape_cy, shape_oy, shape_cz, shape_oz;
1739  mutable Vector dshape_cx, dshape_cy, dshape_cz;
1740 #endif
1741  Array<int> dof_map, dof2nk;
1742 
1743 public:
1744  RT_HexahedronElement(const int p);
1745  virtual void CalcVShape(const IntegrationPoint &ip,
1746  DenseMatrix &shape) const;
1748  DenseMatrix &shape) const
1749  { CalcVShape_RT(Trans, shape); }
1750  virtual void CalcDivShape(const IntegrationPoint &ip,
1751  Vector &divshape) const;
1753  DenseMatrix &I) const
1754  { LocalInterpolation_RT(nk, dof2nk, Trans, I); }
1755  using FiniteElement::Project;
1756  virtual void Project(VectorCoefficient &vc,
1757  ElementTransformation &Trans, Vector &dofs) const
1758  { Project_RT(nk, dof2nk, vc, Trans, dofs); }
1760  DenseMatrix &I) const
1761  { Project_RT(nk, dof2nk, fe, Trans, I); }
1762  virtual void ProjectCurl(const FiniteElement &fe,
1764  DenseMatrix &curl) const
1765  { ProjectCurl_RT(nk, dof2nk, fe, Trans, curl); }
1766 };
1767 
1768 
1770 {
1771  static const double nk[6], c;
1772 
1773 #ifndef MFEM_THREAD_SAFE
1774  mutable Vector shape_x, shape_y, shape_l;
1775  mutable Vector dshape_x, dshape_y, dshape_l;
1776  mutable DenseMatrix u;
1777  mutable Vector divu;
1778 #endif
1779  Array<int> dof2nk;
1780  DenseMatrixInverse Ti;
1781 
1782 public:
1783  RT_TriangleElement(const int p);
1784  virtual void CalcVShape(const IntegrationPoint &ip,
1785  DenseMatrix &shape) const;
1787  DenseMatrix &shape) const
1788  { CalcVShape_RT(Trans, shape); }
1789  virtual void CalcDivShape(const IntegrationPoint &ip,
1790  Vector &divshape) const;
1792  DenseMatrix &I) const
1793  { LocalInterpolation_RT(nk, dof2nk, Trans, I); }
1794  using FiniteElement::Project;
1795  virtual void Project(VectorCoefficient &vc,
1796  ElementTransformation &Trans, Vector &dofs) const
1797  { Project_RT(nk, dof2nk, vc, Trans, dofs); }
1799  DenseMatrix &I) const
1800  { Project_RT(nk, dof2nk, fe, Trans, I); }
1801  virtual void ProjectGrad(const FiniteElement &fe,
1803  DenseMatrix &grad) const
1804  { ProjectGrad_RT(nk, dof2nk, fe, Trans, grad); }
1805 };
1806 
1807 
1809 {
1810  static const double nk[12], c;
1811 
1812 #ifndef MFEM_THREAD_SAFE
1813  mutable Vector shape_x, shape_y, shape_z, shape_l;
1814  mutable Vector dshape_x, dshape_y, dshape_z, dshape_l;
1815  mutable DenseMatrix u;
1816  mutable Vector divu;
1817 #endif
1818  Array<int> dof2nk;
1819  DenseMatrixInverse Ti;
1820 
1821 public:
1822  RT_TetrahedronElement(const int p);
1823  virtual void CalcVShape(const IntegrationPoint &ip,
1824  DenseMatrix &shape) const;
1826  DenseMatrix &shape) const
1827  { CalcVShape_RT(Trans, shape); }
1828  virtual void CalcDivShape(const IntegrationPoint &ip,
1829  Vector &divshape) const;
1831  DenseMatrix &I) const
1832  { LocalInterpolation_RT(nk, dof2nk, Trans, I); }
1833  using FiniteElement::Project;
1834  virtual void Project(VectorCoefficient &vc,
1835  ElementTransformation &Trans, Vector &dofs) const
1836  { Project_RT(nk, dof2nk, vc, Trans, dofs); }
1838  DenseMatrix &I) const
1839  { Project_RT(nk, dof2nk, fe, Trans, I); }
1840  virtual void ProjectCurl(const FiniteElement &fe,
1842  DenseMatrix &curl) const
1843  { ProjectCurl_RT(nk, dof2nk, fe, Trans, curl); }
1844 };
1845 
1846 
1848 {
1849  static const double tk[18];
1850 
1851  Poly_1D::Basis &cbasis1d, &obasis1d;
1852 #ifndef MFEM_THREAD_SAFE
1853  mutable Vector shape_cx, shape_ox, shape_cy, shape_oy, shape_cz, shape_oz;
1854  mutable Vector dshape_cx, dshape_cy, dshape_cz;
1855 #endif
1856  Array<int> dof_map, dof2tk;
1857 
1858 public:
1859  ND_HexahedronElement(const int p);
1860 
1861  virtual void CalcVShape(const IntegrationPoint &ip,
1862  DenseMatrix &shape) const;
1863 
1865  DenseMatrix &shape) const
1866  { CalcVShape_ND(Trans, shape); }
1867 
1868  virtual void CalcCurlShape(const IntegrationPoint &ip,
1869  DenseMatrix &curl_shape) const;
1870 
1872  DenseMatrix &I) const
1873  { LocalInterpolation_ND(tk, dof2tk, Trans, I); }
1874 
1875  using FiniteElement::Project;
1876 
1877  virtual void Project(VectorCoefficient &vc,
1878  ElementTransformation &Trans, Vector &dofs) const
1879  { Project_ND(tk, dof2tk, vc, Trans, dofs); }
1880 
1881  virtual void Project(const FiniteElement &fe,
1883  DenseMatrix &I) const
1884  { Project_ND(tk, dof2tk, fe, Trans, I); }
1885 
1886  virtual void ProjectGrad(const FiniteElement &fe,
1888  DenseMatrix &grad) const
1889  { ProjectGrad_ND(tk, dof2tk, fe, Trans, grad); }
1890 
1891  virtual void ProjectCurl(const FiniteElement &fe,
1893  DenseMatrix &curl) const
1894  { ProjectCurl_ND(tk, dof2tk, fe, Trans, curl); }
1895 };
1896 
1897 
1899 {
1900  static const double tk[8];
1901 
1902  Poly_1D::Basis &cbasis1d, &obasis1d;
1903 #ifndef MFEM_THREAD_SAFE
1904  mutable Vector shape_cx, shape_ox, shape_cy, shape_oy;
1905  mutable Vector dshape_cx, dshape_cy;
1906 #endif
1907  Array<int> dof_map, dof2tk;
1908 
1909 public:
1910  ND_QuadrilateralElement(const int p);
1911  virtual void CalcVShape(const IntegrationPoint &ip,
1912  DenseMatrix &shape) const;
1914  DenseMatrix &shape) const
1915  { CalcVShape_ND(Trans, shape); }
1916  virtual void CalcCurlShape(const IntegrationPoint &ip,
1917  DenseMatrix &curl_shape) const;
1919  DenseMatrix &I) const
1920  { LocalInterpolation_ND(tk, dof2tk, Trans, I); }
1921  using FiniteElement::Project;
1922  virtual void Project(VectorCoefficient &vc,
1923  ElementTransformation &Trans, Vector &dofs) const
1924  { Project_ND(tk, dof2tk, vc, Trans, dofs); }
1925  virtual void Project(const FiniteElement &fe,
1927  DenseMatrix &I) const
1928  { Project_ND(tk, dof2tk, fe, Trans, I); }
1929  virtual void ProjectGrad(const FiniteElement &fe,
1931  DenseMatrix &grad) const
1932  { ProjectGrad_ND(tk, dof2tk, fe, Trans, grad); }
1933 };
1934 
1935 
1937 {
1938  static const double tk[18], c;
1939 
1940 #ifndef MFEM_THREAD_SAFE
1941  mutable Vector shape_x, shape_y, shape_z, shape_l;
1942  mutable Vector dshape_x, dshape_y, dshape_z, dshape_l;
1943  mutable DenseMatrix u;
1944 #endif
1945  Array<int> dof2tk;
1946  DenseMatrixInverse Ti;
1947 
1948 public:
1949  ND_TetrahedronElement(const int p);
1950  virtual void CalcVShape(const IntegrationPoint &ip,
1951  DenseMatrix &shape) const;
1953  DenseMatrix &shape) const
1954  { CalcVShape_ND(Trans, shape); }
1955  virtual void CalcCurlShape(const IntegrationPoint &ip,
1956  DenseMatrix &curl_shape) const;
1958  DenseMatrix &I) const
1959  { LocalInterpolation_ND(tk, dof2tk, Trans, I); }
1960  using FiniteElement::Project;
1961  virtual void Project(VectorCoefficient &vc,
1962  ElementTransformation &Trans, Vector &dofs) const
1963  { Project_ND(tk, dof2tk, vc, Trans, dofs); }
1964  virtual void Project(const FiniteElement &fe,
1966  DenseMatrix &I) const
1967  { Project_ND(tk, dof2tk, fe, Trans, I); }
1968  virtual void ProjectGrad(const FiniteElement &fe,
1970  DenseMatrix &grad) const
1971  { ProjectGrad_ND(tk, dof2tk, fe, Trans, grad); }
1972 
1973  virtual void ProjectCurl(const FiniteElement &fe,
1975  DenseMatrix &curl) const
1976  { ProjectCurl_ND(tk, dof2tk, fe, Trans, curl); }
1977 };
1978 
1980 {
1981  static const double tk[8], c;
1982 
1983 #ifndef MFEM_THREAD_SAFE
1984  mutable Vector shape_x, shape_y, shape_l;
1985  mutable Vector dshape_x, dshape_y, dshape_l;
1986  mutable DenseMatrix u;
1987  mutable Vector curlu;
1988 #endif
1989  Array<int> dof2tk;
1990  DenseMatrixInverse Ti;
1991 
1992 public:
1993  ND_TriangleElement(const int p);
1994  virtual void CalcVShape(const IntegrationPoint &ip,
1995  DenseMatrix &shape) const;
1997  DenseMatrix &shape) const
1998  { CalcVShape_ND(Trans, shape); }
1999  virtual void CalcCurlShape(const IntegrationPoint &ip,
2000  DenseMatrix &curl_shape) const;
2002  DenseMatrix &I) const
2003  { LocalInterpolation_ND(tk, dof2tk, Trans, I); }
2004  using FiniteElement::Project;
2005  virtual void Project(VectorCoefficient &vc,
2006  ElementTransformation &Trans, Vector &dofs) const
2007  { Project_ND(tk, dof2tk, vc, Trans, dofs); }
2008  virtual void Project(const FiniteElement &fe,
2010  DenseMatrix &I) const
2011  { Project_ND(tk, dof2tk, fe, Trans, I); }
2012  virtual void ProjectGrad(const FiniteElement &fe,
2014  DenseMatrix &grad) const
2015  { ProjectGrad_ND(tk, dof2tk, fe, Trans, grad); }
2016 };
2017 
2018 
2020 {
2021  static const double tk[1];
2022 
2023  Poly_1D::Basis &obasis1d;
2024  Array<int> dof2tk;
2025 
2026 public:
2027  ND_SegmentElement(const int p);
2028  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
2029  { obasis1d.Eval(ip.x, shape); }
2030  virtual void CalcVShape(const IntegrationPoint &ip,
2031  DenseMatrix &shape) const;
2033  DenseMatrix &shape) const
2034  { CalcVShape_ND(Trans, shape); }
2035  // virtual void CalcCurlShape(const IntegrationPoint &ip,
2036  // DenseMatrix &curl_shape) const;
2038  DenseMatrix &I) const
2039  { LocalInterpolation_ND(tk, dof2tk, Trans, I); }
2040  using FiniteElement::Project;
2041  virtual void Project(VectorCoefficient &vc,
2042  ElementTransformation &Trans, Vector &dofs) const
2043  { Project_ND(tk, dof2tk, vc, Trans, dofs); }
2044  virtual void Project(const FiniteElement &fe,
2046  DenseMatrix &I) const
2047  { Project_ND(tk, dof2tk, fe, Trans, I); }
2048  virtual void ProjectGrad(const FiniteElement &fe,
2050  DenseMatrix &grad) const
2051  { ProjectGrad_ND(tk, dof2tk, fe, Trans, grad); }
2052 };
2053 
2054 
2056 {
2057 protected:
2059  mutable int *ijk, patch, elem;
2060  mutable Vector weights;
2061 
2062 public:
2063  NURBSFiniteElement(int D, int G, int Do, int O, int F)
2064  : FiniteElement(D, G, Do, O, F)
2065  {
2066  ijk = NULL;
2067  patch = elem = -1;
2068  kv.SetSize(Dim);
2069  weights.SetSize(Dof);
2070  weights = 1.0;
2071  }
2072 
2073  void Reset () const { patch = elem = -1; }
2074  void SetIJK (int *IJK) const { ijk = IJK; }
2075  int GetPatch () const { return patch; }
2076  void SetPatch (int p) const { patch = p; }
2077  int GetElement () const { return elem; }
2078  void SetElement (int e) const { elem = e; }
2079  Array <KnotVector*> &KnotVectors() const { return kv; }
2080  Vector &Weights () const { return weights; }
2081 };
2082 
2084 {
2085 protected:
2086  mutable Vector shape_x;
2087 
2088 public:
2090  : NURBSFiniteElement(1, Geometry::SEGMENT, p + 1, p, FunctionSpace::Qk),
2091  shape_x(p + 1) { }
2092 
2093  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
2094  virtual void CalcDShape(const IntegrationPoint &ip,
2095  DenseMatrix &dshape) const;
2096 };
2097 
2099 {
2100 protected:
2102 
2103 public:
2105  : NURBSFiniteElement(2, Geometry::SQUARE, (p + 1)*(p + 1), p,
2106  FunctionSpace::Qk), u(Dof),
2107  shape_x(p + 1), shape_y(p + 1), dshape_x(p + 1), dshape_y(p + 1) { }
2108 
2109  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
2110  virtual void CalcDShape(const IntegrationPoint &ip,
2111  DenseMatrix &dshape) const;
2112 };
2113 
2115 {
2116 protected:
2118 
2119 public:
2121  : NURBSFiniteElement(3, Geometry::CUBE, (p + 1)*(p + 1)*(p + 1), p,
2122  FunctionSpace::Qk), u(Dof),
2123  shape_x(p + 1), shape_y(p + 1), shape_z(p + 1),
2124  dshape_x(p + 1), dshape_y(p + 1), dshape_z(p + 1) { }
2125 
2126  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
2127  virtual void CalcDShape(const IntegrationPoint &ip,
2128  DenseMatrix &dshape) const;
2129 };
2130 
2131 }
2132 
2133 #endif
Abstract class for Finite Elements.
Definition: fe.hpp:44
RefinedLinear3DFiniteElement()
Construct a quadratic FE on tetrahedron.
Definition: fe.cpp:4073
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1191
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:1189
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:8301
RT_QuadrilateralElement(const int p)
Definition: fe.cpp:9132
DenseMatrix curlshape_J
Definition: fe.hpp:280
const Array< int > & GetDofMap() const
Definition: fe.hpp:1384
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:5872
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.hpp:1756
int GetDim() const
Returns the space dimension for the finite element.
Definition: fe.hpp:82
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:7171
L2_TetrahedronElement(const int p, const int _type=0)
Definition: fe.cpp:8924
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:2714
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1163
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:3651
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:937
static void CalcBernstein(const int p, const double x, double *u, double *d)
Definition: fe.hpp:1300
Class for integration rule.
Definition: intrules.hpp:83
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 void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:8721
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:11040
Linear3DFiniteElement()
Construct a linear FE on tetrahedron.
Definition: fe.cpp:2112
static double CalcDelta(const int p, const double x)
Definition: fe.hpp:1279
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:2896
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1786
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1830
H1_HexahedronElement(const int p)
Definition: fe.cpp:6851
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:8460
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:9903
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:30
RefinedBiLinear2DFiniteElement()
Construct a biquadratic FE on quadrilateral.
Definition: fe.cpp:4305
ND_QuadrilateralElement(const int p)
Definition: fe.cpp:10354
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:1507
L2Pos_TriangleElement(const int p)
Definition: fe.cpp:8873
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.hpp:1929
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:756
H1Pos_SegmentElement(const int p)
Definition: fe.cpp:7109
int GetPatch() const
Definition: fe.hpp:2075
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.hpp:1840
virtual void ProjectDiv(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &div) const
Definition: fe.cpp:116
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:648
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Definition: fe.cpp:4917
virtual void Project(Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:75
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:559
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:2350
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:2129
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:8912
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.hpp:1973
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:2668
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:2948
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:1303
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:7151
FiniteElement(D, G, Do, O, F)
L2Pos_TetrahedronElement(const int p)
Definition: fe.cpp:9075
Array< KnotVector * > kv
Definition: fe.hpp:2058
void CalcVShape_RT(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.cpp:372
void SetSize(int s)
Resize the vector if the new size is different.
Definition: vector.hpp:263
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:9290
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:4979
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:3438
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:858
virtual void Project(Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:177
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:1327
void LocalInterpolation_ND(const double *tk, const Array< int > &d2t, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:687
virtual void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const
Definition: fe.cpp:1685
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:1596
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:2433
static void UniformPoints(const int p, double *x)
Definition: fe.cpp:6150
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
evaluate derivatives of shape function - constant 0
Definition: fe.cpp:2083
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1111
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1864
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:747
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:6679
H1Pos_TetrahedronElement(const int p)
Definition: fe.cpp:7893
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:10988
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:8320
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1715
LagrangeHexFiniteElement(int degree)
Definition: fe.cpp:3703
int GetOrder() const
Returns the order of the finite element.
Definition: fe.hpp:91
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:2446
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:3559
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:3865
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1913
Data type dense matrix using column-major storage.
Definition: densemat.hpp:22
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:3513
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.cpp:259
void Reset() const
Definition: fe.hpp:2073
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.hpp:1795
Class for quadratic FE on triangle.
Definition: fe.hpp:485
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:11002
Basis & OpenBasis(const int p)
Definition: fe.cpp:6548
const double * ClosedPoints(const int p)
Definition: fe.cpp:6523
PositiveFiniteElement(int D, int G, int Do, int O, int F=FunctionSpace::Pk)
Definition: fe.hpp:256
Class for quadratic FE on interval.
Definition: fe.hpp:456
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1710
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:9760
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:5145
void ProjectCurl_ND(const double *tk, const Array< int > &d2t, const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.cpp:509
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:2098
int Space() const
Returns the type of space on each element.
Definition: fe.hpp:94
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:1756
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:69
virtual void Project(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1759
RefinedTriLinear3DFiniteElement()
Construct a biquadratic FE on quadrilateral.
Definition: fe.cpp:4452
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:3690
Quadratic3DFiniteElement()
Construct a quadratic FE on tetrahedron.
Definition: fe.cpp:2168
int GetMapType() const
Definition: fe.hpp:98
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:8174
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:5495
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:5682
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:2545
const Array< int > & GetDofMap() const
Definition: fe.hpp:1362
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:9040
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.hpp:1891
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:1808
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:9105
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:2375
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:2367
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:8244
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:6763
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:961
void ProjectCurl_RT(const double *nk, const Array< int > &d2n, const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.cpp:547
H1Pos_HexahedronElement(const int p)
Definition: fe.cpp:7281
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.cpp:100
static void CalcShape(const int p, const double x, const double y, const double z, double *shape)
Definition: fe.cpp:7995
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:8982
RT_HexahedronElement(const int p)
Definition: fe.cpp:9340
virtual void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const
Definition: fe.cpp:819
int GetElement() const
Definition: fe.hpp:2077
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:7253
Linear2DFiniteElement()
Construct a linear FE on triangle.
Definition: fe.cpp:761
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.hpp:1719
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:2225
Class for refined bi-linear FE on quadrilateral.
Definition: fe.hpp:1040
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:5817
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:4641
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:10885
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:728
L2_SegmentElement(const int p, const int _type=0)
Definition: fe.cpp:8145
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:1828
L2Pos_HexahedronElement(const int p)
Definition: fe.cpp:8648
ND_TriangleElement(const int p)
Definition: fe.cpp:10805
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.hpp:1635
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:9939
static void CalcDShape(const int p, const double x, const double y, const double z, double *dshape_1d, double *dshape)
Definition: fe.cpp:8028
static void CalcBasis(const int p, const double x, double *u)
Definition: fe.hpp:1263
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:783
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
evaluate shape function - constant 1
Definition: fe.cpp:2077
virtual void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const
Definition: fe.cpp:1856
static void CalcBinomTerms(const int p, const double x, const double y, double *u)
Compute the terms in the expansion of the binomial (x + y)^p.
Definition: fe.cpp:6333
virtual void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const
Definition: fe.cpp:63
H1Pos_TriangleElement(const int p)
Definition: fe.cpp:7748
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:3150
Cubic3DFiniteElement()
Construct a cubic FE on tetrahedron.
Definition: fe.cpp:1903
Linear1DFiniteElement()
Construct a linear FE on interval.
Definition: fe.cpp:740
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:4014
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:5805
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:1376
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.cpp:108
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:994
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1752
static void CalcShape(const int p, const double x, const double y, double *shape)
Definition: fe.cpp:7803
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:2399
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:1996
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:3456
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:451
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 void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:8180
Class for refined linear FE on triangle.
Definition: fe.hpp:1007
Class for refined linear FE on tetrahedron.
Definition: fe.hpp:1027
void SetPatch(int p) const
Definition: fe.hpp:2076
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:954
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:965
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:1213
Class for constant FE on triangle.
Definition: fe.hpp:624
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.hpp:2048
ND_SegmentElement(const int p)
Definition: fe.cpp:10964
H1_TriangleElement(const int p)
Definition: fe.cpp:7459
RefinedLinear1DFiniteElement()
Construct a quadratic FE on interval.
Definition: fe.cpp:3903
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:7024
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:7274
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:1480
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:8815
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:7526
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.hpp:1762
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:3930
VectorFiniteElement(int D, int G, int Do, int O, int M, int F=FunctionSpace::Pk)
Definition: fe.hpp:331
void LocalInterpolation_RT(const double *nk, const Array< int > &d2n, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:651
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:3089
static void ChebyshevPoints(const int p, double *x)
Definition: fe.cpp:6301
NURBSFiniteElement(int D, int G, int Do, int O, int F)
Definition: fe.hpp:2063
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:5942
const Array< int > & GetDofMap() const
Definition: fe.hpp:1343
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:3675
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:5783
Class for refined trilinear FE on a hexahedron.
Definition: fe.hpp:1060
ND_TetrahedronElement(const int p)
Definition: fe.cpp:10543
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:8901
Class for quadratic FE on tetrahedron.
Definition: fe.hpp:679
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const =0
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1952
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:734
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1825
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:1771
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1957
const IntegrationRule & GetNodes() const
Definition: fe.hpp:113
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:5032
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:1035
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:8556
DenseMatrix vshape
Definition: fe.hpp:280
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:2203
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.hpp:1725
L2Pos_QuadrilateralElement(const int p)
Definition: fe.cpp:8396
P0SegmentFiniteElement(int Ord=0)
Definition: fe.cpp:2338
static void CalcDShape(const int p, const double x, const double y, double *dshape_1d, double *dshape)
Definition: fe.cpp:7829
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:2691
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:11075
DenseMatrix curlshape
Definition: fe.hpp:280
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.hpp:1961
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:863
Class for linear FE on tetrahedron.
Definition: fe.hpp:654
Class for linear FE on interval.
Definition: fe.hpp:354
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:9501
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:5083
Crouzeix-Raviart finite element on quadrilateral.
Definition: fe.hpp:727
Class for linear FE on triangle.
Definition: fe.hpp:374
Class for quadratic FE on triangle with nodes at the &quot;Gaussian&quot; points.
Definition: fe.hpp:509
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:8793
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:9243
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:2001
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:1635
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:3348
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1918
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:810
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:2138
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:9575
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:7689
virtual void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const
Definition: fe.cpp:1059
void SetElement(int e) const
Definition: fe.hpp:2078
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:1137
void ProjectGrad_ND(const double *tk, const Array< int > &d2t, const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.cpp:630
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:8895
virtual void Project(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1837
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:810
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.hpp:1562
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:11099
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:1151
int GetGeomType() const
Returns the geometry type:
Definition: fe.hpp:85
DenseMatrix m_dshape
Definition: fe.hpp:1466
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Definition: fe.cpp:10489
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:8421
ND_HexahedronElement(const int p)
Definition: fe.cpp:9985
IntegrationRule Nodes
Definition: fe.hpp:48
const Array< int > & GetDofMap() const
Definition: fe.hpp:1422
virtual void Project(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:2044
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:1021
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:8676
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:7131
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:7452
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:3400
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Definition: fe.cpp:10919
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:2560
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:4328
virtual ~FiniteElement()
Definition: fe.hpp:196
H1_QuadrilateralElement(const int p)
Definition: fe.cpp:6707
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:2457
const Array< int > & GetDofMap() const
Definition: fe.hpp:1403
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Overrides the scalar CalcShape function to print an error.
Definition: fe.hpp:2028
void CalcVShape_ND(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.cpp:385
DenseMatrix Jinv
Definition: fe.hpp:279
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1135
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.hpp:2041
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:8191
int GetDof() const
Returns the degrees of freedom in the FE space.
Definition: fe.hpp:88
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:7714
Class for bi-quadratic FE on quadrilateral.
Definition: fe.hpp:525
Class for tri-linear FE on cube.
Definition: fe.hpp:692
Base class Coefficient that may optionally depend on time.
Definition: coefficient.hpp:31
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:905
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:950
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Definition: fe.cpp:5113
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:9098
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:722
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:772
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:4498
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:4195
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:636
static const int * Binom(const int p)
Definition: fe.cpp:6133
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:3250
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:44
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:840
virtual void Project(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1722
static void GaussLobattoPoints(const int p, double *x)
Definition: fe.cpp:6225
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1089
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:9007
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:10201
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:2624
FiniteElement(int D, int G, int Do, int O, int F=FunctionSpace::Pk)
Definition: fe.cpp:22
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:1262
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:2283
Class for bilinear FE on quadrilateral.
Definition: fe.hpp:396
Quad1DFiniteElement()
Construct a quadratic FE on interval.
Definition: fe.cpp:946
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:2032
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:7407
virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const
Definition: fe.cpp:58
Class for linear FE on triangle with nodes at the 3 &quot;Gaussian&quot; points.
Definition: fe.hpp:420
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:894
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1871
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:884
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:801
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Definition: fe.cpp:51
H1Pos_QuadrilateralElement(const int p)
Definition: fe.cpp:7178
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.hpp:2005
Array< KnotVector * > & KnotVectors() const
Definition: fe.hpp:2079
Class for integration point with weight.
Definition: intrules.hpp:25
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:754
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:8844
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:7428
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:8250
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:230
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:2104
Bi-quadratic element on quad with nodes at the 9 Gaussian points.
Definition: fe.hpp:564
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:837
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:8697
void Project_RT(const double *nk, const Array< int > &d2n, VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:404
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:6979
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:929
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:3064
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:10691
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:2824
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:10440
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:391
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.hpp:2012
virtual void ProjectDiv(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &div) const
Definition: fe.cpp:288
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:2572
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:8261
RefinedLinear2DFiniteElement()
Construct a quadratic FE on triangle.
Definition: fe.cpp:3949
void NodalLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I, const NodalFiniteElement &fine_fe) const
Definition: fe.cpp:125
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:5737
void ProjectCurl_2D(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Definition: fe.cpp:158
void SetIJK(int *IJK) const
Definition: fe.hpp:2074
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.hpp:1922
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.hpp:1886
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:3911
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:3639
Array< int > dof_map
Definition: fe.hpp:1468
RT_TetrahedronElement(const int p)
Definition: fe.cpp:9802
L2_TriangleElement(const int p, const int _type=0)
Definition: fe.cpp:8740
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.cpp:5293
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:3844
Vector & Weights() const
Definition: fe.hpp:2080
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:6660
H1_TetrahedronElement(const int p)
Definition: fe.cpp:7578
BiQuad2DFiniteElement()
Construct a biquadratic FE on quadrilateral.
Definition: fe.cpp:1166
virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const
Definition: fe.cpp:2158
Tensor products of 1D FEs (only degree 2 is functional)
Definition: fe.hpp:966
Basis & ClosedBasis(const int p)
Definition: fe.cpp:6571
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:2861
P0TriangleFiniteElement()
Construct P0 triangle finite element.
Definition: fe.cpp:2070
virtual void Project(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:2008
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:2344
virtual void Project(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1798
RT_TriangleElement(const int p)
Definition: fe.cpp:9650
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:6641
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:10979
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.hpp:1877
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Definition: fe.cpp:10275
void Eval(const double x, Vector &u) const
Definition: fe.cpp:6023
L2Pos_SegmentElement(const int p)
Definition: fe.cpp:8223
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:4108
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:5348
Class for cubic FE on tetrahedron.
Definition: fe.hpp:611
Vector data type.
Definition: vector.hpp:33
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:5279
static void CalcBernstein(const int p, const double x, double *u)
Definition: fe.hpp:1298
virtual void Project(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1881
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:7233
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:8532
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:1087
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:9727
Quad2DFiniteElement()
Construct a quadratic FE on triangle.
Definition: fe.cpp:1004
virtual void Project(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1925
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.hpp:1834
Describes the space on each element.
Definition: fe.hpp:27
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:2299
const double * OpenPoints(const int p)
Definition: fe.cpp:6498
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:2037
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:5198
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:8440
static void GaussPoints(const int p, double *x)
Definition: fe.cpp:6165
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Definition: fe.cpp:10733
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:3696
virtual void Project(Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:1405
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:780
H1_SegmentElement(const int p)
Definition: fe.cpp:6618
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:1968
NodalFiniteElement(int D, int G, int Do, int O, int F=FunctionSpace::Pk)
Definition: fe.hpp:215
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:2410
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:6802
L2_HexahedronElement(const int p, const int _type=0)
Definition: fe.cpp:8475
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:710
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:5249
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:3468
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:415
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:4379
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const =0
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:11021
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1791
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.hpp:1968
void Project_ND(const double *tk, const Array< int > &d2t, VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:566
int GetRangeType() const
Definition: fe.hpp:96
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:883
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:2766
static void CalcBasis(const int p, const double x, double *u, double *d)
Definition: fe.hpp:1272
virtual void Project(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const
Definition: fe.hpp:1964
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:3669
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:7000
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.hpp:1801
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:6782
Basis(const int p, const double *nodes, const int _mode=1)
Definition: fe.cpp:5977
static void CalcDBinomTerms(const int p, const double x, const double y, double *d)
Definition: fe.cpp:6397
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Definition: fe.cpp:4862
TriLinear3DFiniteElement()
Construct a tri-linear FE on cube.
Definition: fe.cpp:2247
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:87
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:984
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:8511
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1747
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:3966
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:3429
BiLinear2DFiniteElement()
Construct a bilinear FE on quadrilateral.
Definition: fe.cpp:788
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe.hpp:1996
Poly_1D poly1d
Definition: fe.cpp:6614
Crouzeix-Raviart finite element on triangle.
Definition: fe.hpp:715
L2_QuadrilateralElement(const int p, const int _type=0)
Definition: fe.cpp:8268
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Definition: fe.cpp:5925
Lagrange1DFiniteElement(int degree)
Definition: fe.cpp:3481
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:2509
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Definition: fe.cpp:5618
const Array< int > & GetDofMap() const
Definition: fe.hpp:1324
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:9116
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.cpp:8340
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:850
virtual void ProjectDelta(int vertex, Vector &dofs) const
Definition: fe.hpp:672
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Definition: fe.cpp:7548
void ProjectGrad_RT(const double *nk, const Array< int > &d2n, const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Definition: fe.cpp:482
virtual void Project(Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const
Definition: fe.cpp:320