MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
fe_fixed_order.hpp
Go to the documentation of this file.
1// Copyright (c) 2010-2025, 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
19namespace mfem
20{
21
22/// A 0D point finite element
24{
25public:
26 /// Construct the PointFiniteElement
28
29 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
30
31 void CalcDShape(const IntegrationPoint &ip,
32 DenseMatrix &dshape) const override;
33};
34
35/// A 1D linear element with nodes on the endpoints
37{
38public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
52 DenseMatrix &dshape) const override;
53};
54
55/// A 2D linear element on triangle with nodes at the vertices of the triangle
57{
58public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
72 DenseMatrix &dshape) const override;
73 void ProjectDelta(int vertex, Vector &dofs) const override
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{
80public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
94 DenseMatrix &dshape) const override;
95 void CalcHessian(const IntegrationPoint &ip,
96 DenseMatrix &h) const override;
97 void ProjectDelta(int vertex, Vector &dofs) const override
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{
104public:
105 /// Construct the GaussLinear2DFiniteElement
107 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
108 void CalcDShape(const IntegrationPoint &ip,
109 DenseMatrix &dshape) const override;
110 void ProjectDelta(int vertex, Vector &dofs) const override;
111};
112
113/// A 2D bi-linear element on a square with nodes at the "Gaussian" points
115{
116private:
117 static const real_t p[2];
118
119public:
120 /// Construct the FiniteElement
122 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
123 void CalcDShape(const IntegrationPoint &ip,
124 DenseMatrix &dshape) const override;
125 void ProjectDelta(int vertex, Vector &dofs) const override;
126};
127
128/** @brief A 2D linear element on a square with 3 nodes at the
129 vertices of the lower left triangle */
131{
132public:
133 /// Construct the P1OnQuadFiniteElement
135 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
136 void CalcDShape(const IntegrationPoint &ip,
137 DenseMatrix &dshape) const override;
138 void ProjectDelta(int vertex, Vector &dofs) const override
139 { dofs = 1.0; }
140};
141
142/// A 1D quadratic finite element with uniformly spaced nodes
144{
145public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
159 DenseMatrix &dshape) const override;
160};
161
162/** @brief A 2D quadratic element on triangle with nodes at the
163 vertices and midpoints of the triangle. */
165{
166public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
180 DenseMatrix &dshape) const override;
181
182 void CalcHessian(const IntegrationPoint &ip,
183 DenseMatrix &h) const override;
184 void ProjectDelta(int vertex, Vector &dofs) const override;
185};
186
187/// A quadratic element on triangle with nodes at the "Gaussian" points
189{
190private:
191 static const real_t p[2];
192 DenseMatrix A;
193 mutable DenseMatrix D;
194 mutable Vector pol;
195public:
196 /// Construct the GaussQuad2DFiniteElement
198 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
199 void CalcDShape(const IntegrationPoint &ip,
200 DenseMatrix &dshape) const override;
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{
207public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
221 DenseMatrix &dshape) const override;
222 void ProjectDelta(int vertex, Vector &dofs) const override;
223};
224
225
226/// A 2D bi-quadratic element on a square with nodes at the 9 "Gaussian" points
228{
229public:
230 /// Construct the GaussBiQuad2DFiniteElement
232 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
233 void CalcDShape(const IntegrationPoint &ip,
234 DenseMatrix &dshape) const override;
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{
242public:
243 /// Construct the BiCubic2DFiniteElement
245 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
246 void CalcDShape(const IntegrationPoint &ip,
247 DenseMatrix &dshape) const override;
248
249 /// Compute the Hessian of second order partial derivatives at @a ip.
250 void CalcHessian(const IntegrationPoint &ip,
251 DenseMatrix &h) const override;
252};
253
254/// A 1D cubic element with uniformly spaced nodes
256{
257public:
258 /// Construct the Cubic1DFiniteElement
260
261 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
262
263 void CalcDShape(const IntegrationPoint &ip,
264 DenseMatrix &dshape) const override;
265};
266
267/// A 2D cubic element on a triangle with uniformly spaced nodes
269{
270public:
271 /// Construct the Cubic2DFiniteElement
273
274 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
275
276 void CalcDShape(const IntegrationPoint &ip,
277 DenseMatrix &dshape) const override;
278
279 void CalcHessian(const IntegrationPoint &ip,
280 DenseMatrix &h) const override;
281};
282
283/// A 3D cubic element on a tetrahedron with 20 nodes at the thirds of the
284/// tetrahedron
286{
287public:
288 /// Construct the Cubic3DFiniteElement
290
291 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
292
293 void CalcDShape(const IntegrationPoint &ip,
294 DenseMatrix &dshape) const override;
295};
296
297/// A linear element defined on a triangular prism
299{
300public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
314 DenseMatrix &dshape) const override;
315
316 void ProjectDelta(int vertex, Vector &dofs) const override
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 void GetFaceDofs(int face, int **dofs, int *ndofs) const override;
324};
325
326/// A linear element defined on a square pyramid
328{
329public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
343 DenseMatrix &dshape) const override;
344
345 void ProjectDelta(int vertex, Vector &dofs) const override
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 void GetFaceDofs(int face, int **dofs, int *ndofs) const override;
353};
354
355/// A 2D constant element on a triangle
357{
358public:
359 /// Construct the P0TriangleFiniteElement
361
362 /// evaluate shape function - constant 1
363 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
364
365 /// evaluate derivatives of shape function - constant 0
366 void CalcDShape(const IntegrationPoint &ip,
367 DenseMatrix &dshape) const override;
368 void ProjectDelta(int vertex, Vector &dofs) const override
369 { dofs(0) = 1.0; }
370};
371
372
373/// A 2D constant element on a square
375{
376public:
377 /// Construct the P0QuadFiniteElement
379 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
380 void CalcDShape(const IntegrationPoint &ip,
381 DenseMatrix &dshape) const override;
382 void ProjectDelta(int vertex, Vector &dofs) const override
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{
391public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
405 DenseMatrix &dshape) const override;
406
407 void ProjectDelta(int vertex, Vector &dofs) const override
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 void GetFaceDofs(int face, int **dofs, int *ndofs) const override;
415};
416
417/// A 3D quadratic element on a tetrahedron with uniformly spaced nodes
419{
420public:
421 /// Construct the Quadratic3DFiniteElement
423
424 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
425
426 void CalcDShape(const IntegrationPoint &ip,
427 DenseMatrix &dshape) const override;
428};
429
430/// A 3D tri-linear element on a cube with nodes at the vertices of the cube
432{
433public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
447 DenseMatrix &dshape) const override;
448
449 void ProjectDelta(int vertex, Vector &dofs) const override
450 { dofs = 0.0; dofs(vertex) = 1.0; }
451};
452
453
454/// A 2D Crouzeix-Raviart element on triangle
456{
457public:
458 /// Construct the CrouzeixRaviartFiniteElement
460 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
461 void CalcDShape(const IntegrationPoint &ip,
462 DenseMatrix &dshape) const override;
463 void ProjectDelta(int vertex, Vector &dofs) const override
464 { dofs = 1.0; }
465};
466
467/// A 2D Crouzeix-Raviart finite element on square
469{
470public:
471 /// Construct the CrouzeixRaviartQuadFiniteElement
473 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
474 void CalcDShape(const IntegrationPoint &ip,
475 DenseMatrix &dshape) const override;
476};
477
478
479/// A 1D constant element on a segment
481{
482public:
483 /// Construct the P0SegmentFiniteElement with dummy order @a Ord
484 P0SegmentFiniteElement(int Ord = 0);
485 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
486 void CalcDShape(const IntegrationPoint &ip,
487 DenseMatrix &dshape) const override;
488};
489
490/** @brief A 2D 1st order Raviart-Thomas vector element on a triangle */
492{
493private:
494 static const real_t nk[3][2];
495
496public:
497 /// Construct the RT0TriangleFiniteElement
499
500 void CalcVShape(const IntegrationPoint &ip,
501 DenseMatrix &shape) const override;
502
504 DenseMatrix &shape) const override
505 { CalcVShape_RT(Trans, shape); }
506
507 void CalcDivShape(const IntegrationPoint &ip,
508 Vector &divshape) const override;
509
511 DenseMatrix &I) const override;
512
514
515 void Project(VectorCoefficient &vc,
516 ElementTransformation &Trans, Vector &dofs) const override;
517};
518
519/** @brief A 2D 1st order Raviart-Thomas vector element on a square*/
521{
522private:
523 static const real_t nk[4][2];
524
525public:
526 /// Construct the RT0QuadFiniteElement
528
529 void CalcVShape(const IntegrationPoint &ip,
530 DenseMatrix &shape) const override;
531
533 DenseMatrix &shape) const override
534 { CalcVShape_RT(Trans, shape); }
535
536 void CalcDivShape(const IntegrationPoint &ip,
537 Vector &divshape) const override;
538
540 DenseMatrix &I) const override;
541
543
544 void Project(VectorCoefficient &vc,
545 ElementTransformation &Trans, Vector &dofs) const override;
546};
547
548/** @brief A 2D 2nd order Raviart-Thomas vector element on a triangle */
550{
551private:
552 static const real_t nk[8][2];
553
554public:
555 /// Construct the RT1TriangleFiniteElement
557
558 void CalcVShape(const IntegrationPoint &ip,
559 DenseMatrix &shape) const override;
560
562 DenseMatrix &shape) const override
563 { CalcVShape_RT(Trans, shape); }
564
565 void CalcDivShape(const IntegrationPoint &ip,
566 Vector &divshape) const override;
567
569 DenseMatrix &I) const override;
570
572
573 void Project(VectorCoefficient &vc,
574 ElementTransformation &Trans, Vector &dofs) const override;
575};
576
577/** @brief A 2D 2nd order Raviart-Thomas vector element on a square */
579{
580private:
581 static const real_t nk[12][2];
582
583public:
584 /// Construct the RT1QuadFiniteElement
586
587 void CalcVShape(const IntegrationPoint &ip,
588 DenseMatrix &shape) const override;
589
591 DenseMatrix &shape) const override
592 { CalcVShape_RT(Trans, shape); }
593
594 void CalcDivShape(const IntegrationPoint &ip,
595 Vector &divshape) const override;
596
598 DenseMatrix &I) const override;
599
601
602 void Project(VectorCoefficient &vc,
603 ElementTransformation &Trans, Vector &dofs) const override;
604};
605
606/** @brief A 2D 3rd order Raviart-Thomas vector element on a triangle */
608{
609private:
610 static const real_t M[15][15];
611public:
612 /// Construct the RT2TriangleFiniteElement
614
615 void CalcVShape(const IntegrationPoint &ip,
616 DenseMatrix &shape) const override;
617
619 DenseMatrix &shape) const override
620 { CalcVShape_RT(Trans, shape); }
621
622 void CalcDivShape(const IntegrationPoint &ip,
623 Vector &divshape) const override;
624};
625
626/** @brief A 2D 3rd order Raviart-Thomas vector element on a square */
628{
629private:
630 static const real_t nk[24][2];
631 static const real_t pt[4];
632 static const real_t dpt[3];
633
634public:
635 /// Construct the RT2QuadFiniteElement
637
638 void CalcVShape(const IntegrationPoint &ip,
639 DenseMatrix &shape) const override;
640
642 DenseMatrix &shape) const override
643 { CalcVShape_RT(Trans, shape); }
644
645 void CalcDivShape(const IntegrationPoint &ip,
646 Vector &divshape) const override;
647
649 DenseMatrix &I) const override;
650
652
653 void Project(VectorCoefficient &vc,
654 ElementTransformation &Trans, Vector &dofs) const override;
655};
656
657/// A 1D linear element with nodes at 1/3 and 2/3 (trace of RT1)
659{
660public:
661 /// Construct the P1SegmentFiniteElement
663 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
664 void CalcDShape(const IntegrationPoint &ip,
665 DenseMatrix &dshape) const override;
666};
667
668/// A 1D quadratic element with nodes at the Gaussian points (trace of RT2)
670{
671public:
672 /// Construct the P2SegmentFiniteElement
674 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
675 void CalcDShape(const IntegrationPoint &ip,
676 DenseMatrix &dshape) const override;
677};
678
679/// A 1D element with uniform nodes
681{
682private:
683 Vector rwk;
684#ifndef MFEM_THREAD_SAFE
685 mutable Vector rxxk;
686#endif
687public:
688 /// Construct the Lagrange1DFiniteElement with the provided @a degree
689 Lagrange1DFiniteElement(int degree);
690 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
691 void CalcDShape(const IntegrationPoint &ip,
692 DenseMatrix &dshape) const override;
693};
694
695/// A 3D Crouzeix-Raviart element on the tetrahedron.
697{
698public:
699 /// Construct the P1TetNonConfFiniteElement
701 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
702 void CalcDShape(const IntegrationPoint &ip,
703 DenseMatrix &dshape) const override;
704};
705
706/// A 3D constant element on a tetrahedron
708{
709public:
710 /// Construct the P0TetFiniteElement
712 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
713 void CalcDShape(const IntegrationPoint &ip,
714 DenseMatrix &dshape) const override;
715 void ProjectDelta(int vertex, Vector &dofs) const override
716 { dofs(0) = 1.0; }
717};
718
719/// A 3D constant element on a cube
721{
722public:
723 /// Construct the P0HexFiniteElement
725 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
726 void CalcDShape(const IntegrationPoint &ip,
727 DenseMatrix &dshape) const override;
728 void ProjectDelta(int vertex, Vector &dofs) const override
729 { dofs(0) = 1.0; }
730};
731
732/// A 3D constant element on a wedge
734{
735public:
736 /// Construct the P0WdgFiniteElement
738 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
739 void CalcDShape(const IntegrationPoint &ip,
740 DenseMatrix &dshape) const override;
741 void ProjectDelta(int vertex, Vector &dofs) const override
742 { dofs(0) = 1.0; }
743};
744
745/// A 3D constant element on a pyramid
747{
748public:
749 /// Construct the P0PyrFiniteElement
751 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
752 void CalcDShape(const IntegrationPoint &ip,
753 DenseMatrix &dshape) const override;
754 void ProjectDelta(int vertex, Vector &dofs) const override
755 { dofs(0) = 1.0; }
756};
757
758/** @brief Tensor products of 1D Lagrange1DFiniteElement
759 (only degree 2 is functional) */
761{
762private:
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
771public:
772 /// Construct the LagrangeHexFiniteElement with the provided @a degree
773 LagrangeHexFiniteElement(int degree);
774 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
775 void CalcDShape(const IntegrationPoint &ip,
776 DenseMatrix &dshape) const override;
778};
779
780
781/// A 1D refined linear element
783{
784public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
798 DenseMatrix &dshape) const override;
799};
800
801/// A 2D refined linear element on a triangle
803{
804public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
818 DenseMatrix &dshape) const override;
819};
820
821/// A 2D refined linear element on a tetrahedron
823{
824public:
825 /// Construct the RefinedLinear3DFiniteElement
827
828 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
829
830 void CalcDShape(const IntegrationPoint &ip,
831 DenseMatrix &dshape) const override;
832};
833
834/// A 2D refined bi-linear FE on a square
836{
837public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
851 DenseMatrix &dshape) const override;
852};
853
854/// A 3D refined tri-linear element on a cube
856{
857public:
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 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
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 void CalcDShape(const IntegrationPoint &ip,
871 DenseMatrix &dshape) const override;
872};
873
874
875/// Class for linear FE on wedge
877{
878public:
879 /// Construct a linear FE on wedge
881};
882
883/// Class for quadratic FE on wedge
885{
886public:
887 /// Construct a quadratic FE on wedge
889};
890
891/// Class for cubic FE on wedge
893{
894public:
895 /// Construct a cubic FE on wedge
897};
898
899
900/// A 0th order L2 element on a Wedge
902{
903public:
904 /// Construct the P0WedgeFiniteElement
906};
907
908
909/// A 3D 1st order Nedelec element on a cube
911{
912private:
913 static const real_t tk[12][3];
914
915public:
916 /// Construct the Nedelec1HexFiniteElement
918 void CalcVShape(const IntegrationPoint &ip,
919 DenseMatrix &shape) const override;
921 DenseMatrix &shape) const override
922 { CalcVShape_ND(Trans, shape); }
923 void CalcCurlShape(const IntegrationPoint &ip,
924 DenseMatrix &curl_shape) const override;
926 DenseMatrix &I) const override;
928 void Project(VectorCoefficient &vc,
929 ElementTransformation &Trans, Vector &dofs) const override;
930
931 void ProjectGrad(const FiniteElement &fe,
933 DenseMatrix &grad) const override;
934};
935
936
937/// A 3D 1st order Nedelec element on a tetrahedron
939{
940private:
941 static const real_t tk[6][3];
942
943public:
944 /// Construct the Nedelec1TetFiniteElement
946 void CalcVShape(const IntegrationPoint &ip,
947 DenseMatrix &shape) const override;
949 DenseMatrix &shape) const override
950 { CalcVShape_ND(Trans, shape); }
951 void CalcCurlShape(const IntegrationPoint &ip,
952 DenseMatrix &curl_shape) const override;
954 DenseMatrix &I) const override;
956 void Project(VectorCoefficient &vc,
957 ElementTransformation &Trans, Vector &dofs) const override;
958
959 void ProjectGrad(const FiniteElement &fe,
961 DenseMatrix &grad) const override;
962};
963
964
965/// A 3D 1st order Nedelec element on a wedge
967{
968private:
969 static const real_t tk[9][3];
970
971public:
972 /// Construct the Nedelec1WdgFiniteElement
974 void CalcVShape(const IntegrationPoint &ip,
975 DenseMatrix &shape) const override;
977 DenseMatrix &shape) const override
978 { CalcVShape_ND(Trans, shape); }
979 void CalcCurlShape(const IntegrationPoint &ip,
980 DenseMatrix &curl_shape) const override;
982 DenseMatrix &I) const override;
984 void Project(VectorCoefficient &vc,
985 ElementTransformation &Trans, Vector &dofs) const override;
986
987 void ProjectGrad(const FiniteElement &fe,
989 DenseMatrix &grad) const override;
990};
991
992
993/// A 3D 1st order Nedelec element on a pyramid
995{
996private:
997 static const real_t tk[8][3];
998
999public:
1000 /// Construct the Nedelec1PyrFiniteElement
1002 void CalcVShape(const IntegrationPoint &ip,
1003 DenseMatrix &shape) const override;
1005 DenseMatrix &shape) const override
1006 { CalcVShape_ND(Trans, shape); }
1007 void CalcCurlShape(const IntegrationPoint &ip,
1008 DenseMatrix &curl_shape) const override;
1010 DenseMatrix &I) const override;
1012 void Project(VectorCoefficient &vc,
1013 ElementTransformation &Trans, Vector &dofs) const override;
1014
1015 void ProjectGrad(const FiniteElement &fe,
1016 ElementTransformation &Trans,
1017 DenseMatrix &grad) const override;
1018};
1019
1020
1021/// A 3D 2nd order Nedelec element on a pyramid
1023{
1024private:
1025 static const real_t tk[28][3];
1026
1027public:
1028 /// Construct the Nedelec2PyrFiniteElement
1030 virtual void CalcVShape(const IntegrationPoint &ip,
1031 DenseMatrix &shape) const;
1033 DenseMatrix &shape) const
1034 { CalcVShape_ND(Trans, shape); }
1035 virtual void CalcCurlShape(const IntegrationPoint &ip,
1036 DenseMatrix &curl_shape) const;
1037 virtual void GetLocalInterpolation (ElementTransformation &Trans,
1038 DenseMatrix &I) const;
1040 virtual void Project (VectorCoefficient &vc,
1041 ElementTransformation &Trans, Vector &dofs) const;
1042
1043 virtual void ProjectGrad(const FiniteElement &fe,
1044 ElementTransformation &Trans,
1045 DenseMatrix &grad) const;
1046};
1047
1048
1049/// A 3D 0th order Raviert-Thomas element on a cube
1051{
1052private:
1053 static const real_t nk[6][3];
1054
1055public:
1056 /// Construct the RT0HexFiniteElement
1058
1059 void CalcVShape(const IntegrationPoint &ip,
1060 DenseMatrix &shape) const override;
1061
1063 DenseMatrix &shape) const override
1064 { CalcVShape_RT(Trans, shape); }
1065
1066 void CalcDivShape(const IntegrationPoint &ip,
1067 Vector &divshape) const override;
1068
1070 DenseMatrix &I) const override;
1071
1073
1074 void Project(VectorCoefficient &vc,
1075 ElementTransformation &Trans, Vector &dofs) const override;
1076};
1077
1078
1079/// A 3D 1st order Raviert-Thomas element on a cube
1081{
1082private:
1083 static const real_t nk[36][3];
1084
1085public:
1086 /// Construct the RT1HexFiniteElement
1088
1089 void CalcVShape(const IntegrationPoint &ip,
1090 DenseMatrix &shape) const override;
1091
1093 DenseMatrix &shape) const override
1094 { CalcVShape_RT(Trans, shape); }
1095
1096 void CalcDivShape(const IntegrationPoint &ip,
1097 Vector &divshape) const override;
1098
1100 DenseMatrix &I) const override;
1101
1103
1104 void Project(VectorCoefficient &vc,
1105 ElementTransformation &Trans, Vector &dofs) const override;
1106};
1107
1108
1109/// A 3D 0th order Raviert-Thomas element on a tetrahedron
1111{
1112private:
1113 static const real_t nk[4][3];
1114
1115public:
1116 /// Construct the RT0TetFiniteElement
1118
1119 void CalcVShape(const IntegrationPoint &ip,
1120 DenseMatrix &shape) const override;
1121
1123 DenseMatrix &shape) const override
1124 { CalcVShape_RT(Trans, shape); }
1125
1126 void CalcDivShape(const IntegrationPoint &ip,
1127 Vector &divshape) const override;
1128
1130 DenseMatrix &I) const override;
1131
1133
1134 void Project(VectorCoefficient &vc,
1135 ElementTransformation &Trans, Vector &dofs) const override;
1136};
1137
1138
1139/// A 3D 0th order Raviert-Thomas element on a wedge
1141{
1142private:
1143 static const real_t nk[5][3];
1144
1145public:
1146 /// Construct the RT0WdgFiniteElement
1148
1149 void CalcVShape(const IntegrationPoint &ip,
1150 DenseMatrix &shape) const override;
1151
1153 DenseMatrix &shape) const override
1154 { CalcVShape_RT(Trans, shape); }
1155
1156 void CalcDivShape(const IntegrationPoint &ip,
1157 Vector &divshape) const override;
1158
1160 DenseMatrix &I) const override;
1161
1163
1164 void Project(VectorCoefficient &vc,
1165 ElementTransformation &Trans, Vector &dofs) const override;
1166
1167 void ProjectCurl(const FiniteElement &fe,
1168 ElementTransformation &Trans,
1169 DenseMatrix &curl) const override;
1170};
1171
1172
1173/// A 3D 0th order Raviert-Thomas element on a pyramid
1175{
1176private:
1177 static const real_t nk[5][3];
1178
1179 // If true match RT0TetFiniteElement rather than RT_TetrahedronElement(0)
1180 bool rt0;
1181
1182public:
1183 /// Construct the RT0PyrFiniteElement
1184 RT0PyrFiniteElement(bool rt0tets = true);
1185
1186 void CalcVShape(const IntegrationPoint &ip,
1187 DenseMatrix &shape) const override;
1188
1190 DenseMatrix &shape) const override
1191 { CalcVShape_RT(Trans, shape); }
1192
1193 void CalcDivShape(const IntegrationPoint &ip,
1194 Vector &divshape) const override;
1195
1197 DenseMatrix &I) const override;
1198
1200
1201 void Project(VectorCoefficient &vc,
1202 ElementTransformation &Trans, Vector &dofs) const override;
1203
1204 void ProjectCurl(const FiniteElement &fe,
1205 ElementTransformation &Trans,
1206 DenseMatrix &curl) const override;
1207};
1208
1209
1211{
1212public:
1213 /// Construct the RotTriLinearHexFiniteElement
1215 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
1216 void CalcDShape(const IntegrationPoint &ip,
1217 DenseMatrix &dshape) const override;
1218};
1219
1220
1221} // namespace mfem
1222
1223#endif
A 2D bi-cubic element on a square with uniformly spaces nodes.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
BiCubic2DFiniteElement()
Construct the BiCubic2DFiniteElement.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const override
Compute the Hessian of second order partial derivatives at ip.
Class for cubic FE on wedge.
BiCubic3DFiniteElement()
Construct a cubic FE on wedge.
A 2D bi-linear element on a square with nodes at the vertices of the square.
void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const override
Evaluate the Hessians of all shape functions of a scalar finite element in reference space at the giv...
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
BiLinear2DFiniteElement()
Construct the BiLinear2DFiniteElement.
Class for linear FE on wedge.
BiLinear3DFiniteElement()
Construct a linear FE on wedge.
A 2D bi-quadratic element on a square with uniformly spaced nodes.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
BiQuad2DFiniteElement()
Construct the BiQuad2DFiniteElement.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Class for quadratic FE on wedge.
BiQuadratic3DFiniteElement()
Construct a quadratic FE on wedge.
A 2D Crouzeix-Raviart element on triangle.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
CrouzeixRaviartFiniteElement()
Construct the CrouzeixRaviartFiniteElement.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
A 2D Crouzeix-Raviart finite element on square.
CrouzeixRaviartQuadFiniteElement()
Construct the CrouzeixRaviartQuadFiniteElement.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 1D cubic element with uniformly spaced nodes.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
Cubic1DFiniteElement()
Construct the Cubic1DFiniteElement.
A 2D cubic element on a triangle with uniformly spaced nodes.
Cubic2DFiniteElement()
Construct the Cubic2DFiniteElement.
void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const override
Evaluate the Hessians of all shape functions of a scalar finite element in reference space at the giv...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
Cubic3DFiniteElement()
Construct the Cubic3DFiniteElement.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
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:24
Abstract class for all finite elements.
Definition fe_base.hpp:244
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:126
A 2D bi-linear element on a square with nodes at the "Gaussian" points.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
GaussBiLinear2DFiniteElement()
Construct the FiniteElement.
A 2D bi-quadratic element on a square with nodes at the 9 "Gaussian" points.
GaussBiQuad2DFiniteElement()
Construct the GaussBiQuad2DFiniteElement.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
A linear element on a triangle with nodes at the 3 "Gaussian" points.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
GaussLinear2DFiniteElement()
Construct the GaussLinear2DFiniteElement.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A quadratic element on triangle with nodes at the "Gaussian" points.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
GaussQuad2DFiniteElement()
Construct the GaussQuad2DFiniteElement.
Arbitrary order H1 elements in 3D on a wedge.
Definition fe_h1.hpp:132
Class for integration point with weight.
Definition intrules.hpp:35
Arbitrary order L2 elements in 3D on a wedge.
Definition fe_l2.hpp:167
A 1D element with uniform nodes.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
Lagrange1DFiniteElement(int degree)
Construct the Lagrange1DFiniteElement with the provided degree.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
Tensor products of 1D Lagrange1DFiniteElement (only degree 2 is functional)
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
LagrangeHexFiniteElement(int degree)
Construct the LagrangeHexFiniteElement with the provided degree.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
A 1D linear element with nodes on the endpoints.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Linear1DFiniteElement()
Construct the Linear1DFiniteElement.
A 2D linear element on triangle with nodes at the vertices of the triangle.
Linear2DFiniteElement()
Construct the Linear2DFiniteElement.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
A 3D linear element on a tetrahedron with nodes at the vertices of the tetrahedron.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
virtual function which evaluates the values of all partial derivatives of all shape functions at a gi...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
virtual function which evaluates the values of all shape functions at a given point ip and stores the...
void GetFaceDofs(int face, int **dofs, int *ndofs) const override
Get the dofs associated with the given face. *dofs is set to an internal array of the local dofc on t...
Linear3DFiniteElement()
Construct the Linear3DFiniteElement.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
A linear element defined on a square pyramid.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
virtual function which evaluates the values of all partial derivatives of all shape functions at a gi...
LinearPyramidFiniteElement()
Construct the LinearPyramidFiniteElement.
void GetFaceDofs(int face, int **dofs, int *ndofs) const override
Get the dofs associated with the given face. *dofs is set to an internal array of the local dofc on t...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
virtual function which evaluates the values of all shape functions at a given point ip and stores the...
A linear element defined on a triangular prism.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void GetFaceDofs(int face, int **dofs, int *ndofs) const override
Get the dofs associated with the given face. *dofs is set to an internal array of the local dofc on t...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
virtual function which evaluates the values of all partial derivatives of all shape functions at a gi...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
virtual function which evaluates the values of all shape functions at a given point ip and stores the...
LinearWedgeFiniteElement()
Construct the LinearWedgeFiniteElement.
A 3D 1st order Nedelec element on a cube.
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
Nedelec1HexFiniteElement()
Construct the Nedelec1HexFiniteElement.
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const override
Evaluate the curl of all shape functions of a vector finite element in reference space at the given p...
void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const override
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement....
A 3D 1st order Nedelec element on a pyramid.
Nedelec1PyrFiniteElement()
Construct the Nedelec1PyrFiniteElement.
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const override
Evaluate the curl of all shape functions of a vector finite element in reference space at the given p...
void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const override
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement....
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
A 3D 1st order Nedelec element on a tetrahedron.
Nedelec1TetFiniteElement()
Construct the Nedelec1TetFiniteElement.
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const override
Evaluate the curl of all shape functions of a vector finite element in reference space at the given p...
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const override
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement....
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
A 3D 1st order Nedelec element on a wedge.
void ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const override
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement....
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void CalcCurlShape(const IntegrationPoint &ip, DenseMatrix &curl_shape) const override
Evaluate the curl of all shape functions of a vector finite element in reference space at the given p...
Nedelec1WdgFiniteElement()
Construct the Nedelec1WdgFiniteElement.
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
A 3D 2nd order Nedelec element on a pyramid.
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...
Nedelec2PyrFiniteElement()
Construct the Nedelec2PyrFiniteElement.
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 ProjectGrad(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &grad) const
Compute the discrete gradient matrix from the given FiniteElement onto 'this' FiniteElement....
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 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...
Class for standard nodal finite elements.
Definition fe_base.hpp:721
A 3D constant element on a cube.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
P0HexFiniteElement()
Construct the P0HexFiniteElement.
A 3D constant element on a pyramid.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
P0PyrFiniteElement()
Construct the P0PyrFiniteElement.
A 2D constant element on a square.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
P0QuadFiniteElement()
Construct the P0QuadFiniteElement.
A 1D constant element on a segment.
P0SegmentFiniteElement(int Ord=0)
Construct the P0SegmentFiniteElement with dummy order Ord.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
A 3D constant element on a tetrahedron.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
P0TetFiniteElement()
Construct the P0TetFiniteElement.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 2D constant element on a triangle.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
evaluate derivatives of shape function - constant 0
P0TriangleFiniteElement()
Construct the P0TriangleFiniteElement.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
evaluate shape function - constant 1
A 3D constant element on a wedge.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
P0WdgFiniteElement()
Construct the P0WdgFiniteElement.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 0th order L2 element on a Wedge.
P0WedgeFiniteElement()
Construct the P0WedgeFiniteElement.
A 2D linear element on a square with 3 nodes at the vertices of the lower left triangle.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
P1OnQuadFiniteElement()
Construct the P1OnQuadFiniteElement.
A 1D linear element with nodes at 1/3 and 2/3 (trace of RT1)
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
P1SegmentFiniteElement()
Construct the P1SegmentFiniteElement.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
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.
P1TetNonConfFiniteElement()
Construct the P1TetNonConfFiniteElement.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
A 1D quadratic element with nodes at the Gaussian points (trace of RT2)
P2SegmentFiniteElement()
Construct the P2SegmentFiniteElement.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
A 0D point finite element.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
PointFiniteElement()
Construct the PointFiniteElement.
A 1D quadratic finite element with uniformly spaced nodes.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Quad1DFiniteElement()
Construct the Quad1DFiniteElement.
A 2D quadratic element on triangle with nodes at the vertices and midpoints of the triangle.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Quad2DFiniteElement()
Construct the Quad2DFiniteElement.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcHessian(const IntegrationPoint &ip, DenseMatrix &h) const override
Evaluate the Hessians of all shape functions of a scalar finite element in reference space at the giv...
A 3D quadratic element on a tetrahedron with uniformly spaced nodes.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
Quadratic3DFiniteElement()
Construct the Quadratic3DFiniteElement.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 3D 0th order Raviert-Thomas element on a cube.
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
RT0HexFiniteElement()
Construct the RT0HexFiniteElement.
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
A 3D 0th order Raviert-Thomas element on a pyramid.
void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const override
Compute the discrete curl matrix from the given FiniteElement onto 'this' FiniteElement....
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
RT0PyrFiniteElement(bool rt0tets=true)
Construct the RT0PyrFiniteElement.
A 2D 1st order Raviart-Thomas vector element on a square.
RT0QuadFiniteElement()
Construct the RT0QuadFiniteElement.
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
A 3D 0th order Raviert-Thomas element on a tetrahedron.
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
RT0TetFiniteElement()
Construct the RT0TetFiniteElement.
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
A 2D 1st order Raviart-Thomas vector element on a triangle.
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
RT0TriangleFiniteElement()
Construct the RT0TriangleFiniteElement.
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
A 3D 0th order Raviert-Thomas element on a wedge.
void ProjectCurl(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &curl) const override
Compute the discrete curl matrix from the given FiniteElement onto 'this' FiniteElement....
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
RT0WdgFiniteElement()
Construct the RT0WdgFiniteElement.
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
A 3D 1st order Raviert-Thomas element on a cube.
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
RT1HexFiniteElement()
Construct the RT1HexFiniteElement.
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
A 2D 2nd order Raviart-Thomas vector element on a square.
RT1QuadFiniteElement()
Construct the RT1QuadFiniteElement.
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
A 2D 2nd order Raviart-Thomas vector element on a triangle.
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
RT1TriangleFiniteElement()
Construct the RT1TriangleFiniteElement.
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
A 2D 3rd order Raviart-Thomas vector element on a square.
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
RT2QuadFiniteElement()
Construct the RT2QuadFiniteElement.
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
void GetLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I) const override
Return the local interpolation matrix I (Dof x Dof) where the fine element is the image of the base g...
void Project(VectorCoefficient &vc, ElementTransformation &Trans, Vector &dofs) const override
Given a vector coefficient and a transformation, compute its projection (approximation) in the local ...
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
A 2D 3rd order Raviart-Thomas vector element on a triangle.
void CalcDivShape(const IntegrationPoint &ip, Vector &divshape) const override
Evaluate the divergence of all shape functions of a vector finite element in reference space at the g...
void CalcVShape(const IntegrationPoint &ip, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in reference space at the given...
void CalcVShape(ElementTransformation &Trans, DenseMatrix &shape) const override
Evaluate the values of all shape functions of a vector finite element in physical space at the point ...
RT2TriangleFiniteElement()
Construct the RT2TriangleFiniteElement.
A 2D refined bi-linear FE on a square.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
RefinedBiLinear2DFiniteElement()
Construct the RefinedBiLinear2DFiniteElement.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
A 1D refined linear element.
RefinedLinear1DFiniteElement()
Construct the RefinedLinear1DFiniteElement.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
A 2D refined linear element on a triangle.
RefinedLinear2DFiniteElement()
Construct the RefinedLinear2DFiniteElement.
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
A 2D refined linear element on a tetrahedron.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
RefinedLinear3DFiniteElement()
Construct the RefinedLinear3DFiniteElement.
A 3D refined tri-linear element on a cube.
RefinedTriLinear3DFiniteElement()
Construct the RefinedTriLinear3DFiniteElement.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
RotTriLinearHexFiniteElement()
Construct the RotTriLinearHexFiniteElement.
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Evaluate the values of all shape functions of a scalar finite element in reference space at the given...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
Evaluate the gradients of all shape functions of a scalar finite element in reference space at the gi...
A 3D tri-linear element on a cube with nodes at the vertices of the cube.
TriLinear3DFiniteElement()
Construct the TriLinear3DFiniteElement.
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
void CalcDShape(const IntegrationPoint &ip, DenseMatrix &dshape) const override
void CalcShape(const IntegrationPoint &ip, Vector &shape) const override
Base class for vector Coefficients that optionally depend on time and space.
Intermediate class for finite elements whose basis functions return vector values.
Definition fe_base.hpp:813
void CalcVShape_ND(ElementTransformation &Trans, DenseMatrix &shape) const
Definition fe_base.cpp:1071
void CalcVShape_RT(ElementTransformation &Trans, DenseMatrix &shape) const
Definition fe_base.cpp:1059
Vector data type.
Definition vector.hpp:82
float real_t
Definition config.hpp:43