v0.8.20
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 23 of file TriPolynomialBase.cpp.

23 {}

◆ ~TriPolynomialBase()

TriPolynomialBase::~TriPolynomialBase ( )

Definition at line 24 of file TriPolynomialBase.cpp.

24 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 456 of file TriPolynomialBase.cpp.

457  {
459 
461  CHKERR ctx_ptr->query_interface(IDD_TRI_BASE_FUNCTION, &iface);
462  cTx = reinterpret_cast<EntPolynomialBaseCtx *>(iface);
463 
464  int nb_gauss_pts = pts.size2();
465  if (!nb_gauss_pts) {
467  }
468 
469  if (pts.size1() < 2) {
470  SETERRQ(
471  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
472  "Wrong dimension of pts, should be at least 3 rows with coordinates");
473  }
474 
475  const FieldApproximationBase base = cTx->bAse;
477  if (cTx->copyNodeBase == LASTBASE) {
478  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 3, false);
480  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
481  &pts(0, 0), &pts(1, 0), nb_gauss_pts);
482  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(3, 2, false);
484  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
485  } else {
486  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
487  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
488  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(base) =
489  data.dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(cTx->copyNodeBase);
490  }
491  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
492  (unsigned int)nb_gauss_pts) {
493  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
494  "Base functions or nodes has wrong number of integration points "
495  "for base %s",
496  ApproximationBaseNames[base]);
497  }
498 
499  if(1) {
500  // In linear geometry derivatives are constant,
501  // this in expense of efficiency makes implementation
502  // consistent between vertices and other types of entities
503  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(3, 2, false);
505  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin());
506  MatrixDouble diffN(nb_gauss_pts, 6);
507  for (int gg = 0; gg != nb_gauss_pts; ++gg) {
508  for (int nn = 0; nn != 3; ++nn) {
509  for (int dd = 0; dd != 2; ++dd) {
510  diffN(gg, nn * 2 + dd) =
511  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(nn, dd);
512  }
513  }
514  }
515  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(
516  diffN.size1(), diffN.size2(), false);
517  data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().swap(diffN.data());
518  }
519 
520  switch (cTx->sPace) {
521  case H1:
522  CHKERR getValueH1(pts);
523  break;
524  case HDIV:
525  CHKERR getValueHdiv(pts);
526  break;
527  case HCURL:
528  CHKERR getValueHcurl(pts);
529  break;
530  case L2:
531  CHKERR getValueL2(pts);
532  break;
533  default:
534  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
535  }
536 
538 }
field with continuous normal traction
Definition: definitions.h:171
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...
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:77
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:475
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
MoFEMErrorCode getValueH1(MatrixDouble &pts)
const FieldApproximationBase copyNodeBase
static const char *const ApproximationBaseNames[]
Definition: definitions.h:156
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:141
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:170
#define CHKERR
Inline error check.
Definition: definitions.h:594
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:405
continuous field
Definition: definitions.h:169
field with C-1 continuity
Definition: definitions.h:172

◆ getValueH1()

MoFEMErrorCode TriPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 41 of file TriPolynomialBase.cpp.

41  {
43 
45  const FieldApproximationBase base = cTx->bAse;
46  if (cTx->basePolynomialsType0 == NULL)
47  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
48  "Polynomial type not set");
49  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
50  double *diffL, const int dim) =
52 
53  int nb_gauss_pts = pts.size2();
54 
55  if (data.spacesOnEntities[MBEDGE].test(H1)) {
56  // edges
57  if (data.dataOnEntities[MBEDGE].size() != 3) {
58  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
59  }
60  int sense[3], order[3];
61  double *H1edgeN[3], *diffH1edgeN[3];
62  for (int ee = 0; ee < 3; ee++) {
63  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
64  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
65  "data inconsistency");
66  }
67  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
68  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
69  int nb_dofs = NBEDGE_H1(data.dataOnEntities[MBEDGE][ee].getDataOrder());
70  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
71  false);
72  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,
73  2 * nb_dofs, false);
74  H1edgeN[ee] = &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
75  diffH1edgeN[ee] =
76  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
77  }
79  sense, order,
80  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
81  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
82  H1edgeN, diffH1edgeN, nb_gauss_pts, base_polynomials);
83  }
84 
85  if (data.spacesOnEntities[MBTRI].test(H1)) {
86  // face
87  if (data.dataOnEntities[MBTRI].size() != 1) {
88  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
89  }
90  int nb_dofs = NBFACETRI_H1(data.dataOnEntities[MBTRI][0].getDataOrder());
91  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, nb_dofs,
92  false);
93  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
94  2 * nb_dofs, false);
95  const int face_nodes[] = {0, 1, 2};
97  face_nodes, data.dataOnEntities[MBTRI][0].getDataOrder(),
98  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
99  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
100  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
101  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
102  nb_gauss_pts, base_polynomials);
103  }
104 
106 }
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:475
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:141
DataForcesAndSourcesCore & dAta
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
continuous field
Definition: definitions.h:169

◆ getValueHcurl()

MoFEMErrorCode TriPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 437 of file TriPolynomialBase.cpp.

437  {
439 
440  switch (cTx->bAse) {
444  break;
447  break;
448  default:
449  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
450  }
451 
453 }
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:475
EntPolynomialBaseCtx * cTx
MoFEMErrorCode getValueHcurlDemkowiczBase(MatrixDouble &pts)
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:143
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ getValueHcurlAinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlAinsworthBase ( MatrixDouble pts)
private

Definition at line 258 of file TriPolynomialBase.cpp.

258  {
260 
262  const FieldApproximationBase base = cTx->bAse;
263  if (data.dataOnEntities[MBTRI].size() != 1) {
264  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
265  }
266  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
267  double *diffL, const int dim) =
269 
270  int nb_gauss_pts = pts.size2();
271 
272  // Calculation H-curl on triangle faces
273  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
274  if (data.dataOnEntities[MBEDGE].size() != 3) {
275  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
276  }
277  int sense[3], order[3];
278  double *hcurl_edge_n[3];
279  double *diff_hcurl_edge_n[3];
280  for (int ee = 0; ee < 3; ee++) {
281  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
282  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
283  "data inconsistency");
284  }
285  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
286  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
287  int nb_dofs = NBEDGE_AINSWORTH_HCURL(
288  data.dataOnEntities[MBEDGE][ee].getDataOrder());
289  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
290  3 * nb_dofs, false);
291  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
292  nb_gauss_pts, 2 * 3 * nb_dofs, false);
293  hcurl_edge_n[ee] =
294  &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
295  diff_hcurl_edge_n[ee] =
296  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
297  }
299  sense, order,
300  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
301  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
302  hcurl_edge_n, diff_hcurl_edge_n, nb_gauss_pts, base_polynomials);
303  } else {
304  for (int ee = 0; ee < 3; ee++) {
305  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, 0, false);
306  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts, 0,
307  false);
308  }
309  }
310 
311  if (data.spacesOnEntities[MBTRI].test(HCURL)) {
312 
313  // cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << endl;
314  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
315  //
316  // face
317  if (data.dataOnEntities[MBTRI].size() != 1) {
318  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
319  }
320  int order = data.dataOnEntities[MBTRI][0].getDataOrder();
321  int nb_dofs = NBFACETRI_AINSWORTH_HCURL(order);
322  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
323  false);
324  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
325  3 * 2 * nb_dofs, false);
326  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
327  int face_nodes[] = {0, 1, 2};
329  face_nodes, order,
330  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
331  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
332  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
333  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
334  nb_gauss_pts, base_polynomials);
335  // cerr << data.dataOnEntities[MBTRI][0].getN(base) << endl;
336  } else {
337  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
338  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
339  }
340 
342 }
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:1251
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:233
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:141
DataForcesAndSourcesCore & dAta
#define NBEDGE_AINSWORTH_HCURL(P)
field with continuous tangents
Definition: definitions.h:170
#define CHKERR
Inline error check.
Definition: definitions.h:594
#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:405

◆ getValueHcurlDemkowiczBase()

MoFEMErrorCode TriPolynomialBase::getValueHcurlDemkowiczBase ( MatrixDouble pts)
private

Definition at line 345 of file TriPolynomialBase.cpp.

345  {
347 
349  const FieldApproximationBase base = cTx->bAse;
350 
351  int nb_gauss_pts = pts.size2();
352 
353  // Calculation H-curl on triangle faces
354  if (data.spacesOnEntities[MBEDGE].test(HCURL)) {
355  if (data.dataOnEntities[MBEDGE].size() != 3) {
356  SETERRQ1(
357  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
358  "wrong number of data structures on edges, should be three but is %d",
359  data.dataOnEntities[MBEDGE].size());
360  }
361  int sense[3], order[3];
362  double *hcurl_edge_n[3];
363  double *diff_hcurl_edge_n[3];
364  for (int ee = 0; ee != 3; ++ee) {
365  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
366  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
367  "orientation (sense) of edge is not set");
368  }
369  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
370  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
371  int nb_dofs = NBEDGE_DEMKOWICZ_HCURL(
372  data.dataOnEntities[MBEDGE][ee].getDataOrder());
373  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,
374  3 * nb_dofs, false);
375  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
376  nb_gauss_pts, 2 * 3 * nb_dofs, false);
377  hcurl_edge_n[ee] =
378  &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
379  diff_hcurl_edge_n[ee] =
380  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
381  }
383  sense, order,
384  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
385  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
386  hcurl_edge_n, diff_hcurl_edge_n, nb_gauss_pts);
387  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
388  // cerr << data.dataOnEntities[MBEDGE][0].getDiffN(base) << endl;
389  // cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << endl;
390  // cerr << data.dataOnEntities[MBEDGE][0].getN(base) << endl;
391  } else {
392  // No DOFs on faces, resize base function matrices, indicating that no
393  // dofs on them.
394  for (int ee = 0; ee != 3; ++ee) {
395  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, 0, false);
396  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts, 0,
397  false);
398  }
399  }
400 
401  if (data.spacesOnEntities[MBTRI].test(HCURL)) {
402 
403  // cerr << data.dataOnEntities[MBVERTEX][0].getN(base) << endl;
404  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
405  //
406  // face
407  if (data.dataOnEntities[MBTRI].size() != 1) {
408  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
409  "No data struture to keep base functions on face");
410  }
411  int order = data.dataOnEntities[MBTRI][0].getDataOrder();
412  int nb_dofs = NBFACETRI_DEMKOWICZ_HCURL(order);
413  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 3 * nb_dofs,
414  false);
415  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts,
416  3 * 2 * nb_dofs, false);
417  // cerr << data.dataOnEntities[MBVERTEX][0].getDiffN(base) << endl;
418  int face_nodes[] = {0, 1, 2};
420  face_nodes, order,
421  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
422  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
423  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
424  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
425  nb_gauss_pts);
426  // cerr << data.dataOnEntities[MBTRI][0].getN(base) << endl;
427  } else {
428  // No DOFs on faces, resize base function matrices, indicating that no
429  // dofs on them.
430  data.dataOnEntities[MBTRI][0].getN(base).resize(nb_gauss_pts, 0, false);
431  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(nb_gauss_pts, 0, false);
432  }
433 
435 }
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:2765
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:141
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:170
#define CHKERR
Inline error check.
Definition: definitions.h:594
#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:405
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:2237

◆ getValueHdiv()

MoFEMErrorCode TriPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 241 of file TriPolynomialBase.cpp.

241  {
243 
244  switch (cTx->bAse) {
247  return getValueHdivAinsworthBase(pts);
249  return getValueHdivDemkowiczBase(pts);
250  default:
251  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
252  }
253 
255 }
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
EntPolynomialBaseCtx * cTx
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:143
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:405

◆ getValueHdivAinsworthBase()

MoFEMErrorCode TriPolynomialBase::getValueHdivAinsworthBase ( MatrixDouble pts)
private

Definition at line 140 of file TriPolynomialBase.cpp.

140  {
142 
144  const FieldApproximationBase base = cTx->bAse;
145  if (cTx->basePolynomialsType0 == NULL)
146  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
147  "Polynomial type not set");
148  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
149  double *diffL, const int dim) =
151 
152  int nb_gauss_pts = pts.size2();
153 
154  double *PHI_f_e[3];
155  double *PHI_f;
156 
157  N_face_edge.resize(1, 3, false);
158  N_face_bubble.resize(1, false);
159  int face_order = data.dataOnEntities[MBTRI][0].getDataOrder();
160  // three edges on face
161  for (int ee = 0; ee < 3; ee++) {
162  N_face_edge(0, ee).resize(
163  nb_gauss_pts, 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(face_order), false);
164  PHI_f_e[ee] = &((N_face_edge(0, ee))(0, 0));
165  }
166  N_face_bubble[0].resize(nb_gauss_pts,
167  3 * NBFACETRI_AINSWORTH_FACE_HDIV(face_order), false);
168  PHI_f = &*(N_face_bubble[0].data().begin());
169 
170  int face_nodes[3] = {0, 1, 2};
172  face_nodes, face_order,
173  &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0), NULL, PHI_f_e, NULL,
174  nb_gauss_pts, 3, base_polynomials);
176  face_nodes, face_order,
177  &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0), NULL, PHI_f, NULL,
178  nb_gauss_pts, 3, base_polynomials);
179 
180  // set shape functions into data structure
181  if (data.dataOnEntities[MBTRI].size() != 1) {
182  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
183  }
184  data.dataOnEntities[MBTRI][0].getN(base).resize(
185  nb_gauss_pts, 3 * NBFACETRI_AINSWORTH_HDIV(face_order), false);
186  int col = 0;
187  for (int oo = 0; oo < face_order; oo++) {
188  for (int ee = 0; ee < 3; ee++) {
189  for (int dd = 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(oo);
190  dd < 3 * NBFACETRI_AINSWORTH_EDGE_HDIV(oo + 1); dd++, col++) {
191  for (int gg = 0; gg < nb_gauss_pts; gg++) {
192  data.dataOnEntities[MBTRI][0].getN(base)(gg, col) =
193  N_face_edge(0, ee)(gg, dd);
194  }
195  }
196  }
197  for (int dd = 3 * NBFACETRI_AINSWORTH_FACE_HDIV(oo);
198  dd < 3 * NBFACETRI_AINSWORTH_FACE_HDIV(oo + 1); dd++, col++) {
199  for (int gg = 0; gg < nb_gauss_pts; gg++) {
200  data.dataOnEntities[MBTRI][0].getN(base)(gg, col) =
201  N_face_bubble[0](gg, dd);
202  }
203  }
204  }
205 
207 }
#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:475
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:141
DataForcesAndSourcesCore & dAta
ublas::matrix< MatrixDouble > N_face_edge
#define CHKERR
Inline error check.
Definition: definitions.h:594
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:37
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:160
#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:405
#define NBFACETRI_AINSWORTH_HDIV(P)

◆ getValueHdivDemkowiczBase()

MoFEMErrorCode TriPolynomialBase::getValueHdivDemkowiczBase ( MatrixDouble pts)
private

Definition at line 209 of file TriPolynomialBase.cpp.

209  {
211 
213  // set shape functions into data structure
214  if (data.dataOnEntities[MBTRI].size() != 1) {
215  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
216  }
217  const FieldApproximationBase base = cTx->bAse;
218  if (base != DEMKOWICZ_JACOBI_BASE) {
219  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
220  "This should be used only with DEMKOWICZ_JACOBI_BASE "
221  "but base is %s",
222  ApproximationBaseNames[base]);
223  }
224  int order = data.dataOnEntities[MBTRI][0].getDataOrder();
225  int nb_gauss_pts = pts.size2();
226  data.dataOnEntities[MBTRI][0].getN(base).resize(
227  nb_gauss_pts, 3 * NBFACETRI_DEMKOWICZ_HDIV(order), false);
228  double *phi_f = &*data.dataOnEntities[MBTRI][0].getN(base).data().begin();
229  if (NBFACETRI_DEMKOWICZ_HDIV(order) == 0)
231  int face_nodes[3] = {0, 1, 2};
233  face_nodes, order,
234  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
235  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(), phi_f,
236  NULL, nb_gauss_pts, 3);
237 
239 }
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:617
const FieldApproximationBase bAse
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
static const char *const ApproximationBaseNames[]
Definition: definitions.h:156
EntPolynomialBaseCtx * cTx
#define NBFACETRI_DEMKOWICZ_HDIV(P)
FieldApproximationBase
approximation base
Definition: definitions.h:141
DataForcesAndSourcesCore & dAta
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405

◆ getValueL2()

MoFEMErrorCode TriPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 108 of file TriPolynomialBase.cpp.

108  {
110 
112  const FieldApproximationBase base = cTx->bAse;
113  if (cTx->basePolynomialsType0 == NULL)
114  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
115  "Polynomial type not set");
116  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
117  double *diffL, const int dim) =
119 
120  int nb_gauss_pts = pts.size2();
121 
122  data.dataOnEntities[MBTRI][0].getN(base).resize(
123  nb_gauss_pts, NBFACETRI_L2(data.dataOnEntities[MBTRI][0].getDataOrder()),
124  false);
125  data.dataOnEntities[MBTRI][0].getDiffN(base).resize(
126  nb_gauss_pts,
127  2 * NBFACETRI_L2(data.dataOnEntities[MBTRI][0].getDataOrder()), false);
128 
130  data.dataOnEntities[MBTRI][0].getDataOrder(),
131  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin(),
132  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin(),
133  &*data.dataOnEntities[MBTRI][0].getN(base).data().begin(),
134  &*data.dataOnEntities[MBTRI][0].getDiffN(base).data().begin(),
135  nb_gauss_pts, base_polynomials);
136 
138 }
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:475
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
EntPolynomialBaseCtx * cTx
FieldApproximationBase
approximation base
Definition: definitions.h:141
DataForcesAndSourcesCore & dAta
#define NBFACETRI_L2(P)
Number of base functions on triangle for L2 space.
#define CHKERR
Inline error check.
Definition: definitions.h:594
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
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 27 of file TriPolynomialBase.cpp.

28  {
30  *iface = NULL;
31  if (uuid == IDD_TET_BASE_FUNCTION) {
32  *iface = const_cast<TriPolynomialBase *>(this);
34  } else {
35  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
36  }
39 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
static const MOFEMuuid IDD_TET_BASE_FUNCTION
#define CHKERR
Inline error check.
Definition: definitions.h:594
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:405

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: