v0.13.0
Public Member Functions | Private Member Functions | Private Attributes | List of all members
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. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
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

- Static Public Member Functions inherited from MoFEM::UnknownInterface
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...
 

Detailed Description

Calculate base functions on triangle.

Definition at line 33 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 446 of file QuadPolynomialBase.cpp.

447  {
449 
450  cTx = ctx_ptr->getInterface<EntPolynomialBaseCtx>();
451 
452  int nb_gauss_pts = pts.size2();
453  if (!nb_gauss_pts)
455 
456  if (pts.size1() < 2)
457  SETERRQ(
458  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
459  "Wrong dimension of pts, should be at least 3 rows with coordinates");
460 
461  const FieldApproximationBase base = cTx->bAse;
462  EntitiesFieldData &data = cTx->dAta;
463  if (cTx->copyNodeBase != NOBASE)
464  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
465  "Shape base has to be on NOBASE", ApproximationBaseNames[base]);
466 
467  auto &base_shape = data.dataOnEntities[MBVERTEX][0].getN(cTx->copyNodeBase);
468  auto &diff_base =
469  data.dataOnEntities[MBVERTEX][0].getDiffN(cTx->copyNodeBase);
470 
471  if (base_shape.size1() != pts.size2())
472  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
473  "Number of base functions integration points not equal number of "
474  "set integration point");
475  if (base_shape.size2() != 4)
476  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
477  "Number of shape functions should be four");
478  if (diff_base.size1() != pts.size2())
479  SETERRQ2(
480  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
481  "Number of diff base functions integration points not equal number of "
482  "set integration point %d != %d",
483  diff_base.size1(), pts.size2());
484  if (diff_base.size2() != 8)
485  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
486  "Number of shape functions should be four");
487 
488  switch (cTx->sPace) {
489  case H1:
490  CHKERR getValueH1(pts);
491  break;
492  case L2:
493  CHKERR getValueL2(pts);
494  break;
495  case HCURL:
496  CHKERR getValueHcurl(pts);
497  break;
498  case HDIV:
499  CHKERR getValueHdiv(pts);
500  break;
501  default:
502  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
503  }
504 
506 }
FieldApproximationBase
approximation base
Definition: definitions.h:71
@ NOBASE
Definition: definitions.h:72
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
@ L2
field with C-1 continuity
Definition: definitions.h:101
@ H1
continuous field
Definition: definitions.h:98
@ HCURL
field with continuous tangents
Definition: definitions.h:99
@ HDIV
field with continuous normal traction
Definition: definitions.h:100
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:45
static const char *const ApproximationBaseNames[]
Definition: definitions.h:85
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
Class used to pass element data to calculate base functions on tet,triangle,edge.
const FieldApproximationBase copyNodeBase
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 refernce to pointer of interface.

◆ getValueH1()

MoFEMErrorCode QuadPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 29 of file QuadPolynomialBase.cpp.

29  {
31 
32  const FieldApproximationBase base = cTx->bAse;
33  EntitiesFieldData &data = cTx->dAta;
34  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
35  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
36  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(base) =
37  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(cTx->copyNodeBase);
38 
39  switch (cTx->bAse) {
43  break;
46  break;
48  default:
49  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
50  }
51 
53 }
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:73
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:75
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:79
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:77
MoFEMErrorCode getValueH1AinsworthBase(MatrixDouble &pts)
MoFEMErrorCode getValueH1DemkowiczBase(MatrixDouble &pts)

◆ getValueH1AinsworthBase()

MoFEMErrorCode QuadPolynomialBase::getValueH1AinsworthBase ( MatrixDouble pts)
private

Definition at line 55 of file QuadPolynomialBase.cpp.

55  {
57 
58  EntitiesFieldData &data = cTx->dAta;
59  const FieldApproximationBase base = cTx->bAse;
60  const auto copy_base = cTx->copyNodeBase;
61 
62  if (cTx->basePolynomialsType0 == NULL)
63  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
64  "Polynomial type not set");
65  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
66  double *diffL, const int dim) =
68 
69  auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
70  auto &copy_diff_base_fun =
71  data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
72 
73  int nb_gauss_pts = pts.size2();
74  auto &vert_dat = data.dataOnEntities[MBVERTEX][0];
75 
76  if (data.spacesOnEntities[MBEDGE].test(H1)) {
77  // edges
78  if (data.dataOnEntities[MBEDGE].size() != 4)
79  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
80  "should be four edges on quad");
81 
82  int sense[4], order[4];
83  double *H1edgeN[4], *diffH1edgeN[4];
84  for (int ee = 0; ee != 4; ++ee) {
85  auto &ent_dat = data.dataOnEntities[MBEDGE][ee];
86  if (ent_dat.getSense() == 0)
87  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "sense not set");
88 
89  sense[ee] = ent_dat.getSense();
90  order[ee] = ent_dat.getOrder();
91  int nb_dofs = NBEDGE_H1(ent_dat.getOrder());
92  ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
93  ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
94  H1edgeN[ee] = &*ent_dat.getN(base).data().begin();
95  diffH1edgeN[ee] = &*ent_dat.getDiffN(base).data().begin();
96  }
98  sense, order, &*copy_base_fun.data().begin(),
99  &*copy_diff_base_fun.data().begin(), H1edgeN, diffH1edgeN, nb_gauss_pts,
100  base_polynomials);
101  }
102 
103  if (data.spacesOnEntities[MBQUAD].test(H1)) {
104 
105  // face
106  if (data.dataOnEntities[MBQUAD].size() != 1)
107  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
108  "should be one quad to store bubble base on quad");
109 
110  auto &ent_dat = data.dataOnEntities[MBQUAD][0];
111  int nb_dofs = NBFACEQUAD_H1(ent_dat.getOrder());
112  ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
113  ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
114  int face_nodes[] = {0, 1, 2, 3};
116  face_nodes, ent_dat.getOrder(),
117  &*vert_dat.getN(base).data().begin(),
118  &*vert_dat.getDiffN(base).data().begin(),
119  &*ent_dat.getN(base).data().begin(),
120  &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts,
121  base_polynomials);
122  }
123 
125 }
static Index< 'p', 3 > p
const int dim
#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:969
#define NBEDGE_H1(P)
Numer 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:1101
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 127 of file QuadPolynomialBase.cpp.

127  {
129 
130  EntitiesFieldData &data = cTx->dAta;
131  const FieldApproximationBase base = cTx->bAse;
132  const auto copy_base = cTx->copyNodeBase;
133 
134  int nb_gauss_pts = pts.size2();
135  auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
136  auto &copy_diff_base_fun =
137  data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
138 
139  if (data.spacesOnEntities[MBEDGE].test(H1)) {
140  // edges
141  if (data.dataOnEntities[MBEDGE].size() != 4)
142  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
143  "should be four edges on quad");
144 
145  int sense[4], order[4];
146  double *H1edgeN[4], *diffH1edgeN[4];
147  for (int ee = 0; ee != 4; ++ee) {
148  auto &ent_dat = data.dataOnEntities[MBEDGE][ee];
149  if (ent_dat.getSense() == 0)
150  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "sense not set");
151 
152  sense[ee] = ent_dat.getSense();
153  order[ee] = ent_dat.getOrder();
154  int nb_dofs = NBEDGE_H1(ent_dat.getOrder());
155  ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
156  ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
157  H1edgeN[ee] = &*ent_dat.getN(base).data().begin();
158  diffH1edgeN[ee] = &*ent_dat.getDiffN(base).data().begin();
159  }
161  sense, order, &*copy_base_fun.data().begin(),
162  &*copy_diff_base_fun.data().begin(), H1edgeN, diffH1edgeN,
163  nb_gauss_pts);
164  }
165 
166  if (data.spacesOnEntities[MBQUAD].test(H1)) {
167 
168  // face
169  if (data.dataOnEntities[MBQUAD].size() != 1)
170  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
171  "should be one quad to store bubble base on quad");
172 
173  auto &ent_dat = data.dataOnEntities[MBQUAD][0];
174  int nb_dofs = NBFACEQUAD_H1(ent_dat.getOrder());
175  int p = ent_dat.getOrder();
176  int order[2] = {p, p};
177  ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
178  ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
179 
180  int face_nodes[] = {0, 1, 2, 3};
181  if (nb_dofs) {
183  face_nodes, order, &*copy_base_fun.data().begin(),
184  &*copy_diff_base_fun.data().begin(),
185  &*ent_dat.getN(base).data().begin(),
186  &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts);
187  }
188  }
189 
191 }
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 242 of file QuadPolynomialBase.cpp.

242  {
244 
245  switch (cTx->bAse) {
248  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
249  "Ainsworth not implemented on quad for Hcurl base");
250  break;
253  break;
255  default:
256  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
257  }
258 
260 }
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode QuadPolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 263 of file QuadPolynomialBase.cpp.

263  {
265 
266  EntitiesFieldData &data = cTx->dAta;
267  const FieldApproximationBase base = cTx->bAse;
268  const auto copy_base = cTx->copyNodeBase;
269 
270  int nb_gauss_pts = pts.size2();
271  auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
272  auto &copy_diff_base_fun =
273  data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
274 
275  // Calculation H-curl on quad edges
276  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
277 
278  if (data.dataOnEntities[MBEDGE].size() != 4)
279  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
280  "wrong number of edges on quad, should be 4 but is %d",
281  data.dataOnEntities[MBEDGE].size());
282 
283  int sense[4], order[4];
284  double *hcurl_edge_n[4];
285  double *diff_hcurl_edge_n[4];
286 
287  for (int ee = 0; ee != 4; ++ee) {
288 
289  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0)
290  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
291  "orientation (sense) of edge is not set");
292 
293  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
294  order[ee] = data.dataOnEntities[MBEDGE][ee].getOrder();
295  int nb_dofs = NBEDGE_DEMKOWICZ_HCURL(
296  data.dataOnEntities[MBEDGE][ee].getOrder());
297 
298  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
299  3 * nb_dofs, false);
300  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
301  nb_gauss_pts, 3 * 2 * nb_dofs, false);
302 
303  hcurl_edge_n[ee] =
304  &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
305  diff_hcurl_edge_n[ee] =
306  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
307  }
309  sense, order, &*copy_base_fun.data().begin(),
310  &*copy_diff_base_fun.data().begin(), hcurl_edge_n, diff_hcurl_edge_n,
311  nb_gauss_pts);
312  }
313 
314  if (data.spacesOnEntities[MBQUAD].test(HCURL)) {
315 
316  // face
317  if (data.dataOnEntities[MBQUAD].size() != 1)
318  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
319  "No data struture to keep base functions on face");
320 
321  int p = data.dataOnEntities[MBQUAD][0].getOrder();
322  const int nb_dofs_family = NBFACEQUAD_DEMKOWICZ_FAMILY_HCURL(p, p);
323  if (nb_dofs_family) {
324  faceFamily.resize(2, 3 * nb_dofs_family * nb_gauss_pts, false);
325  diffFaceFamily.resize(2, 6 * nb_dofs_family * nb_gauss_pts, false);
326 
327  int order[2] = {p, p};
328  double *face_family_ptr[] = {&faceFamily(0, 0), &faceFamily(1, 0)};
329  double *diff_face_family_ptr[] = {&diffFaceFamily(0, 0),
330  &diffFaceFamily(1, 0)};
331  int face_nodes[] = {0, 1, 2, 3};
333  face_nodes, order, &*copy_base_fun.data().begin(),
334  &*copy_diff_base_fun.data().begin(), face_family_ptr,
335  diff_face_family_ptr, nb_gauss_pts);
336  }
337 
338  // put family back
339 
340  const int nb_dofs = NBFACEQUAD_DEMKOWICZ_HCURL(p);
341  auto &face_n = data.dataOnEntities[MBQUAD][0].getN(base);
342  auto &diff_face_n = data.dataOnEntities[MBQUAD][0].getDiffN(base);
343  face_n.resize(nb_gauss_pts, 3 * nb_dofs, false);
344  diff_face_n.resize(nb_gauss_pts, 3 * 2 * nb_dofs, false);
345 
346  if (nb_dofs) {
347  double *ptr_f0 = &faceFamily(0, 0);
348  double *ptr_f1 = &faceFamily(1, 0);
349  double *ptr = &face_n(0, 0);
350  for (int n = 0; n != faceFamily.size2() / 3; ++n) {
351  for (int j = 0; j != 3; ++j) {
352  *ptr = *ptr_f0;
353  ++ptr;
354  ++ptr_f0;
355  }
356  for (int j = 0; j != 3; ++j) {
357  *ptr = *ptr_f1;
358  ++ptr;
359  ++ptr_f1;
360  }
361  }
362 
363  double *diff_ptr_f0 = &diffFaceFamily(0, 0);
364  double *diff_ptr_f1 = &diffFaceFamily(1, 0);
365  double *diff_ptr = &diff_face_n(0, 0);
366  for (int n = 0; n != diffFaceFamily.size2() / 6; ++n) {
367  for (int j = 0; j != 6; ++j) {
368  *diff_ptr = *diff_ptr_f0;
369  ++diff_ptr;
370  ++diff_ptr_f0;
371  }
372  for (int j = 0; j != 6; ++j) {
373  *diff_ptr = *diff_ptr_f1;
374  ++diff_ptr;
375  ++diff_ptr_f1;
376  }
377  }
378  }
379  }
381 }
static Index< 'n', 3 > n
#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)
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 383 of file QuadPolynomialBase.cpp.

383  {
385 
386  switch (cTx->bAse) {
389  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
390  "Ainsworth not implemented on quad for Hdiv base");
391  break;
394  break;
396  default:
397  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
398  }
399 
401 }
MoFEMErrorCode getValueHdivDemkowiczBase(MatrixDouble &pts)

◆ getValueHdivDemkowiczBase()

MoFEMErrorCode QuadPolynomialBase::getValueHdivDemkowiczBase ( MatrixDouble pts)
private

Definition at line 404 of file QuadPolynomialBase.cpp.

404  {
406 
407  EntitiesFieldData &data = cTx->dAta;
408  const FieldApproximationBase base = cTx->bAse;
409  const auto copy_base = cTx->copyNodeBase;
410  int nb_gauss_pts = pts.size2();
411 
412  auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
413  auto &copy_diff_base_fun =
414  data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
415 
416  if (data.spacesOnEntities[MBQUAD].test(HDIV)) {
417 
418  // face
419  if (data.dataOnEntities[MBQUAD].size() != 1)
420  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
421  "No data struture to keep base functions on face");
422 
423  int p = data.dataOnEntities[MBQUAD][0].getOrder();
424  const int nb_dofs = NBFACEQUAD_DEMKOWICZ_HDIV(p);
425  auto &face_n = data.dataOnEntities[MBQUAD][0].getN(base);
426  auto &diff_face_n = data.dataOnEntities[MBQUAD][0].getDiffN(base);
427  face_n.resize(nb_gauss_pts, 3 * nb_dofs, false);
428  diff_face_n.resize(nb_gauss_pts, 6 * nb_dofs, false);
429 
430  if (nb_dofs) {
431 
432  std::array<int, 2> order = {p, p};
433  std::array<int, 6> face_nodes = {0, 1, 2, 3};
435  face_nodes.data(), order.data(), &*copy_base_fun.data().begin(),
436  &*copy_diff_base_fun.data().begin(), &*face_n.data().begin(),
437  &*diff_face_n.data().begin(), nb_gauss_pts);
438 
439  }
440  }
441 
443 }
#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 193 of file QuadPolynomialBase.cpp.

193  {
195 
196  switch (cTx->bAse) {
199  break;
201  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
202  "Ainsworth not implemented on quad for L2 base");
203  break;
206  break;
208  default:
209  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
210  }
211 
213 }
MoFEMErrorCode getValueL2DemkowiczBase(MatrixDouble &pts)

◆ getValueL2DemkowiczBase()

MoFEMErrorCode QuadPolynomialBase::getValueL2DemkowiczBase ( MatrixDouble pts)
private

Definition at line 215 of file QuadPolynomialBase.cpp.

215  {
217 
218  EntitiesFieldData &data = cTx->dAta;
219  const FieldApproximationBase base = cTx->bAse;
220  const auto copy_base = cTx->copyNodeBase;
221 
222  int nb_gauss_pts = pts.size2();
223  auto &copy_base_fun = data.dataOnEntities[MBVERTEX][0].getN(copy_base);
224  auto &copy_diff_base_fun =
225  data.dataOnEntities[MBVERTEX][0].getDiffN(copy_base);
226 
227  auto &ent_dat = data.dataOnEntities[MBQUAD][0];
228  int p = ent_dat.getOrder();
229  int order[2] = {p, p};
230  int nb_dofs = NBFACEQUAD_L2(p);
231  ent_dat.getN(base).resize(nb_gauss_pts, nb_dofs, false);
232  ent_dat.getDiffN(base).resize(nb_gauss_pts, 2 * nb_dofs, false);
233 
235  order, &*copy_base_fun.data().begin(),
236  &*copy_diff_base_fun.data().begin(), &*ent_dat.getN(base).data().begin(),
237  &*ent_dat.getDiffN(base).data().begin(), nb_gauss_pts);
238 
240 }
#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 23 of file QuadPolynomialBase.cpp.

24  {
25  *iface = const_cast<QuadPolynomialBase *>(this);
26  return 0;
27 }
Calculate base functions on triangle.

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* MoFEM::QuadPolynomialBase::cTx
private

Definition at line 44 of file QuadPolynomialBase.hpp.

◆ diffFaceFamily

MatrixDouble MoFEM::QuadPolynomialBase::diffFaceFamily
private

Definition at line 58 of file QuadPolynomialBase.hpp.

◆ faceFamily

MatrixDouble MoFEM::QuadPolynomialBase::faceFamily
private

Definition at line 57 of file QuadPolynomialBase.hpp.


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