v0.8.23
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 71 of file FatPrismPolynomialBase.hpp.

Constructor & Destructor Documentation

◆ FatPrismPolynomialBase()

FatPrismPolynomialBase::FatPrismPolynomialBase ( )

Definition at line 76 of file FatPrismPolynomialBase.cpp.

76 {}

◆ ~FatPrismPolynomialBase()

FatPrismPolynomialBase::~FatPrismPolynomialBase ( )

Definition at line 75 of file FatPrismPolynomialBase.cpp.

75 {}

Member Function Documentation

◆ getValue()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 79 of file FatPrismPolynomialBase.cpp.

80  {
81 
83 
85  ierr = ctx_ptr->query_interface(IDD_FATPRISM_BASE_FUNCTION, &iface);
86  CHKERRG(ierr);
87  cTx = reinterpret_cast<FatPrismPolynomialBaseCtx *>(iface);
88  if (!cTx->fePtr) {
89  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
90  "Pointer to element should be given "
91  "when EntPolynomialBaseCtx is constructed "
92  "(use different constructor)");
93  }
94 
95  int nb_gauss_pts = pts.size2();
96  if (!nb_gauss_pts) {
98  }
99 
100  if (pts.size1() < 1) {
101  SETERRQ(
102  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
103  "Wrong dimension of pts, should be at least 3 rows with coordinates");
104  }
105 
106  const FieldApproximationBase base = cTx->bAse;
108 
109  if (cTx->copyNodeBase == LASTBASE) {
110  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
111  "It is assumed that base for vertices is calculated");
112  } else {
113  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
114  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
115  }
116  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 6, false);
117  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 12,
118  false);
119  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
120  (unsigned int)nb_gauss_pts) {
121  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
122  "Base functions or nodes has wrong number of integration points "
123  "for base %s",
124  ApproximationBaseNames[base]);
125  }
126  if (cTx->gaussPtsTrianglesOnly.size2() *
127  cTx->gaussPtsThroughThickness.size2() !=
128  pts.size2()) {
129  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
130  }
131 
132  switch (cTx->sPace) {
133  case H1:
135  CHKERRG(ierr);
137  CHKERRG(ierr);
138  ierr = getValueH1(pts);
139  CHKERRG(ierr);
140  break;
141  case HDIV:
142  ierr = getValueHdiv(pts);
143  CHKERRG(ierr);
144  break;
145  case HCURL:
146  ierr = getValueHcurl(pts);
147  CHKERRG(ierr);
148  break;
149  case L2:
150  ierr = getValueL2(pts);
151  CHKERRG(ierr);
152  break;
153  default:
154  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
155  }
156 
158 }
FatPrismPolynomialBaseCtx * cTx
field with continuous normal traction
Definition: definitions.h:173
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:501
const FieldApproximationBase bAse
base class for all interface classes
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
const FieldApproximationBase copyNodeBase
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
static const char *const ApproximationBaseNames[]
Definition: definitions.h:158
FieldApproximationBase
approximation base
Definition: definitions.h:143
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
DataForcesAndSourcesCore & dAta
field with continuous tangents
Definition: definitions.h:172
const NumeredEntFiniteElement * fePtr
static const MOFEMuuid IDD_FATPRISM_BASE_FUNCTION
MoFEMErrorCode getValueH1(MatrixDouble &pts)
continuous field
Definition: definitions.h:171
MoFEMErrorCode getValueL2(MatrixDouble &pts)
field with C-1 continuity
Definition: definitions.h:174

◆ getValueH1()

MoFEMErrorCode FatPrismPolynomialBase::getValueH1 ( MatrixDouble pts)
private

Definition at line 225 of file FatPrismPolynomialBase.cpp.

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

◆ getValueH1ThroughThickness()

MoFEMErrorCode FatPrismPolynomialBase::getValueH1ThroughThickness ( )
private

Definition at line 178 of file FatPrismPolynomialBase.cpp.

178  {
179 
181 
182  // DataForcesAndSourcesCore& data = cTx->dAta;
183  const FieldApproximationBase base = cTx->bAse;
184  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
185  double *diffL, const int dim) =
187 
188  int nb_gauss_pts_through_thickness = cTx->gaussPtsThroughThickness.size2();
189 
190  // Calculate Legendre approx. on edges
191  for (unsigned int ee = 3; ee <= 5; ee++) {
192  int sense = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getSense();
193  int order =
194  cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDataOrder() - 2;
195  cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getN(base).resize(
196  nb_gauss_pts_through_thickness, order < 0 ? 0 : 1 + order, false);
197  cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(
198  nb_gauss_pts_through_thickness, order < 0 ? 0 : 1 + order, false);
199  if (order < 0)
200  continue;
201  double diff_s = 2.; // s = s(xi), ds/dxi = 2., because change of basis
202  for (int gg = 0; gg < nb_gauss_pts_through_thickness; gg++) {
203  double s =
204  2 * cTx->gaussPtsThroughThickness(0, gg) - 1; // makes form -1..1
205  if (!sense) {
206  s *= -1;
207  diff_s *= -1;
208  }
209  // calculate Legendre polynomials at integration points on edges thorough
210  // thickness
211  ierr = base_polynomials(
212  order, s, &diff_s,
213  &cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getN(base)(gg,
214  0),
215  &cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee].getDiffN(base)(
216  gg, 0),
217  1);
218  CHKERRG(ierr);
219  }
220  }
221 
223 }
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:501
DataForcesAndSourcesCore & dataTroughThickness
const FieldApproximationBase bAse
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
boost::ptr_vector< EntData > dataOnEntities[MBMAXTYPE]
FieldApproximationBase
approximation base
Definition: definitions.h:143
const int order
approximation order

◆ getValueH1TrianglesOnly()

MoFEMErrorCode FatPrismPolynomialBase::getValueH1TrianglesOnly ( )
private

Definition at line 160 of file FatPrismPolynomialBase.cpp.

160  {
161 
163 
164  const FieldApproximationBase base = cTx->bAse;
165  // MoFEMErrorCode (*base_polynomials)(
166  // int p,double s,double *diff_s,double *L,double *diffL,const int dim
167  // ) = cTx->basePolynomialsType0;
168 
171  boost::shared_ptr<BaseFunctionCtx>(new FlatPrismPolynomialBaseCtx(
172  cTx->dataTrianglesOnly, cTx->mOab, cTx->fePtr, H1, base, NOBASE)));
173  CHKERRG(ierr);
174 
176 }
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:501
const FieldApproximationBase bAse
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (t...
DataForcesAndSourcesCore & dataTrianglesOnly
FieldApproximationBase
approximation base
Definition: definitions.h:143
Class used to pass element data to calculate base functions on flat prism.
const NumeredEntFiniteElement * fePtr
continuous field
Definition: definitions.h:171

◆ getValueHcurl()

MoFEMErrorCode FatPrismPolynomialBase::getValueHcurl ( MatrixDouble pts)
private

Definition at line 523 of file FatPrismPolynomialBase.cpp.

523  {
524  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
525 }

◆ getValueHdiv()

MoFEMErrorCode FatPrismPolynomialBase::getValueHdiv ( MatrixDouble pts)
private

Definition at line 519 of file FatPrismPolynomialBase.cpp.

519  {
520  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
521 }

◆ getValueL2()

MoFEMErrorCode FatPrismPolynomialBase::getValueL2 ( MatrixDouble pts)
private

Definition at line 513 of file FatPrismPolynomialBase.cpp.

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

◆ query_interface()

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

Reimplemented from MoFEM::BaseFunction.

Definition at line 60 of file FatPrismPolynomialBase.cpp.

61  {
63  *iface = NULL;
64  if (uuid == IDD_FATPRISM_BASE_FUNCTION) {
65  *iface = const_cast<FatPrismPolynomialBase *>(this);
67  } else {
68  SETERRQ(PETSC_COMM_WORLD, MOFEM_DATA_INCONSISTENCY, "wrong interference");
69  }
70  ierr = BaseFunction::query_interface(uuid, iface);
71  CHKERRG(ierr);
73 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:544
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
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 83 of file FatPrismPolynomialBase.hpp.

◆ diffN

MatrixDouble MoFEM::FatPrismPolynomialBase::diffN
private

Definition at line 99 of file FatPrismPolynomialBase.hpp.

◆ N

MatrixDouble MoFEM::FatPrismPolynomialBase::N
private

Definition at line 98 of file FatPrismPolynomialBase.hpp.


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