MFEM  v4.5.1
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
fe_fixed_order.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2022, Lawrence Livermore National Security, LLC. Produced
2 // at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3 // LICENSE and NOTICE for details. LLNL-CODE-806117.
4 //
5 // This file is part of the MFEM library. For more information and source code
6 // availability visit https://mfem.org.
7 //
8 // MFEM is free software; you can redistribute it and/or modify it under the
9 // terms of the BSD-3 license. We welcome feedback and contributions, see file
10 // CONTRIBUTING.md for details.
11 
12 #ifndef MFEM_FE_FIXED_ORDER
13 #define MFEM_FE_FIXED_ORDER
14 
15 #include "fe_base.hpp"
16 #include "fe_h1.hpp"
17 #include "fe_l2.hpp"
18 
19 namespace mfem
20 {
21 
22 /// A 0D point finite element
24 {
25 public:
26  /// Construct the PointFiniteElement
28 
29  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
30 
31  virtual void CalcDShape(const IntegrationPoint &ip,
32  DenseMatrix &dshape) const;
33 };
34 
35 /// A 1D linear element with nodes on the endpoints
37 {
38 public:
39  /// Construct the Linear1DFiniteElement
41 
42  /** virtual function which evaluates the values of all
43  shape functions at a given point ip and stores
44  them in the vector shape of dimension Dof (2) */
45  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
46 
47  /** virtual function which evaluates the derivatives of all
48  shape functions at a given point ip and stores them in
49  the matrix dshape (Dof x Dim) (2 x 1) so that each row
50  contains the derivative of one shape function */
51  virtual void CalcDShape(const IntegrationPoint &ip,
52  DenseMatrix &dshape) const;
53 };
54 
55 /// A 2D linear element on triangle with nodes at the vertices of the triangle
57 {
58 public:
59  /// Construct the Linear2DFiniteElement
61 
62  /** virtual function which evaluates the values of all
63  shape functions at a given point ip and stores
64  them in the vector shape of dimension Dof (3) */
65  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
66 
67  /** virtual function which evaluates the values of all
68  partial derivatives of all shape functions at a given
69  point ip and stores them in the matrix dshape (Dof x Dim) (3 x 2)
70  so that each row contains the derivatives of one shape function */
71  virtual void CalcDShape(const IntegrationPoint &ip,
72  DenseMatrix &dshape) const;
73  virtual void ProjectDelta(int vertex, Vector &dofs) const
74  { dofs = 0.0; dofs(vertex) = 1.0; }
75 };
76 
77 /// A 2D bi-linear element on a square with nodes at the vertices of the square
79 {
80 public:
81  /// Construct the BiLinear2DFiniteElement
83 
84  /** virtual function which evaluates the values of all
85  shape functions at a given point ip and stores
86  them in the vector shape of dimension Dof (4) */
87  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
88 
89  /** virtual function which evaluates the values of all
90  partial derivatives of all shape functions at a given
91  point ip and stores them in the matrix dshape (Dof x Dim) (4 x 2)
92  so that each row contains the derivatives of one shape function */
93  virtual void CalcDShape(const IntegrationPoint &ip,
94  DenseMatrix &dshape) const;
95  virtual void CalcHessian (const IntegrationPoint &ip,
96  DenseMatrix &h) const;
97  virtual void ProjectDelta(int vertex, Vector &dofs) const
98  { dofs = 0.0; dofs(vertex) = 1.0; } // { dofs = 1.0; }
99 };
100 
101 /// A linear element on a triangle with nodes at the 3 "Gaussian" points
103 {
104 public:
105  /// Construct the GaussLinear2DFiniteElement
107  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
108  virtual void CalcDShape(const IntegrationPoint &ip,
109  DenseMatrix &dshape) const;
110  virtual void ProjectDelta(int vertex, Vector &dofs) const;
111 };
112 
113 /// A 2D bi-linear element on a square with nodes at the "Gaussian" points
115 {
116 private:
117  static const double p[2];
118 
119 public:
120  /// Construct the FiniteElement
122  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
123  virtual void CalcDShape(const IntegrationPoint &ip,
124  DenseMatrix &dshape) const;
125  virtual void ProjectDelta(int vertex, Vector &dofs) const;
126 };
127 
128 /** @brief A 2D linear element on a square with 3 nodes at the
129  vertices of the lower left triangle */
131 {
132 public:
133  /// Construct the P1OnQuadFiniteElement
135  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
136  virtual void CalcDShape(const IntegrationPoint &ip,
137  DenseMatrix &dshape) const;
138  virtual void ProjectDelta(int vertex, Vector &dofs) const
139  { dofs = 1.0; }
140 };
141 
142 /// A 1D quadratic finite element with uniformly spaced nodes
144 {
145 public:
146  /// Construct the Quad1DFiniteElement
148 
149  /** virtual function which evaluates the values of all
150  shape functions at a given point ip and stores
151  them in the vector shape of dimension Dof (3) */
152  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
153 
154  /** virtual function which evaluates the derivatives of all
155  shape functions at a given point ip and stores them in
156  the matrix dshape (Dof x Dim) (3 x 1) so that each row
157  contains the derivative of one shape function */
158  virtual void CalcDShape(const IntegrationPoint &ip,
159  DenseMatrix &dshape) const;
160 };
161 
162 /** @brief A 2D quadratic element on triangle with nodes at the
163  vertices and midpoints of the triangle. */
165 {
166 public:
167  /// Construct the Quad2DFiniteElement
169 
170  /** virtual function which evaluates the values of all
171  shape functions at a given point ip and stores
172  them in the vector shape of dimension Dof (6) */
173  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
174 
175  /** virtual function which evaluates the values of all
176  partial derivatives of all shape functions at a given
177  point ip and stores them in the matrix dshape (Dof x Dim) (6 x 2)
178  so that each row contains the derivatives of one shape function */
179  virtual void CalcDShape(const IntegrationPoint &ip,
180  DenseMatrix &dshape) const;
181 
182  virtual void CalcHessian (const IntegrationPoint &ip,
183  DenseMatrix &h) const;
184  virtual void ProjectDelta(int vertex, Vector &dofs) const;
185 };
186 
187 /// A quadratic element on triangle with nodes at the "Gaussian" points
189 {
190 private:
191  static const double p[2];
192  DenseMatrix A;
193  mutable DenseMatrix D;
194  mutable Vector pol;
195 public:
196  /// Construct the GaussQuad2DFiniteElement
198  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
199  virtual void CalcDShape(const IntegrationPoint &ip,
200  DenseMatrix &dshape) const;
201  // virtual void ProjectDelta(int vertex, Vector &dofs) const;
202 };
203 
204 /// A 2D bi-quadratic element on a square with uniformly spaced nodes
206 {
207 public:
208  /// Construct the BiQuad2DFiniteElement
210 
211  /** virtual function which evaluates the values of all
212  shape functions at a given point ip and stores
213  them in the vector shape of dimension Dof (9) */
214  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
215 
216  /** virtual function which evaluates the values of all
217  partial derivatives of all shape functions at a given
218  point ip and stores them in the matrix dshape (Dof x Dim) (9 x 2)
219  so that each row contains the derivatives of one shape function */
220  virtual void CalcDShape(const IntegrationPoint &ip,
221  DenseMatrix &dshape) const;
222  virtual void ProjectDelta(int vertex, Vector &dofs) const;
223 };
224 
225 
226 /// A 2D bi-quadratic element on a square with nodes at the 9 "Gaussian" points
228 {
229 public:
230  /// Construct the GaussBiQuad2DFiniteElement
232  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
233  virtual void CalcDShape(const IntegrationPoint &ip,
234  DenseMatrix &dshape) const;
235  // virtual void ProjectDelta(int vertex, Vector &dofs) const { dofs = 1.; }
236 };
237 
238 
239 /// A 2D bi-cubic element on a square with uniformly spaces nodes
241 {
242 public:
243  /// Construct the BiCubic2DFiniteElement
245  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
246  virtual void CalcDShape(const IntegrationPoint &ip,
247  DenseMatrix &dshape) const;
248 
249  /// Compute the Hessian of second order partial derivatives at @a ip.
250  virtual void CalcHessian (const IntegrationPoint &ip,
251  DenseMatrix &h) const;
252 };
253 
254 /// A 1D cubic element with uniformly spaced nodes
256 {
257 public:
258  /// Construct the Cubic1DFiniteElement
260 
261  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
262 
263  virtual void CalcDShape(const IntegrationPoint &ip,
264  DenseMatrix &dshape) const;
265 };
266 
267 /// A 2D cubic element on a triangle with uniformly spaced nodes
269 {
270 public:
271  /// Construct the Cubic2DFiniteElement
273 
274  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
275 
276  virtual void CalcDShape(const IntegrationPoint &ip,
277  DenseMatrix &dshape) const;
278 
279  virtual void CalcHessian (const IntegrationPoint &ip,
280  DenseMatrix &h) const;
281 };
282 
283 /// A 3D cubic element on a tetrahedron with 20 nodes at the thirds of the
284 /// tetrahedron
286 {
287 public:
288  /// Construct the Cubic3DFiniteElement
290 
291  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
292 
293  virtual void CalcDShape(const IntegrationPoint &ip,
294  DenseMatrix &dshape) const;
295 };
296 
297 /// A linear element defined on a triangular prism
299 {
300 public:
301  /// Construct the LinearWedgeFiniteElement
303 
304  /** @brief virtual function which evaluates the values of all
305  shape functions at a given point ip and stores
306  them in the vector shape of dimension Dof (4) */
307  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
308 
309  /** @brief virtual function which evaluates the values of all
310  partial derivatives of all shape functions at a given
311  point ip and stores them in the matrix dshape (Dof x Dim) (4 x 3)
312  so that each row contains the derivatives of one shape function */
313  virtual void CalcDShape(const IntegrationPoint &ip,
314  DenseMatrix &dshape) const;
315 
316  virtual void ProjectDelta(int vertex, Vector &dofs) const
317  { dofs = 0.0; dofs(vertex) = 1.0; }
318 
319  /** @brief Get the dofs associated with the given @a face.
320  @a *dofs is set to an internal array of the local dofc on the
321  face, while *ndofs is set to the number of dofs on that face.
322  */
323  virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const;
324 };
325 
326 /// A linear element defined on a square pyramid
328 {
329 public:
330  /// Construct the LinearPyramidFiniteElement
332 
333  /** @brief virtual function which evaluates the values of all
334  shape functions at a given point ip and stores
335  them in the vector shape of dimension Dof (4) */
336  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
337 
338  /** @brief virtual function which evaluates the values of all
339  partial derivatives of all shape functions at a given
340  point ip and stores them in the matrix dshape (Dof x Dim) (4 x 3)
341  so that each row contains the derivatives of one shape function */
342  virtual void CalcDShape(const IntegrationPoint &ip,
343  DenseMatrix &dshape) const;
344 
345  virtual void ProjectDelta(int vertex, Vector &dofs) const
346  { dofs = 0.0; dofs(vertex) = 1.0; }
347 
348  /** @brief Get the dofs associated with the given @a face.
349  @a *dofs is set to an internal array of the local dofc on the
350  face, while *ndofs is set to the number of dofs on that face.
351  */
352  virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const;
353 };
354 
355 /// A 2D constant element on a triangle
357 {
358 public:
359  /// Construct the P0TriangleFiniteElement
361 
362  /// evaluate shape function - constant 1
363  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
364 
365  /// evaluate derivatives of shape function - constant 0
366  virtual void CalcDShape(const IntegrationPoint &ip,
367  DenseMatrix &dshape) const;
368  virtual void ProjectDelta(int vertex, Vector &dofs) const
369  { dofs(0) = 1.0; }
370 };
371 
372 
373 /// A 2D constant element on a square
375 {
376 public:
377  /// Construct the P0QuadFiniteElement
379  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
380  virtual void CalcDShape(const IntegrationPoint &ip,
381  DenseMatrix &dshape) const;
382  virtual void ProjectDelta(int vertex, Vector &dofs) const
383  { dofs(0) = 1.0; }
384 };
385 
386 
387 /** @brief A 3D linear element on a tetrahedron with nodes at the
388  vertices of the tetrahedron */
390 {
391 public:
392  /// Construct the Linear3DFiniteElement
394 
395  /** @brief virtual function which evaluates the values of all
396  shape functions at a given point ip and stores
397  them in the vector shape of dimension Dof (4) */
398  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
399 
400  /** @brief virtual function which evaluates the values of all
401  partial derivatives of all shape functions at a given
402  point ip and stores them in the matrix dshape (Dof x Dim) (4 x 3)
403  so that each row contains the derivatives of one shape function */
404  virtual void CalcDShape(const IntegrationPoint &ip,
405  DenseMatrix &dshape) const;
406 
407  virtual void ProjectDelta(int vertex, Vector &dofs) const
408  { dofs = 0.0; dofs(vertex) = 1.0; }
409 
410  /** @brief Get the dofs associated with the given @a face.
411  @a *dofs is set to an internal array of the local dofc on the
412  face, while *ndofs is set to the number of dofs on that face.
413  */
414  virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const;
415 };
416 
417 /// A 3D quadratic element on a tetrahedron with uniformly spaced nodes
419 {
420 public:
421  /// Construct the Quadratic3DFiniteElement
423 
424  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
425 
426  virtual void CalcDShape(const IntegrationPoint &ip,
427  DenseMatrix &dshape) const;
428 };
429 
430 /// A 3D tri-linear element on a cube with nodes at the vertices of the cube
432 {
433 public:
434  /// Construct the TriLinear3DFiniteElement
436 
437  /** virtual function which evaluates the values of all
438  shape functions at a given point ip and stores
439  them in the vector shape of dimension Dof (8) */
440  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
441 
442  /** virtual function which evaluates the values of all
443  partial derivatives of all shape functions at a given
444  point ip and stores them in the matrix dshape (Dof x Dim) (8 x 3)
445  so that each row contains the derivatives of one shape function */
446  virtual void CalcDShape(const IntegrationPoint &ip,
447  DenseMatrix &dshape) const;
448 
449  virtual void ProjectDelta(int vertex, Vector &dofs) const
450  { dofs = 0.0; dofs(vertex) = 1.0; }
451 };
452 
453 
454 /// A 2D Crouzeix-Raviart element on triangle
456 {
457 public:
458  /// Construct the CrouzeixRaviartFiniteElement
460  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
461  virtual void CalcDShape(const IntegrationPoint &ip,
462  DenseMatrix &dshape) const;
463  virtual void ProjectDelta(int vertex, Vector &dofs) const
464  { dofs = 1.0; }
465 };
466 
467 /// A 2D Crouzeix-Raviart finite element on square
469 {
470 public:
471  /// Construct the CrouzeixRaviartQuadFiniteElement
473  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
474  virtual void CalcDShape(const IntegrationPoint &ip,
475  DenseMatrix &dshape) const;
476 };
477 
478 
479 /// A 1D constant element on a segment
481 {
482 public:
483  /// Construct the P0SegmentFiniteElement with dummy order @a Ord
484  P0SegmentFiniteElement(int Ord = 0);
485  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
486  virtual void CalcDShape(const IntegrationPoint &ip,
487  DenseMatrix &dshape) const;
488 };
489 
490 /** @brief A 2D 1st order Raviart-Thomas vector element on a triangle */
492 {
493 private:
494  static const double nk[3][2];
495 
496 public:
497  /// Construct the RT0TriangleFiniteElement
499 
500  virtual void CalcVShape(const IntegrationPoint &ip,
501  DenseMatrix &shape) const;
502 
504  DenseMatrix &shape) const
505  { CalcVShape_RT(Trans, shape); }
506 
507  virtual void CalcDivShape(const IntegrationPoint &ip,
508  Vector &divshape) const;
509 
511  DenseMatrix &I) const;
512 
514 
515  virtual void Project (VectorCoefficient &vc,
516  ElementTransformation &Trans, Vector &dofs) const;
517 };
518 
519 /** @brief A 2D 1st order Raviart-Thomas vector element on a square*/
521 {
522 private:
523  static const double nk[4][2];
524 
525 public:
526  /// Construct the RT0QuadFiniteElement
528 
529  virtual void CalcVShape(const IntegrationPoint &ip,
530  DenseMatrix &shape) const;
531 
533  DenseMatrix &shape) const
534  { CalcVShape_RT(Trans, shape); }
535 
536  virtual void CalcDivShape(const IntegrationPoint &ip,
537  Vector &divshape) const;
538 
540  DenseMatrix &I) const;
541 
543 
544  virtual void Project (VectorCoefficient &vc,
545  ElementTransformation &Trans, Vector &dofs) const;
546 };
547 
548 /** @brief A 2D 2nd order Raviart-Thomas vector element on a triangle */
550 {
551 private:
552  static const double nk[8][2];
553 
554 public:
555  /// Construct the RT1TriangleFiniteElement
557 
558  virtual void CalcVShape(const IntegrationPoint &ip,
559  DenseMatrix &shape) const;
560 
562  DenseMatrix &shape) const
563  { CalcVShape_RT(Trans, shape); }
564 
565  virtual void CalcDivShape(const IntegrationPoint &ip,
566  Vector &divshape) const;
567 
569  DenseMatrix &I) const;
570 
572 
573  virtual void Project (VectorCoefficient &vc,
574  ElementTransformation &Trans, Vector &dofs) const;
575 };
576 
577 /** @brief A 2D 2nd order Raviart-Thomas vector element on a square */
579 {
580 private:
581  static const double nk[12][2];
582 
583 public:
584  /// Construct the RT1QuadFiniteElement
586 
587  virtual void CalcVShape(const IntegrationPoint &ip,
588  DenseMatrix &shape) const;
589 
591  DenseMatrix &shape) const
592  { CalcVShape_RT(Trans, shape); }
593 
594  virtual void CalcDivShape(const IntegrationPoint &ip,
595  Vector &divshape) const;
596 
598  DenseMatrix &I) const;
599 
601 
602  virtual void Project (VectorCoefficient &vc,
603  ElementTransformation &Trans, Vector &dofs) const;
604 };
605 
606 /** @brief A 2D 3rd order Raviart-Thomas vector element on a triangle */
608 {
609 private:
610  static const double M[15][15];
611 public:
612  /// Construct the RT2TriangleFiniteElement
614 
615  virtual void CalcVShape(const IntegrationPoint &ip,
616  DenseMatrix &shape) const;
617 
619  DenseMatrix &shape) const
620  { CalcVShape_RT(Trans, shape); }
621 
622  virtual void CalcDivShape(const IntegrationPoint &ip,
623  Vector &divshape) const;
624 };
625 
626 /** @brief A 2D 3rd order Raviart-Thomas vector element on a square */
628 {
629 private:
630  static const double nk[24][2];
631  static const double pt[4];
632  static const double dpt[3];
633 
634 public:
635  /// Construct the RT2QuadFiniteElement
637 
638  virtual void CalcVShape(const IntegrationPoint &ip,
639  DenseMatrix &shape) const;
640 
642  DenseMatrix &shape) const
643  { CalcVShape_RT(Trans, shape); }
644 
645  virtual void CalcDivShape(const IntegrationPoint &ip,
646  Vector &divshape) const;
647 
649  DenseMatrix &I) const;
650 
652 
653  virtual void Project (VectorCoefficient &vc,
654  ElementTransformation &Trans, Vector &dofs) const;
655 };
656 
657 /// A 1D linear element with nodes at 1/3 and 2/3 (trace of RT1)
659 {
660 public:
661  /// Construct the P1SegmentFiniteElement
663  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
664  virtual void CalcDShape(const IntegrationPoint &ip,
665  DenseMatrix &dshape) const;
666 };
667 
668 /// A 1D quadratic element with nodes at the Gaussian points (trace of RT2)
670 {
671 public:
672  /// Construct the P2SegmentFiniteElement
674  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
675  virtual void CalcDShape(const IntegrationPoint &ip,
676  DenseMatrix &dshape) const;
677 };
678 
679 /// A 1D element with uniform nodes
681 {
682 private:
683  Vector rwk;
684 #ifndef MFEM_THREAD_SAFE
685  mutable Vector rxxk;
686 #endif
687 public:
688  /// Construct the Lagrange1DFiniteElement with the provided @a degree
689  Lagrange1DFiniteElement (int degree);
690  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
691  virtual void CalcDShape(const IntegrationPoint &ip,
692  DenseMatrix &dshape) const;
693 };
694 
695 /// A 3D Crouzeix-Raviart element on the tetrahedron.
697 {
698 public:
699  /// Construct the P1TetNonConfFiniteElement
701  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
702  virtual void CalcDShape(const IntegrationPoint &ip,
703  DenseMatrix &dshape) const;
704 };
705 
706 /// A 3D constant element on a tetrahedron
708 {
709 public:
710  /// Construct the P0TetFiniteElement
712  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
713  virtual void CalcDShape(const IntegrationPoint &ip,
714  DenseMatrix &dshape) const;
715  virtual void ProjectDelta(int vertex, Vector &dofs) const
716  { dofs(0) = 1.0; }
717 };
718 
719 /// A 3D constant element on a cube
721 {
722 public:
723  /// Construct the P0HexFiniteElement
725  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
726  virtual void CalcDShape(const IntegrationPoint &ip,
727  DenseMatrix &dshape) const;
728  virtual void ProjectDelta(int vertex, Vector &dofs) const
729  { dofs(0) = 1.0; }
730 };
731 
732 /// A 3D constant element on a wedge
734 {
735 public:
736  /// Construct the P0WdgFiniteElement
738  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
739  virtual void CalcDShape(const IntegrationPoint &ip,
740  DenseMatrix &dshape) const;
741  virtual void ProjectDelta(int vertex, Vector &dofs) const
742  { dofs(0) = 1.0; }
743 };
744 
745 /// A 3D constant element on a pyramid
747 {
748 public:
749  /// Construct the P0PyrFiniteElement
751  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
752  virtual void CalcDShape(const IntegrationPoint &ip,
753  DenseMatrix &dshape) const;
754  virtual void ProjectDelta(int vertex, Vector &dofs) const
755  { dofs(0) = 1.0; }
756 };
757 
758 /** @brief Tensor products of 1D Lagrange1DFiniteElement
759  (only degree 2 is functional) */
761 {
762 private:
764  int dof1d;
765  int *I, *J, *K;
766 #ifndef MFEM_THREAD_SAFE
767  mutable Vector shape1dx, shape1dy, shape1dz;
768  mutable DenseMatrix dshape1dx, dshape1dy, dshape1dz;
769 #endif
770 
771 public:
772  /// Construct the LagrangeHexFiniteElement with the provided @a degree
773  LagrangeHexFiniteElement (int degree);
774  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
775  virtual void CalcDShape(const IntegrationPoint &ip,
776  DenseMatrix &dshape) const;
778 };
779 
780 
781 /// A 1D refined linear element
783 {
784 public:
785  /// Construct the RefinedLinear1DFiniteElement
787 
788  /** virtual function which evaluates the values of all
789  shape functions at a given point ip and stores
790  them in the vector shape of dimension Dof (3) */
791  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
792 
793  /** virtual function which evaluates the derivatives of all
794  shape functions at a given point ip and stores them in
795  the matrix dshape (Dof x Dim) (3 x 1) so that each row
796  contains the derivative of one shape function */
797  virtual void CalcDShape(const IntegrationPoint &ip,
798  DenseMatrix &dshape) const;
799 };
800 
801 /// A 2D refined linear element on a triangle
803 {
804 public:
805  /// Construct the RefinedLinear2DFiniteElement
807 
808  /** virtual function which evaluates the values of all
809  shape functions at a given point ip and stores
810  them in the vector shape of dimension Dof (6) */
811  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
812 
813  /** virtual function which evaluates the values of all
814  partial derivatives of all shape functions at a given
815  point ip and stores them in the matrix dshape (Dof x Dim) (6 x 2)
816  so that each row contains the derivatives of one shape function */
817  virtual void CalcDShape(const IntegrationPoint &ip,
818  DenseMatrix &dshape) const;
819 };
820 
821 /// A 2D refined linear element on a tetrahedron
823 {
824 public:
825  /// Construct the RefinedLinear3DFiniteElement
827 
828  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
829 
830  virtual void CalcDShape(const IntegrationPoint &ip,
831  DenseMatrix &dshape) const;
832 };
833 
834 /// A 2D refined bi-linear FE on a square
836 {
837 public:
838  /// Construct the RefinedBiLinear2DFiniteElement
840 
841  /** virtual function which evaluates the values of all
842  shape functions at a given point ip and stores
843  them in the vector shape of dimension Dof (9) */
844  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
845 
846  /** virtual function which evaluates the values of all
847  partial derivatives of all shape functions at a given
848  point ip and stores them in the matrix dshape (Dof x Dim) (9 x 2)
849  so that each row contains the derivatives of one shape function */
850  virtual void CalcDShape(const IntegrationPoint &ip,
851  DenseMatrix &dshape) const;
852 };
853 
854 /// A 3D refined tri-linear element on a cube
856 {
857 public:
858  /// Construct the RefinedTriLinear3DFiniteElement
860 
861  /** virtual function which evaluates the values of all
862  shape functions at a given point ip and stores
863  them in the vector shape of dimension Dof (9) */
864  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
865 
866  /** virtual function which evaluates the values of all
867  partial derivatives of all shape functions at a given
868  point ip and stores them in the matrix dshape (Dof x Dim) (9 x 2)
869  so that each row contains the derivatives of one shape function */
870  virtual void CalcDShape(const IntegrationPoint &ip,
871  DenseMatrix &dshape) const;
872 };
873 
874 
875 /// Class for linear FE on wedge
877 {
878 public:
879  /// Construct a linear FE on wedge
881 };
882 
883 /// Class for quadratic FE on wedge
885 {
886 public:
887  /// Construct a quadratic FE on wedge
889 };
890 
891 /// Class for cubic FE on wedge
893 {
894 public:
895  /// Construct a cubic FE on wedge
897 };
898 
899 
900 /// A 0th order L2 element on a Wedge
902 {
903 public:
904  /// Construct the P0WedgeFiniteElement
906 };
907 
908 
909 /// A 3D 1st order Nedelec element on a cube
911 {
912 private:
913  static const double tk[12][3];
914 
915 public:
916  /// Construct the Nedelec1HexFiniteElement
918  virtual void CalcVShape(const IntegrationPoint &ip,
919  DenseMatrix &shape) const;
921  DenseMatrix &shape) const
922  { CalcVShape_ND(Trans, shape); }
923  virtual void CalcCurlShape(const IntegrationPoint &ip,
924  DenseMatrix &curl_shape) const;
926  DenseMatrix &I) const;
928  virtual void Project (VectorCoefficient &vc,
929  ElementTransformation &Trans, Vector &dofs) const;
930 
931  virtual void ProjectGrad(const FiniteElement &fe,
933  DenseMatrix &grad) const;
934 };
935 
936 
937 /// A 3D 1st order Nedelec element on a tetrahedron
939 {
940 private:
941  static const double tk[6][3];
942 
943 public:
944  /// Construct the Nedelec1TetFiniteElement
946  virtual void CalcVShape(const IntegrationPoint &ip,
947  DenseMatrix &shape) const;
949  DenseMatrix &shape) const
950  { CalcVShape_ND(Trans, shape); }
951  virtual void CalcCurlShape(const IntegrationPoint &ip,
952  DenseMatrix &curl_shape) const;
954  DenseMatrix &I) const;
956  virtual void Project (VectorCoefficient &vc,
957  ElementTransformation &Trans, Vector &dofs) const;
958 
959  virtual void ProjectGrad(const FiniteElement &fe,
961  DenseMatrix &grad) const;
962 };
963 
964 
965 /// A 3D 1st order Nedelec element on a wedge
967 {
968 private:
969  static const double tk[9][3];
970 
971 public:
972  /// Construct the Nedelec1WdgFiniteElement
974  virtual void CalcVShape(const IntegrationPoint &ip,
975  DenseMatrix &shape) const;
977  DenseMatrix &shape) const
978  { CalcVShape_ND(Trans, shape); }
979  virtual void CalcCurlShape(const IntegrationPoint &ip,
980  DenseMatrix &curl_shape) const;
982  DenseMatrix &I) const;
984  virtual void Project (VectorCoefficient &vc,
985  ElementTransformation &Trans, Vector &dofs) const;
986 
987  virtual void ProjectGrad(const FiniteElement &fe,
989  DenseMatrix &grad) const;
990 };
991 
992 
993 /// A 3D 1st order Nedelec element on a pyramid
995 {
996 private:
997  static const double tk[8][3];
998 
999 public:
1000  /// Construct the Nedelec1PyrFiniteElement
1002  virtual void CalcVShape(const IntegrationPoint &ip,
1003  DenseMatrix &shape) const;
1005  DenseMatrix &shape) const
1006  { CalcVShape_ND(Trans, shape); }
1007  virtual void CalcCurlShape(const IntegrationPoint &ip,
1008  DenseMatrix &curl_shape) const;
1010  DenseMatrix &I) const;
1011  using FiniteElement::Project;
1012  virtual void Project (VectorCoefficient &vc,
1013  ElementTransformation &Trans, Vector &dofs) const;
1014 
1015  virtual void ProjectGrad(const FiniteElement &fe,
1017  DenseMatrix &grad) const;
1018 };
1019 
1020 
1021 /// A 3D 0th order Raviert-Thomas element on a cube
1023 {
1024 private:
1025  static const double nk[6][3];
1026 
1027 public:
1028  /// Construct the RT0HexFiniteElement
1030 
1031  virtual void CalcVShape(const IntegrationPoint &ip,
1032  DenseMatrix &shape) const;
1033 
1035  DenseMatrix &shape) const
1036  { CalcVShape_RT(Trans, shape); }
1037 
1038  virtual void CalcDivShape(const IntegrationPoint &ip,
1039  Vector &divshape) const;
1040 
1042  DenseMatrix &I) const;
1043 
1044  using FiniteElement::Project;
1045 
1046  virtual void Project (VectorCoefficient &vc,
1047  ElementTransformation &Trans, Vector &dofs) const;
1048 };
1049 
1050 
1051 /// A 3D 1st order Raviert-Thomas element on a cube
1053 {
1054 private:
1055  static const double nk[36][3];
1056 
1057 public:
1058  /// Construct the RT1HexFiniteElement
1060 
1061  virtual void CalcVShape(const IntegrationPoint &ip,
1062  DenseMatrix &shape) const;
1063 
1065  DenseMatrix &shape) const
1066  { CalcVShape_RT(Trans, shape); }
1067 
1068  virtual void CalcDivShape(const IntegrationPoint &ip,
1069  Vector &divshape) const;
1070 
1072  DenseMatrix &I) const;
1073 
1074  using FiniteElement::Project;
1075 
1076  virtual void Project (VectorCoefficient &vc,
1077  ElementTransformation &Trans, Vector &dofs) const;
1078 };
1079 
1080 
1081 /// A 3D 0th order Raviert-Thomas element on a tetrahedron
1083 {
1084 private:
1085  static const double nk[4][3];
1086 
1087 public:
1088  /// Construct the RT0TetFiniteElement
1090 
1091  virtual void CalcVShape(const IntegrationPoint &ip,
1092  DenseMatrix &shape) const;
1093 
1095  DenseMatrix &shape) const
1096  { CalcVShape_RT(Trans, shape); }
1097 
1098  virtual void CalcDivShape(const IntegrationPoint &ip,
1099  Vector &divshape) const;
1100 
1102  DenseMatrix &I) const;
1103 
1104  using FiniteElement::Project;
1105 
1106  virtual void Project (VectorCoefficient &vc,
1107  ElementTransformation &Trans, Vector &dofs) const;
1108 };
1109 
1110 
1111 /// A 3D 0th order Raviert-Thomas element on a wedge
1113 {
1114 private:
1115  static const double nk[5][3];
1116 
1117 public:
1118  /// Construct the RT0WdgFiniteElement
1120 
1121  virtual void CalcVShape(const IntegrationPoint &ip,
1122  DenseMatrix &shape) const;
1123 
1125  DenseMatrix &shape) const
1126  { CalcVShape_RT(Trans, shape); }
1127 
1128  virtual void CalcDivShape(const IntegrationPoint &ip,
1129  Vector &divshape) const;
1130 
1132  DenseMatrix &I) const;
1133 
1134  using FiniteElement::Project;
1135 
1136  virtual void Project (VectorCoefficient &vc,
1137  ElementTransformation &Trans, Vector &dofs) const;
1138 
1139  virtual void ProjectCurl(const FiniteElement &fe,
1141  DenseMatrix &curl) const;
1142 };
1143 
1144 
1145 /// A 3D 0th order Raviert-Thomas element on a pyramid
1147 {
1148 private:
1149  static const double nk[5][3];
1150 
1151  // If true match RT0TetFiniteElement rather than RT_TetrahedronElement(0)
1152  bool rt0;
1153 
1154 public:
1155  /// Construct the RT0PyrFiniteElement
1156  RT0PyrFiniteElement(bool rt0tets = true);
1157 
1158  virtual void CalcVShape(const IntegrationPoint &ip,
1159  DenseMatrix &shape) const;
1160 
1162  DenseMatrix &shape) const
1163  { CalcVShape_RT(Trans, shape); }
1164 
1165  virtual void CalcDivShape(const IntegrationPoint &ip,
1166  Vector &divshape) const;
1167 
1169  DenseMatrix &I) const;
1170 
1171  using FiniteElement::Project;
1172 
1173  virtual void Project (VectorCoefficient &vc,
1174  ElementTransformation &Trans, Vector &dofs) const;
1175 
1176  virtual void ProjectCurl(const FiniteElement &fe,
1178  DenseMatrix &curl) const;
1179 };
1180 
1181 
1183 {
1184 public:
1185  /// Construct the RotTriLinearHexFiniteElement
1187  virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
1188  virtual void CalcDShape(const IntegrationPoint &ip,
1189  DenseMatrix &dshape) const;
1190 };
1191 
1192 
1193 } // namespace mfem
1194 
1195 #endif
Abstract class for all finite elements.
Definition: fe_base.hpp:235
A 2D 2nd order Raviart-Thomas vector element on a triangle.
RefinedLinear3DFiniteElement()
Construct the RefinedLinear3DFiniteElement.
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Arbitrary order L2 elements in 3D on a wedge.
Definition: fe_l2.hpp:140
RT0TriangleFiniteElement()
Construct the RT0TriangleFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual function which evaluates the values of all shape functions at a given point ip and stores the...
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
RT2TriangleFiniteElement()
Construct the RT2TriangleFiniteElement.
GaussQuad2DFiniteElement()
Construct the GaussQuad2DFiniteElement.
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 2D 1st order Raviart-Thomas vector element on a square.
A 1D constant element on a segment.
A 1D linear element with nodes at 1/3 and 2/3 (trace of RT1)
A 1D quadratic element with nodes at the Gaussian points (trace of RT2)
CrouzeixRaviartQuadFiniteElement()
Construct the CrouzeixRaviartQuadFiniteElement.
P2SegmentFiniteElement()
Construct the P2SegmentFiniteElement.
Linear3DFiniteElement()
Construct the Linear3DFiniteElement.
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
RT1TriangleFiniteElement()
Construct the RT1TriangleFiniteElement.
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
RefinedBiLinear2DFiniteElement()
Construct the RefinedBiLinear2DFiniteElement.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
P0WedgeFiniteElement()
Construct the P0WedgeFiniteElement.
A 3D 1st order Nedelec element on a tetrahedron.
Base class for vector Coefficients that optionally depend on time and space.
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Evaluate the curl of all shape functions of a vector finite element in reference space at the given p...
virtual void Project(Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const
Given a coefficient and a transformation, compute its projection (approximation) in the local finite ...
Definition: fe_base.cpp:125
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the matrix depends on it.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual function which evaluates the values of all shape functions at a given point ip and stores the...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Compute the discrete curl matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the matrix depends on it.
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
A 3D 1st order Nedelec element on a pyramid.
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the matrix depends on it.
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
PointFiniteElement()
Construct the PointFiniteElement.
P0WdgFiniteElement()
Construct the P0WdgFiniteElement.
void CalcVShape_RT(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe_base.cpp:903
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 2D 3rd order Raviart-Thomas vector element on a square.
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const
Compute the Hessian of second order partial derivatives at ip.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
RT0WdgFiniteElement()
Construct the RT0WdgFiniteElement.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
evaluate derivatives of shape function - constant 0
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
BiQuadratic3DFiniteElement()
Construct a quadratic FE on wedge.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
LagrangeHexFiniteElement(int degree)
Construct the LagrangeHexFiniteElement with the provided degree.
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
Data type dense matrix using column-major storage.
Definition: densemat.hpp:23
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
GaussBiQuad2DFiniteElement()
Construct the GaussBiQuad2DFiniteElement.
A 2D quadratic element on triangle with nodes at the vertices and midpoints of the triangle...
A 0D point finite element.
A 3D 0th order Raviert-Thomas element on a cube.
A 1D quadratic finite element with uniformly spaced nodes.
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
virtual function which evaluates the values of all partial derivatives of all shape functions at a gi...
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the matrix depends on it.
RefinedTriLinear3DFiniteElement()
Construct the RefinedTriLinear3DFiniteElement.
P0QuadFiniteElement()
Construct the P0QuadFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
Quadratic3DFiniteElement()
Construct the Quadratic3DFiniteElement.
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
A 3D constant element on a tetrahedron.
Intermediate class for finite elements whose basis functions return vector values.
Definition: fe_base.hpp:786
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 3D Crouzeix-Raviart element on the tetrahedron.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
RT1HexFiniteElement()
Construct the RT1HexFiniteElement.
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
A 2D bi-cubic element on a square with uniformly spaces nodes.
Class for cubic FE on wedge.
RT0QuadFiniteElement()
Construct the RT0QuadFiniteElement.
RT1QuadFiniteElement()
Construct the RT1QuadFiniteElement.
virtual void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const
Evaluate the Hessians of all shape functions of a scalar finite element in reference space at the giv...
Linear2DFiniteElement()
Construct the Linear2DFiniteElement.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 2D refined bi-linear FE on a square.
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Cubic2DFiniteElement()
Construct the Cubic2DFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
P0PyrFiniteElement()
Construct the P0PyrFiniteElement.
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
evaluate shape function - constant 1
virtual void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const
Evaluate the Hessians of all shape functions of a scalar finite element in reference space at the giv...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
Cubic3DFiniteElement()
Construct the Cubic3DFiniteElement.
Linear1DFiniteElement()
Construct the Linear1DFiniteElement.
A 1D element with uniform nodes.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
RT0HexFiniteElement()
Construct the RT0HexFiniteElement.
Class for standard nodal finite elements.
Definition: fe_base.hpp:706
CrouzeixRaviartFiniteElement()
Construct the CrouzeixRaviartFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Nedelec1PyrFiniteElement()
Construct the Nedelec1PyrFiniteElement.
virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const
Get the dofs associated with the given face. *dofs is set to an internal array of the local dofc on t...
A 2D bi-linear element on a square with nodes at the "Gaussian" points.
A 1D refined linear element.
A 2D refined linear element on a triangle.
A 2D refined linear element on a tetrahedron.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
A 2D constant element on a triangle.
RefinedLinear1DFiniteElement()
Construct the RefinedLinear1DFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
A linear element defined on a triangular prism.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
A 3D 0th order Raviert-Thomas element on a wedge.
Cubic1DFiniteElement()
Construct the Cubic1DFiniteElement.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
A 3D refined tri-linear element on a cube.
Nedelec1TetFiniteElement()
Construct the Nedelec1TetFiniteElement.
RT2QuadFiniteElement()
Construct the RT2QuadFiniteElement.
A 3D quadratic element on a tetrahedron with uniformly spaced nodes.
virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const
Get the dofs associated with the given face. *dofs is set to an internal array of the local dofc on t...
GaussLinear2DFiniteElement()
Construct the GaussLinear2DFiniteElement.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the matrix depends on it.
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
virtual void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const
Compute the discrete curl matrix from the given FiniteElement onto 'this' FiniteElement. The ElementTransformation is included to support cases when the matrix depends on it.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
P1TetNonConfFiniteElement()
Construct the P1TetNonConfFiniteElement.
P0SegmentFiniteElement(int Ord=0)
Construct the P0SegmentFiniteElement with dummy order Ord.
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
A 3D linear element on a tetrahedron with nodes at the vertices of the tetrahedron.
Nedelec1HexFiniteElement()
Construct the Nedelec1HexFiniteElement.
A 1D linear element with nodes on the endpoints.
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
A 2D Crouzeix-Raviart finite element on square.
A 2D linear element on triangle with nodes at the vertices of the triangle.
A quadratic element on triangle with nodes at the "Gaussian" points.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
virtual function which evaluates the values of all partial derivatives of all shape functions at a gi...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const
Evaluate the Hessians of all shape functions of a scalar finite element in reference space at the giv...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
RT0TetFiniteElement()
Construct the RT0TetFiniteElement.
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
LinearPyramidFiniteElement()
Construct the LinearPyramidFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
P0HexFiniteElement()
Construct the P0HexFiniteElement.
A 2D 2nd order Raviart-Thomas vector element on a square.
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
P1SegmentFiniteElement()
Construct the P1SegmentFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
BiCubic2DFiniteElement()
Construct the BiCubic2DFiniteElement.
void CalcVShape_ND(ElementTransformation &Trans, DenseMatrix &shape) const
Definition: fe_base.cpp:915
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
A 3D 1st order Raviert-Thomas element on a cube.
BiCubic3DFiniteElement()
Construct a cubic FE on wedge.
A 2D bi-quadratic element on a square with uniformly spaced nodes.
P0TetFiniteElement()
Construct the P0TetFiniteElement.
A 3D tri-linear element on a cube with nodes at the vertices of the cube.
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
A 3D constant element on a wedge.
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Evaluate the curl of all shape functions of a vector finite element in reference space at the given p...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
A 3D constant element on a cube.
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Class for quadratic FE on wedge.
A 2D bi-linear element on a square with nodes at the vertices of the square.
Class for linear FE on wedge.
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Evaluate the curl of all shape functions of a vector finite element in reference space at the given p...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
Quad1DFiniteElement()
Construct the Quad1DFiniteElement.
Nedelec1WdgFiniteElement()
Construct the Nedelec1WdgFiniteElement.
A linear element on a triangle with nodes at the 3 "Gaussian" points.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
virtual function which evaluates the values of all partial derivatives of all shape functions at a gi...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Class for integration point with weight.
Definition: intrules.hpp:25
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 2D bi-quadratic element on a square with nodes at the 9 "Gaussian" points.
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
RT0PyrFiniteElement(bool rt0tets=true)
Construct the RT0PyrFiniteElement.
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
A 3D 1st order Nedelec element on a wedge.
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
RefinedLinear2DFiniteElement()
Construct the RefinedLinear2DFiniteElement.
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
BiQuad2DFiniteElement()
Construct the BiQuad2DFiniteElement.
virtual void GetFaceDofs(int face, int **dofs, int *ndofs) const
Get the dofs associated with the given face. *dofs is set to an internal array of the local dofc on t...
Tensor products of 1D Lagrange1DFiniteElement (only degree 2 is functional)
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
P0TriangleFiniteElement()
Construct the P0TriangleFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
A 3D constant element on a pyramid.
A 2D constant element on a square.
A linear element defined on a square pyramid.
BiLinear3DFiniteElement()
Construct a linear FE on wedge.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
Vector data type.
Definition: vector.hpp:60
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
A 2D cubic element on a triangle with uniformly spaced nodes.
A 1D cubic element with uniformly spaced nodes.
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const
Evaluate the curl of all shape functions of a vector finite element in reference space at the given p...
Quad2DFiniteElement()
Construct the Quad2DFiniteElement.
A 2D 1st order Raviart-Thomas vector element on a triangle.
A 2D linear element on a square with 3 nodes at the vertices of the lower left triangle.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
A 3D 0th order Raviert-Thomas element on a pyramid.
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
GaussBiLinear2DFiniteElement()
Construct the FiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual function which evaluates the values of all shape functions at a given point ip and stores the...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
Arbitrary order H1 elements in 3D on a wedge.
Definition: fe_h1.hpp:130
A 0th order L2 element on a Wedge.
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
A 3D 0th order Raviert-Thomas element on a tetrahedron.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
P1OnQuadFiniteElement()
Construct the P1OnQuadFiniteElement.
A 2D 3rd order Raviart-Thomas vector element on a triangle.
RotTriLinearHexFiniteElement()
Construct the RotTriLinearHexFiniteElement.
LinearWedgeFiniteElement()
Construct the LinearWedgeFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
virtual void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
virtual void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
virtual void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
TriLinear3DFiniteElement()
Construct the TriLinear3DFiniteElement.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
BiLinear2DFiniteElement()
Construct the BiLinear2DFiniteElement.
A 3D 1st order Nedelec element on a cube.
A 2D Crouzeix-Raviart element on triangle.
virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
Lagrange1DFiniteElement(int degree)
Construct the Lagrange1DFiniteElement with the provided degree.
virtual void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
virtual void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
virtual void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...
virtual void ProjectDelta(int vertex, Vector &dofs) const
Project a delta function centered on the given vertex in the local finite dimensional space represent...