|
| v0.14.0
|
Calculate base functions on tetrahedral.
More...
#include <src/approximation/EdgePolynomialBase.hpp>
|
using | DofsSideMap = multi_index_container< DofsSideMapData, indexed_by< ordered_non_unique< tag< TypeSide_mi_tag >, composite_key< DofsSideMapData, member< DofsSideMapData, EntityType, &DofsSideMapData::type >, member< DofsSideMapData, int, &DofsSideMapData::side > >>, ordered_unique< tag< EntDofIdx_mi_tag >, member< DofsSideMapData, int, &DofsSideMapData::dof > > > > |
| Map entity stype and side to element/entity dof index. More...
|
|
static MoFEMErrorCode | getLibVersion (Version &version) |
| Get library version. More...
|
|
static MoFEMErrorCode | getFileVersion (moab::Interface &moab, Version &version) |
| Get database major version. More...
|
|
static MoFEMErrorCode | setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD)) |
| Get database major version. More...
|
|
static MoFEMErrorCode | getInterfaceVersion (Version &version) |
| Get database major version. More...
|
|
Calculate base functions on tetrahedral.
Definition at line 19 of file EdgePolynomialBase.hpp.
◆ EdgePolynomialBase()
MoFEM::EdgePolynomialBase::EdgePolynomialBase |
( |
| ) |
|
|
default |
◆ ~EdgePolynomialBase()
virtual MoFEM::EdgePolynomialBase::~EdgePolynomialBase |
( |
| ) |
|
|
virtualdefault |
◆ getValue()
Reimplemented from MoFEM::BaseFunction.
Definition at line 15 of file EdgePolynomialBase.cpp.
21 int nb_gauss_pts = pts.size2();
28 "Wrong dimension of pts, should be at least 3 rows with coordinates");
35 data.
dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 2,
39 &pts(0, 0), nb_gauss_pts);
45 (
unsigned int)nb_gauss_pts)
47 "Base functions or nodes has wrong number of integration points "
51 data.
dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 2 * 1,
53 for (
auto gg = 0; gg != nb_gauss_pts; ++gg)
◆ getValueH1()
◆ getValueH1AinsworthBase()
Definition at line 196 of file EdgePolynomialBase.cpp.
201 PetscErrorCode (*base_polynomials)(
int p,
double s,
double *diff_s,
double *
L,
202 double *
diffL,
const int dim) =
205 int nb_gauss_pts = pts.size2();
207 const int side_number = 0;
224 for (
int gg = 0; gg != nb_gauss_pts; gg++) {
226 double s = 2 * pts(0, gg) - 1;
234 &*
L.data().begin(), &*
diffL.data().begin(), 1);
236 for (
unsigned int pp = 0;
◆ getValueH1BernsteinBezierBase()
Definition at line 101 of file EdgePolynomialBase.cpp.
106 const int nb_gauss_pts = pts.size2();
109 auto &ptr = data.getBBAlphaIndicesSharedPtr(
field_name);
123 auto &ptr = data.getBBDiffNSharedPtr(
field_name);
130 auto &get_diff_n = get_diff_base(data.
dataOnEntities[MBVERTEX][0]);
131 get_n.resize(nb_gauss_pts, 2,
false);
132 get_diff_n.resize(nb_gauss_pts, 2,
false);
137 (
unsigned int)nb_gauss_pts)
139 "Base functions or nodes has wrong number of integration points "
145 vertex_alpha.resize(2, 2,
false);
146 vertex_alpha(0, 0) = data.
dataOnEntities[MBVERTEX][0].getBBNodeOrder()[0];
147 vertex_alpha(0, 1) = 0;
148 vertex_alpha(1, 0) = 0;
149 vertex_alpha(1, 1) = data.
dataOnEntities[MBVERTEX][0].getBBNodeOrder()[1];
152 1, nb_gauss_pts, vertex_alpha.size1(), &vertex_alpha(0, 0), &
lambda(0, 0),
154 std::array<double, 2>
f = {
155 boost::math::factorial<double>(
157 boost::math::factorial<double>(
160 for (
int g = 0;
g != nb_gauss_pts; ++
g)
161 for (
int n = 0;
n != 2; ++
n) {
163 get_diff_n(
g,
n) *=
f[
n];
169 "Wrong size ent of ent data");
176 get_n.resize(nb_gauss_pts, nb_dofs,
false);
177 get_diff_n.resize(nb_gauss_pts, nb_dofs,
false);
181 edge_alpha.resize(nb_dofs, 2);
184 order, nb_gauss_pts, edge_alpha.size1(), &edge_alpha(0, 0),
189 get_n.resize(nb_gauss_pts, 0,
false);
190 get_diff_n.resize(nb_gauss_pts, 0,
false);
◆ getValueH1DemkowiczBase()
Definition at line 260 of file EdgePolynomialBase.cpp.
266 int nb_gauss_pts = pts.size2();
268 const int side_number = 0;
279 double diff_shape[] = {-1, 1};
283 for (
int gg = 0; gg != nb_gauss_pts; gg++) {
284 const double mu0 = 1.0 - pts(0, gg);
285 const double mu1 = pts(0, gg);
291 order, &*shape.data().begin(), diff_shape,
292 &*data.
dataOnEntities[MBEDGE][side_number].getN(base).data().begin(),
◆ getValueHcurl()
◆ getValueHcurlAinsworthBase()
Definition at line 379 of file EdgePolynomialBase.cpp.
385 PetscErrorCode (*base_polynomials)(
int p,
double s,
double *diff_s,
double *
L,
386 double *
diffL,
const int dim) =
389 int nb_gauss_pts = pts.size2();
399 data.
dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
401 data.
dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
407 &*data.
dataOnEntities[MBEDGE][0].getN(base).data().begin(),
nullptr,
408 nb_gauss_pts, base_polynomials);
411 data.
dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0,
false);
412 data.
dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
◆ getValueHcurlDemkowiczBase()
Definition at line 420 of file EdgePolynomialBase.cpp.
426 int nb_gauss_pts = pts.size2();
431 "No data structure to store base functions");
437 data.
dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
439 data.
dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
450 data.
dataOnEntities[MBEDGE][0].getN(base).resize(nb_gauss_pts, 0,
false);
451 data.
dataOnEntities[MBEDGE][0].getDiffN(base).resize(nb_gauss_pts, 0,
◆ getValueHdiv()
◆ getValueL2()
◆ getValueL2AinsworthBase()
Definition at line 321 of file EdgePolynomialBase.cpp.
327 PetscErrorCode (*base_polynomials)(
int p,
double s,
double *diff_s,
double *
L,
328 double *
diffL,
const int dim) =
334 int nb_gauss_pts = pts.size2();
336 constexpr
int side_number = 0;
345 &*data.
dataOnEntities[MBEDGE][side_number].getN(base).data().begin();
347 &*data.
dataOnEntities[MBEDGE][side_number].getDiffN(base).data().begin();
354 for (
int gg = 0; gg != nb_gauss_pts; gg++) {
355 double mu = 2 * pts(0, gg) - 1;
359 fun_n[qd_shift +
n] =
l[
n];
360 diff_fun_n[qd_shift +
n] = diff_l[
n];
◆ getValueL2DemkowiczBase()
◆ query_interface()
◆ cTx
◆ diffL
The documentation for this struct was generated from the following files:
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
MoFEMErrorCode Hcurl_Ainsworth_EdgeBaseFunctions_MBTET_ON_EDGE(int sense, int p, double *N, double *diffN, double *edgeN, double *diff_edgeN, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Edge based H-curl base functions on edge.
#define NBEDGE_H1(P)
Number of base function on edge for H1 space.
Class used to pass element data to calculate base functions on tet,triangle,edge.
MoFEMErrorCode getValueL2(MatrixDouble &pts)
MoFEMErrorCode getValueH1AinsworthBase(MatrixDouble &pts)
const static char *const ApproximationBaseNames[]
@ L2
field with C-1 continuity
UBlasMatrix< double > MatrixDouble
MoFEMErrorCode Hcurl_Demkowicz_EdgeBaseFunctions_MBEDGE(int sense, int p, double *n, double *diff_n, double *phi, double *diff_phi, int nb_integration_pts)
Edge based H-curl base functions on edge.
MoFEMErrorCode getValueH1(MatrixDouble &pts)
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
PetscErrorCode Lobatto_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Lobatto base functions .
const FieldApproximationBase copyNodeBase
#define CHKERR
Inline error check.
static MoFEMErrorCode baseFunctionsEdge(const int N, const int gdim, const int n_alpha, const int *alpha, const double *lambda, const double *grad_lambda, double *base, double *grad_base)
MoFEMErrorCode getValueH1DemkowiczBase(MatrixDouble &pts)
#define NBEDGE_L2(P)
Number of base functions on edge from L2 space.
UBlasMatrix< int > MatrixInt
MoFEMErrorCode getValueH1BernsteinBezierBase(MatrixDouble &pts)
MoFEMErrorCode getValueL2DemkowiczBase(MatrixDouble &pts)
#define NBEDGE_DEMKOWICZ_HCURL(P)
const FieldApproximationBase bAse
constexpr auto field_name
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getValueL2AinsworthBase(MatrixDouble &pts)
const double v
phase velocity of light in medium (cm/ns)
Calculate base functions on tetrahedral.
#define NBEDGE_AINSWORTH_HCURL(P)
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
@ HCURL
field with continuous tangents
@ MOFEM_DATA_INCONSISTENCY
const std::string fieldName
FieldApproximationBase
approximation base
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
data structure for finite element entity
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ HDIV
field with continuous normal traction
static MoFEMErrorCode generateIndicesEdgeEdge(const int N, int *alpha)
MoFEMErrorCode H1_BubbleShapeFunctions_ONSEGMENT(int p, double *N, double *diffN, double *bubbleN, double *diff_bubbleN, int nb_integration_pts)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
FTensor::Index< 'l', 3 > l
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)