v0.8.16
Public Member Functions | Private Member Functions | Private Attributes | List of all members
MoFEM::TriPolynomialBase Struct Reference

Calculate base functions on triangle. More...

#include <src/approximation/TriPolynomialBase.hpp>

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

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 TriPolynomialBase ()
 
 ~TriPolynomialBase ()
 
MoFEMErrorCode getValue (MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::BaseFunction
 BaseFunction ()
 
 ~BaseFunction ()
 
virtual MoFEMErrorCode getValue (MatrixDouble &pts_x, MatrixDouble &pts_t, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of 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 ()
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 

Private Member Functions

MoFEMErrorCode getValueH1 (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2 (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdiv (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurl (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdivAinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurlAinsworthBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdivDemkowiczBase (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurlDemkowiczBase (MatrixDouble &pts)
 

Private Attributes

EntPolynomialBaseCtxcTx
 
ublas::matrix< MatrixDoubleN_face_edge
 
ublas::vector< MatrixDoubleN_face_bubble
 
ublas::matrix< MatrixDoublediffN_face_edge
 
ublas::vector< MatrixDoublediffN_face_bubble
 

Additional Inherited Members

- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Detailed Description

Calculate base functions on triangle.

Definition at line 33 of file TriPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ TriPolynomialBase()

TriPolynomialBase::TriPolynomialBase ( )

Definition at line 51 of file TriPolynomialBase.cpp.

51 {}

◆ ~TriPolynomialBase()

TriPolynomialBase::~TriPolynomialBase ( )

Definition at line 52 of file TriPolynomialBase.cpp.

52 {}

Member Function Documentation

◆ getValue()

MoFEMErrorCode TriPolynomialBase::getValue ( MatrixDouble pts,
boost::shared_ptr< BaseFunctionCtx ctx_ptr 
)
virtual

Reimplemented from MoFEM::BaseFunction.

Definition at line 484 of file TriPolynomialBase.cpp.

485  {
487 
489  CHKERR ctx_ptr->query_interface(IDD_TRI_BASE_FUNCTION, &iface);
490  cTx = reinterpret_cast<EntPolynomialBaseCtx *>(iface);
491 
492  int nb_gauss_pts = pts.size2();
493  if (!nb_gauss_pts) {
495  }
496 
497  if (pts.size1() < 2) {
498  SETERRQ(
499  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
500  "Wrong dimension of pts, should be at least 3 rows with coordinates");
501  }
502 
503  const FieldApproximationBase base = cTx->bAse;
505  if (cTx->copyNodeBase == LASTBASE) {
506  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 3, false);
508  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
509  &pts(0, 0), &pts(1, 0), nb_gauss_pts);
510  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(3, 2, false);
512  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
513  } else {
514  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
515  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
516  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(base) =
517  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(cTx->copyNodeBase);
518  }
519  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
520  (unsigned int)nb_gauss_pts) {
521  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
522  "Base functions or nodes has wrong number of integration points "
523  "for base %s",
524  ApproximationBaseNames[base]);
525  }
526 
527  if(1) {
528  // In linear geometry derivatives are constant,
529  // this in expense of efficiency makes implementation
530  // consistent between vertices and other types of entities
531  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(3, 2, false);
533  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
534  MatrixDouble diffN(nb_gauss_pts, 6);
535  for (int gg = 0; gg != nb_gauss_pts; ++gg) {
536  for (int nn = 0; nn != 3; ++nn) {
537  for (int dd = 0; dd != 2; ++dd) {
538  diffN(gg, nn * 2 + dd) =
539  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(nn, dd);
540  }
541  }
542  }
543  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(
544  diffN.size1(), diffN.size2(), false);
545  data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().swap(diffN.data());
546  }
547 
548  switch (cTx->sPace) {
549  case H1:
550  CHKERR getValueH1(pts);
551  break;
552  case HDIV:
553  CHKERR getValueHdiv(pts);
554  break;
555  case HCURL:
556  CHKERR getValueHcurl(pts);
557  break;
558  case L2:
559  CHKERR getValueL2(pts);
560  break;
561  default:
562  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
563  }
564 
566 }
field with continuous normal traction
Definition: definitions.h:170
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
data structure for finite element entityIt keeps that about indices of degrees of freedom...
Class used to pass element data to calculate base functions on tet,triangle,edge. ...
const FieldApproximationBase bAse
base class for all interface classes
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
MoFEMErrorCode getValueH1(MatrixDouble &pts)
const FieldApproximationBase copyNodeBase
static const char *const ApproximationBaseNames[]
Definition: definitions.h:155
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getValueL2(MatrixDouble &pts)
PetscErrorCode ShapeDiffMBTRI(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:186
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Common.hpp:212
field with continuous tangents
Definition: definitions.h:169
#define CHKERR
Inline error check.
Definition: definitions.h:578
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition: fem_tools.c:175
static const MOFEMuuid IDD_TRI_BASE_FUNCTION
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
continuous field
Definition: definitions.h:168
field with C-1 continuity
Definition: definitions.h:171

◆ getValueH1()

MoFEMErrorCode TriPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 69 of file TriPolynomialBase.cpp.

69  {
71 
73  const FieldApproximationBase base = cTx->bAse;
74  if (cTx->basePolynomialsType0 == NULL)
75  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
76  "Polynomial type not set");
77  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
78  double *diffL, const int dim) =
80 
81  int nb_gauss_pts = pts.size2();
82 
83  if (data.spacesOnEntities[MBEDGE].test(H1)) {
84  // edges
85  if (data.dataOnEntities[MBEDGE].size() != 3) {
86  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
87  }
88  int sense[3], order[3];
89  double *H1edgeN[3], *diffH1edgeN[3];
90  for (int ee = 0; ee < 3; ee++) {
91  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
92  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
93  "data inconsistency");
94  }
95  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
96  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
97  int nb_dofs = NBEDGE_H1(data.dataOnEntities[MBEDGE][ee].getDataOrder());
98  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
99  false);
100  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,
101  2 * nb_dofs, false);
102  H1edgeN[ee] = &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
103  diffH1edgeN[ee] =
104  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
105  }
107  sense, order,
108  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
109  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
110  H1edgeN, diffH1edgeN, nb_gauss_pts, base_polynomials);
111  }
112 
113  if (data.spacesOnEntities[MBTRI].test(H1)) {
114  // face
115  if (data.dataOnEntities[MBTRI].size() != 1) {
116  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
117  }
118  int nb_dofs = NBFACETRI_H1(data.dataOnEntities[MBTRI][0].getDataOrder());
119  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, nb_dofs,
120  false);
121  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
122  2 * nb_dofs, false);
123  const int face_nodes[] = {0, 1, 2};
125  face_nodes, data.dataOnEntities[MBTRI][0].getDataOrder(),
126  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
127  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
128  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
129  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
130  nb_gauss_pts, base_polynomials);
131  }
132 
134 }
std::bitset< LASTSPACE > spacesOnEntities[MBMAXTYPE]
spaces on entity types
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
data structure for finite element entityIt keeps that about indices of degrees of freedom...
PetscErrorCode H1_EdgeShapeFunctions_MBTRI(int *sense, int *p, double *N, double *diffN, double *edgeN[3], double *diff_edgeN[3], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
H1_EdgeShapeFunctions_MBTRI.
Definition: h1.c:27
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
PetscErrorCode H1_FaceShapeFunctions_MBTRI(const int *face_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:162
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
continuous field
Definition: definitions.h:168

◆ getValueHcurl()

MoFEMErrorCode TriPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 465 of file TriPolynomialBase.cpp.

465  {
467 
468  switch (cTx->bAse) {
472  break;
475  break;
476  default:
477  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
478  }
479 
481 }
MoFEMErrorCode getValueHcurlAinsworthBase(MatrixDouble &pts)
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:142
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

◆ getValueHcurlAinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlAinsworthBase ( MatrixDouble pts)
private

Definition at line 286 of file TriPolynomialBase.cpp.

286  {
288 
290  const FieldApproximationBase base = cTx->bAse;
291  if (data.dataOnEntities[MBTRI].size() != 1) {
292  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
293  }
294  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
295  double *diffL, const int dim) =
297 
298  int nb_gauss_pts = pts.size2();
299 
300  // Calculation H-curl on triangle faces
301  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
302  if (data.dataOnEntities[MBEDGE].size() != 3) {
303  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
304  }
305  int sense[3], order[3];
306  double *hcurl_edge_n[3];
307  double *diff_hcurl_edge_n[3];
308  for (int ee = 0; ee < 3; ee++) {
309  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
310  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
311  "data inconsistency");
312  }
313  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
314  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
315  int nb_dofs = NBEDGE_AINSWORTH_HCURL(
316  data.dataOnEntities[MBEDGE][ee].getDataOrder());
317  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
318  3 * nb_dofs, false);
319  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
320  nb_gauss_pts, 2 * 3 * nb_dofs, false);
321  hcurl_edge_n[ee] =
322  &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
323  diff_hcurl_edge_n[ee] =
324  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
325  }
327  sense, order,
328  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
329  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
330  hcurl_edge_n, diff_hcurl_edge_n, nb_gauss_pts, base_polynomials);
331  } else {
332  for (int ee = 0; ee < 3; ee++) {
333  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, 0, false);
334  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts, 0,
335  false);
336  }
337  }
338 
339  if (data.spacesOnEntities[MBTRI].test(HCURL)) {
340 
341  // cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << endl;
342  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
343  //
344  // face
345  if (data.dataOnEntities[MBTRI].size() != 1) {
346  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
347  }
348  int order = data.dataOnEntities[MBTRI][0].getDataOrder();
349  int nb_dofs = NBFACETRI_AINSWORTH_HCURL(order);
350  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
351  false);
352  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
353  3 * 2 * nb_dofs, false);
354  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
355  int face_nodes[] = {0, 1, 2};
357  face_nodes, order,
358  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
359  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
360  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
361  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
362  nb_gauss_pts, base_polynomials);
363  // cerr << data.dataOnEntities[MBTRI][0].getN(base) << endl;
364  } else {
365  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
366  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
367  }
368 
370 }
std::bitset< LASTSPACE > spacesOnEntities[MBMAXTYPE]
spaces on entity types
MoFEMErrorCode Hcurl_Ainsworth_FaceFunctions_MBTET_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int nb_integration_pts, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face H-curl functions.
Definition: Hcurl.cpp:1258
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
data structure for finite element entityIt keeps that about indices of degrees of freedom...
MoFEMErrorCode Hcurl_Ainsworth_EdgeBaseFunctions_MBTET_ON_FACE(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 face.
Definition: Hcurl.cpp:240
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
#define NBEDGE_AINSWORTH_HCURL(P)
field with continuous tangents
Definition: definitions.h:169
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define NBFACETRI_AINSWORTH_HCURL(P)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 373 of file TriPolynomialBase.cpp.

373  {
375 
377  const FieldApproximationBase base = cTx->bAse;
378 
379  int nb_gauss_pts = pts.size2();
380 
381  // Calculation H-curl on triangle faces
382  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
383  if (data.dataOnEntities[MBEDGE].size() != 3) {
384  SETERRQ1(
385  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
386  "wrong number of data structures on edges, should be three but is %d",
387  data.dataOnEntities[MBEDGE].size());
388  }
389  int sense[3], order[3];
390  double *hcurl_edge_n[3];
391  double *diff_hcurl_edge_n[3];
392  for (int ee = 0; ee != 3; ++ee) {
393  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
394  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
395  "orientation (sense) of edge is not set");
396  }
397  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
398  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
399  int nb_dofs = NBEDGE_DEMKOWICZ_HCURL(
400  data.dataOnEntities[MBEDGE][ee].getDataOrder());
401  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
402  3 * nb_dofs, false);
403  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
404  nb_gauss_pts, 2 * 3 * nb_dofs, false);
405  hcurl_edge_n[ee] =
406  &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
407  diff_hcurl_edge_n[ee] =
408  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
409  }
411  sense, order,
412  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
413  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
414  hcurl_edge_n, diff_hcurl_edge_n, nb_gauss_pts);
415  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
416  // cerr << data.dataOnEntities[MBEDGE][0].getDiffN(base) << endl;
417  // cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << endl;
418  // cerr << data.dataOnEntities[MBEDGE][0].getN(base) << endl;
419  } else {
420  // No DOFs on faces, resize base function matrices, indicating that no
421  // dofs on them.
422  for (int ee = 0; ee != 3; ++ee) {
423  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, 0, false);
424  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts, 0,
425  false);
426  }
427  }
428 
429  if (data.spacesOnEntities[MBTRI].test(HCURL)) {
430 
431  // cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << endl;
432  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
433  //
434  // face
435  if (data.dataOnEntities[MBTRI].size() != 1) {
436  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
437  "No data struture to keep base functions on face");
438  }
439  int order = data.dataOnEntities[MBTRI][0].getDataOrder();
440  int nb_dofs = NBFACETRI_DEMKOWICZ_HCURL(order);
441  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
442  false);
443  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
444  3 * 2 * nb_dofs, false);
445  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
446  int face_nodes[] = {0, 1, 2};
448  face_nodes, order,
449  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
450  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
451  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
452  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
453  nb_gauss_pts);
454  // cerr << data.dataOnEntities[MBTRI][0].getN(base) << endl;
455  } else {
456  // No DOFs on faces, resize base function matrices, indicating that no
457  // dofs on them.
458  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
459  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
460  }
461 
463 }
std::bitset< LASTSPACE > spacesOnEntities[MBMAXTYPE]
spaces on entity types
data structure for finite element entityIt keeps that about indices of degrees of freedom...
MoFEMErrorCode Hcurl_Demkowicz_FaceBaseFunctions_MBTRI(int *faces_nodes, int p, double *n, double *diff_n, double *phi, double *diff_phi, int nb_integration_pts)
Face base interior function.
Definition: Hcurl.cpp:2771
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:169
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define NBFACETRI_DEMKOWICZ_HCURL(P)
#define NBEDGE_DEMKOWICZ_HCURL(P)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
MoFEMErrorCode Hcurl_Demkowicz_EdgeBaseFunctions_MBTRI(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 teriangle.
Definition: Hcurl.cpp:2243

◆ getValueHdiv()

MoFEMErrorCode TriPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 269 of file TriPolynomialBase.cpp.

269  {
271 
272  switch (cTx->bAse) {
275  return getValueHdivAinsworthBase(pts);
277  return getValueHdivDemkowiczBase(pts);
278  default:
279  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
280  }
281 
283 }
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
EntPolynomialBaseCtx * cTx
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:142
MoFEMErrorCode getValueHdivDemkowiczBase(MatrixDouble &pts)
MoFEMErrorCode getValueHdivAinsworthBase(MatrixDouble &pts)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

◆ getValueHdivAinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueHdivAinsworthBase ( MatrixDouble pts)
private

Definition at line 168 of file TriPolynomialBase.cpp.

168  {
170 
172  const FieldApproximationBase base = cTx->bAse;
173  if (cTx->basePolynomialsType0 == NULL)
174  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
175  "Polynomial type not set");
176  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
177  double *diffL, const int dim) =
179 
180  int nb_gauss_pts = pts.size2();
181 
182  double *PHI_f_e[3];
183  double *PHI_f;
184 
185  N_face_edge.resize(1, 3, false);
186  N_face_bubble.resize(1, false);
187  int face_order = data.dataOnEntities[MBTRI][0].getDataOrder();
188  // three edges on face
189  for (int ee = 0; ee < 3; ee++) {
190  N_face_edge(0, ee).resize(
191  nb_gauss_pts, 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(face_order), false);
192  PHI_f_e[ee] = &((N_face_edge(0, ee))(0, 0));
193  }
194  N_face_bubble[0].resize(nb_gauss_pts,
195  3 * NBFACETRI_AINSWORTH_FACE_HDIV(face_order), false);
196  PHI_f = &*(N_face_bubble[0].data().begin());
197 
198  int face_nodes[3] = {0, 1, 2};
200  face_nodes, face_order,
201  &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0), NULL, PHI_f_e, NULL,
202  nb_gauss_pts, 3, base_polynomials);
204  face_nodes, face_order,
205  &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0), NULL, PHI_f, NULL,
206  nb_gauss_pts, 3, base_polynomials);
207 
208  // set shape functions into data structure
209  if (data.dataOnEntities[MBTRI].size() != 1) {
210  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
211  }
212  data.dataOnEntities[MBTRI][0].getN(base).resize(
213  nb_gauss_pts, 3 * NBFACETRI_AINSWORTH_HDIV(face_order), false);
214  int col = 0;
215  for (int oo = 0; oo < face_order; oo++) {
216  for (int ee = 0; ee < 3; ee++) {
217  for (int dd = 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(oo);
218  dd < 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(oo + 1); dd++, col++) {
219  for (int gg = 0; gg < nb_gauss_pts; gg++) {
220  data.dataOnEntities[MBTRI][0].getN(base)(gg, col) =
221  N_face_edge(0, ee)(gg, dd);
222  }
223  }
224  }
225  for (int dd = 3 * NBFACETRI_AINSWORTH_FACE_HDIV(oo);
226  dd < 3 * NBFACETRI_AINSWORTH_FACE_HDIV(oo + 1); dd++, col++) {
227  for (int gg = 0; gg < nb_gauss_pts; gg++) {
228  data.dataOnEntities[MBTRI][0].getN(base)(gg, col) =
229  N_face_bubble[0](gg, dd);
230  }
231  }
232  }
233 
235 }
#define NBFACETRI_AINSWORTH_EDGE_HDIV(P)
ublas::vector< MatrixDouble > N_face_bubble
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T *> &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
data structure for finite element entityIt keeps that about indices of degrees of freedom...
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
ublas::matrix< MatrixDouble > N_face_edge
#define CHKERR
Inline error check.
Definition: definitions.h:578
MoFEMErrorCode Hdiv_Ainsworth_EdgeFaceShapeFunctions_MBTET_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f_e[3], double *diff_phi_f_e[3], int gdim, int nb, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Hdiv base functions, Edge-based face functions by Ainsworth .
Definition: Hdiv.cpp:43
MoFEMErrorCode Hdiv_Ainsworth_FaceBubbleShapeFunctions_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int gdim, int nb, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Face bubble functions by Ainsworth .
Definition: Hdiv.cpp:164
#define NBFACETRI_AINSWORTH_FACE_HDIV(P)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
#define NBFACETRI_AINSWORTH_HDIV(P)

◆ getValueHdivDemkowiczBase()

MoFEMErrorCode TriPolynomialBase::getValueHdivDemkowiczBase ( MatrixDouble pts)
private

Definition at line 237 of file TriPolynomialBase.cpp.

237  {
239 
241  // set shape functions into data structure
242  if (data.dataOnEntities[MBTRI].size() != 1) {
243  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
244  }
245  const FieldApproximationBase base = cTx->bAse;
246  if (base != DEMKOWICZ_JACOBI_BASE) {
247  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
248  "This should be used only with DEMKOWICZ_JACOBI_BASE "
249  "but base is %s",
250  ApproximationBaseNames[base]);
251  }
252  int order = data.dataOnEntities[MBTRI][0].getDataOrder();
253  int nb_gauss_pts = pts.size2();
254  data.dataOnEntities[MBTRI][0].getN(base).resize(
255  nb_gauss_pts, 3 * NBFACETRI_DEMKOWICZ_HDIV(order), false);
256  double *phi_f = &*data.dataOnEntities[MBTRI][0].getN(base).data().begin();
257  if (NBFACETRI_DEMKOWICZ_HDIV(order) == 0)
259  int face_nodes[3] = {0, 1, 2};
261  face_nodes, order,
262  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
263  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(), phi_f,
264  NULL, nb_gauss_pts, 3);
265 
267 }
data structure for finite element entityIt keeps that about indices of degrees of freedom...
MoFEMErrorCode Hdiv_Demkowicz_Face_MBTET_ON_FACE(int *faces_nodes, int p, double *N, double *diffN, double *phi_f, double *diff_phi_f, int gdim, int nb)
Definition: Hdiv.cpp:621
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
static const char *const ApproximationBaseNames[]
Definition: definitions.h:155
EntPolynomialBaseCtx * cTx
#define NBFACETRI_DEMKOWICZ_HDIV(P)
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

◆ getValueL2()

MoFEMErrorCode TriPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 136 of file TriPolynomialBase.cpp.

136  {
138 
140  const FieldApproximationBase base = cTx->bAse;
141  if (cTx->basePolynomialsType0 == NULL)
142  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
143  "Polynomial type not set");
144  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
145  double *diffL, const int dim) =
147 
148  int nb_gauss_pts = pts.size2();
149 
150  data.dataOnEntities[MBTRI][0].getN(base).resize(
151  nb_gauss_pts, NBFACETRI_L2(data.dataOnEntities[MBTRI][0].getDataOrder()),
152  false);
153  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(
154  nb_gauss_pts,
155  2 * NBFACETRI_L2(data.dataOnEntities[MBTRI][0].getDataOrder()), false);
156 
158  data.dataOnEntities[MBTRI][0].getDataOrder(),
159  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
160  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
161  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
162  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
163  nb_gauss_pts, base_polynomials);
164 
166 }
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
data structure for finite element entityIt keeps that about indices of degrees of freedom...
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
#define NBFACETRI_L2(P)
Number of base functions on triangle for L2 space.
#define CHKERR
Inline error check.
Definition: definitions.h:578
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
PetscErrorCode L2_Ainsworth_ShapeFunctions_MBTRI(int p, double *N, double *diffN, double *L2N, double *diff_L2N, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Get base functions on triangle for L2 space.
Definition: l2.c:29

◆ query_interface()

MoFEMErrorCode TriPolynomialBase::query_interface ( const MOFEMuuid uuid,
MoFEM::UnknownInterface **  iface 
) const
virtual

Reimplemented from MoFEM::BaseFunction.

Definition at line 55 of file TriPolynomialBase.cpp.

56  {
58  *iface = NULL;
59  if (uuid == IDD_TET_BASE_FUNCTION) {
60  *iface = const_cast<TriPolynomialBase *>(this);
62  } else {
63  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
64  }
67 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:459
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
static const MOFEMuuid IDD_TET_BASE_FUNCTION
#define CHKERR
Inline error check.
Definition: definitions.h:578
Calculate base functions on triangle.
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, MoFEM::UnknownInterface **iface) const
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403

Member Data Documentation

◆ cTx

EntPolynomialBaseCtx* MoFEM::TriPolynomialBase::cTx
private

Definition at line 45 of file TriPolynomialBase.hpp.

◆ diffN_face_bubble

ublas::vector<MatrixDouble> MoFEM::TriPolynomialBase::diffN_face_bubble
private

Definition at line 54 of file TriPolynomialBase.hpp.

◆ diffN_face_edge

ublas::matrix<MatrixDouble> MoFEM::TriPolynomialBase::diffN_face_edge
private

Definition at line 53 of file TriPolynomialBase.hpp.

◆ N_face_bubble

ublas::vector<MatrixDouble> MoFEM::TriPolynomialBase::N_face_bubble
private

Definition at line 52 of file TriPolynomialBase.hpp.

◆ N_face_edge

ublas::matrix<MatrixDouble> MoFEM::TriPolynomialBase::N_face_edge
private

Definition at line 51 of file TriPolynomialBase.hpp.


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