MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
fe_pos.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_POS
13#define MFEM_FE_POS
14
15#include "fe_base.hpp"
16#include "fe_pyramid.hpp"
17
18namespace mfem
19{
20
21/** @brief Class for finite elements utilizing the
22 always positive Bernstein basis. */
24{
25public:
26 /** @brief Construct PositiveFiniteElement with given
27 @param D Reference space dimension
28 @param G Geometry type (of type Geometry::Type)
29 @param Do Number of degrees of freedom in the FiniteElement
30 @param O Order/degree of the FiniteElement
31 @param F FunctionSpace type of the FiniteElement
32 */
33 PositiveFiniteElement(int D, Geometry::Type G, int Do, int O,
34 int F = FunctionSpace::Pk) :
35 ScalarFiniteElement(D, G, Do, O, F)
36 { }
37
39 DenseMatrix &I) const override
40 { ScalarLocalInterpolation(Trans, I, *this); }
41
43 DenseMatrix &R) const override
44 { ScalarLocalL2Restriction(Trans, R, *this); }
45
48 DenseMatrix &I) const override
49 { CheckScalarFE(fe).ScalarLocalInterpolation(Trans, I, *this); }
50
52
53 // Low-order monotone "projection" (actually it is not a projection): the
54 // dofs are set to be the Coefficient values at the nodes.
55 void Project(Coefficient &coeff,
56 ElementTransformation &Trans, Vector &dofs) const override;
57
58 void Project (VectorCoefficient &vc,
59 ElementTransformation &Trans, Vector &dofs) const override;
60
61 void Project(const FiniteElement &fe, ElementTransformation &Trans,
62 DenseMatrix &I) const override;
63};
64
65
68{
69public:
70 PositiveTensorFiniteElement(const int dims, const int p,
71 const DofMapType dmtype);
72
74 DofToQuad::Mode mode) const override
75 {
76 return (mode == DofToQuad::FULL) ?
78 GetTensorDofToQuad(*this, ir, mode, basis1d, true, dof2quad_array);
79 }
80
81 void GetFaceMap(const int face_id, Array<int> &face_map) const override;
82};
83
84
85/// A 2D positive bi-quadratic element on a square utilizing the 2nd order
86/// Bernstein basis
88{
89public:
90 /// Construct the BiQuadPos2DFiniteElement
92 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
93 void CalcDShape(const IntegrationPoint &ip,
94 DenseMatrix &dshape) const override;
96 DenseMatrix &I) const override;
98 void Project(Coefficient &coeff, ElementTransformation &Trans,
99 Vector &dofs) const override;
101 Vector &dofs) const override;
102 void ProjectDelta(int vertex, Vector &dofs) const override
103 { dofs = 0.; dofs(vertex) = 1.; }
104};
105
106
107/// A 1D quadratic positive element utilizing the 2nd order Bernstein basis
109{
110public:
111 /// Construct the QuadPos1DFiniteElement
113 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
114 void CalcDShape(const IntegrationPoint &ip,
115 DenseMatrix &dshape) const override;
116};
117
118
119/// Arbitrary order H1 elements in 1D utilizing the Bernstein basis
121{
122private:
123#ifndef MFEM_THREAD_SAFE
124 // This is to share scratch space between invocations, which helps speed
125 // things up, but with OpenMP, we need one copy per thread. Right now, we
126 // solve this by allocating this space within each function call every time
127 // we call it. Alternatively, we should do some sort thread private thing.
128 // Brunner, Jan 2014
129 mutable Vector shape_x, dshape_x;
130#endif
131
132public:
133 /// Construct the H1Pos_SegmentElement of order @a p
134 H1Pos_SegmentElement(const int p);
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};
140
141
142/// Arbitrary order H1 elements in 2D utilizing the Bernstein basis on a square
144{
145private:
146#ifndef MFEM_THREAD_SAFE
147 // See comment in H1Pos_SegmentElement
148 mutable Vector shape_x, shape_y, dshape_x, dshape_y;
149#endif
150
151public:
152 /// Construct the H1Pos_QuadrilateralElement of order @a p
153 H1Pos_QuadrilateralElement(const int p);
154 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
155 void CalcDShape(const IntegrationPoint &ip,
156 DenseMatrix &dshape) const override;
157 void ProjectDelta(int vertex, Vector &dofs) const override;
158};
159
160
161/// Arbitrary order H1 elements in 3D utilizing the Bernstein basis on a cube
163{
164private:
165#ifndef MFEM_THREAD_SAFE
166 // See comment in H1Pos_SegmentElement.
167 mutable Vector shape_x, shape_y, shape_z, dshape_x, dshape_y, dshape_z;
168#endif
169
170public:
171 /// Construct the H1Pos_HexahedronElement of order @a p
172 H1Pos_HexahedronElement(const int p);
173 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
174 void CalcDShape(const IntegrationPoint &ip,
175 DenseMatrix &dshape) const override;
176 void ProjectDelta(int vertex, Vector &dofs) const override;
177};
178
179
180/// Arbitrary order H1 elements in 2D utilizing the Bernstein basis on a triangle
182{
183protected:
184#ifndef MFEM_THREAD_SAFE
187#endif
189
190public:
191 /// Construct the H1Pos_TriangleElement of order @a p
192 H1Pos_TriangleElement(const int p);
193
194 // The size of shape is (p+1)(p+2)/2 (dof).
195 static void CalcShape(const int p, const real_t x, const real_t y,
196 real_t *shape);
197
198 // The size of dshape_1d is p+1; the size of dshape is (dof x dim).
199 static void CalcDShape(const int p, const real_t x, const real_t y,
200 real_t *dshape_1d, real_t *dshape);
201
202 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
203 void CalcDShape(const IntegrationPoint &ip,
204 DenseMatrix &dshape) const override;
205};
206
207
208/// Arbitrary order H1 elements in 3D utilizing the Bernstein basis on a
209/// tetrahedron
211{
212protected:
213#ifndef MFEM_THREAD_SAFE
216#endif
218
219public:
220 /// Construct the H1Pos_TetrahedronElement of order @a p
221 H1Pos_TetrahedronElement(const int p);
222
223 // The size of shape is (p+1)(p+2)(p+3)/6 (dof).
224 static void CalcShape(const int p, const real_t x, const real_t y,
225 const real_t z, real_t *shape);
226
227 // The size of dshape_1d is p+1; the size of dshape is (dof x dim).
228 static void CalcDShape(const int p, const real_t x, const real_t y,
229 const real_t z, real_t *dshape_1d, real_t *dshape);
230
231 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
232 void CalcDShape(const IntegrationPoint &ip,
233 DenseMatrix &dshape) const override;
234};
235
236
237/// Arbitrary order H1 elements in 3D utilizing the Bernstein basis on a wedge
239{
240protected:
241#ifndef MFEM_THREAD_SAFE
244#endif
246
249
250public:
251 /// Construct the H1Pos_WedgeElement of order @a p
252 H1Pos_WedgeElement(const int p);
253
254 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
255 void CalcDShape(const IntegrationPoint &ip,
256 DenseMatrix &dshape) const override;
257};
258
259
260/// Arbitrary order H1 elements in 3D utilizing the Bernstein basis on a pyramid
261///
262/// The pyramid affine-related coordinates $\lambda_i$ for $i=1,\ldots,5$ can
263/// be used to define a positive H1 basis by noting that $\lambda_i \ge 0$
264/// inside the pyramid for all $i$ and that $\sum_{i=1}^5\lambda_i=1$. This
265/// leads to $1 = (\sum_{i=1}^5\lambda_i)^p$. The terms of this product,
266/// expanded as a polynomial in the $\lambda_i$, can be used as a Bernstein
267/// basis of order $p$ on a pyramid.
269{
270protected:
271 const int nterms;
272#ifndef MFEM_THREAD_SAFE
276#endif
277 std::map<int,int> dof_map;
278
279 struct Index
280 {
281 Index() = default;
282 int operator()(int i1, int i2, int i3, int i4, int i5)
283 {
284 const int p = i1 + i2 + i3 + i4 + i5;
285 const int min24 = std::min(i2,i4);
286 i1 += min24;
287 i2 -= min24;
288 i3 += min24;
289 i4 -= min24;
290 return i2 + i3 * (p - i4 - i5) - i4 * i5 * (p + 2)
291 - ((i3 - 3) * i3) / 2 + i4 * ((p + 1) * (p + 2)) / 2
292 + (i4 * i5 * (i4 + i5)) / 2 + ((i4 - 1) * i4 * (i4 +1)) / 6
293 - (p + 2) * ((i4 - 1) * i4) / 2
294 + (i5 * (5 + 2 * p - i5) * (i5 * i5 - i5 * (5 + 2 * p)
295 + 2 * (5 + 5 * p + p * p))) / 24;
296 }
297 };
298
299public:
300 /// Construct the H1Pos_PyramidElement of order @a p
301 H1Pos_PyramidElement(const int p);
302
303 // The size of shape is (p+1)(p+2)(p+3)(p+4)/24.
304 // The size of shape_1d should be at least p+1.
305 static void CalcShape(const int p, const real_t x, const real_t y,
306 const real_t z, real_t *shape_1d, real_t *shape);
307
308 // The size of dshape is (p+1)(p+2)(p+3)(p+4)/24 by 3.
309 // The size of dshape_1d should be at least p+1.
310 static void CalcDShape(const int p, const real_t x, const real_t y,
311 const real_t z, real_t *dshape_1d, real_t *dshape);
312
313 virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
314 virtual void CalcDShape(const IntegrationPoint &ip,
315 DenseMatrix &dshape) const;
316
317 // Returns (p+1)(p+2)(p+3)(p+4)/24 which is the size of the temporary arrays
318 // needed above
319 int GetNumTerms() const { return nterms; }
320};
321
322
323/// Arbitrary order L2 elements in 1D utilizing the Bernstein basis on a segment
325{
326private:
327#ifndef MFEM_THREAD_SAFE
328 mutable Vector shape_x, dshape_x;
329#endif
330
331public:
332 /// Construct the L2Pos_SegmentElement of order @a p
333 L2Pos_SegmentElement(const int p);
334 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
335 void CalcDShape(const IntegrationPoint &ip,
336 DenseMatrix &dshape) const override;
337 void ProjectDelta(int vertex, Vector &dofs) const override;
338};
339
340
341/// Arbitrary order L2 elements in 2D utilizing the Bernstein basis on a square
343{
344private:
345#ifndef MFEM_THREAD_SAFE
346 mutable Vector shape_x, shape_y, dshape_x, dshape_y;
347#endif
348
349public:
350 /// Construct the L2Pos_QuadrilateralElement of order @a p
351 L2Pos_QuadrilateralElement(const int p);
352 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
353 void CalcDShape(const IntegrationPoint &ip,
354 DenseMatrix &dshape) const override;
355 void ProjectDelta(int vertex, Vector &dofs) const override;
356};
357
358
359/// Arbitrary order L2 elements in 3D utilizing the Bernstein basis on a cube
361{
362private:
363#ifndef MFEM_THREAD_SAFE
364 mutable Vector shape_x, shape_y, shape_z, dshape_x, dshape_y, dshape_z;
365#endif
366
367public:
368 /// Construct the L2Pos_HexahedronElement of order @a p
369 L2Pos_HexahedronElement(const int p);
370 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
371 void CalcDShape(const IntegrationPoint &ip,
372 DenseMatrix &dshape) const override;
373 void ProjectDelta(int vertex, Vector &dofs) const override;
374};
375
376
377/// Arbitrary order L2 elements in 2D utilizing the Bernstein basis on a triangle
379{
380private:
381#ifndef MFEM_THREAD_SAFE
382 mutable Vector dshape_1d;
383#endif
384
385public:
386 /// Construct the L2Pos_TriangleElement of order @a p
387 L2Pos_TriangleElement(const int p);
388 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
389 void CalcDShape(const IntegrationPoint &ip,
390 DenseMatrix &dshape) const override;
391 void ProjectDelta(int vertex, Vector &dofs) const override;
392};
393
394
395/// Arbitrary order L2 elements in 3D utilizing the Bernstein basis on a
396/// tetrahedron
398{
399private:
400#ifndef MFEM_THREAD_SAFE
401 mutable Vector dshape_1d;
402#endif
403
404public:
405 /// Construct the L2Pos_TetrahedronElement of order @a p
406 L2Pos_TetrahedronElement(const int p);
407 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
408 void CalcDShape(const IntegrationPoint &ip,
409 DenseMatrix &dshape) const override;
410 void ProjectDelta(int vertex, Vector &dofs) const override;
411};
412
413
414/// Arbitrary order L2 elements in 3D utilizing the Bernstein basis on a wedge
416{
417protected:
418#ifndef MFEM_THREAD_SAFE
421#endif
423
426
427public:
428 /// Construct the L2Pos_WedgeElement of order @a p
429 L2Pos_WedgeElement(const int p);
430
431 void CalcShape(const IntegrationPoint &ip, Vector &shape) const override;
432 void CalcDShape(const IntegrationPoint &ip,
433 DenseMatrix &dshape) const override;
434};
435
436/// Arbitrary order L2 elements in 3D utilizing the Bernstein basis on a pyramid
438{
439protected:
440 const int nterms;
441#ifndef MFEM_THREAD_SAFE
445#endif
446 std::map<int,int> dof_map;
447
448 struct Index
449 {
450 Index() = default;
451 int operator()(int i1, int i2, int i3, int i4, int i5)
452 {
453 const int p = i1 + i2 + i3 + i4 + i5;
454 const int min24 = std::min(i2,i4);
455 i1 += min24;
456 i2 -= min24;
457 i3 += min24;
458 i4 -= min24;
459 return i2 + i3 * (p - i4 - i5) - i4 * i5 * (p + 2)
460 - ((i3 - 3) * i3) / 2 + i4 * ((p + 1) * (p + 2)) / 2
461 + (i4 * i5 * (i4 + i5)) / 2 + ((i4 - 1) * i4 * (i4 +1)) / 6
462 - (p + 2) * ((i4 - 1) * i4) / 2
463 + (i5 * (5 + 2 * p - i5) * (i5 * i5 - i5 * (5 + 2 * p)
464 + 2 * (5 + 5 * p + p * p))) / 24;
465 }
466 };
467
468 // Returns (p+1)(p+2)(p+3)(p+4)/24 which is the size of the temporary arrays
469 // needed below
470 int GetNumTerms() const { return nterms; }
471
472 // The size of shape is (p+1)(p+2)(p+3)(p+4)/24.
473 // The size of shape_1d should be at least p+1.
474 static void CalcShape(const int p, const real_t x, const real_t y,
475 const real_t z, real_t *shape_1d, real_t *shape);
476
477 // The size of dshape is (p+1)(p+2)(p+3)(p+4)/24 by 3.
478 // The size of dshape_1d should be at least p+1.
479 static void CalcDShape(const int p, const real_t x, const real_t y,
480 const real_t z, real_t *dshape_1d, real_t *dshape);
481
482public:
483 /// Construct the L2Pos_PyramidElement of order @a p
484 L2Pos_PyramidElement(const int p);
485
486 virtual void CalcShape(const IntegrationPoint &ip, Vector &shape) const;
487 virtual void CalcDShape(const IntegrationPoint &ip,
488 DenseMatrix &dshape) const;
489};
490
491} // namespace mfem
492
493#endif
BiQuadPos2DFiniteElement()
Construct the BiQuadPos2DFiniteElement.
Definition fe_pos.cpp:95
void Project(Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const override
Given a coefficient and a transformation, compute its projection (approximation) in the local finite ...
Definition fe_pos.cpp:220
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...
Definition fe_pos.cpp:118
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...
Definition fe_pos.cpp:191
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Definition fe_pos.hpp:102
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...
Definition fe_pos.cpp:142
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Data type dense matrix using column-major storage.
Definition densemat.hpp:24
Structure representing the matrices/tensors needed to evaluate (in reference space) the values,...
Definition fe_base.hpp:141
Mode
Type of data stored in the arrays B, Bt, G, and Gt.
Definition fe_base.hpp:154
@ FULL
Full multidimensional representation which does not use tensor product structure. The ordering of the...
Definition fe_base.hpp:158
Abstract class for all finite elements.
Definition fe_base.hpp:244
virtual const DofToQuad & GetDofToQuad(const IntegrationRule &ir, DofToQuad::Mode mode) const
Return a DofToQuad structure corresponding to the given IntegrationRule using the given DofToQuad::Mo...
Definition fe_base.cpp:365
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
Array< DofToQuad * > dof2quad_array
Container for all DofToQuad objects created by the FiniteElement.
Definition fe_base.hpp:263
@ Pk
Polynomials of order k.
Definition fe_base.hpp:230
Arbitrary order H1 elements in 3D utilizing the Bernstein basis on a cube.
Definition fe_pos.hpp:163
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...
Definition fe_pos.cpp:475
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...
Definition fe_pos.cpp:454
H1Pos_HexahedronElement(const int p)
Construct the H1Pos_HexahedronElement of order p.
Definition fe_pos.cpp:432
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Definition fe_pos.cpp:499
H1Pos_PyramidElement(const int p)
Construct the H1Pos_PyramidElement of order p.
Definition fe_pos.cpp:1045
static void CalcShape(const int p, const real_t x, const real_t y, const real_t z, real_t *shape_1d, real_t *shape)
Definition fe_pos.cpp:1193
static void CalcDShape(const int p, const real_t x, const real_t y, const real_t z, real_t *dshape_1d, real_t *dshape)
Definition fe_pos.cpp:1246
std::map< int, int > dof_map
Definition fe_pos.hpp:277
Arbitrary order H1 elements in 2D utilizing the Bernstein basis on a square.
Definition fe_pos.hpp:144
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...
Definition fe_pos.cpp:384
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Definition fe_pos.cpp:425
H1Pos_QuadrilateralElement(const int p)
Construct the H1Pos_QuadrilateralElement of order p.
Definition fe_pos.cpp:364
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...
Definition fe_pos.cpp:404
Arbitrary order H1 elements in 1D utilizing the Bernstein basis.
Definition fe_pos.hpp:121
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...
Definition fe_pos.cpp:337
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...
Definition fe_pos.cpp:317
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Definition fe_pos.cpp:357
H1Pos_SegmentElement(const int p)
Construct the H1Pos_SegmentElement of order p.
Definition fe_pos.cpp:299
H1Pos_TetrahedronElement(const int p)
Construct the H1Pos_TetrahedronElement of order p.
Definition fe_pos.cpp:651
static void CalcShape(const int p, const real_t x, const real_t y, const real_t z, real_t *shape)
Definition fe_pos.cpp:753
static void CalcDShape(const int p, const real_t x, const real_t y, const real_t z, real_t *dshape_1d, real_t *dshape)
Definition fe_pos.cpp:786
Arbitrary order H1 elements in 2D utilizing the Bernstein basis on a triangle.
Definition fe_pos.hpp:182
H1Pos_TriangleElement(const int p)
Construct the H1Pos_TriangleElement of order p.
Definition fe_pos.cpp:506
static void CalcDShape(const int p, const real_t x, const real_t y, real_t *dshape_1d, real_t *dshape)
Definition fe_pos.cpp:587
static void CalcShape(const int p, const real_t x, const real_t y, real_t *shape)
Definition fe_pos.cpp:561
Arbitrary order H1 elements in 3D utilizing the Bernstein basis on a wedge.
Definition fe_pos.hpp:239
H1Pos_TriangleElement TriangleFE
Definition fe_pos.hpp:247
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...
Definition fe_pos.cpp:1020
H1Pos_SegmentElement SegmentFE
Definition fe_pos.hpp:248
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...
Definition fe_pos.cpp:1001
H1Pos_WedgeElement(const int p)
Construct the H1Pos_WedgeElement of order p.
Definition fe_pos.cpp:903
Class for integration point with weight.
Definition intrules.hpp:35
Class for an integration rule - an Array of IntegrationPoint.
Definition intrules.hpp:100
Arbitrary order L2 elements in 3D utilizing the Bernstein basis on a cube.
Definition fe_pos.hpp:361
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...
Definition fe_pos.cpp:1624
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Definition fe_pos.cpp:1669
L2Pos_HexahedronElement(const int p)
Construct the L2Pos_HexahedronElement of order p.
Definition fe_pos.cpp:1597
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...
Definition fe_pos.cpp:1645
Arbitrary order L2 elements in 3D utilizing the Bernstein basis on a pyramid.
Definition fe_pos.hpp:438
static void CalcDShape(const int p, const real_t x, const real_t y, const real_t z, real_t *dshape_1d, real_t *dshape)
Definition fe_pos.cpp:1975
std::map< int, int > dof_map
Definition fe_pos.hpp:446
static void CalcShape(const int p, const real_t x, const real_t y, const real_t z, real_t *shape_1d, real_t *shape)
Definition fe_pos.cpp:1922
L2Pos_PyramidElement(const int p)
Construct the L2Pos_PyramidElement of order p.
Definition fe_pos.cpp:1880
Arbitrary order L2 elements in 2D utilizing the Bernstein basis on a square.
Definition fe_pos.hpp:343
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...
Definition fe_pos.cpp:1543
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Definition fe_pos.cpp:1582
L2Pos_QuadrilateralElement(const int p)
Construct the L2Pos_QuadrilateralElement of order p.
Definition fe_pos.cpp:1519
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...
Definition fe_pos.cpp:1562
Arbitrary order L2 elements in 1D utilizing the Bernstein basis on a segment.
Definition fe_pos.hpp:325
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...
Definition fe_pos.cpp:1495
L2Pos_SegmentElement(const int p)
Construct the L2Pos_SegmentElement of order p.
Definition fe_pos.cpp:1474
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...
Definition fe_pos.cpp:1501
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Definition fe_pos.cpp:1512
L2Pos_TetrahedronElement(const int p)
Construct the L2Pos_TetrahedronElement of order p.
Definition fe_pos.cpp:1739
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Definition fe_pos.cpp:1780
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...
Definition fe_pos.cpp:1762
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...
Definition fe_pos.cpp:1769
Arbitrary order L2 elements in 2D utilizing the Bernstein basis on a triangle.
Definition fe_pos.hpp:379
L2Pos_TriangleElement(const int p)
Construct the L2Pos_TriangleElement of order p.
Definition fe_pos.cpp:1688
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...
Definition fe_pos.cpp:1710
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...
Definition fe_pos.cpp:1716
void ProjectDelta(int vertex, Vector &dofs) const override
Project a delta function centered on the given vertex in the local finite dimensional space represent...
Definition fe_pos.cpp:1727
Arbitrary order L2 elements in 3D utilizing the Bernstein basis on a wedge.
Definition fe_pos.hpp:416
L2Pos_TriangleElement TriangleFE
Definition fe_pos.hpp:424
L2Pos_WedgeElement(const int p)
Construct the L2Pos_WedgeElement of order p.
Definition fe_pos.cpp:1793
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...
Definition fe_pos.cpp:1855
L2Pos_SegmentElement SegmentFE
Definition fe_pos.hpp:425
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...
Definition fe_pos.cpp:1836
Class for finite elements utilizing the always positive Bernstein basis.
Definition fe_pos.hpp:24
void GetLocalRestriction(ElementTransformation &Trans, DenseMatrix &R) const override
Return a local restriction matrix R (Dof x Dof) mapping fine dofs to coarse dofs.
Definition fe_pos.hpp:42
PositiveFiniteElement(int D, Geometry::Type G, int Do, int O, int F=FunctionSpace::Pk)
Construct PositiveFiniteElement with given.
Definition fe_pos.hpp:33
void GetTransferMatrix(const FiniteElement &fe, ElementTransformation &Trans, DenseMatrix &I) const override
Return interpolation matrix, I, which maps dofs from a coarse element, fe, to the fine dofs on this f...
Definition fe_pos.hpp:46
void Project(Coefficient &coeff, ElementTransformation &Trans, Vector &dofs) const override
Given a coefficient and a transformation, compute its projection (approximation) in the local finite ...
Definition fe_pos.cpp:25
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...
Definition fe_pos.hpp:38
PositiveTensorFiniteElement(const int dims, const int p, const DofMapType dmtype)
Definition fe_pos.cpp:81
void GetFaceMap(const int face_id, Array< int > &face_map) const override
Return the mapping from lexicographic face DOFs to lexicographic element DOFs for the given local fac...
Definition fe_pos.cpp:88
const DofToQuad & GetDofToQuad(const IntegrationRule &ir, DofToQuad::Mode mode) const override
Return a DofToQuad structure corresponding to the given IntegrationRule using the given DofToQuad::Mo...
Definition fe_pos.hpp:73
A 1D quadratic positive element utilizing the 2nd order Bernstein basis.
Definition fe_pos.hpp:109
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...
Definition fe_pos.cpp:278
QuadPos1DFiniteElement()
Construct the QuadPos1DFiniteElement.
Definition fe_pos.cpp:270
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...
Definition fe_pos.cpp:288
Class for finite elements with basis functions that return scalar values.
Definition fe_base.hpp:662
void ScalarLocalInterpolation(ElementTransformation &Trans, DenseMatrix &I, const ScalarFiniteElement &fine_fe) const
Get matrix I "Interpolation" defined through local L2-projection in the space defined by the fine_fe.
Definition fe_base.cpp:561
void ScalarLocalL2Restriction(ElementTransformation &Trans, DenseMatrix &R, const ScalarFiniteElement &coarse_fe) const
Get restriction matrix R defined through local L2-projection in the space defined by the coarse_fe.
Definition fe_base.cpp:602
static const ScalarFiniteElement & CheckScalarFE(const FiniteElement &fe)
Definition fe_base.hpp:664
Poly_1D::Basis & basis1d
Definition fe_base.hpp:1251
static const DofToQuad & GetTensorDofToQuad(const FiniteElement &fe, const IntegrationRule &ir, DofToQuad::Mode mode, const Poly_1D::Basis &basis, bool closed, Array< DofToQuad * > &dof2quad_array)
Definition fe_base.cpp:2599
Base class for vector Coefficients that optionally depend on time and space.
Vector data type.
Definition vector.hpp:82
float real_t
Definition config.hpp:43
real_t p(const Vector &x, real_t t)
int operator()(int i1, int i2, int i3, int i4, int i5)
Definition fe_pos.hpp:282
int operator()(int i1, int i2, int i3, int i4, int i5)
Definition fe_pos.hpp:451