v0.15.0
Loading...
Searching...
No Matches
MoFEM::QuadPolynomialBase Struct Reference

Calculate base functions on triangle. More...

#include "src/approximation/QuadPolynomialBase.hpp"

Inheritance diagram for MoFEM::QuadPolynomialBase:
[legend]
Collaboration diagram for MoFEM::QuadPolynomialBase:
[legend]

Public Member Functions

MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 QuadPolynomialBase ()=default
 
 ~QuadPolynomialBase ()=default
 
MoFEMErrorCode getValue (MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::BaseFunction
virtual MoFEMErrorCode getValue (MatrixDouble &pts_x, MatrixDouble &pts_t, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::BaseFunctionUnknownInterface
virtual ~BaseFunctionUnknownInterface ()=default
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface.
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface reference to pointer of interface.
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface.
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface.
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface.
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface.
 
virtual ~UnknownInterface ()=default
 

Private Member Functions

MoFEMErrorCode getValueH1 (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2 (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurl (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdiv (MatrixDouble &pts)
 
MoFEMErrorCode getValueH1AinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueH1DemkowiczBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2DemkowiczBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurlDemkowiczBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdivDemkowiczBase (MatrixDouble &pts)
 

Private Attributes

EntPolynomialBaseCtxcTx
 
MatrixDouble faceFamily
 
MatrixDouble diffFaceFamily
 

Additional Inherited Members

- Public Types inherited from MoFEM::BaseFunction
using DofsSideMap
 Map entity stype and side to element/entity dof index.
 
- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version.
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version.
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version.
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version.
 

Detailed Description

Calculate base functions on triangle.

Definition at line 21 of file QuadPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ QuadPolynomialBase()

MoFEM::QuadPolynomialBase::QuadPolynomialBase ( )
default

◆ ~QuadPolynomialBase()

MoFEM::QuadPolynomialBase::~QuadPolynomialBase ( )
default

Member Function Documentation

◆ getValue()

MoFEMErrorCode QuadPolynomialBase::getValue ( MatrixDouble & pts,
boost::shared_ptr< BaseFunctionCtx > ctx_ptr )
virtual

Reimplemented from MoFEM::BaseFunction.

Definition at line 430 of file QuadPolynomialBase.cpp.

431 {
433
435
436 int nb_gauss_pts = pts.size2();
437 if (!nb_gauss_pts)
439
440 if (pts.size1() < 2)
441 SETERRQ(
442 PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
443 "Wrong dimension of pts, should be at least 3 rows with coordinates");
444
445 const FieldApproximationBase base = cTx->bAse;
446 EntitiesFieldData &data = cTx->dAta;
447 if (cTx->copyNodeBase != NOBASE)
448 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
449 "Shape base has to be on NOBASE <%s>",
451
452 auto &base_shape = data.dataOnEntities[MBVERTEX][0].getN(cTx->copyNodeBase);
453 auto &diff_base =
454 data.dataOnEntities[MBVERTEX][0].getDiffN(cTx->copyNodeBase);
455
456 if (base_shape.size1() != pts.size2())
457 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
458 "Number of base functions integration points not equal number of "
459 "set integration point");
460 if (base_shape.size2() != 4)
461 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
462 "Number of shape functions should be four");
463 if (diff_base.size1() != pts.size2())
464 SETERRQ(
465 PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
466 "Number of diff base functions integration points not equal number of "
467 "set integration point %ld != %ld",
468 diff_base.size1(), pts.size2());
469 if (diff_base.size2() != 8)
470 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
471 "Number of shape functions should be four");
472
473 switch (cTx->spaceContinuity) {
474 case CONTINUOUS:
475 switch (cTx->sPace) {
476 case H1:
477 CHKERR getValueH1(pts);
478 break;
479 case L2:
480 CHKERR getValueL2(pts);
481 break;
482 case HCURL:
484 break;
485 case HDIV:
486 CHKERR getValueHdiv(pts);
487 break;
488 default:
489 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
490 }
491
492 break;
493
494 default:
495 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Unknown continuity");
496 }
497
499}
FieldApproximationBase
approximation base
Definition definitions.h:58
@ NOBASE
Definition definitions.h:59
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
@ L2
field with C-1 continuity
Definition definitions.h:88
@ H1
continuous field
Definition definitions.h:85
@ HCURL
field with continuous tangents
Definition definitions.h:86
@ HDIV
field with continuous normal traction
Definition definitions.h:87
@ CONTINUOUS
Regular field.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ MOFEM_DATA_INCONSISTENCY
Definition definitions.h:31
@ MOFEM_NOT_IMPLEMENTED
Definition definitions.h:32
static const char *const ApproximationBaseNames[]
Definition definitions.h:72
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
Class used to pass element data to calculate base functions on tet,triangle,edge.
const FieldApproximationBase copyNodeBase
const FieldContinuity spaceContinuity
const FieldApproximationBase bAse
data structure for finite element entity
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
MoFEMErrorCode getValueH1(MatrixDouble &pts)
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
MoFEMErrorCode getValueL2(MatrixDouble &pts)
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.

◆ getValueH1()

MoFEMErrorCode QuadPolynomialBase::getValueH1 ( MatrixDouble & pts)
private

Definition at line 15 of file QuadPolynomialBase.cpp.

15 {
17
18 const FieldApproximationBase base = cTx->bAse;
19 EntitiesFieldData &data = cTx->dAta;
20 data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
21 data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
22 data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(base) =
23 data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(cTx->copyNodeBase);
24
25 switch (cTx->bAse) {
29 break;
32 break;
34 default:
35 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
36 }
37
39}
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base nme:nme847.
Definition definitions.h:60
@ AINSWORTH_LOBATTO_BASE
Definition definitions.h:62
@ DEMKOWICZ_JACOBI_BASE
Definition definitions.h:66
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition definitions.h:64
MoFEMErrorCode getValueH1AinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueH1DemkowiczBase(MatrixDouble &pts)

◆ getValueH1AinsworthBase()

MoFEMErrorCode QuadPolynomialBase::getValueH1AinsworthBase ( MatrixDouble & pts)
private

Definition at line 41 of file QuadPolynomialBase.cpp.

41 {
43
44 EntitiesFieldData &data = cTx->dAta;
45 const FieldApproximationBase base = cTx->bAse;
46 const auto copy_base = cTx->copyNodeBase;
47
48 if (cTx->basePolynomialsType0 == NULL)
49 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
50 "Polynomial type not set");
51 PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
52 double *diffL, const int dim) =
54
55 auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
56 auto &copy_diff_base_fun =
57 data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
58
59 int nb_gauss_pts = pts.size2();
60 auto &vert_dat = data.dataOnEntities[MBVERTEX][0];
61
62 if (data.spacesOnEntities[MBEDGE].test(H1)) {
63 // edges
64 if (data.dataOnEntities[MBEDGE].size() != 4)
65 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
66 "should be four edges on quad");
67
68 int sense[4], order[4];
69 double *H1edgeN[4], *diffH1edgeN[4];
70 for (int ee = 0; ee != 4; ++ee) {
71 auto &ent_dat = data.dataOnEntities[MBEDGE][ee];
72 if (ent_dat.getSense() == 0)
73 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "sense not set");
74
75 sense[ee] = ent_dat.getSense();
76 order[ee] = ent_dat.getOrder();
77 int nb_dofs = NBEDGE_H1(ent_dat.getOrder());
78 ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
79 ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
80 H1edgeN[ee] = &*ent_dat.getN(base).data().begin();
81 diffH1edgeN[ee] = &*ent_dat.getDiffN(base).data().begin();
82 }
84 sense, order, &*copy_base_fun.data().begin(),
85 &*copy_diff_base_fun.data().begin(), H1edgeN, diffH1edgeN, nb_gauss_pts,
86 base_polynomials);
87 }
88
89 if (data.spacesOnEntities[MBQUAD].test(H1)) {
90
91 // face
92 if (data.dataOnEntities[MBQUAD].size() != 1)
93 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
94 "should be one quad to store bubble base on quad");
95
96 auto &ent_dat = data.dataOnEntities[MBQUAD][0];
97 int nb_dofs = NBFACEQUAD_H1(ent_dat.getOrder());
98 ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
99 ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
100 int face_nodes[] = {0, 1, 2, 3};
102 face_nodes, ent_dat.getOrder(), &*vert_dat.getN(base).data().begin(),
103 &*vert_dat.getDiffN(base).data().begin(),
104 &*ent_dat.getN(base).data().begin(),
105 &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts,
106 base_polynomials);
107 }
108
110}
constexpr int order
#define NBFACEQUAD_H1(P)
Number of base functions on quad for H1 space.
PetscErrorCode H1_QuadShapeFunctions_MBQUAD(int *faces_nodes, int p, double *N, double *diffN, double *faceN, double *diff_faceN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition h1.c:959
#define NBEDGE_H1(P)
Number of base function on edge for H1 space.
PetscErrorCode H1_EdgeShapeFunctions_MBQUAD(int *sense, int *p, double *N, double *diffN, double *edgeN[4], double *diff_edgeN[4], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition h1.c:1091
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types

◆ getValueH1DemkowiczBase()

MoFEMErrorCode QuadPolynomialBase::getValueH1DemkowiczBase ( MatrixDouble & pts)
private

Definition at line 112 of file QuadPolynomialBase.cpp.

112 {
114
115 EntitiesFieldData &data = cTx->dAta;
116 const FieldApproximationBase base = cTx->bAse;
117 const auto copy_base = cTx->copyNodeBase;
118
119 int nb_gauss_pts = pts.size2();
120 auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
121 auto &copy_diff_base_fun =
122 data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
123
124 if (data.spacesOnEntities[MBEDGE].test(H1)) {
125 // edges
126 if (data.dataOnEntities[MBEDGE].size() != 4)
127 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
128 "should be four edges on quad");
129
130 int sense[4], order[4];
131 double *H1edgeN[4], *diffH1edgeN[4];
132 for (int ee = 0; ee != 4; ++ee) {
133 auto &ent_dat = data.dataOnEntities[MBEDGE][ee];
134 if (ent_dat.getSense() == 0)
135 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "sense not set");
136
137 sense[ee] = ent_dat.getSense();
138 order[ee] = ent_dat.getOrder();
139 int nb_dofs = NBEDGE_H1(ent_dat.getOrder());
140 ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
141 ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
142 H1edgeN[ee] = &*ent_dat.getN(base).data().begin();
143 diffH1edgeN[ee] = &*ent_dat.getDiffN(base).data().begin();
144 }
146 sense, order, &*copy_base_fun.data().begin(),
147 &*copy_diff_base_fun.data().begin(), H1edgeN, diffH1edgeN,
148 nb_gauss_pts);
149 }
150
151 if (data.spacesOnEntities[MBQUAD].test(H1)) {
152
153 // face
154 if (data.dataOnEntities[MBQUAD].size() != 1)
155 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
156 "should be one quad to store bubble base on quad");
157
158 auto &ent_dat = data.dataOnEntities[MBQUAD][0];
159 int nb_dofs = NBFACEQUAD_H1(ent_dat.getOrder());
160 int p = ent_dat.getOrder();
161 int order[2] = {p, p};
162 ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
163 ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
164
165 int face_nodes[] = {0, 1, 2, 3};
166 if (nb_dofs) {
168 face_nodes, order, &*copy_base_fun.data().begin(),
169 &*copy_diff_base_fun.data().begin(),
170 &*ent_dat.getN(base).data().begin(),
171 &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts);
172 }
173 }
174
176}
MoFEMErrorCode H1_FaceShapeFunctions_ONQUAD(int *face_nodes, int *p, double *N, double *diffN, double *faceN, double *diff_faceN, int nb_integration_pts)
H1 Face bubble functions on Quad.
MoFEMErrorCode H1_EdgeShapeFunctions_ONQUAD(int *sense, int *p, double *N, double *diffN, double *edgeN[4], double *edgeDiffN[4], int nb_integration_pts)
H1 Edge base functions on Quad.

◆ getValueHcurl()

MoFEMErrorCode QuadPolynomialBase::getValueHcurl ( MatrixDouble & pts)
private

Definition at line 227 of file QuadPolynomialBase.cpp.

227 {
229
230 switch (cTx->bAse) {
233 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
234 "Ainsworth not implemented on quad for Hcurl base");
235 break;
238 break;
240 default:
241 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
242 }
243
245}
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode QuadPolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble & pts)
private

Definition at line 248 of file QuadPolynomialBase.cpp.

248 {
250
251 EntitiesFieldData &data = cTx->dAta;
252 const FieldApproximationBase base = cTx->bAse;
253 const auto copy_base = cTx->copyNodeBase;
254
255 int nb_gauss_pts = pts.size2();
256 auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
257 auto &copy_diff_base_fun =
258 data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
259
260 // Calculation H-curl on quad edges
261 if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
262
263 if (data.dataOnEntities[MBEDGE].size() != 4)
264 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
265 "wrong number of edges on quad, should be 4 but is %ld",
266 data.dataOnEntities[MBEDGE].size());
267
268 int sense[4], order[4];
269 double *hcurl_edge_n[4];
270 double *diff_hcurl_edge_n[4];
271
272 for (int ee = 0; ee != 4; ++ee) {
273
274 if (data.dataOnEntities[MBEDGE][ee].getSense() == 0)
275 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
276 "orientation (sense) of edge is not set");
277
278 sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
279 order[ee] = data.dataOnEntities[MBEDGE][ee].getOrder();
280 int nb_dofs =
281 NBEDGE_DEMKOWICZ_HCURL(data.dataOnEntities[MBEDGE][ee].getOrder());
282
283 data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
284 3 * nb_dofs, false);
285 data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
286 nb_gauss_pts, 3 * 2 * nb_dofs, false);
287
288 hcurl_edge_n[ee] =
289 &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
290 diff_hcurl_edge_n[ee] =
291 &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
292 }
294 sense, order, &*copy_base_fun.data().begin(),
295 &*copy_diff_base_fun.data().begin(), hcurl_edge_n, diff_hcurl_edge_n,
296 nb_gauss_pts);
297 }
298
299 if (data.spacesOnEntities[MBQUAD].test(HCURL)) {
300
301 // face
302 if (data.dataOnEntities[MBQUAD].size() != 1)
303 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
304 "No data struture to keep base functions on face");
305
306 int p = data.dataOnEntities[MBQUAD][0].getOrder();
307 const int nb_dofs_family = NBFACEQUAD_DEMKOWICZ_FAMILY_HCURL(p, p);
308 if (nb_dofs_family) {
309 faceFamily.resize(2, 3 * nb_dofs_family * nb_gauss_pts, false);
310 diffFaceFamily.resize(2, 6 * nb_dofs_family * nb_gauss_pts, false);
311
312 int order[2] = {p, p};
313 double *face_family_ptr[] = {&faceFamily(0, 0), &faceFamily(1, 0)};
314 double *diff_face_family_ptr[] = {&diffFaceFamily(0, 0),
315 &diffFaceFamily(1, 0)};
316 int face_nodes[] = {0, 1, 2, 3};
318 face_nodes, order, &*copy_base_fun.data().begin(),
319 &*copy_diff_base_fun.data().begin(), face_family_ptr,
320 diff_face_family_ptr, nb_gauss_pts);
321 }
322
323 // put family back
324
325 const int nb_dofs = NBFACEQUAD_DEMKOWICZ_HCURL(p);
326 auto &face_n = data.dataOnEntities[MBQUAD][0].getN(base);
327 auto &diff_face_n = data.dataOnEntities[MBQUAD][0].getDiffN(base);
328 face_n.resize(nb_gauss_pts, 3 * nb_dofs, false);
329 diff_face_n.resize(nb_gauss_pts, 3 * 2 * nb_dofs, false);
330
331 if (nb_dofs) {
332 double *ptr_f0 = &faceFamily(0, 0);
333 double *ptr_f1 = &faceFamily(1, 0);
334 double *ptr = &face_n(0, 0);
335 for (int n = 0; n != faceFamily.size2() / 3; ++n) {
336 for (int j = 0; j != 3; ++j) {
337 *ptr = *ptr_f0;
338 ++ptr;
339 ++ptr_f0;
340 }
341 for (int j = 0; j != 3; ++j) {
342 *ptr = *ptr_f1;
343 ++ptr;
344 ++ptr_f1;
345 }
346 }
347
348 double *diff_ptr_f0 = &diffFaceFamily(0, 0);
349 double *diff_ptr_f1 = &diffFaceFamily(1, 0);
350 double *diff_ptr = &diff_face_n(0, 0);
351 for (int n = 0; n != diffFaceFamily.size2() / 6; ++n) {
352 for (int j = 0; j != 6; ++j) {
353 *diff_ptr = *diff_ptr_f0;
354 ++diff_ptr;
355 ++diff_ptr_f0;
356 }
357 for (int j = 0; j != 6; ++j) {
358 *diff_ptr = *diff_ptr_f1;
359 ++diff_ptr;
360 ++diff_ptr_f1;
361 }
362 }
363 }
364 }
366}
#define NBFACEQUAD_DEMKOWICZ_FAMILY_HCURL(P, Q)
Number of base functions on quad for Hcurl space.
#define NBEDGE_DEMKOWICZ_HCURL(P)
#define NBFACEQUAD_DEMKOWICZ_HCURL(P)
const double n
refractive index of diffusive medium
FTensor::Index< 'j', 3 > j
MoFEMErrorCode Hcurl_EdgeShapeFunctions_ONQUAD(int *sense, int *p, double *N, double *diffN, double *edgeN[4], double *curl_edgeN[4], int nb_integration_pts)
MoFEMErrorCode Hcurl_FaceShapeFunctions_ONQUAD(int *face_nodes, int *p, double *N, double *diffN, double *faceN[], double *diff_faceN[], int nb_integration_pts)

◆ getValueHdiv()

MoFEMErrorCode QuadPolynomialBase::getValueHdiv ( MatrixDouble & pts)
private

Definition at line 368 of file QuadPolynomialBase.cpp.

368 {
370
371 switch (cTx->bAse) {
374 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
375 "Ainsworth not implemented on quad for Hdiv base");
376 break;
379 break;
381 default:
382 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
383 }
384
386}
MoFEMErrorCode getValueHdivDemkowiczBase(MatrixDouble &pts)

◆ getValueHdivDemkowiczBase()

MoFEMErrorCode QuadPolynomialBase::getValueHdivDemkowiczBase ( MatrixDouble & pts)
private

Definition at line 389 of file QuadPolynomialBase.cpp.

389 {
391
392 EntitiesFieldData &data = cTx->dAta;
393 const FieldApproximationBase base = cTx->bAse;
394 const auto copy_base = cTx->copyNodeBase;
395 int nb_gauss_pts = pts.size2();
396
397 auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
398 auto &copy_diff_base_fun =
399 data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
400
401 if (data.spacesOnEntities[MBQUAD].test(HDIV)) {
402
403 // face
404 if (data.dataOnEntities[MBQUAD].size() != 1)
405 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
406 "No data struture to keep base functions on face");
407
408 int p = data.dataOnEntities[MBQUAD][0].getOrder();
409 const int nb_dofs = NBFACEQUAD_DEMKOWICZ_HDIV(p);
410 auto &face_n = data.dataOnEntities[MBQUAD][0].getN(base);
411 auto &diff_face_n = data.dataOnEntities[MBQUAD][0].getDiffN(base);
412 face_n.resize(nb_gauss_pts, 3 * nb_dofs, false);
413 diff_face_n.resize(nb_gauss_pts, 6 * nb_dofs, false);
414
415 if (nb_dofs) {
416
417 std::array<int, 2> order = {p, p};
418 std::array<int, 6> face_nodes = {0, 1, 2, 3};
420 face_nodes.data(), order.data(), &*copy_base_fun.data().begin(),
421 &*copy_diff_base_fun.data().begin(), &*face_n.data().begin(),
422 &*diff_face_n.data().begin(), nb_gauss_pts);
423 }
424 }
425
427}
#define NBFACEQUAD_DEMKOWICZ_HDIV(P)
MoFEMErrorCode Hdiv_FaceShapeFunctions_ONQUAD(int *face_nodes, int *p, double *N, double *diffN, double *faceN, double *diff_faceN, int nb_integration_pts)

◆ getValueL2()

MoFEMErrorCode QuadPolynomialBase::getValueL2 ( MatrixDouble & pts)
private

Definition at line 178 of file QuadPolynomialBase.cpp.

178 {
180
181 switch (cTx->bAse) {
184 break;
186 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
187 "Ainsworth not implemented on quad for L2 base");
188 break;
191 break;
193 default:
194 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
195 }
196
198}
MoFEMErrorCode getValueL2DemkowiczBase(MatrixDouble &pts)

◆ getValueL2DemkowiczBase()

MoFEMErrorCode QuadPolynomialBase::getValueL2DemkowiczBase ( MatrixDouble & pts)
private

Definition at line 200 of file QuadPolynomialBase.cpp.

200 {
202
203 EntitiesFieldData &data = cTx->dAta;
204 const FieldApproximationBase base = cTx->bAse;
205 const auto copy_base = cTx->copyNodeBase;
206
207 int nb_gauss_pts = pts.size2();
208 auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
209 auto &copy_diff_base_fun =
210 data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
211
212 auto &ent_dat = data.dataOnEntities[MBQUAD][0];
213 int p = ent_dat.getOrder();
214 int order[2] = {p, p};
215 int nb_dofs = NBFACEQUAD_L2(p);
216 ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
217 ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
218
220 order, &*copy_base_fun.data().begin(),
221 &*copy_diff_base_fun.data().begin(), &*ent_dat.getN(base).data().begin(),
222 &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts);
223
225}
#define NBFACEQUAD_L2(P)
Number of base functions on quad for L2 space.
MoFEMErrorCode L2_FaceShapeFunctions_ONQUAD(int *p, double *N, double *diffN, double *face_buble, double *diff_face_bubble, int nb_integration_pts)
L2 Face base functions on Quad.

◆ query_interface()

MoFEMErrorCode QuadPolynomialBase::query_interface ( boost::typeindex::type_index type_index,
UnknownInterface ** iface ) const
virtual

Reimplemented from MoFEM::BaseFunction.

Definition at line 9 of file QuadPolynomialBase.cpp.

10 {
11 *iface = const_cast<QuadPolynomialBase *>(this);
12 return 0;
13}
Calculate base functions on triangle.

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* MoFEM::QuadPolynomialBase::cTx
private

Definition at line 32 of file QuadPolynomialBase.hpp.

◆ diffFaceFamily

MatrixDouble MoFEM::QuadPolynomialBase::diffFaceFamily
private

Definition at line 46 of file QuadPolynomialBase.hpp.

◆ faceFamily

MatrixDouble MoFEM::QuadPolynomialBase::faceFamily
private

Definition at line 45 of file QuadPolynomialBase.hpp.


The documentation for this struct was generated from the following files: