v0.8.19
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 102 of file FatPrismPolynomialBase.cpp.

102 {}

◆ ~FatPrismPolynomialBase()

FatPrismPolynomialBase::~FatPrismPolynomialBase ( )

Definition at line 101 of file FatPrismPolynomialBase.cpp.

101 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 105 of file FatPrismPolynomialBase.cpp.

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

◆ getValueH1()

MoFEMErrorCode FatPrismPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 251 of file FatPrismPolynomialBase.cpp.

251  {
252 
254 
256  const FieldApproximationBase base = cTx->bAse;
257  // MoFEMErrorCode (*base_polynomials)(
258  // int p,double s,double *diff_s,double *L,double *diffL,const int dim
259  // ) = cTx->basePolynomialsType0;
260 
261  int nb_gauss_pts = pts.size2();
262  int nb_gauss_pts_on_faces = cTx->gaussPtsTrianglesOnly.size2();
263  int nb_gauss_pts_through_thickness = cTx->gaussPtsThroughThickness.size2();
264 
265  try {
266  // nodes
267  // linear for xi,eta and zeta
268  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 6, false);
269  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 18);
270  noalias(data.dataOnEntities[MBVERTEX][0].getN(base)) =
271  data.dataOnEntities[MBVERTEX][0].getN(NOBASE);
272  noalias(data.dataOnEntities[MBVERTEX][0].getDiffN(base)) =
273  data.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
274 
275  // edges on triangles
276  for (int ee = 0; ee < 9; ee++) {
277  if (ee >= 3 && ee <= 5) {
278  // through thickness ho approximation
279  // linear xi,eta, ho terms for zeta
280  int order =
281  cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDataOrder();
282  int nb_dofs = NBEDGE_H1(order);
283  if ((unsigned int)nb_dofs !=
285  .getN(base)
286  .size2()) {
287  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
288  "nb_dofs != nb_dofs %d != %d", nb_dofs,
290  .getN(base)
291  .size2());
292  }
293  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
294  false);
295  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
296  nb_gauss_pts, 3 * nb_dofs, false);
297  if (nb_dofs == 0)
298  continue;
299  const int prism_edge_map[9][2] = {{0, 1}, {1, 2}, {2, 0},
300  {0, 3}, {1, 4}, {2, 5},
301  {3, 4}, {4, 5}, {5, 3}};
302  int gg = 0;
303  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
304  double tri_n =
305  cTx->dataTrianglesOnly.dataOnEntities[MBVERTEX][0].getN(base)(
306  ggf, prism_edge_map[ee][0]);
307  double dksi_tri_n[2];
308  for (int kk = 0; kk < 2; kk++) {
309  dksi_tri_n[kk] =
310  cTx->dataTrianglesOnly.dataOnEntities[MBVERTEX][0].getDiffN(
311  base)(ggf, 2 * prism_edge_map[ee][0] + kk);
312  }
313  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
314 
315  double zeta = cTx->gaussPtsThroughThickness(0, ggt);
316  double n0 = N_MBEDGE0(zeta);
317  double n1 = N_MBEDGE1(zeta);
318  double n0n1 = n0 * n1;
319 
320  for (int dd = 0; dd < nb_dofs; dd++) {
321 
322  double l =
323  cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getN(
324  base)(ggt, dd);
325  double diff_l =
326  cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDiffN(
327  base)(ggt, dd);
328 
329  double edge_m = n0n1 * l;
330  double dzeta_edge_m =
331  (diffN_MBEDGE0 * n1 + n0 * diffN_MBEDGE1) * l + n0n1 * diff_l;
332  data.dataOnEntities[MBEDGE][ee].getN(base)(gg, dd) =
333  tri_n * edge_m;
334  for (int kk = 0; kk < 2; kk++) {
335  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(
336  gg, 3 * dd + kk) = dksi_tri_n[kk] * edge_m;
337  }
338  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 2) =
339  tri_n * dzeta_edge_m;
340  }
341  }
342  }
343  } else {
344  // on triangles ho approximation
345  // ho terms on edges, linear zeta
346  int nb_dofs = cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee]
347  .getN(base)
348  .size2();
349  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
350  false);
351  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
352  nb_gauss_pts, 3 * nb_dofs, false);
353  for (int dd = 0; dd < nb_dofs; dd++) {
354  int gg = 0;
355  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
356  double tri_n =
357  cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getN(base)(
358  ggf, dd);
359  double dksi_tri_n =
360  cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getDiffN(
361  base)(ggf, 2 * dd + 0);
362  double deta_tri_n =
363  cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getDiffN(
364  base)(ggf, 2 * dd + 1);
365  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness;
366  ggt++, gg++) {
367  double zeta = cTx->gaussPtsThroughThickness(0, ggt);
368  double dzeta, edge_shape;
369  if (ee < 3) {
370  dzeta = diffN_MBEDGE0;
371  edge_shape = N_MBEDGE0(zeta);
372  } else {
373  dzeta = diffN_MBEDGE1;
374  edge_shape = N_MBEDGE1(zeta);
375  }
376  data.dataOnEntities[MBEDGE][ee].getN(base)(gg, dd) =
377  tri_n * edge_shape;
378  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 0) =
379  dksi_tri_n * edge_shape;
380  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 1) =
381  deta_tri_n * edge_shape;
382  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 2) =
383  tri_n * dzeta;
384  }
385  }
386  }
387  }
388  }
389  // triangles
390  // ho on triangles, linear zeta
391  for (int ff = 3; ff <= 4; ff++) {
392  int nb_dofs;
393  nb_dofs =
394  cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getN(base).size2();
395  data.dataOnEntities[MBTRI][ff].getN(base).resize(nb_gauss_pts, nb_dofs,
396  false);
397  data.dataOnEntities[MBTRI][ff].getDiffN(base).resize(nb_gauss_pts,
398  3 * nb_dofs, false);
399  for (int dd = 0; dd < nb_dofs; dd++) {
400  int gg = 0;
401  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
402  double tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getN(
403  base)(ggf, dd);
404  double dksi_tri_n[2];
405  for (int kk = 0; kk < 2; kk++) {
406  dksi_tri_n[kk] =
407  cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getDiffN(base)(
408  ggf, 2 * dd + kk);
409  }
410  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
411  double zeta = cTx->gaussPtsThroughThickness(0, ggt);
412  double dzeta, edge_shape;
413  if (ff == 3) {
414  dzeta = diffN_MBEDGE0;
415  edge_shape = N_MBEDGE0(zeta);
416  } else {
417  dzeta = diffN_MBEDGE1;
418  edge_shape = N_MBEDGE1(zeta);
419  }
420  data.dataOnEntities[MBTRI][ff].getN(base)(gg, dd) =
421  tri_n * edge_shape;
422  for (int kk = 0; kk < 2; kk++) {
423  data.dataOnEntities[MBTRI][ff].getDiffN(base)(gg, 3 * dd + kk) =
424  dksi_tri_n[kk] * edge_shape;
425  }
426  data.dataOnEntities[MBTRI][ff].getDiffN(base)(gg, 3 * dd + 2) =
427  tri_n * dzeta;
428  }
429  }
430  }
431  }
432  // quads
433  // higher order edges and zeta
434  {
435 
436  int quads_nodes[3 * 4];
437  int quad_order[3] = {0, 0, 0};
438  double *quad_n[3], *diff_quad_n[3];
439  SideNumber_multiIndex &side_table =
441  SideNumber_multiIndex::nth_index<1>::type::iterator siit;
442  siit = side_table.get<1>().lower_bound(boost::make_tuple(MBQUAD, 0));
443  SideNumber_multiIndex::nth_index<1>::type::iterator hi_siit;
444  hi_siit = side_table.get<1>().upper_bound(boost::make_tuple(MBQUAD, 3));
445  EntityHandle ent = cTx->fePtr->getEnt();
446  int num_nodes;
447  const EntityHandle *conn;
448  rval = cTx->mOab.get_connectivity(ent, conn, num_nodes, true);
450  // std::cerr << "\n\n" << std::endl;
451  // const int quad_nodes[3][4] = { {0,1,4,3}, {1,2,5,4}, {0,2,5,3} };
452  for (; siit != hi_siit; siit++) {
453  // std::cerr << "sn " << siit->side_number << std::endl;
454  int num_nodes_quad;
455  const EntityHandle *conn_quad;
456  EntityHandle quad = siit->get()->ent;
457  rval =
458  cTx->mOab.get_connectivity(quad, conn_quad, num_nodes_quad, true);
460  for (int nn = 0; nn < num_nodes_quad; nn++) {
461  quads_nodes[4 * siit->get()->side_number + nn] =
462  std::distance(conn, std::find(conn, conn + 6, conn_quad[nn]));
463  // std::cerr
464  // << "quad " << quad
465  // << " side number " << siit->side_number
466  // << " " << quads_nodes[4*siit->side_number+nn]
467  // << " " << conn[quads_nodes[4*siit->side_number+nn]]
468  // << " " << conn_quad[nn]
469  // << std::endl;
470  }
471  int order = data.dataOnEntities[MBQUAD][siit->get()->side_number]
472  .getDataOrder();
473  quad_order[(int)siit->get()->side_number] = order;
474  data.dataOnEntities[MBQUAD][siit->get()->side_number].getN(base).resize(
475  nb_gauss_pts, NBFACEQUAD_H1(order), false);
476  data.dataOnEntities[MBQUAD][siit->get()->side_number]
477  .getDiffN(base)
478  .resize(nb_gauss_pts, 3 * NBFACEQUAD_H1(order), false);
479  if (data.dataOnEntities[MBQUAD][siit->get()->side_number]
480  .getN(base)
481  .size2() > 0) {
482  quad_n[(int)siit->get()->side_number] =
483  &*data.dataOnEntities[MBQUAD][siit->get()->side_number]
484  .getN(base)
485  .data()
486  .begin();
487  diff_quad_n[(int)siit->get()->side_number] =
488  &*data.dataOnEntities[MBQUAD][siit->get()->side_number]
489  .getDiffN(base)
490  .data()
491  .begin();
492  } else {
493  quad_n[(int)siit->get()->side_number] = NULL;
494  diff_quad_n[(int)siit->get()->side_number] = NULL;
495  }
496  }
497  if (quad_order[0] > 0 || quad_order[1] > 0 || quad_order[2] > 0) {
498  double *vertex_n =
499  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin();
500  double *diff_vertex_n =
501  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin();
503  quads_nodes, quad_order, vertex_n, diff_vertex_n, quad_n,
504  diff_quad_n, nb_gauss_pts, Legendre_polynomials);
505  CHKERRG(ierr);
506  }
507  }
508  // prism
509  {
510  int order = data.dataOnEntities[MBPRISM][0].getDataOrder();
511  double *vertex_n = &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0);
512  double *diff_vertex_n =
513  &data.dataOnEntities[MBVERTEX][0].getDiffN(base)(0, 0);
514  data.dataOnEntities[MBPRISM][0].getN(base).resize(
515  nb_gauss_pts, NBVOLUMEPRISM_H1(order), false);
516  data.dataOnEntities[MBPRISM][0].getDiffN(base).resize(
517  nb_gauss_pts, 3 * NBVOLUMEPRISM_H1(order), false);
518  if (NBVOLUMEPRISM_H1(order) > 0) {
520  order, vertex_n, diff_vertex_n,
521  &data.dataOnEntities[MBPRISM][0].getN(base)(0, 0),
522  &data.dataOnEntities[MBPRISM][0].getDiffN(base)(0, 0), nb_gauss_pts,
524  CHKERRG(ierr);
525  }
526  }
527  } catch (MoFEMException const &e) {
528  SETERRQ(PETSC_COMM_SELF, e.errorCode, e.errorMessage);
529  } catch (std::exception &ex) {
530  std::ostringstream ss;
531  ss << "thorw in method: " << ex.what() << " at line " << __LINE__
532  << " in file " << __FILE__;
533  SETERRQ(PETSC_COMM_SELF, MOFEM_STD_EXCEPTION_THROW, ss.str().c_str());
534  }
535 
537 }
#define CHKERRQ_MOAB(a)
check error code of MoAB function
Definition: definitions.h:513
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:499
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:542
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:506
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:141
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 204 of file FatPrismPolynomialBase.cpp.

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

◆ getValueH1TrianglesOnly()

MoFEMErrorCode FatPrismPolynomialBase::getValueH1TrianglesOnly ( )
private

Definition at line 186 of file FatPrismPolynomialBase.cpp.

186  {
187 
189 
190  const FieldApproximationBase base = cTx->bAse;
191  // MoFEMErrorCode (*base_polynomials)(
192  // int p,double s,double *diff_s,double *L,double *diffL,const int dim
193  // ) = cTx->basePolynomialsType0;
194 
197  boost::shared_ptr<BaseFunctionCtx>(new FlatPrismPolynomialBaseCtx(
198  cTx->dataTrianglesOnly, cTx->mOab, cTx->fePtr, H1, base, NOBASE)));
199  CHKERRG(ierr);
200 
202 }
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:499
const FieldApproximationBase bAse
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
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:141
Class used to pass element data to calculate base functions on flat prism.
const NumeredEntFiniteElement * fePtr
continuous field
Definition: definitions.h:169

◆ getValueHcurl()

MoFEMErrorCode FatPrismPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 549 of file FatPrismPolynomialBase.cpp.

549  {
550  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
551 }

◆ getValueHdiv()

MoFEMErrorCode FatPrismPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 545 of file FatPrismPolynomialBase.cpp.

545  {
546  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
547 }

◆ getValueL2()

MoFEMErrorCode FatPrismPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 539 of file FatPrismPolynomialBase.cpp.

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

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 86 of file FatPrismPolynomialBase.cpp.

87  {
89  *iface = NULL;
90  if (uuid == IDD_FATPRISM_BASE_FUNCTION) {
91  *iface = const_cast<FatPrismPolynomialBase *>(this);
93  } else {
94  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
95  }
96  ierr = BaseFunction::query_interface(uuid, iface);
97  CHKERRG(ierr);
99 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:499
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:542
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
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 81 of file FatPrismPolynomialBase.hpp.

◆ diffN

MatrixDouble MoFEM::FatPrismPolynomialBase::diffN
private

Definition at line 97 of file FatPrismPolynomialBase.hpp.

◆ N

MatrixDouble MoFEM::FatPrismPolynomialBase::N
private

Definition at line 96 of file FatPrismPolynomialBase.hpp.


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