v0.8.13
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 59 of file FlatPrismPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ FlatPrismPolynomialBase()

FlatPrismPolynomialBase::FlatPrismPolynomialBase ( )

Definition at line 97 of file FlatPrismPolynomialBase.cpp.

97 {}

◆ ~FlatPrismPolynomialBase()

FlatPrismPolynomialBase::~FlatPrismPolynomialBase ( )

Definition at line 96 of file FlatPrismPolynomialBase.cpp.

96 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 99 of file FlatPrismPolynomialBase.cpp.

102  {
103 
105 
107  ierr = ctx_ptr->query_interface(IDD_FLATPRISM_BASE_FUNCTION,&iface); CHKERRG(ierr);
108  cTx = reinterpret_cast<FlatPrismPolynomialBaseCtx*>(iface);
109  if(!cTx->fePtr) {
110  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,
111  "Pointer to element should be given "
112  "when EntPolynomialBaseCtx is constructed "
113  "(use different constructor)"
114  );
115  }
116 
117  int nb_gauss_pts = pts.size2();
118  if(!nb_gauss_pts) {
120  }
121 
122  if(pts.size1()<1) {
123  SETERRQ(
124  PETSC_COMM_SELF,
126  "Wrong dimension of pts, should be at least 3 rows with coordinates"
127  );
128  }
129 
130  const FieldApproximationBase base = cTx->bAse;
132 
133  if(cTx->copyNodeBase==LASTBASE) {
134  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"Not implemented");
135  } else {
136  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) = data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
137  }
138  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts,6,false);
139  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts,12,false);
140  if(data.dataOnEntities[MBVERTEX][0].getN(base).size1()!=(unsigned int)nb_gauss_pts) {
141  SETERRQ1(
142  PETSC_COMM_SELF,
144  "Base functions or nodes has wrong number of integration points for base %s",
146  );
147  }
148 
149  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts,6,false);
150  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts,12,false);
151  N.resize(nb_gauss_pts,3,false);
152  diffN.resize(3,2,false);
153  ierr = ShapeMBTRI(&*N.data().begin(),&pts(0,0),&pts(1,0),nb_gauss_pts); CHKERRG(ierr);
154  ierr = ShapeDiffMBTRI(&*diffN.data().begin()); CHKERRG(ierr);
155 
156  // This is needed to have proper order of nodes on faces
157 
158  rval = cTx->mOab.get_connectivity(cTx->fePtr->getEnt(),connPrism,numNodes,true); CHKERRQ_MOAB(rval);
160  SideNumber_multiIndex::nth_index<1>::type::iterator siit3 = side_table.get<1>().find(boost::make_tuple(MBTRI,3));
161  SideNumber_multiIndex::nth_index<1>::type::iterator siit4 = side_table.get<1>().find(boost::make_tuple(MBTRI,4));
162  if(siit3==side_table.get<1>().end()) SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"data inconsistency");
163  if(siit4==side_table.get<1>().end()) SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"data inconsistency");
164  rval = cTx->mOab.get_connectivity(siit3->get()->ent,connFace3,numNodes,true); CHKERRQ_MOAB(rval);
165  rval = cTx->mOab.get_connectivity(siit4->get()->ent,connFace4,numNodes,true); CHKERRQ_MOAB(rval);
166 
167  for(int nn = 0;nn<3;nn++) {
168  faceNodes[0][nn] = std::distance(connPrism,std::find(connPrism,connPrism+3,connFace3[nn]));
169  faceNodes[1][nn] = std::distance(connPrism+3,std::find(connPrism+3,connPrism+6,connFace4[nn]));
170  for(int gg = 0;gg<nb_gauss_pts;gg++) {
171  double val = N(gg,nn);
172  double val_x = diffN(nn,0);
173  double val_y = diffN(nn,1);
174  data.dataOnEntities[MBVERTEX][0].getN(base)(gg,nn) = val;
175  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg,2*nn+0) = val_x;
176  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg,2*nn+1) = val_y;
177  data.dataOnEntities[MBVERTEX][0].getN(base)(gg,3+nn) = val;
178  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg,6+2*nn+0) = val_x;
179  data.dataOnEntities[MBVERTEX][0].getDiffN(base)(gg,6+2*nn+1) = val_y;
180  }
181  }
182  // for(int nn = 0;nn<3;nn++) {
183  // if(faceNodes[0][nn]!=faceNodes[1][nn]) {
184  // SETERRQ(
185  // PETSC_COMM_SELF,
186  // MOFEM_DATA_INCONSISTENCY,
187  // "Node order different on both faces"
188  // );
189  // }
190  // }
191 
192  switch (cTx->sPace) {
193  case H1:
194  ierr = getValueH1(pts); CHKERRG(ierr);
195  break;
196  case HDIV:
197  ierr = getValueHdiv(pts); CHKERRG(ierr);
198  break;
199  case HCURL:
200  ierr = getValueHcurl(pts); CHKERRG(ierr);
201  break;
202  case L2:
203  ierr = getValueL2(pts); CHKERRG(ierr);
204  break;
205  default:
206  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"Not yet implemented");
207  }
208 
210 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:533
field with continuous normal traction
Definition: definitions.h:181
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:519
const FieldApproximationBase bAse
base class for all interface classes
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
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:526
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
const FieldApproximationBase copyNodeBase
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
static const char *const ApproximationBaseNames[]
Definition: definitions.h:155
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
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
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
Class used to pass element data to calculate base functions on flat prism.
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:180
const NumeredEntFiniteElement * fePtr
SideNumber_multiIndex & getSideNumberTable() const
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:179
EntityHandle getEnt() const
Get element entity.
MoFEMErrorCode getValueL2(MatrixDouble &pts)
MoFEMErrorCode getValueH1(MatrixDouble &pts)
field with C-1 continuity
Definition: definitions.h:182

◆ getValueH1()

MoFEMErrorCode FlatPrismPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 212 of file FlatPrismPolynomialBase.cpp.

212  {
213 
215 
217  const FieldApproximationBase base = cTx->bAse;
218  PetscErrorCode (*base_polynomials)(
219  int p,double s,double *diff_s,double *L,double *diffL,const int dim
221 
222  int nb_gauss_pts = pts.size2();
223 
224  //edges
225  if(data.dataOnEntities[MBEDGE].size()!=9) {
226  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"data inconsistency");
227  }
228  int valid_edges[] = { 1,1,1, 0,0,0, 1,1,1 };
229  int sense[9],order[9];
230  double *H1edgeN[9],*diffH1edgeN[9];
231  if((data.spacesOnEntities[MBEDGE]).test(H1)) {
232  for(int ee = 0;ee<9;ee++) {
233  if(!valid_edges[ee]) continue;
234  if(data.dataOnEntities[MBEDGE][ee].getSense() == 0) {
235  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"data inconsistency");
236  }
237  sense[ee] = data.dataOnEntities[MBEDGE][ee].getSense();
238  order[ee] = data.dataOnEntities[MBEDGE][ee].getDataOrder();
239  int nb_dofs = NBEDGE_H1(data.dataOnEntities[MBEDGE][ee].getDataOrder());
240  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,nb_dofs,false);
241  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,2*nb_dofs,false);
242  H1edgeN[ee] = &*data.dataOnEntities[MBEDGE][ee].getN(base).data().begin();
243  diffH1edgeN[ee] = &*data.dataOnEntities[MBEDGE][ee].getDiffN(base).data().begin();
244  }
245  //shape functions on face 3
247  &sense[0],
248  &order[0],
249  &*N.data().begin(),
250  &*diffN.data().begin(),
251  &H1edgeN[0],
252  &diffH1edgeN[0],
253  nb_gauss_pts,base_polynomials
254  ); CHKERRG(ierr);
255  //shape functions on face 4
257  &sense[6],
258  &order[6],
259  &*N.data().begin(),
260  &*diffN.data().begin(),
261  &H1edgeN[6],
262  &diffH1edgeN[6],
263  nb_gauss_pts,base_polynomials
264  ); CHKERRG(ierr);
265  }
266 
267  //face
268  if(data.dataOnEntities[MBTRI].size()!=5) {
269  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"data inconsistency");
270  }
271  if((data.spacesOnEntities[MBTRI]).test(H1)) {
272  for(int ff = 3;ff<=4;ff++) {
273  int nb_dofs = NBFACETRI_H1(data.dataOnEntities[MBTRI][ff].getDataOrder());
274  data.dataOnEntities[MBTRI][ff].getN(base).resize(nb_gauss_pts,nb_dofs,false);
275  data.dataOnEntities[MBTRI][ff].getDiffN(base).resize(nb_gauss_pts,2*nb_dofs,false);
277  faceNodes[ff-3],
278  data.dataOnEntities[MBTRI][ff].getDataOrder(),
279  &*N.data().begin(),
280  &*diffN.data().begin(),
281  &*data.dataOnEntities[MBTRI][ff].getN(base).data().begin(),
282  &*data.dataOnEntities[MBTRI][ff].getDiffN(base).data().begin(),
283  nb_gauss_pts,base_polynomials
284  ); CHKERRG(ierr);
285  }
286  }
287 
289 }
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:519
const FieldApproximationBase bAse
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
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:526
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
FieldApproximationBase
approximation base
Definition: definitions.h:140
DataForcesAndSourcesCore & dAta
#define NBFACETRI_H1(P)
Number of base function on triangle for H1 space.
FlatPrismPolynomialBaseCtx * cTx
continuous field
Definition: definitions.h:179

◆ getValueHcurl()

MoFEMErrorCode FlatPrismPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 301 of file FlatPrismPolynomialBase.cpp.

301  {
302  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"Not yet implemented");
303 }

◆ getValueHdiv()

MoFEMErrorCode FlatPrismPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 297 of file FlatPrismPolynomialBase.cpp.

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

◆ getValueL2()

MoFEMErrorCode FlatPrismPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 291 of file FlatPrismPolynomialBase.cpp.

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

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 80 of file FlatPrismPolynomialBase.cpp.

82  {
83 
85  *iface = NULL;
86  if(uuid == IDD_FLATPRISM_BASE_FUNCTION) {
87  *iface = const_cast<FlatPrismPolynomialBase*>(this);
89  } else {
90  SETERRQ(PETSC_COMM_WORLD,MOFEM_DATA_INCONSISTENCY,"wrong interference");
91  }
94 }
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:519
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:562
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:526
Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (t...
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, MoFEM::UnknownInterface **iface) const

Member Data Documentation

◆ connFace3

const EntityHandle* MoFEM::FlatPrismPolynomialBase::connFace3
private

Definition at line 84 of file FlatPrismPolynomialBase.hpp.

◆ connFace4

const EntityHandle* MoFEM::FlatPrismPolynomialBase::connFace4
private

Definition at line 85 of file FlatPrismPolynomialBase.hpp.

◆ connPrism

const EntityHandle* MoFEM::FlatPrismPolynomialBase::connPrism
private

Definition at line 83 of file FlatPrismPolynomialBase.hpp.

◆ cTx

FlatPrismPolynomialBaseCtx* MoFEM::FlatPrismPolynomialBase::cTx
private

Definition at line 72 of file FlatPrismPolynomialBase.hpp.

◆ diffN

MatrixDouble MoFEM::FlatPrismPolynomialBase::diffN
private

Definition at line 88 of file FlatPrismPolynomialBase.hpp.

◆ faceNodes

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

Definition at line 86 of file FlatPrismPolynomialBase.hpp.

◆ N

MatrixDouble MoFEM::FlatPrismPolynomialBase::N
private

Definition at line 87 of file FlatPrismPolynomialBase.hpp.

◆ numNodes

int MoFEM::FlatPrismPolynomialBase::numNodes
private

Definition at line 82 of file FlatPrismPolynomialBase.hpp.


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