v0.14.0
FatPrismPolynomialBase.cpp
Go to the documentation of this file.
1 /** \file FatPrismPolynomialBase.cpp
2 \brief Implementation of Ainsworth-Cole H1 base on edge
3 
4 \todo Prism element can be integrated exploiting tonsorial product. Current
5 implementation do not take that opportunity. That can be viewed as a bug.
6 
7 */
8 
9 
10 
11 using namespace MoFEM;
12 
14  boost::typeindex::type_index type_index, UnknownInterface **iface) const {
15  *iface = const_cast<FatPrismPolynomialBaseCtx *>(this);
16  return 0;
17 }
18 
20  EntitiesFieldData &data, EntitiesFieldData &data_triangles_only,
21  EntitiesFieldData &data_trough_thickness,
22  MatrixDouble &gauss_pts_triangles_only,
23  MatrixDouble &gauss_pts_through_thickness, moab::Interface &moab,
24  const NumeredEntFiniteElement *fe_ptr, const FieldSpace space,
25  const FieldApproximationBase base,
26  const FieldApproximationBase copy_node_base)
27  : EntPolynomialBaseCtx(data, space, CONTINUOUS, base, copy_node_base),
28  dataTrianglesOnly(data_triangles_only),
29  dataTroughThickness(data_trough_thickness),
30  gaussPtsTrianglesOnly(gauss_pts_triangles_only),
31  gaussPtsThroughThickness(gauss_pts_through_thickness), mOab(moab),
32  fePtr(fe_ptr) {
33 
34  ierr = setBase();
35  CHKERRABORT(PETSC_COMM_WORLD, ierr);
36 }
38 
40 FatPrismPolynomialBase::query_interface(boost::typeindex::type_index type_index,
41  UnknownInterface **iface) const {
42  *iface = const_cast<FatPrismPolynomialBase *>(this);
43  return 0;
44 }
45 
48 
51  boost::shared_ptr<BaseFunctionCtx> ctx_ptr) {
52 
54 
56 
57  if (!cTx->fePtr)
58  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
59  "Pointer to element should be given "
60  "when EntPolynomialBaseCtx is constructed "
61  "(use different constructor)");
62 
63  int nb_gauss_pts = pts.size2();
64  if (!nb_gauss_pts)
66 
67  if (pts.size1() < 3)
68  SETERRQ(
69  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
70  "Wrong dimension of pts, should be at least 3 rows with coordinates");
71 
72  const FieldApproximationBase base = cTx->bAse;
73  EntitiesFieldData &data = cTx->dAta;
74 
75  if (cTx->copyNodeBase == LASTBASE)
76  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
77  "It is assumed that base for vertices is calculated");
78  else
79  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
80  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
81 
82  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 6, false);
83  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 12,
84  false);
85  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
86  (unsigned int)nb_gauss_pts)
87  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
88  "Base functions or nodes has wrong number of integration points "
89  "for base %s",
91 
92  if (cTx->gaussPtsTrianglesOnly.size2() *
93  cTx->gaussPtsThroughThickness.size2() !=
94  pts.size2())
95  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
96 
97  switch (cTx->sPace) {
98  case H1:
101  CHKERR getValueH1(pts);
102  break;
103  case HDIV:
104  CHKERR getValueHdiv(pts);
105  break;
106  case HCURL:
107  CHKERR getValueHcurl(pts);
108  break;
109  case L2:
110  CHKERR getValueL2(pts);
111  break;
112  default:
113  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
114  }
115 
117 }
118 
121  const FieldApproximationBase base = cTx->bAse;
124  boost::shared_ptr<BaseFunctionCtx>(new FlatPrismPolynomialBaseCtx(
125  cTx->dataTrianglesOnly, cTx->mOab, cTx->fePtr, H1, base, NOBASE)));
127 }
128 
131 
132  // EntitiesFieldData& data = cTx->dAta;
133  const FieldApproximationBase base = cTx->bAse;
134  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
135  double *diffL, const int dim) =
137 
138  int nb_gauss_pts_through_thickness = cTx->gaussPtsThroughThickness.size2();
139  // Calculate Legendre approx. on edges
140  for (unsigned int ee = 3; ee <= 5; ee++) {
141  auto &ent_data = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee];
142  int sense = ent_data.getSense();
143  int order = ent_data.getOrder();
144  int nb_dofs = NBEDGE_H1(order);
145  ent_data.getN(base).resize(nb_gauss_pts_through_thickness, nb_dofs, false);
146  ent_data.getDiffN(base).resize(nb_gauss_pts_through_thickness, nb_dofs,
147  false);
148 
149  if (nb_dofs > 0) {
150  double diff_s = 1.; // s = s(xi), ds/dxi = 2., because change of basis
151  for (int gg = 0; gg < nb_gauss_pts_through_thickness; gg++) {
152  double s =
153  2 * cTx->gaussPtsThroughThickness(0, gg) - 1; // makes form -1..1
154  if (sense == -1) {
155  s *= -2;
156  diff_s *= -2;
157  }
158  // calculate Legendre polynomials at integration points on edges
159  // thorough thickness
160  CHKERR base_polynomials(order - 2, s, &diff_s,
161  &ent_data.getN(base)(gg, 0),
162  &ent_data.getDiffN(base)(gg, 0), 1);
163  }
164  }
165  }
166 
168 }
169 
172 
173  EntitiesFieldData &data = cTx->dAta;
174  const FieldApproximationBase base = cTx->bAse;
175 
176  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
177  double *diffL, const int dim) =
179 
180  int nb_gauss_pts = pts.size2();
181  int nb_gauss_pts_on_faces = cTx->gaussPtsTrianglesOnly.size2();
182  int nb_gauss_pts_through_thickness = cTx->gaussPtsThroughThickness.size2();
183 
184  // nodes
185  // linear for xi, eta and zeta
186  auto &vert_dat = data.dataOnEntities[MBVERTEX][0];
187  vert_dat.getN(base).resize(nb_gauss_pts, 6, false);
188  vert_dat.getDiffN(base).resize(nb_gauss_pts, 18);
189  noalias(vert_dat.getN(base)) = data.dataOnEntities[MBVERTEX][0].getN(NOBASE);
190  noalias(vert_dat.getDiffN(base)) =
191  data.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
192 
193  auto &vert_n = vert_dat.getN(base);
194  auto &vert_diff_n = vert_dat.getDiffN(base);
195 
196  constexpr int prism_edge_map[9][2] = {{0, 1}, {1, 2}, {2, 0}, {0, 3}, {1, 4},
197  {2, 5}, {3, 4}, {4, 5}, {5, 3}};
198 
199  auto edge_through_thickness = [&](const int ee) {
201  // through thickness ho approximation
202  // linear xi,eta, ho terms for zeta
203 
204  auto &thickness_ent = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee];
205  auto &prism_ent = data.dataOnEntities[MBEDGE][ee];
206 
207  int order = thickness_ent.getOrder();
208  int nb_dofs = NBEDGE_H1(order);
209  if ((unsigned int)nb_dofs != thickness_ent.getN(base).size2())
210  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
211  "nb_dofs != nb_dofs %d != %d", nb_dofs,
212  thickness_ent.getN(base).size2());
213  prism_ent.getN(base).resize(nb_gauss_pts, nb_dofs, false);
214  prism_ent.getDiffN(base).resize(nb_gauss_pts, 3 * nb_dofs, false);
215  if (nb_dofs == 0)
217 
218 
219  int gg = 0;
220  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
221 
222  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
223  double extrude = vert_n(gg, prism_edge_map[ee][0]) +
224  vert_n(gg, prism_edge_map[ee][1]);
225 
226  double diff_extrude[3];
227  for (auto d : {0, 1, 2})
228  diff_extrude[d] = vert_diff_n(gg, 3 * prism_edge_map[ee][0] + d) +
229  vert_diff_n(gg, 3 * prism_edge_map[ee][1] + d);
230 
231  double n0 = vert_n(gg, 0) + vert_n(gg, 1) + vert_n(gg, 2);
232  double n1 = vert_n(gg, 3) + vert_n(gg, 4) + vert_n(gg, 5);
233 
234  // Calculate base through thickness directly from shape fuctions. I
235  // leave that bit of the code, could be useful for integration on the
236  // skeleton.
237 
238  // double l[order + 1], diff_l[3 * (order + 1)];
239  // double ksi = n1 - n0;
240  // double diff_ksi[3];
241  // for (auto d : {0, 1, 2})
242  // diff_ksi[d] =
243  // vert_diff_n(gg, 3 * 3 + d) + vert_diff_n(gg, 3 * 4 + d) +
244  // vert_diff_n(gg, 3 * 5 + d) - vert_diff_n(gg, 3 * 0 + d) -
245  // vert_diff_n(gg, 3 * 1 + d) - vert_diff_n(gg, 3 * 2 + d);
246  // if(sense == -1) {
247  // ksi *= -1;
248  // for (auto d : {0, 1, 2})
249  // diff_ksi[d] *= -1;
250  // }
251  // CHKERR base_polynomials(order - 2, ksi, diff_ksi, l, diff_l, 3);
252 
253  double *l = &(thickness_ent.getN(base)(ggt, 0));
254  double *diff_l_1d = &(thickness_ent.getDiffN(base)(ggt, 0));
255 
256  double bubble = n0 * n1;
257  double diff_bubble[3];
258  for (auto d : {0, 1, 2})
259  diff_bubble[d] =
260  n0 * (vert_diff_n(gg, 3 * 3 + d) + vert_diff_n(gg, 3 * 4 + d) +
261  vert_diff_n(gg, 3 * 5 + d)) +
262 
263  n1 * (vert_diff_n(gg, 3 * 0 + d) + vert_diff_n(gg, 3 * 1 + d) +
264  vert_diff_n(gg, 3 * 2 + d));
265 
266  for (int dd = 0; dd != nb_dofs; dd++) {
267 
268  prism_ent.getN(base)(gg, dd) = extrude * bubble * l[dd];
269  for (auto d : {0, 1, 2})
270  prism_ent.getDiffN(base)(gg, 3 * dd + d) =
271  diff_extrude[d] * bubble * l[dd] +
272  extrude * diff_bubble[d] * l[dd];
273 
274  prism_ent.getDiffN(base)(gg, 3 * dd + 2) +=
275  extrude * bubble * 2 * diff_l_1d[dd];
276 
277  }
278  }
279  }
281  };
282 
283  auto edge_on_the_triangle = [&](const int ee) {
285  // on triangles ho approximation
286  // ho terms on edges, linear zeta
287  int nb_dofs =
288  cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getN(base).size2();
289  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
290  false);
291  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,
292  3 * nb_dofs, false);
293  for (int dd = 0; dd < nb_dofs; dd++) {
294  int gg = 0;
295  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
296  double tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getN(
297  base)(ggf, dd);
298  double dksi_tri_n =
299  cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getDiffN(base)(
300  ggf, 2 * dd + 0);
301  double deta_tri_n =
302  cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getDiffN(base)(
303  ggf, 2 * dd + 1);
304  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
305  double zeta = cTx->gaussPtsThroughThickness(0, ggt);
306  double dzeta, edge_shape;
307  if (ee < 3) {
308  dzeta = diffN_MBEDGE0;
309  edge_shape = N_MBEDGE0(zeta);
310  } else {
311  dzeta = diffN_MBEDGE1;
312  edge_shape = N_MBEDGE1(zeta);
313  }
314  data.dataOnEntities[MBEDGE][ee].getN(base)(gg, dd) =
315  tri_n * edge_shape;
316  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 0) =
317  dksi_tri_n * edge_shape;
318  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 1) =
319  deta_tri_n * edge_shape;
320  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 2) =
321  tri_n * dzeta;
322  }
323  }
324  }
326  };
327 
328  auto trinagle_through_thickness = [&](const int ff) {
330  int nb_dofs;
331  nb_dofs =
332  cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getN(base).size2();
333  data.dataOnEntities[MBTRI][ff].getN(base).resize(nb_gauss_pts, nb_dofs,
334  false);
335  data.dataOnEntities[MBTRI][ff].getDiffN(base).resize(nb_gauss_pts,
336  3 * nb_dofs, false);
337  for (int dd = 0; dd < nb_dofs; dd++) {
338  int gg = 0;
339  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
340  double tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getN(
341  base)(ggf, dd);
342  double dksi_tri_n[2];
343  for (int kk = 0; kk < 2; kk++) {
344  dksi_tri_n[kk] =
345  cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getDiffN(base)(
346  ggf, 2 * dd + kk);
347  }
348  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
349  double zeta = cTx->gaussPtsThroughThickness(0, ggt);
350  double dzeta, edge_shape;
351  if (ff == 3) {
352  dzeta = diffN_MBEDGE0;
353  edge_shape = N_MBEDGE0(zeta);
354  } else {
355  dzeta = diffN_MBEDGE1;
356  edge_shape = N_MBEDGE1(zeta);
357  }
358  data.dataOnEntities[MBTRI][ff].getN(base)(gg, dd) =
359  tri_n * edge_shape;
360  for (auto kk : {0, 1})
361  data.dataOnEntities[MBTRI][ff].getDiffN(base)(gg, 3 * dd + kk) =
362  dksi_tri_n[kk] * edge_shape;
363  data.dataOnEntities[MBTRI][ff].getDiffN(base)(gg, 3 * dd + 2) =
364  tri_n * dzeta;
365  }
366  }
367  }
369  };
370 
371  auto quads_base = [&]() {
373  int quads_nodes[3 * 4];
374  int quad_order[3] = {0, 0, 0};
375  double *quad_n[3], *diff_quad_n[3];
376  SideNumber_multiIndex &side_table =
378  SideNumber_multiIndex::nth_index<1>::type::iterator siit;
379  siit = side_table.get<1>().lower_bound(boost::make_tuple(MBQUAD, 0));
380  SideNumber_multiIndex::nth_index<1>::type::iterator hi_siit;
381  hi_siit = side_table.get<1>().upper_bound(boost::make_tuple(MBQUAD, 3));
382  if (std::distance(siit, hi_siit) != 3)
383  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
384  "Expected three quads on prisms");
385  EntityHandle ent = cTx->fePtr->getEnt();
386  int num_nodes;
387  const EntityHandle *conn;
388  CHKERR cTx->mOab.get_connectivity(ent, conn, num_nodes, true);
389  for (; siit != hi_siit; ++siit) {
390  int side = siit->get()->side_number;
391  if(side < 0 && side > 2)
392  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
393  "Side in range 0,1,2 expected");
394  int num_nodes_quad;
395  const EntityHandle *conn_quad;
396  EntityHandle quad = siit->get()->ent;
397  CHKERR cTx->mOab.get_connectivity(quad, conn_quad, num_nodes_quad, true);
398  for (int nn = 0; nn < num_nodes_quad; nn++) {
399  quads_nodes[4 * side + nn] =
400  std::distance(conn, std::find(conn, conn + 6, conn_quad[nn]));
401  }
402  int order = data.dataOnEntities[MBQUAD][side].getOrder();
403  quad_order[side] = order;
404  data.dataOnEntities[MBQUAD][side].getN(base).resize(
405  nb_gauss_pts, NBFACEQUAD_H1(order), false);
406  data.dataOnEntities[MBQUAD][side].getDiffN(base).resize(
407  nb_gauss_pts, 3 * NBFACEQUAD_H1(order), false);
408  if (data.dataOnEntities[MBQUAD][side].getN(base).size2() > 0) {
409  quad_n[side] =
410  &*data.dataOnEntities[MBQUAD][side].getN(base).data().begin();
411  diff_quad_n[side] =
412  &*data.dataOnEntities[MBQUAD][side].getDiffN(base).data().begin();
413  } else {
414  quad_n[side] = NULL;
415  diff_quad_n[side] = NULL;
416  }
417  }
418  if (quad_order[0] > 0 || quad_order[1] > 0 || quad_order[2] > 0) {
419  double *vertex_n =
420  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin();
421  double *diff_vertex_n =
422  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin();
423  CHKERR H1_QuadShapeFunctions_MBPRISM(quads_nodes, quad_order, vertex_n,
424  diff_vertex_n, quad_n, diff_quad_n,
425  nb_gauss_pts, base_polynomials);
426  }
428  };
429 
430  auto prim_base = [&]() {
432  int order = data.dataOnEntities[MBPRISM][0].getOrder();
433  double *vertex_n = &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0);
434  double *diff_vertex_n =
435  &data.dataOnEntities[MBVERTEX][0].getDiffN(base)(0, 0);
436  data.dataOnEntities[MBPRISM][0].getN(base).resize(
437  nb_gauss_pts, NBVOLUMEPRISM_H1(order), false);
438  data.dataOnEntities[MBPRISM][0].getDiffN(base).resize(
439  nb_gauss_pts, 3 * NBVOLUMEPRISM_H1(order), false);
440  if (NBVOLUMEPRISM_H1(order) > 0) {
442  order, vertex_n, diff_vertex_n,
443  &data.dataOnEntities[MBPRISM][0].getN(base)(0, 0),
444  &data.dataOnEntities[MBPRISM][0].getDiffN(base)(0, 0), nb_gauss_pts,
445  base_polynomials);
446  }
448  };
449 
450  // edges on triangles
451  int ee = 0;
452  for (; ee <= 2; ++ee)
453  CHKERR edge_on_the_triangle(ee);
454  for (; ee <= 5; ++ee)
455  CHKERR edge_through_thickness(ee);
456  for (; ee <= 8; ++ee)
457  CHKERR edge_on_the_triangle(ee);
458 
459  // triangles
460  // ho on triangles, linear zeta
461  for (int ff = 3; ff <= 4; ++ff)
462  CHKERR trinagle_through_thickness(ff);
463 
464  // quads
465  // higher order edges and zeta
466  CHKERR quads_base();
467 
468  // prism
469  CHKERR prim_base();
470 
472 }
473 
476  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
478 }
479 
481  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
482 }
483 
485  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
486 }
UBlasMatrix< double >
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
MoFEM::UnknownInterface::getInterface
MoFEMErrorCode getInterface(IFACE *&iface) const
Get interface reference to pointer of interface.
Definition: UnknownInterface.hpp:93
diffN_MBEDGE0
#define diffN_MBEDGE0
derivative of edge shape function
Definition: fem_tools.h:107
MoFEM::FatPrismPolynomialBaseCtx::gaussPtsTrianglesOnly
MatrixDouble & gaussPtsTrianglesOnly
Definition: FatPrismPolynomialBase.hpp:31
MoFEM::FatPrismPolynomialBaseCtx::query_interface
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: FatPrismPolynomialBase.cpp:13
H1
@ H1
continuous field
Definition: definitions.h:85
NBEDGE_H1
#define NBEDGE_H1(P)
Numer of base function on edge for H1 space.
Definition: h1_hdiv_hcurl_l2.h:55
MoFEM::FatPrismPolynomialBase::getValue
MoFEMErrorCode getValue(MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
Definition: FatPrismPolynomialBase.cpp:50
LASTBASE
@ LASTBASE
Definition: definitions.h:69
EntityHandle
MoFEM::FlatPrismPolynomialBase::getValue
MoFEMErrorCode getValue(MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
Definition: FlatPrismPolynomialBase.cpp:37
MoFEM::EntPolynomialBaseCtx
Class used to pass element data to calculate base functions on tet,triangle,edge.
Definition: EntPolynomialBaseCtx.hpp:22
MoFEM::interface_EntFiniteElement::getSideNumberTable
SideNumber_multiIndex & getSideNumberTable() const
Definition: FEMultiIndices.hpp:705
NOBASE
@ NOBASE
Definition: definitions.h:59
ApproximationBaseNames
const static char *const ApproximationBaseNames[]
Definition: definitions.h:72
L2
@ L2
field with C-1 continuity
Definition: definitions.h:88
MoFEM::FatPrismPolynomialBaseCtx::dataTroughThickness
EntitiesFieldData & dataTroughThickness
Definition: FatPrismPolynomialBase.hpp:29
MoFEM::FatPrismPolynomialBase::getValueH1ThroughThickness
MoFEMErrorCode getValueH1ThroughThickness()
Definition: FatPrismPolynomialBase.cpp:129
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::EntPolynomialBaseCtx::dAta
EntitiesFieldData & dAta
Definition: EntPolynomialBaseCtx.hpp:36
MoFEM::FatPrismPolynomialBaseCtx::FatPrismPolynomialBaseCtx
FatPrismPolynomialBaseCtx(EntitiesFieldData &data, EntitiesFieldData &data_triangles_only, EntitiesFieldData &data_trough_thickness, MatrixDouble &gauss_pts_triangles_only, MatrixDouble &gauss_pts_through_thickness, moab::Interface &moab, const NumeredEntFiniteElement *fe_ptr, const FieldSpace space, const FieldApproximationBase base, const FieldApproximationBase copy_node_base=LASTBASE)
Definition: FatPrismPolynomialBase.cpp:19
MoFEM::interface_RefEntity::getEnt
EntityHandle getEnt() const
Get the entity handle.
Definition: RefEntsMultiIndices.hpp:603
MoFEM::FatPrismPolynomialBase::getValueHdiv
MoFEMErrorCode getValueHdiv(MatrixDouble &pts)
Definition: FatPrismPolynomialBase.cpp:480
N_MBEDGE0
#define N_MBEDGE0(x)
edge shape function
Definition: fem_tools.h:105
MoFEM::FatPrismPolynomialBase
Calculate base functions on tetrahedral.
Definition: FatPrismPolynomialBase.hpp:56
zeta
double zeta
Viscous hardening.
Definition: plastic.cpp:126
MoFEM::FatPrismPolynomialBaseCtx::fePtr
const NumeredEntFiniteElement * fePtr
Definition: FatPrismPolynomialBase.hpp:35
MoFEM::FatPrismPolynomialBase::getValueHcurl
MoFEMErrorCode getValueHcurl(MatrixDouble &pts)
Definition: FatPrismPolynomialBase.cpp:484
MoFEM::FatPrismPolynomialBase::getValueL2
MoFEMErrorCode getValueL2(MatrixDouble &pts)
Definition: FatPrismPolynomialBase.cpp:474
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:2002
FieldSpace
FieldSpace
approximation spaces
Definition: definitions.h:82
MoFEM::EntPolynomialBaseCtx::copyNodeBase
const FieldApproximationBase copyNodeBase
Definition: EntPolynomialBaseCtx.hpp:41
MoFEM::FatPrismPolynomialBase::FatPrismPolynomialBase
FatPrismPolynomialBase()
Definition: FatPrismPolynomialBase.cpp:47
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:548
H1_QuadShapeFunctions_MBPRISM
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:642
MoFEM::FatPrismPolynomialBase::~FatPrismPolynomialBase
~FatPrismPolynomialBase()
Definition: FatPrismPolynomialBase.cpp:46
N_MBEDGE1
#define N_MBEDGE1(x)
edge shape function
Definition: fem_tools.h:106
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
H1_VolumeShapeFunctions_MBPRISM
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:790
SideNumber_multiIndex
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< ordered_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, signed char, &SideNumber::side_number > > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
Definition: RefEntsMultiIndices.hpp:101
MoFEM::FlatPrismPolynomialBase
Calculate base functions on tetrahedral.
Definition: FlatPrismPolynomialBase.hpp:44
MoFEM::EntPolynomialBaseCtx::sPace
const FieldSpace sPace
Definition: EntPolynomialBaseCtx.hpp:37
MoFEM::FatPrismPolynomialBaseCtx
Class used to pass element data to calculate base functions on fat prism.
Definition: FatPrismPolynomialBase.hpp:23
MoFEM::FatPrismPolynomialBase::getValueH1
MoFEMErrorCode getValueH1(MatrixDouble &pts)
Definition: FatPrismPolynomialBase.cpp:170
MoFEM::L
VectorDouble L
Definition: Projection10NodeCoordsOnField.cpp:124
MoFEM::FatPrismPolynomialBaseCtx::~FatPrismPolynomialBaseCtx
~FatPrismPolynomialBaseCtx()
Definition: FatPrismPolynomialBase.cpp:37
MoFEM::FatPrismPolynomialBase::getValueH1TrianglesOnly
MoFEMErrorCode getValueH1TrianglesOnly()
Definition: FatPrismPolynomialBase.cpp:119
NBFACEQUAD_H1
#define NBFACEQUAD_H1(P)
Number of base functions on quad for H1 space.
Definition: h1_hdiv_hcurl_l2.h:65
MoFEM::EntPolynomialBaseCtx::bAse
const FieldApproximationBase bAse
Definition: EntPolynomialBaseCtx.hpp:39
diffN_MBEDGE1
#define diffN_MBEDGE1
derivative of edge shape function
Definition: fem_tools.h:108
MoFEM::UnknownInterface
base class for all interface classes
Definition: UnknownInterface.hpp:34
MoFEM::FatPrismPolynomialBaseCtx::mOab
moab::Interface & mOab
Definition: FatPrismPolynomialBase.hpp:34
FTensor::dd
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
MoFEM::FatPrismPolynomialBase::cTx
FatPrismPolynomialBaseCtx * cTx
Definition: FatPrismPolynomialBase.hpp:68
MoFEM::Exceptions::ierr
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
MoFEM::FatPrismPolynomialBase::query_interface
MoFEMErrorCode query_interface(boost::typeindex::type_index type_index, UnknownInterface **iface) const
Definition: FatPrismPolynomialBase.cpp:40
HCURL
@ HCURL
field with continuous tangents
Definition: definitions.h:86
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::FlatPrismPolynomialBaseCtx
Class used to pass element data to calculate base functions on flat prism.
Definition: FlatPrismPolynomialBase.hpp:21
sdf_hertz_2d_axisymm_plane.d
float d
Definition: sdf_hertz_2d_axisymm_plane.py:4
FieldApproximationBase
FieldApproximationBase
approximation base
Definition: definitions.h:58
MoFEM::EntitiesFieldData::dataOnEntities
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
Definition: EntitiesFieldData.hpp:57
MoFEM::NumeredEntFiniteElement
Partitioned (Indexed) Finite Element in Problem.
Definition: FEMultiIndices.hpp:728
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
NBVOLUMEPRISM_H1
#define NBVOLUMEPRISM_H1(P)
Number of base functions on prism for H1 space.
Definition: h1_hdiv_hcurl_l2.h:80
MoFEM::EntitiesFieldData
data structure for finite element entity
Definition: EntitiesFieldData.hpp:40
CONTINUOUS
@ CONTINUOUS
Regular field.
Definition: definitions.h:100
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
HDIV
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEM::FatPrismPolynomialBaseCtx::gaussPtsThroughThickness
MatrixDouble & gaussPtsThroughThickness
Definition: FatPrismPolynomialBase.hpp:32
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
MoFEM::FatPrismPolynomialBaseCtx::dataTrianglesOnly
EntitiesFieldData & dataTrianglesOnly
Definition: FatPrismPolynomialBase.hpp:28
MoFEM::EntPolynomialBaseCtx::basePolynomialsType0
PetscErrorCode(* basePolynomialsType0)(int p, double s, double *diff_s, double *L, double *diffL, const int dim)
Definition: EntPolynomialBaseCtx.hpp:27
MoFEM::EntPolynomialBaseCtx::setBase
MoFEMErrorCode setBase()
Definition: EntPolynomialBaseCtx.cpp:38