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,
21  EntitiesFieldData &data_triangles_only,
22  EntitiesFieldData &data_trough_thickness,
23  MatrixDouble &gauss_pts_triangles_only,
24  MatrixDouble &gauss_pts_through_thickness, moab::Interface &moab,
25  const NumeredEntFiniteElement *fe_ptr, const FieldSpace space,
26  const FieldApproximationBase base,
27  const FieldApproximationBase copy_node_base)
28  : EntPolynomialBaseCtx(data, space, base, copy_node_base),
29  dataTrianglesOnly(data_triangles_only),
30  dataTroughThickness(data_trough_thickness),
31  gaussPtsTrianglesOnly(gauss_pts_triangles_only),
32  gaussPtsThroughThickness(gauss_pts_through_thickness), mOab(moab),
33  fePtr(fe_ptr) {
34 
35  ierr = setBase();
36  CHKERRABORT(PETSC_COMM_WORLD, ierr);
37 }
39 
41 FatPrismPolynomialBase::query_interface(boost::typeindex::type_index type_index,
42  UnknownInterface **iface) const {
43  *iface = const_cast<FatPrismPolynomialBase *>(this);
44  return 0;
45 }
46 
49 
52  boost::shared_ptr<BaseFunctionCtx> ctx_ptr) {
53 
55 
57 
58  if (!cTx->fePtr)
59  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
60  "Pointer to element should be given "
61  "when EntPolynomialBaseCtx is constructed "
62  "(use different constructor)");
63 
64  int nb_gauss_pts = pts.size2();
65  if (!nb_gauss_pts)
67 
68  if (pts.size1() < 3)
69  SETERRQ(
70  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
71  "Wrong dimension of pts, should be at least 3 rows with coordinates");
72 
73  const FieldApproximationBase base = cTx->bAse;
74  EntitiesFieldData &data = cTx->dAta;
75 
76  if (cTx->copyNodeBase == LASTBASE)
77  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
78  "It is assumed that base for vertices is calculated");
79  else
80  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(base) =
81  data.dataOnEntities[MBVERTEX][0].getNSharedPtr(cTx->copyNodeBase);
82 
83  data.dataOnEntities[MBVERTEX][0].getN(base).resize(nb_gauss_pts, 6, false);
84  data.dataOnEntities[MBVERTEX][0].getDiffN(base).resize(nb_gauss_pts, 12,
85  false);
86  if (data.dataOnEntities[MBVERTEX][0].getN(base).size1() !=
87  (unsigned int)nb_gauss_pts)
88  SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
89  "Base functions or nodes has wrong number of integration points "
90  "for base %s",
92 
93  if (cTx->gaussPtsTrianglesOnly.size2() *
94  cTx->gaussPtsThroughThickness.size2() !=
95  pts.size2())
96  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
97 
98  switch (cTx->sPace) {
99  case H1:
102  CHKERR getValueH1(pts);
103  break;
104  case HDIV:
105  CHKERR getValueHdiv(pts);
106  break;
107  case HCURL:
108  CHKERR getValueHcurl(pts);
109  break;
110  case L2:
111  CHKERR getValueL2(pts);
112  break;
113  default:
114  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
115  }
116 
118 }
119 
122  const FieldApproximationBase base = cTx->bAse;
125  boost::shared_ptr<BaseFunctionCtx>(new FlatPrismPolynomialBaseCtx(
126  cTx->dataTrianglesOnly, cTx->mOab, cTx->fePtr, H1, base, NOBASE)));
128 }
129 
132 
133  // EntitiesFieldData& data = cTx->dAta;
134  const FieldApproximationBase base = cTx->bAse;
135  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
136  double *diffL, const int dim) =
138 
139  int nb_gauss_pts_through_thickness = cTx->gaussPtsThroughThickness.size2();
140  // Calculate Legendre approx. on edges
141  for (unsigned int ee = 3; ee <= 5; ee++) {
142  auto &ent_data = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee];
143  int sense = ent_data.getSense();
144  int order = ent_data.getOrder();
145  int nb_dofs = NBEDGE_H1(order);
146  ent_data.getN(base).resize(nb_gauss_pts_through_thickness, nb_dofs, false);
147  ent_data.getDiffN(base).resize(nb_gauss_pts_through_thickness, nb_dofs,
148  false);
149 
150  if (nb_dofs > 0) {
151  double diff_s = 1.; // s = s(xi), ds/dxi = 2., because change of basis
152  for (int gg = 0; gg < nb_gauss_pts_through_thickness; gg++) {
153  double s =
154  2 * cTx->gaussPtsThroughThickness(0, gg) - 1; // makes form -1..1
155  if (sense == -1) {
156  s *= -2;
157  diff_s *= -2;
158  }
159  // calculate Legendre polynomials at integration points on edges
160  // thorough thickness
161  CHKERR base_polynomials(order - 2, s, &diff_s,
162  &ent_data.getN(base)(gg, 0),
163  &ent_data.getDiffN(base)(gg, 0), 1);
164  }
165  }
166  }
167 
169 }
170 
173 
174  EntitiesFieldData &data = cTx->dAta;
175  const FieldApproximationBase base = cTx->bAse;
176 
177  PetscErrorCode (*base_polynomials)(int p, double s, double *diff_s, double *L,
178  double *diffL, const int dim) =
180 
181  int nb_gauss_pts = pts.size2();
182  int nb_gauss_pts_on_faces = cTx->gaussPtsTrianglesOnly.size2();
183  int nb_gauss_pts_through_thickness = cTx->gaussPtsThroughThickness.size2();
184 
185  // nodes
186  // linear for xi, eta and zeta
187  auto &vert_dat = data.dataOnEntities[MBVERTEX][0];
188  vert_dat.getN(base).resize(nb_gauss_pts, 6, false);
189  vert_dat.getDiffN(base).resize(nb_gauss_pts, 18);
190  noalias(vert_dat.getN(base)) = data.dataOnEntities[MBVERTEX][0].getN(NOBASE);
191  noalias(vert_dat.getDiffN(base)) =
192  data.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
193 
194  auto &vert_n = vert_dat.getN(base);
195  auto &vert_diff_n = vert_dat.getDiffN(base);
196 
197  constexpr int prism_edge_map[9][2] = {{0, 1}, {1, 2}, {2, 0}, {0, 3}, {1, 4},
198  {2, 5}, {3, 4}, {4, 5}, {5, 3}};
199 
200  auto edge_through_thickness = [&](const int ee) {
202  // through thickness ho approximation
203  // linear xi,eta, ho terms for zeta
204 
205  auto &thickness_ent = cTx->dataTroughThickness.dataOnEntities[MBEDGE][ee];
206  auto &prism_ent = data.dataOnEntities[MBEDGE][ee];
207 
208  int order = thickness_ent.getOrder();
209  int nb_dofs = NBEDGE_H1(order);
210  if ((unsigned int)nb_dofs != thickness_ent.getN(base).size2())
211  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
212  "nb_dofs != nb_dofs %d != %d", nb_dofs,
213  thickness_ent.getN(base).size2());
214  prism_ent.getN(base).resize(nb_gauss_pts, nb_dofs, false);
215  prism_ent.getDiffN(base).resize(nb_gauss_pts, 3 * nb_dofs, false);
216  if (nb_dofs == 0)
218 
219 
220  int gg = 0;
221  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
222 
223  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
224  double extrude = vert_n(gg, prism_edge_map[ee][0]) +
225  vert_n(gg, prism_edge_map[ee][1]);
226 
227  double diff_extrude[3];
228  for (auto d : {0, 1, 2})
229  diff_extrude[d] = vert_diff_n(gg, 3 * prism_edge_map[ee][0] + d) +
230  vert_diff_n(gg, 3 * prism_edge_map[ee][1] + d);
231 
232  double n0 = vert_n(gg, 0) + vert_n(gg, 1) + vert_n(gg, 2);
233  double n1 = vert_n(gg, 3) + vert_n(gg, 4) + vert_n(gg, 5);
234 
235  // Calculate base through thickness directly from shape fuctions. I
236  // leave that bit of the code, could be useful for integration on the
237  // skeleton.
238 
239  // double l[order + 1], diff_l[3 * (order + 1)];
240  // double ksi = n1 - n0;
241  // double diff_ksi[3];
242  // for (auto d : {0, 1, 2})
243  // diff_ksi[d] =
244  // vert_diff_n(gg, 3 * 3 + d) + vert_diff_n(gg, 3 * 4 + d) +
245  // vert_diff_n(gg, 3 * 5 + d) - vert_diff_n(gg, 3 * 0 + d) -
246  // vert_diff_n(gg, 3 * 1 + d) - vert_diff_n(gg, 3 * 2 + d);
247  // if(sense == -1) {
248  // ksi *= -1;
249  // for (auto d : {0, 1, 2})
250  // diff_ksi[d] *= -1;
251  // }
252  // CHKERR base_polynomials(order - 2, ksi, diff_ksi, l, diff_l, 3);
253 
254  double *l = &(thickness_ent.getN(base)(ggt, 0));
255  double *diff_l_1d = &(thickness_ent.getDiffN(base)(ggt, 0));
256 
257  double bubble = n0 * n1;
258  double diff_bubble[3];
259  for (auto d : {0, 1, 2})
260  diff_bubble[d] =
261  n0 * (vert_diff_n(gg, 3 * 3 + d) + vert_diff_n(gg, 3 * 4 + d) +
262  vert_diff_n(gg, 3 * 5 + d)) +
263 
264  n1 * (vert_diff_n(gg, 3 * 0 + d) + vert_diff_n(gg, 3 * 1 + d) +
265  vert_diff_n(gg, 3 * 2 + d));
266 
267  for (int dd = 0; dd != nb_dofs; dd++) {
268 
269  prism_ent.getN(base)(gg, dd) = extrude * bubble * l[dd];
270  for (auto d : {0, 1, 2})
271  prism_ent.getDiffN(base)(gg, 3 * dd + d) =
272  diff_extrude[d] * bubble * l[dd] +
273  extrude * diff_bubble[d] * l[dd];
274 
275  prism_ent.getDiffN(base)(gg, 3 * dd + 2) +=
276  extrude * bubble * 2 * diff_l_1d[dd];
277 
278  }
279  }
280  }
282  };
283 
284  auto edge_on_the_triangle = [&](const int ee) {
286  // on triangles ho approximation
287  // ho terms on edges, linear zeta
288  int nb_dofs =
289  cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getN(base).size2();
290  data.dataOnEntities[MBEDGE][ee].getN(base).resize(nb_gauss_pts, nb_dofs,
291  false);
292  data.dataOnEntities[MBEDGE][ee].getDiffN(base).resize(nb_gauss_pts,
293  3 * nb_dofs, false);
294  for (int dd = 0; dd < nb_dofs; dd++) {
295  int gg = 0;
296  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
297  double tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getN(
298  base)(ggf, dd);
299  double dksi_tri_n =
300  cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getDiffN(base)(
301  ggf, 2 * dd + 0);
302  double deta_tri_n =
303  cTx->dataTrianglesOnly.dataOnEntities[MBEDGE][ee].getDiffN(base)(
304  ggf, 2 * dd + 1);
305  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
306  double zeta = cTx->gaussPtsThroughThickness(0, ggt);
307  double dzeta, edge_shape;
308  if (ee < 3) {
309  dzeta = diffN_MBEDGE0;
310  edge_shape = N_MBEDGE0(zeta);
311  } else {
312  dzeta = diffN_MBEDGE1;
313  edge_shape = N_MBEDGE1(zeta);
314  }
315  data.dataOnEntities[MBEDGE][ee].getN(base)(gg, dd) =
316  tri_n * edge_shape;
317  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 0) =
318  dksi_tri_n * edge_shape;
319  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 1) =
320  deta_tri_n * edge_shape;
321  data.dataOnEntities[MBEDGE][ee].getDiffN(base)(gg, 3 * dd + 2) =
322  tri_n * dzeta;
323  }
324  }
325  }
327  };
328 
329  auto trinagle_through_thickness = [&](const int ff) {
331  int nb_dofs;
332  nb_dofs =
333  cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getN(base).size2();
334  data.dataOnEntities[MBTRI][ff].getN(base).resize(nb_gauss_pts, nb_dofs,
335  false);
336  data.dataOnEntities[MBTRI][ff].getDiffN(base).resize(nb_gauss_pts,
337  3 * nb_dofs, false);
338  for (int dd = 0; dd < nb_dofs; dd++) {
339  int gg = 0;
340  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
341  double tri_n = cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getN(
342  base)(ggf, dd);
343  double dksi_tri_n[2];
344  for (int kk = 0; kk < 2; kk++) {
345  dksi_tri_n[kk] =
346  cTx->dataTrianglesOnly.dataOnEntities[MBTRI][ff].getDiffN(base)(
347  ggf, 2 * dd + kk);
348  }
349  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
350  double zeta = cTx->gaussPtsThroughThickness(0, ggt);
351  double dzeta, edge_shape;
352  if (ff == 3) {
353  dzeta = diffN_MBEDGE0;
354  edge_shape = N_MBEDGE0(zeta);
355  } else {
356  dzeta = diffN_MBEDGE1;
357  edge_shape = N_MBEDGE1(zeta);
358  }
359  data.dataOnEntities[MBTRI][ff].getN(base)(gg, dd) =
360  tri_n * edge_shape;
361  for (auto kk : {0, 1})
362  data.dataOnEntities[MBTRI][ff].getDiffN(base)(gg, 3 * dd + kk) =
363  dksi_tri_n[kk] * edge_shape;
364  data.dataOnEntities[MBTRI][ff].getDiffN(base)(gg, 3 * dd + 2) =
365  tri_n * dzeta;
366  }
367  }
368  }
370  };
371 
372  auto quads_base = [&]() {
374  int quads_nodes[3 * 4];
375  int quad_order[3] = {0, 0, 0};
376  double *quad_n[3], *diff_quad_n[3];
377  SideNumber_multiIndex &side_table =
379  SideNumber_multiIndex::nth_index<1>::type::iterator siit;
380  siit = side_table.get<1>().lower_bound(boost::make_tuple(MBQUAD, 0));
381  SideNumber_multiIndex::nth_index<1>::type::iterator hi_siit;
382  hi_siit = side_table.get<1>().upper_bound(boost::make_tuple(MBQUAD, 3));
383  if (std::distance(siit, hi_siit) != 3)
384  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
385  "Expected three quads on prisms");
386  EntityHandle ent = cTx->fePtr->getEnt();
387  int num_nodes;
388  const EntityHandle *conn;
389  CHKERR cTx->mOab.get_connectivity(ent, conn, num_nodes, true);
390  for (; siit != hi_siit; ++siit) {
391  int side = siit->get()->side_number;
392  if(side < 0 && side > 2)
393  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
394  "Side in range 0,1,2 expected");
395  int num_nodes_quad;
396  const EntityHandle *conn_quad;
397  EntityHandle quad = siit->get()->ent;
398  CHKERR cTx->mOab.get_connectivity(quad, conn_quad, num_nodes_quad, true);
399  for (int nn = 0; nn < num_nodes_quad; nn++) {
400  quads_nodes[4 * side + nn] =
401  std::distance(conn, std::find(conn, conn + 6, conn_quad[nn]));
402  }
403  int order = data.dataOnEntities[MBQUAD][side].getOrder();
404  quad_order[side] = order;
405  data.dataOnEntities[MBQUAD][side].getN(base).resize(
406  nb_gauss_pts, NBFACEQUAD_H1(order), false);
407  data.dataOnEntities[MBQUAD][side].getDiffN(base).resize(
408  nb_gauss_pts, 3 * NBFACEQUAD_H1(order), false);
409  if (data.dataOnEntities[MBQUAD][side].getN(base).size2() > 0) {
410  quad_n[side] =
411  &*data.dataOnEntities[MBQUAD][side].getN(base).data().begin();
412  diff_quad_n[side] =
413  &*data.dataOnEntities[MBQUAD][side].getDiffN(base).data().begin();
414  } else {
415  quad_n[side] = NULL;
416  diff_quad_n[side] = NULL;
417  }
418  }
419  if (quad_order[0] > 0 || quad_order[1] > 0 || quad_order[2] > 0) {
420  double *vertex_n =
421  &*data.dataOnEntities[MBVERTEX][0].getN(base).data().begin();
422  double *diff_vertex_n =
423  &*data.dataOnEntities[MBVERTEX][0].getDiffN(base).data().begin();
424  CHKERR H1_QuadShapeFunctions_MBPRISM(quads_nodes, quad_order, vertex_n,
425  diff_vertex_n, quad_n, diff_quad_n,
426  nb_gauss_pts, base_polynomials);
427  }
429  };
430 
431  auto prim_base = [&]() {
433  int order = data.dataOnEntities[MBPRISM][0].getOrder();
434  double *vertex_n = &data.dataOnEntities[MBVERTEX][0].getN(base)(0, 0);
435  double *diff_vertex_n =
436  &data.dataOnEntities[MBVERTEX][0].getDiffN(base)(0, 0);
437  data.dataOnEntities[MBPRISM][0].getN(base).resize(
438  nb_gauss_pts, NBVOLUMEPRISM_H1(order), false);
439  data.dataOnEntities[MBPRISM][0].getDiffN(base).resize(
440  nb_gauss_pts, 3 * NBVOLUMEPRISM_H1(order), false);
441  if (NBVOLUMEPRISM_H1(order) > 0) {
443  order, vertex_n, diff_vertex_n,
444  &data.dataOnEntities[MBPRISM][0].getN(base)(0, 0),
445  &data.dataOnEntities[MBPRISM][0].getDiffN(base)(0, 0), nb_gauss_pts,
446  base_polynomials);
447  }
449  };
450 
451  // edges on triangles
452  int ee = 0;
453  for (; ee <= 2; ++ee)
454  CHKERR edge_on_the_triangle(ee);
455  for (; ee <= 5; ++ee)
456  CHKERR edge_through_thickness(ee);
457  for (; ee <= 8; ++ee)
458  CHKERR edge_on_the_triangle(ee);
459 
460  // triangles
461  // ho on triangles, linear zeta
462  for (int ff = 3; ff <= 4; ++ff)
463  CHKERR trinagle_through_thickness(ff);
464 
465  // quads
466  // higher order edges and zeta
467  CHKERR quads_base();
468 
469  // prism
470  CHKERR prim_base();
471 
473 }
474 
477  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
479 }
480 
482  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
483 }
484 
486  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
487 }
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
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:51
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:130
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::Types::MatrixDouble
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
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:481
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:485
MoFEM::FatPrismPolynomialBase::getValueL2
MoFEMErrorCode getValueL2(MatrixDouble &pts)
Definition: FatPrismPolynomialBase.cpp:475
order
constexpr int order
Definition: dg_projection.cpp:18
MoFEM::Interface
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1975
FieldSpace
FieldSpace
approximation spaces
Definition: definitions.h:82
MoFEM::EntPolynomialBaseCtx::copyNodeBase
const FieldApproximationBase copyNodeBase
Definition: EntPolynomialBaseCtx.hpp:40
MoFEM::FatPrismPolynomialBase::FatPrismPolynomialBase
FatPrismPolynomialBase()
Definition: FatPrismPolynomialBase.cpp:48
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
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:47
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:171
MoFEM::L
VectorDouble L
Definition: Projection10NodeCoordsOnField.cpp:124
MoFEM::FatPrismPolynomialBaseCtx::~FatPrismPolynomialBaseCtx
~FatPrismPolynomialBaseCtx()
Definition: FatPrismPolynomialBase.cpp:38
MoFEM::FatPrismPolynomialBase::getValueH1TrianglesOnly
MoFEMErrorCode getValueH1TrianglesOnly()
Definition: FatPrismPolynomialBase.cpp:120
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:38
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:41
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:56
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:440
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
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
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:346
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:36