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

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

#include <src/approximation/FatPrismPolynomialBase.hpp>

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

Public Member Functions

MoFEMErrorCode query_interface (const MOFEMuuid &uuid, MoFEM::UnknownInterface **iface) const
 
 FatPrismPolynomialBase ()
 
 ~FatPrismPolynomialBase ()
 
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 getValueH1TrianglesOnly ()
 
MoFEMErrorCode getValueH1ThroughThickness ()
 
MoFEMErrorCode getValueH1 (MatrixDouble &pts)
 
MoFEMErrorCode getValueL2 (MatrixDouble &pts)
 
MoFEMErrorCode getValueHdiv (MatrixDouble &pts)
 
MoFEMErrorCode getValueHcurl (MatrixDouble &pts)
 

Private Attributes

FatPrismPolynomialBaseCtxcTx
 
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 69 of file FatPrismPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ FatPrismPolynomialBase()

FatPrismPolynomialBase::FatPrismPolynomialBase ( )

Definition at line 105 of file FatPrismPolynomialBase.cpp.

105 {}

◆ ~FatPrismPolynomialBase()

FatPrismPolynomialBase::~FatPrismPolynomialBase ( )

Definition at line 104 of file FatPrismPolynomialBase.cpp.

104 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 107 of file FatPrismPolynomialBase.cpp.

110  {
111 
113 
115  ierr = ctx_ptr->query_interface(IDD_FATPRISM_BASE_FUNCTION,&iface); CHKERRG(ierr);
116  cTx = reinterpret_cast<FatPrismPolynomialBaseCtx*>(iface);
117  if(!cTx->fePtr) {
118  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,
119  "Pointer to element should be given "
120  "when EntPolynomialBaseCtx is constructed "
121  "(use different constructor)"
122  );
123  }
124 
125  int nb_gauss_pts = pts.size2();
126  if(!nb_gauss_pts) {
128  }
129 
130  if(pts.size1()<1) {
131  SETERRQ(
132  PETSC_COMM_SELF,
134  "Wrong dimension of pts, should be at least 3 rows with coordinates"
135  );
136  }
137 
138  const FieldApproximationBase base = cTx->bAse;
140 
141  if(cTx->copyNodeBase==LASTBASE) {
142  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"It is assumed that base for vertices is calculated");
143  } else {
144  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) = data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
145  }
146  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts,6,false);
147  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts,12,false);
148  if(data.dataOnEntities[MBVERTEX][0].getN(base).size1()!=(unsigned int)nb_gauss_pts) {
149  SETERRQ1(
150  PETSC_COMM_SELF,
152  "Base functions or nodes has wrong number of integration points for base %s",
154  );
155  }
156  if(
158  pts.size2()
159  ) {
160  SETERRQ(PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"data inconsistency");
161  }
162 
163  switch (cTx->sPace) {
164  case H1:
167  ierr = getValueH1(pts); CHKERRG(ierr);
168  break;
169  case HDIV:
170  ierr = getValueHdiv(pts); CHKERRG(ierr);
171  break;
172  case HCURL:
173  ierr = getValueHcurl(pts); CHKERRG(ierr);
174  break;
175  case L2:
176  ierr = getValueL2(pts); CHKERRG(ierr);
177  break;
178  default:
179  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"Not yet implemented");
180  }
181 
183 }
FatPrismPolynomialBaseCtx * cTx
Class used to pass element data to calculate base functions on fat prism.
field with continuous normal traction
Definition: definitions.h:170
data structure for finite element entityIt keeps that about indices of degrees of freedom...
#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
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:526
#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 getValueHdiv(MatrixDouble &pts)
static const char *const ApproximationBaseNames[]
Definition: definitions.h:155
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
FieldApproximationBase
approximation base
Definition: definitions.h:140
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:169
const NumeredEntFiniteElement * fePtr
static const MOFEMuuid IDD_FATPRISM_BASE_FUNCTION
MoFEMErrorCode getValueH1(MatrixDouble &pts)
continuous field
Definition: definitions.h:168
MoFEMErrorCode getValueL2(MatrixDouble &pts)
field with C-1 continuity
Definition: definitions.h:171

◆ getValueH1()

MoFEMErrorCode FatPrismPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 249 of file FatPrismPolynomialBase.cpp.

249  {
250 
252 
254  const FieldApproximationBase base = cTx->bAse;
255  // MoFEMErrorCode (*base_polynomials)(
256  // int p,double s,double *diff_s,double *L,double *diffL,const int dim
257  // ) = cTx->basePolynomialsType0;
258 
259  int nb_gauss_pts = pts.size2();
260  int nb_gauss_pts_on_faces = cTx->gaussPtsTrianglesOnly.size2();
261  int nb_gauss_pts_through_thickness = cTx->gaussPtsThroughThickness.size2();
262 
263  try {
264  // nodes
265  // linear for xi,eta and zeta
266  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts,6,false);
267  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts,18);
268  noalias(data.dataOnEntities[MBVERTEX][0].getN(base)) = data.dataOnEntities[MBVERTEX][0].getN(NOBASE);
269  noalias(data.dataOnEntities[MBVERTEX][0].getDiffN(base)) = data.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
270 
271  // edges on triangles
272  for(int ee = 0;ee<9;ee++) {
273  if(ee>=3&&ee<=5) {
274  // through thickness ho approximation
275  // linear xi,eta, ho terms for zeta
276  int order = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDataOrder();
277  int nb_dofs = NBEDGE_H1(order);
278  if((unsigned int)nb_dofs!=cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getN(base).size2()) {
279  SETERRQ2(
280  PETSC_COMM_SELF,MOFEM_DATA_INCONSISTENCY,"nb_dofs != nb_dofs %d != %d",
281  nb_dofs,cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getN(base).size2()
282  );
283  }
284  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,nb_dofs,false);
285  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,3*nb_dofs,false);
286  if(nb_dofs == 0) continue;
287  const int prism_edge_map[9][2] = {
288  {0,1}, {1,2}, {2,0}, {0,3}, {1,4}, {2,5}, {3,4}, {4,5}, {5,3}
289  };
290  int gg = 0;
291  for(int ggf = 0;ggf<nb_gauss_pts_on_faces;ggf++) {
292  double tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBVERTEX][0].getN(base)(ggf,prism_edge_map[ee][0]);
293  double dksi_tri_n[2];
294  for(int kk = 0;kk<2;kk++) {
295  dksi_tri_n[kk] = cTx->dataTrianglesOnly.dataOnEntities[MBVERTEX][0].getDiffN(base)(
296  ggf,2*prism_edge_map[ee][0]+kk
297  );
298  }
299  for(int ggt = 0;ggt<nb_gauss_pts_through_thickness;ggt++,gg++) {
300 
301  double zeta = cTx->gaussPtsThroughThickness(0,ggt);
302  double n0 = N_MBEDGE0(zeta);
303  double n1 = N_MBEDGE1(zeta);
304  double n0n1 = n0*n1;
305 
306  for(int dd = 0;dd<nb_dofs;dd++) {
307 
308  double l = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getN(base)(ggt,dd);
309  double diff_l = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDiffN(base)(ggt,dd);
310 
311  double edge_m = n0n1*l;
312  double dzeta_edge_m = (diffN_MBEDGE0*n1+n0*diffN_MBEDGE1)*l + n0n1*diff_l;
313  data.dataOnEntities[MBEDGE][ee].getN(base)(gg,dd) = tri_n*edge_m;
314  for(int kk = 0;kk<2;kk++) {
315  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg,3*dd+kk) = dksi_tri_n[kk]*edge_m;
316  }
317  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg,3*dd+2) = tri_n*dzeta_edge_m;
318 
319  }
320  }
321  }
322  } else {
323  // on triangles ho approximation
324  // ho terms on edges, linear zeta
325  int nb_dofs = cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getN(base).size2();
326  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts,nb_dofs,false);
327  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,3*nb_dofs,false);
328  for(int dd = 0;dd<nb_dofs;dd++) {
329  int gg = 0;
330  for(int ggf = 0;ggf<nb_gauss_pts_on_faces;ggf++) {
331  double tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getN(base)(ggf,dd);
332  double dksi_tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getDiffN(base)(ggf,2*dd+0);
333  double deta_tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getDiffN(base)(ggf,2*dd+1);
334  for(int ggt = 0;ggt<nb_gauss_pts_through_thickness;ggt++,gg++) {
335  double zeta = cTx->gaussPtsThroughThickness(0,ggt);
336  double dzeta,edge_shape;
337  if(ee<3) {
338  dzeta = diffN_MBEDGE0;
339  edge_shape = N_MBEDGE0(zeta);
340  } else {
341  dzeta = diffN_MBEDGE1;
342  edge_shape = N_MBEDGE1(zeta);
343  }
344  data.dataOnEntities[MBEDGE][ee].getN(base)(gg,dd) = tri_n*edge_shape;
345  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg,3*dd+0) = dksi_tri_n*edge_shape;
346  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg,3*dd+1) = deta_tri_n*edge_shape;
347  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg,3*dd+2) = tri_n*dzeta;
348  }
349  }
350  }
351  }
352  }
353  // triangles
354  // ho on triangles, linear zeta
355  for(int ff = 3;ff<=4;ff++) {
356  int nb_dofs;
357  nb_dofs = cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getN(base).size2();
358  data.dataOnEntities[MBTRI][ff].getN(base).resize(nb_gauss_pts,nb_dofs,false);
359  data.dataOnEntities[MBTRI][ff].getDiffN(base).resize(nb_gauss_pts,3*nb_dofs,false);
360  for(int dd = 0;dd<nb_dofs;dd++) {
361  int gg = 0;
362  for(int ggf = 0;ggf<nb_gauss_pts_on_faces;ggf++) {
363  double tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getN(base)(ggf,dd);
364  double dksi_tri_n[2];
365  for(int kk = 0;kk<2;kk++) {
366  dksi_tri_n[kk] = cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getDiffN(base)(ggf,2*dd+kk);
367  }
368  for(int ggt = 0;ggt<nb_gauss_pts_through_thickness;ggt++,gg++) {
369  double zeta = cTx->gaussPtsThroughThickness(0,ggt);
370  double dzeta,edge_shape;
371  if(ff == 3) {
372  dzeta = diffN_MBEDGE0;
373  edge_shape = N_MBEDGE0(zeta);
374  } else {
375  dzeta = diffN_MBEDGE1;
376  edge_shape = N_MBEDGE1(zeta);
377  }
378  data.dataOnEntities[MBTRI][ff].getN(base)(gg,dd) = tri_n*edge_shape;
379  for(int kk = 0;kk<2;kk++) {
380  data.dataOnEntities[MBTRI][ff].getDiffN(base)(gg,3*dd+kk) = dksi_tri_n[kk]*edge_shape;
381  }
382  data.dataOnEntities[MBTRI][ff].getDiffN(base)(gg,3*dd+2) = tri_n*dzeta;
383  }
384  }
385  }
386  }
387  // quads
388  // higher order edges and zeta
389  {
390 
391  int quads_nodes[3*4];
392  int quad_order[3] = { 0, 0, 0};
393  double *quad_n[3],*diff_quad_n[3];
395  SideNumber_multiIndex::nth_index<1>::type::iterator siit;
396  siit = side_table.get<1>().lower_bound(boost::make_tuple(MBQUAD,0));
397  SideNumber_multiIndex::nth_index<1>::type::iterator hi_siit;
398  hi_siit = side_table.get<1>().upper_bound(boost::make_tuple(MBQUAD,3));
399  EntityHandle ent = cTx->fePtr->getEnt();
400  int num_nodes;
401  const EntityHandle *conn;
402  rval = cTx->mOab.get_connectivity(ent,conn,num_nodes,true); CHKERRQ_MOAB(rval);
403  // std::cerr << "\n\n" << std::endl;
404  // const int quad_nodes[3][4] = { {0,1,4,3}, {1,2,5,4}, {0,2,5,3} };
405  for(;siit!=hi_siit;siit++) {
406  // std::cerr << "sn " << siit->side_number << std::endl;
407  int num_nodes_quad;
408  const EntityHandle *conn_quad;
409  EntityHandle quad = siit->get()->ent;
410  rval = cTx->mOab.get_connectivity(
411  quad,conn_quad,num_nodes_quad,true
412  ); CHKERRQ_MOAB(rval);
413  for(int nn = 0;nn<num_nodes_quad;nn++) {
414  quads_nodes[4*siit->get()->side_number+nn] = std::distance(conn,std::find(conn,conn+6,conn_quad[nn]));
415  // std::cerr
416  // << "quad " << quad
417  // << " side number " << siit->side_number
418  // << " " << quads_nodes[4*siit->side_number+nn]
419  // << " " << conn[quads_nodes[4*siit->side_number+nn]]
420  // << " " << conn_quad[nn]
421  // << std::endl;
422  }
423  int order = data.dataOnEntities[MBQUAD][siit->get()->side_number].getDataOrder();
424  quad_order[(int)siit->get()->side_number] = order;
425  data.dataOnEntities[MBQUAD][siit->get()->side_number].getN(base).resize(nb_gauss_pts,NBFACEQUAD_H1(order),false);
426  data.dataOnEntities[MBQUAD][siit->get()->side_number].getDiffN(base).resize(nb_gauss_pts,3*NBFACEQUAD_H1(order),false);
427  if(data.dataOnEntities[MBQUAD][siit->get()->side_number].getN(base).size2()>0) {
428  quad_n[(int)siit->get()->side_number] = &*data.dataOnEntities[MBQUAD][siit->get()->side_number].getN(base).data().begin();
429  diff_quad_n[(int)siit->get()->side_number] = &*data.dataOnEntities[MBQUAD][siit->get()->side_number].getDiffN(base).data().begin();
430  } else {
431  quad_n[(int)siit->get()->side_number] = NULL;
432  diff_quad_n[(int)siit->get()->side_number] = NULL;
433  }
434  }
435  if(quad_order[0]>0||quad_order[1]>0||quad_order[2]>0) {
436  double *vertex_n = &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin();
437  double *diff_vertex_n = &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin();
439  quads_nodes,quad_order,vertex_n,diff_vertex_n,quad_n,diff_quad_n,nb_gauss_pts,Legendre_polynomials
440  ); CHKERRG(ierr);
441  }
442  }
443  // prism
444  {
445  int order = data.dataOnEntities[MBPRISM][0].getDataOrder();
446  double *vertex_n = &data.dataOnEntities[MBVERTEX][0].getN(base)(0,0);
447  double *diff_vertex_n = &data.dataOnEntities[MBVERTEX][0].getDiffN(base)(0,0);
448  data.dataOnEntities[MBPRISM][0].getN(base).resize(nb_gauss_pts,NBVOLUMEPRISM_H1(order),false);
449  data.dataOnEntities[MBPRISM][0].getDiffN(base).resize(nb_gauss_pts,3*NBVOLUMEPRISM_H1(order),false);
450  if(NBVOLUMEPRISM_H1(order)>0) {
452  order,
453  vertex_n,
454  diff_vertex_n,
455  &data.dataOnEntities[MBPRISM][0].getN(base)(0,0),
456  &data.dataOnEntities[MBPRISM][0].getDiffN(base)(0,0),
457  nb_gauss_pts,
459  ); CHKERRG(ierr);
460  }
461  }
462  } catch (MoFEMException const &e) {
463  SETERRQ(PETSC_COMM_SELF,e.errorCode,e.errorMessage);
464  } catch (std::exception& ex) {
465  std::ostringstream ss;
466  ss << "thorw in method: " << ex.what() << " at line " << __LINE__ << " in file " << __FILE__;
467  SETERRQ(PETSC_COMM_SELF,MOFEM_STD_EXCEPTION_THROW,ss.str().c_str());
468  }
469 
471 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:497
PetscErrorCode Legendre_polynomials(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Calculate Legendre approximation basis.
FatPrismPolynomialBaseCtx * cTx
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
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
data structure for finite element entityIt keeps that about indices of degrees of freedom...
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
DataForcesAndSourcesCore & dataTroughThickness
const FieldApproximationBase bAse
#define N_MBEDGE1(x)
edge shape function
Definition: fem_tools.h:80
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:526
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.
Exception to catch.
Definition: Common.hpp:26
#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]
#define diffN_MBEDGE1
derivative of edge shape function
Definition: fem_tools.h:82
char errorMessage[255]
Definition: Common.hpp:28
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
DataForcesAndSourcesCore & dataTrianglesOnly
#define N_MBEDGE0(x)
edge shape function
Definition: fem_tools.h:79
FieldApproximationBase
approximation base
Definition: definitions.h:140
static MoFEMErrorCodeGeneric< moab::ErrorCode > rval
Definition: Common.hpp:78
PetscErrorCode H1_QuadShapeFunctions_MBPRISM(int *faces_nodes, int *p, double *N, double *diffN, double *faceN[], double *diff_faceN[], int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:586
#define NBFACEQUAD_H1(P)
Number of base functions on quad for H1 space.
DataForcesAndSourcesCore & dAta
PetscErrorCode H1_VolumeShapeFunctions_MBPRISM(int p, double *N, double *diffN, double *volumeN, double *diff_volumeN, int GDIM, PetscErrorCode(*base_polynomials)(int p, double s, double *diff_s, double *L, double *diffL, const int dim))
Definition: h1.c:696
const int errorCode
Definition: Common.hpp:27
#define diffN_MBEDGE0
derivative of edge shape function
Definition: fem_tools.h:81
SideNumber_multiIndex & getSideNumberTable() const
const NumeredEntFiniteElement * fePtr
#define NBVOLUMEPRISM_H1(P)
Number of base functions on prism for H1 space.
EntityHandle getEnt() const
Get element entity.

◆ getValueH1ThroughThickness()

MoFEMErrorCode FatPrismPolynomialBase::getValueH1ThroughThickness ( )
private

Definition at line 206 of file FatPrismPolynomialBase.cpp.

206  {
207 
209 
210  // DataForcesAndSourcesCore& data = cTx->dAta;
211  const FieldApproximationBase base = cTx->bAse;
212  PetscErrorCode (*base_polynomials)(
213  int p,double s,double *diff_s,double *L,double *diffL,const int dim
215 
216  int nb_gauss_pts_through_thickness = cTx->gaussPtsThroughThickness.size2();
217 
218  // Calculate Legendre approx. on edges
219  for(unsigned int ee = 3;ee<=5;ee++) {
220  int sense = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getSense();
221  int order = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDataOrder()-2;
222  cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getN(base).resize(
223  nb_gauss_pts_through_thickness,order<0?0:1+order,false
224  );
225  cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
226  nb_gauss_pts_through_thickness,order<0?0:1+order,false
227  );
228  if(order<0) continue;
229  double diff_s = 2.; // s = s(xi), ds/dxi = 2., because change of basis
230  for(int gg = 0;gg<nb_gauss_pts_through_thickness;gg++) {
231  double s = 2*cTx->gaussPtsThroughThickness(0,gg)-1; // makes form -1..1
232  if(!sense) {
233  s *= -1;
234  diff_s *= -1;
235  }
236  // calculate Legendre polynomials at integration points on edges thorough thickness
237  ierr = base_polynomials(
238  order,s,&diff_s,
239  &cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getN(base)(gg,0),
240  &cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg,0),
241  1
242  ); CHKERRG(ierr);
243  }
244  }
245 
247 }
FatPrismPolynomialBaseCtx * cTx
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
DataForcesAndSourcesCore & dataTroughThickness
const FieldApproximationBase bAse
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:526
#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 MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
FieldApproximationBase
approximation base
Definition: definitions.h:140

◆ getValueH1TrianglesOnly()

MoFEMErrorCode FatPrismPolynomialBase::getValueH1TrianglesOnly ( )
private

Definition at line 185 of file FatPrismPolynomialBase.cpp.

185  {
186 
188 
189  const FieldApproximationBase base = cTx->bAse;
190  // MoFEMErrorCode (*base_polynomials)(
191  // int p,double s,double *diff_s,double *L,double *diffL,const int dim
192  // ) = cTx->basePolynomialsType0;
193 
196  boost::shared_ptr<BaseFunctionCtx>(
199  )
200  )
201  ); CHKERRG(ierr);
202 
204 }
FatPrismPolynomialBaseCtx * cTx
MoFEMErrorCode getValue(MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
#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
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:526
#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...
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
DataForcesAndSourcesCore & dataTrianglesOnly
FieldApproximationBase
approximation base
Definition: definitions.h:140
Class used to pass element data to calculate base functions on flat prism.
const NumeredEntFiniteElement * fePtr
continuous field
Definition: definitions.h:168

◆ getValueHcurl()

MoFEMErrorCode FatPrismPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 483 of file FatPrismPolynomialBase.cpp.

483  {
484  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"Not yet implemented");
485 }

◆ getValueHdiv()

MoFEMErrorCode FatPrismPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 479 of file FatPrismPolynomialBase.cpp.

479  {
480  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"Not yet implemented");
481 }

◆ getValueL2()

MoFEMErrorCode FatPrismPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 473 of file FatPrismPolynomialBase.cpp.

473  {
475  SETERRQ(PETSC_COMM_SELF,MOFEM_NOT_IMPLEMENTED,"Not yet implemented");
477 }
#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 FatPrismPolynomialBase::query_interface ( const MOFEMuuid uuid,
MoFEM::UnknownInterface **  iface 
) const
virtual

Reimplemented from MoFEM::BaseFunction.

Definition at line 89 of file FatPrismPolynomialBase.cpp.

91  {
93  *iface = NULL;
94  if(uuid == IDD_FATPRISM_BASE_FUNCTION) {
95  *iface = const_cast<FatPrismPolynomialBase*>(this);
97  } else {
98  SETERRQ(PETSC_COMM_WORLD,MOFEM_DATA_INCONSISTENCY,"wrong interference");
99  }
102 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:483
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:526
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:490
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Common.hpp:80
Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (t...
static const MOFEMuuid IDD_FATPRISM_BASE_FUNCTION
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, MoFEM::UnknownInterface **iface) const

Member Data Documentation

◆ cTx

FatPrismPolynomialBaseCtx* MoFEM::FatPrismPolynomialBase::cTx
private

Definition at line 82 of file FatPrismPolynomialBase.hpp.

◆ diffN

MatrixDouble MoFEM::FatPrismPolynomialBase::diffN
private

Definition at line 98 of file FatPrismPolynomialBase.hpp.

◆ N

MatrixDouble MoFEM::FatPrismPolynomialBase::N
private

Definition at line 97 of file FatPrismPolynomialBase.hpp.


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