![]() |
MFEM v4.8.0
Finite element discretization library
|
#include <fe_pyramid.hpp>
Public Member Functions | |
FuentesPyramid ()=default | |
void | phi_E (int p, Vector s, const DenseMatrix &grad_s, Vector &u, DenseMatrix &grad_u) const |
void | phi_Q (int p, Vector s, Vector t, DenseMatrix &u) const |
void | phi_Q (int p, Vector s, const DenseMatrix &grad_s, Vector t, const DenseMatrix &grad_t, DenseMatrix &u, DenseTensor &grad_u) const |
void | phi_T (int p, Vector nu, DenseMatrix &u) const |
void | phi_T (int p, Vector nu, const DenseMatrix &grad_nu, DenseMatrix &u, DenseTensor &grad_u) const |
void | E_E (int p, Vector s, Vector sds, DenseMatrix &u) const |
void | E_E (int p, Vector s, const DenseMatrix &grad_s, DenseMatrix &u, DenseMatrix &curl_u) const |
void | E_Q (int p, Vector s, Vector ds, Vector t, DenseTensor &u) const |
void | E_Q (int p, Vector s, const DenseMatrix &grad_s, Vector t, const DenseMatrix &grad_t, DenseTensor &u, DenseTensor &curl_u) const |
void | E_T (int p, Vector s, Vector sds, DenseTensor &u) const |
void | E_T (int p, Vector s, const DenseMatrix &grad_s, DenseTensor &u, DenseTensor &curl_u) const |
void | V_Q (int p, Vector s, Vector ds, Vector t, Vector dt, DenseTensor &u) const |
void | V_T (int p, Vector s, Vector sdsxds, DenseTensor &u) const |
void | V_T (int p, Vector s, Vector sdsxds, real_t dsdsxds, DenseTensor &u, DenseMatrix &du) const |
void | VT_T (int p, Vector s, Vector sds, Vector sdsxds, real_t mu, Vector grad_mu, DenseTensor &u) const |
void | VT_T (int p, Vector s, Vector sds, Vector sdsxds, Vector grad_s2, real_t mu, Vector grad_mu, DenseTensor &u, DenseMatrix &du) const |
void | V_L (int p, Vector sx, const DenseMatrix &grad_sx, Vector sy, const DenseMatrix &grad_sy, real_t t, Vector grad_t, DenseTensor &u) const |
void | V_R (int p, Vector s, const DenseMatrix &grad_s, real_t mu, Vector grad_mu, real_t t, Vector grad_t, DenseMatrix &u) const |
Static Protected Attributes | |
static constexpr real_t | one = 1.0 |
static constexpr real_t | zero = 0.0 |
static constexpr real_t | apex_tol = 1e-8 |
Base class for arbitrary order basis functions on pyramid-shaped elements
This base class provides a common class to store temporary vectors, matrices, and tensors computed by various functions defined on pyramid-shaped elements.
The function names defined here are chosen to reflect, as closely as possible, those used in the paper "Orientation embedded high order shape functions for the exact sequence elements of all shapes" by Federico Fuentes, Brendan Keith, Leszek Demkowicz, and Sriram Nagaraj, see https://doi.org/10.1016/j.camwa.2015.04.027.
Definition at line 38 of file fe_pyramid.hpp.
|
default |
|
inlinestatic |
u must be at least p + 1 in length
Definition at line 390 of file fe_pyramid.hpp.
|
static |
Definition at line 698 of file fe_pyramid.cpp.
|
static |
Definition at line 708 of file fe_pyramid.cpp.
|
static |
Definition at line 717 of file fe_pyramid.cpp.
|
static |
u must be at least p + 1 in length
Definition at line 485 of file fe_pyramid.cpp.
|
static |
Definition at line 493 of file fe_pyramid.cpp.
|
static |
Definition at line 503 of file fe_pyramid.cpp.
|
static |
Definition at line 512 of file fe_pyramid.cpp.
|
inlinestatic |
u must be at least p + 1 in length
Definition at line 373 of file fe_pyramid.hpp.
|
static |
Definition at line 666 of file fe_pyramid.cpp.
|
static |
Definition at line 676 of file fe_pyramid.cpp.
|
static |
Definition at line 685 of file fe_pyramid.cpp.
|
static |
u must be at least p + 1 in length
Definition at line 447 of file fe_pyramid.cpp.
|
static |
Definition at line 454 of file fe_pyramid.cpp.
|
static |
Definition at line 464 of file fe_pyramid.cpp.
|
static |
Definition at line 472 of file fe_pyramid.cpp.
|
static |
Integrated Jacobi Polynomials.
These are the integrals of the shifted and scaled Jacobi polynomials provided above and defined as:
These polynomials are computed as:
With
u must be at least p + 1 in length
Definition at line 525 of file fe_pyramid.cpp.
|
static |
Definition at line 559 of file fe_pyramid.cpp.
|
inlinestatic |
Definition at line 360 of file fe_pyramid.hpp.
|
inlinestatic |
Definition at line 364 of file fe_pyramid.hpp.
Integrated Legendre Polynomials.
These are the integrals of the shifted and scaled Legendre polynomials provided above and defined as:
These polynomials are computed as:
u must be at least p + 1 in length
Definition at line 369 of file fe_pyramid.cpp.
|
static |
Definition at line 395 of file fe_pyramid.cpp.
Definition at line 427 of file fe_pyramid.cpp.
|
static |
Definition at line 435 of file fe_pyramid.cpp.
|
static |
Shifted and Scaled Jacobi Polynomials.
Implements a scaled and shifted set of Jacobi polynomials
where alpha
u must be at least p + 1 in length
Definition at line 591 of file fe_pyramid.cpp.
|
static |
Definition at line 613 of file fe_pyramid.cpp.
|
static |
Definition at line 645 of file fe_pyramid.cpp.
|
static |
Definition at line 654 of file fe_pyramid.cpp.
Shifted and Scaled Legendre Polynomials.
Implements a scaled and shifted set of Legendre polynomials
where t >= 0.0, x
u must be at least p + 1 in length
Definition at line 293 of file fe_pyramid.cpp.
|
static |
Definition at line 313 of file fe_pyramid.cpp.
Definition at line 351 of file fe_pyramid.cpp.
|
static |
Definition at line 359 of file fe_pyramid.cpp.
|
inlinestatic |
Definition at line 87 of file fe_pyramid.hpp.
Divergences of the above "normal" vector functions divided by 3.
Definition at line 160 of file fe_pyramid.cpp.
Definition at line 175 of file fe_pyramid.cpp.
Definition at line 163 of file fe_pyramid.cpp.
Definition at line 166 of file fe_pyramid.cpp.
Definition at line 172 of file fe_pyramid.cpp.
Definition at line 169 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::E_E | ( | int | p, |
Vector | s, | ||
const DenseMatrix & | grad_s, | ||
DenseMatrix & | u, | ||
DenseMatrix & | curl_u ) const |
Definition at line 983 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::E_E | ( | int | p, |
Vector | s, | ||
Vector | sds, | ||
DenseMatrix & | u ) const |
This is a vector-valued function associated with an edge of a pyramid
The vector s contains two coordinate values and ds is related to the gradient of these coordinates with respect to the reference coordinates i.e. sds = s0 grad s1 - s1 grad s0
Definition at line 960 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::E_Q | ( | int | p, |
Vector | s, | ||
const DenseMatrix & | grad_s, | ||
Vector | t, | ||
const DenseMatrix & | grad_t, | ||
DenseTensor & | u, | ||
DenseTensor & | curl_u ) const |
Definition at line 1066 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::E_Q | ( | int | p, |
Vector | s, | ||
Vector | ds, | ||
Vector | t, | ||
DenseTensor & | u ) const |
Definition at line 1027 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::E_T | ( | int | p, |
Vector | s, | ||
const DenseMatrix & | grad_s, | ||
DenseTensor & | u, | ||
DenseTensor & | curl_u ) const |
Definition at line 1171 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::E_T | ( | int | p, |
Vector | s, | ||
Vector | sds, | ||
DenseTensor & | u ) const |
Definition at line 1137 of file fe_pyramid.cpp.
Gradients of the "Affine" Coordinates.
Definition at line 21 of file fe_pyramid.cpp.
|
static |
Gradients of the above two component vectors.
Definition at line 54 of file fe_pyramid.cpp.
Definition at line 28 of file fe_pyramid.cpp.
|
static |
Definition at line 62 of file fe_pyramid.cpp.
Definition at line 35 of file fe_pyramid.cpp.
|
static |
Definition at line 70 of file fe_pyramid.cpp.
Definition at line 42 of file fe_pyramid.cpp.
|
static |
Definition at line 78 of file fe_pyramid.cpp.
Definition at line 49 of file fe_pyramid.cpp.
Definition at line 172 of file fe_pyramid.hpp.
Definition at line 186 of file fe_pyramid.cpp.
|
static |
Definition at line 178 of file fe_pyramid.cpp.
|
static |
Definition at line 200 of file fe_pyramid.cpp.
Definition at line 174 of file fe_pyramid.hpp.
Definition at line 193 of file fe_pyramid.cpp.
Definition at line 216 of file fe_pyramid.cpp.
|
static |
Definition at line 233 of file fe_pyramid.cpp.
|
static |
Definition at line 241 of file fe_pyramid.cpp.
Definition at line 222 of file fe_pyramid.cpp.
|
static |
Definition at line 250 of file fe_pyramid.cpp.
Definition at line 228 of file fe_pyramid.cpp.
Pyramid "Affine" Coordinates.
Definition at line 90 of file fe_pyramid.hpp.
Three component vectors associated with triangular faces.
Definition at line 131 of file fe_pyramid.hpp.
Vector functions related to the normals to the triangular faces
Computes
Definition at line 118 of file fe_pyramid.cpp.
Definition at line 141 of file fe_pyramid.hpp.
Definition at line 153 of file fe_pyramid.cpp.
Two component vectors associated with edges touching the apex.
Definition at line 109 of file fe_pyramid.hpp.
Computes
Definition at line 86 of file fe_pyramid.cpp.
Definition at line 92 of file fe_pyramid.hpp.
Definition at line 133 of file fe_pyramid.hpp.
Definition at line 125 of file fe_pyramid.cpp.
Definition at line 111 of file fe_pyramid.hpp.
Definition at line 94 of file fe_pyramid.cpp.
Definition at line 94 of file fe_pyramid.hpp.
Definition at line 135 of file fe_pyramid.hpp.
Definition at line 132 of file fe_pyramid.cpp.
Definition at line 113 of file fe_pyramid.hpp.
Definition at line 102 of file fe_pyramid.cpp.
Definition at line 96 of file fe_pyramid.hpp.
Definition at line 139 of file fe_pyramid.hpp.
Definition at line 146 of file fe_pyramid.cpp.
Definition at line 137 of file fe_pyramid.hpp.
Definition at line 139 of file fe_pyramid.cpp.
Definition at line 115 of file fe_pyramid.hpp.
Definition at line 110 of file fe_pyramid.cpp.
Definition at line 98 of file fe_pyramid.hpp.
Definition at line 167 of file fe_pyramid.hpp.
|
inlinestatic |
Definition at line 182 of file fe_pyramid.hpp.
Definition at line 177 of file fe_pyramid.hpp.
Definition at line 190 of file fe_pyramid.hpp.
Definition at line 208 of file fe_pyramid.cpp.
Definition at line 169 of file fe_pyramid.hpp.
|
inlinestatic |
Definition at line 184 of file fe_pyramid.hpp.
Definition at line 196 of file fe_pyramid.hpp.
Definition at line 205 of file fe_pyramid.hpp.
Definition at line 209 of file fe_pyramid.hpp.
Definition at line 275 of file fe_pyramid.cpp.
Definition at line 259 of file fe_pyramid.cpp.
Definition at line 198 of file fe_pyramid.hpp.
Definition at line 207 of file fe_pyramid.hpp.
Definition at line 211 of file fe_pyramid.hpp.
Definition at line 267 of file fe_pyramid.cpp.
Definition at line 199 of file fe_pyramid.hpp.
u must be at least p + 1 in length
Definition at line 730 of file fe_pyramid.cpp.
|
static |
Definition at line 737 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::phi_E | ( | int | p, |
Vector | s, | ||
const DenseMatrix & | grad_s, | ||
Vector & | u, | ||
DenseMatrix & | grad_u ) const |
grad_s must be 2x3
Definition at line 762 of file fe_pyramid.cpp.
Definition at line 744 of file fe_pyramid.cpp.
|
static |
Definition at line 751 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::phi_Q | ( | int | p, |
Vector | s, | ||
const DenseMatrix & | grad_s, | ||
Vector | t, | ||
const DenseMatrix & | grad_t, | ||
DenseMatrix & | u, | ||
DenseTensor & | grad_u ) const |
Definition at line 814 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::phi_Q | ( | int | p, |
Vector | s, | ||
Vector | t, | ||
DenseMatrix & | u ) const |
u must be at least (p+1)x(p+1) in size
Definition at line 786 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::phi_T | ( | int | p, |
Vector | nu, | ||
const DenseMatrix & | grad_nu, | ||
DenseMatrix & | u, | ||
DenseTensor & | grad_u ) const |
Definition at line 900 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::phi_T | ( | int | p, |
Vector | nu, | ||
DenseMatrix & | u ) const |
Definition at line 868 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::V_L | ( | int | p, |
Vector | sx, | ||
const DenseMatrix & | grad_sx, | ||
Vector | sy, | ||
const DenseMatrix & | grad_sy, | ||
real_t | t, | ||
Vector | grad_t, | ||
DenseTensor & | u ) const |
This implements
u must be at least (p+1)x(p+1)x3
Definition at line 1469 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::V_Q | ( | int | p, |
Vector | s, | ||
Vector | ds, | ||
Vector | t, | ||
Vector | dt, | ||
DenseTensor & | u ) const |
This is a vector-valued function associated with the quadrilateral face of a pyramid
The vectors s and t contain pairs of coordinate values and ds and dt are related to derivatives of these coordinates:
ds = s0 grad s1 - s1 grad s0
dt = t0 grad t1 - t1 grad t0
Definition at line 1238 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::V_R | ( | int | p, |
Vector | s, | ||
const DenseMatrix & | grad_s, | ||
real_t | mu, | ||
Vector | grad_mu, | ||
real_t | t, | ||
Vector | grad_t, | ||
DenseMatrix & | u ) const |
This implements
u must be at least (p+1)x3
Definition at line 1538 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::V_T | ( | int | p, |
Vector | s, | ||
Vector | sdsxds, | ||
DenseTensor & | u ) const |
This is a vector-valued function associated with the triangular faces of a pyramid
The vector s contains three coordinate values and sdsxds is related to derivatives of these coordinates with respect to the reference coordinates:
sdsxds = s0 grad s1 x grad s2 + s1 grad s2 x grad s0 + s2 grad s0 x grad s1
Definition at line 1272 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::V_T | ( | int | p, |
Vector | s, | ||
Vector | sdsxds, | ||
real_t | dsdsxds, | ||
DenseTensor & | u, | ||
DenseMatrix & | du ) const |
This computes V_T as above and its divergence
The vector s contains three coordinate values and sdsxds is related to derivatives of these coordinates with respect to the reference coordinates:
sdsxds = s0 grad s1 x grad s2 + s1 grad s2 x grad s0 + s2 grad s0 x grad s1
The scalar dsdsxds is the divergence of sdsxds:
dsdsxds = grad s0 dot (grad s1 x grad s2)
Definition at line 1306 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::VT_T | ( | int | p, |
Vector | s, | ||
Vector | sds, | ||
Vector | sdsxds, | ||
real_t | mu, | ||
Vector | grad_mu, | ||
DenseTensor & | u ) const |
Definition at line 1345 of file fe_pyramid.cpp.
void mfem::FuentesPyramid::VT_T | ( | int | p, |
Vector | s, | ||
Vector | sds, | ||
Vector | sdsxds, | ||
Vector | grad_s2, | ||
real_t | mu, | ||
Vector | grad_mu, | ||
DenseTensor & | u, | ||
DenseMatrix & | du ) const |
Definition at line 1398 of file fe_pyramid.cpp.
|
staticconstexprprotected |
Definition at line 82 of file fe_pyramid.hpp.
|
staticconstexprprotected |
Definition at line 80 of file fe_pyramid.hpp.
|
staticconstexprprotected |
Definition at line 81 of file fe_pyramid.hpp.