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

Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (that not perfect) More...

#include <src/approximation/FlatPrismPolynomialBase.hpp>

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

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 FlatPrismPolynomialBase ()
 
 ~FlatPrismPolynomialBase ()
 
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)
 

Private Attributes

FlatPrismPolynomialBaseCtxcTx
 
int numNodes
 
const EntityHandleconnPrism
 
const EntityHandleconnFace3
 
const EntityHandleconnFace4
 
int faceNodes [2][3]
 
MatrixDouble N
 
MatrixDouble diffN
 

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 tetrahedral

FIXME: Need moab and mofem finite element structure to work (that not perfect)

Definition at line 60 of file FlatPrismPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ FlatPrismPolynomialBase()

FlatPrismPolynomialBase::FlatPrismPolynomialBase ( )

Definition at line 90 of file FlatPrismPolynomialBase.cpp.

90 {}

◆ ~FlatPrismPolynomialBase()

FlatPrismPolynomialBase::~FlatPrismPolynomialBase ( )

Definition at line 89 of file FlatPrismPolynomialBase.cpp.

89 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 93 of file FlatPrismPolynomialBase.cpp.

94  {
95 
97 
99  CHKERR ctx_ptr->query_interface(IDD_FLATPRISM_BASE_FUNCTION, &iface);
100  cTx = reinterpret_cast<FlatPrismPolynomialBaseCtx *>(iface);
101  if (!cTx->fePtr) {
102  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
103  "Pointer to element should be given "
104  "when EntPolynomialBaseCtx is constructed "
105  "(use different constructor)");
106  }
107 
108  int nb_gauss_pts = pts.size2();
109  if (!nb_gauss_pts) {
111  }
112 
113  if (pts.size1() < 1) {
114  SETERRQ(
115  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
116  "Wrong dimension of pts, should be at least 3 rows with coordinates");
117  }
118 
119  const FieldApproximationBase base = cTx->bAse;
121 
122  if (cTx->copyNodeBase == LASTBASE) {
123  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
124  } else {
125  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
126  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
127  }
128  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 6, false);
129  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 12,
130  false);
131  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
132  (unsigned int)nb_gauss_pts) {
133  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
134  "Base functions or nodes has wrong number of integration points "
135  "for base %s",
136  ApproximationBaseNames[base]);
137  }
138 
139  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 6, false);
140  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 12,
141  false);
142  N.resize(nb_gauss_pts, 3, false);
143  diffN.resize(3, 2, false);
144  CHKERR ShapeMBTRI(&*N.data().begin(), &pts(0, 0), &pts(1, 0), nb_gauss_pts);
145  CHKERR ShapeDiffMBTRI(&*diffN.data().begin());
146 
147  // This is needed to have proper order of nodes on faces
148 
149  CHKERR cTx->mOab.get_connectivity(cTx->fePtr->getEnt(), connPrism, numNodes,
150  true);
151  SideNumber_multiIndex &side_table =
153  SideNumber_multiIndex::nth_index<1>::type::iterator siit3 =
154  side_table.get<1>().find(boost::make_tuple(MBTRI, 3));
155  SideNumber_multiIndex::nth_index<1>::type::iterator siit4 =
156  side_table.get<1>().find(boost::make_tuple(MBTRI, 4));
157  if (siit3 == side_table.get<1>().end())
158  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
159  if (siit4 == side_table.get<1>().end())
160  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
161  CHKERR cTx->mOab.get_connectivity(siit3->get()->ent, connFace3, numNodes,
162  true);
163  CHKERR cTx->mOab.get_connectivity(siit4->get()->ent, connFace4, numNodes,
164  true);
165 
166  for (int nn = 0; nn < 3; nn++) {
167  faceNodes[0][nn] = std::distance(
168  connPrism, std::find(connPrism, connPrism + 3, connFace3[nn]));
169  faceNodes[1][nn] = std::distance(
170  connPrism + 3, std::find(connPrism + 3, connPrism + 6, connFace4[nn]));
171  for (int gg = 0; gg < nb_gauss_pts; gg++) {
172  double val = N(gg, nn);
173  double val_x = diffN(nn, 0);
174  double val_y = diffN(nn, 1);
175  data.dataOnEntities[MBVERTEX][0].getN(base)(gg, nn) = val;
176  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg, 2 * nn + 0) = val_x;
177  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg, 2 * nn + 1) = val_y;
178  data.dataOnEntities[MBVERTEX][0].getN(base)(gg, 3 + nn) = val;
179  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg, 6 + 2 * nn + 0) =
180  val_x;
181  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg, 6 + 2 * nn + 1) =
182  val_y;
183  }
184  }
185  // for(int nn = 0;nn<3;nn++) {
186  // if(faceNodes[0][nn]!=faceNodes[1][nn]) {
187  // SETERRQ(
188  // PETSC_COMM_SELF,
189  // MOFEM_DATA_INCONSISTENCY,
190  // "Node order different on both faces"
191  // );
192  // }
193  // }
194 
195  switch (cTx->sPace) {
196  case H1:
197  CHKERR getValueH1(pts);
198  break;
199  case HDIV:
200  CHKERR getValueHdiv(pts);
201  break;
202  case HCURL:
203  CHKERR getValueHcurl(pts);
204  break;
205  case L2:
206  CHKERR getValueL2(pts);
207  break;
208  default:
209  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
210  }
211 
213 }
field with continuous normal traction
Definition: definitions.h:170
data structure for finite element entityIt keeps that about indices of degrees of freedom...
static const MOFEMuuid IDD_FLATPRISM_BASE_FUNCTION
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
const FieldApproximationBase bAse
base class for all interface classes
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< hashed_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, char, &SideNumber::side_number > > >, ordered_non_unique< const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
#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]
const FieldApproximationBase copyNodeBase
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
static const char *const ApproximationBaseNames[]
Definition: definitions.h:155
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
PetscErrorCode ShapeDiffMBTRI(double *diffN)
calculate derivatives of shape functions
Definition: fem_tools.c:186
FieldApproximationBase
approximation base
Definition: definitions.h:140
Class used to pass element data to calculate base functions on flat prism.
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:169
const NumeredEntFiniteElement * fePtr
SideNumber_multiIndex & getSideNumberTable() const
#define CHKERR
Inline error check.
Definition: definitions.h:578
FlatPrismPolynomialBaseCtx * cTx
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition: fem_tools.c:175
continuous field
Definition: definitions.h:168
EntityHandle getEnt() const
Get element entity.
MoFEMErrorCode getValueL2(MatrixDouble &pts)
MoFEMErrorCode getValueH1(MatrixDouble &pts)
field with C-1 continuity
Definition: definitions.h:171

◆ getValueH1()

MoFEMErrorCode FlatPrismPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 215 of file FlatPrismPolynomialBase.cpp.

215  {
216 
218 
220  const FieldApproximationBase base = cTx->bAse;
221  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
222  double *diffL, const int dim) =
224 
225  int nb_gauss_pts = pts.size2();
226 
227  // edges
228  if (data.dataOnEntities[MBEDGE].size() != 9) {
229  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
230  }
231  int valid_edges[] = {1, 1, 1, 0, 0, 0, 1, 1, 1};
232  int sense[9], order[9];
233  double *H1edgeN[9], *diffH1edgeN[9];
234  if ((data.spacesOnEntities[MBEDGE]).test(H1)) {
235  for (int ee = 0; ee < 9; ee++) {
236  if (!valid_edges[ee])
237  continue;
238  if (data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
239  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
240  "data inconsistency");
241  }
242  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
243  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
244  int nb_dofs = NBEDGE_H1(data.dataOnEntities[MBEDGE][ee].getDataOrder());
245  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
246  false);
247  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,
248  2 * nb_dofs, false);
249  H1edgeN[ee] = &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
250  diffH1edgeN[ee] =
251  &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
252  }
253  // shape functions on face 3
255  &sense[0], &order[0], &*N.data().begin(), &*diffN.data().begin(),
256  &H1edgeN[0], &diffH1edgeN[0], nb_gauss_pts, base_polynomials);
257  // shape functions on face 4
259  &sense[6], &order[6], &*N.data().begin(), &*diffN.data().begin(),
260  &H1edgeN[6], &diffH1edgeN[6], nb_gauss_pts, base_polynomials);
261  }
262 
263  // face
264  if (data.dataOnEntities[MBTRI].size() != 5) {
265  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
266  }
267  if ((data.spacesOnEntities[MBTRI]).test(H1)) {
268  for (int ff = 3; ff <= 4; ff++) {
269  int nb_dofs = NBFACETRI_H1(data.dataOnEntities[MBTRI][ff].getDataOrder());
270  data.dataOnEntities[MBTRI][ff].getN(base).resize(nb_gauss_pts, nb_dofs,
271  false);
272  data.dataOnEntities[MBTRI][ff].getDiffN(base).resize(nb_gauss_pts,
273  2 * nb_dofs, false);
275  faceNodes[ff - 3], data.dataOnEntities[MBTRI][ff].getDataOrder(),
276  &*N.data().begin(), &*diffN.data().begin(),
277  &*data.dataOnEntities[MBTRI][ff].getN(base).data().begin(),
278  &*data.dataOnEntities[MBTRI][ff].getDiffN(base).data().begin(),
279  nb_gauss_pts, base_polynomials);
280  }
281  }
282 
284 }
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
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
const FieldApproximationBase bAse
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
#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]
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
FlatPrismPolynomialBaseCtx * cTx
continuous field
Definition: definitions.h:168

◆ getValueHcurl()

MoFEMErrorCode FlatPrismPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 296 of file FlatPrismPolynomialBase.cpp.

296  {
297  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
298 }

◆ getValueHdiv()

MoFEMErrorCode FlatPrismPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 292 of file FlatPrismPolynomialBase.cpp.

292  {
293  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
294 }

◆ getValueL2()

MoFEMErrorCode FlatPrismPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 286 of file FlatPrismPolynomialBase.cpp.

286  {
288  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
290 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 74 of file FlatPrismPolynomialBase.cpp.

75  {
76 
78  *iface = NULL;
79  if (uuid == IDD_FLATPRISM_BASE_FUNCTION) {
80  *iface = const_cast<FlatPrismPolynomialBase *>(this);
82  } else {
83  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
84  }
87 }
static const MOFEMuuid IDD_FLATPRISM_BASE_FUNCTION
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (t...
#define CHKERR
Inline error check.
Definition: definitions.h:578
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, MoFEM::UnknownInterface **iface) const

Member Data Documentation

◆ connFace3

const EntityHandle* MoFEM::FlatPrismPolynomialBase::connFace3
private

Definition at line 85 of file FlatPrismPolynomialBase.hpp.

◆ connFace4

const EntityHandle* MoFEM::FlatPrismPolynomialBase::connFace4
private

Definition at line 86 of file FlatPrismPolynomialBase.hpp.

◆ connPrism

const EntityHandle* MoFEM::FlatPrismPolynomialBase::connPrism
private

Definition at line 84 of file FlatPrismPolynomialBase.hpp.

◆ cTx

FlatPrismPolynomialBaseCtx* MoFEM::FlatPrismPolynomialBase::cTx
private

Definition at line 73 of file FlatPrismPolynomialBase.hpp.

◆ diffN

MatrixDouble MoFEM::FlatPrismPolynomialBase::diffN
private

Definition at line 89 of file FlatPrismPolynomialBase.hpp.

◆ faceNodes

int MoFEM::FlatPrismPolynomialBase::faceNodes[2][3]
private

Definition at line 87 of file FlatPrismPolynomialBase.hpp.

◆ N

MatrixDouble MoFEM::FlatPrismPolynomialBase::N
private

Definition at line 88 of file FlatPrismPolynomialBase.hpp.

◆ numNodes

int MoFEM::FlatPrismPolynomialBase::numNodes
private

Definition at line 83 of file FlatPrismPolynomialBase.hpp.


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