v0.9.0
FatPrismElementForcesAndSourcesCore.cpp
Go to the documentation of this file.
1 /** \file FatPrismElementForcesAndSourcesCore.cpp
2 
3 \brief Implementation of fat prism element
4 
5 */
6 
7 /* This file is part of MoFEM.
8  * MoFEM is free software: you can redistribute it and/or modify it under
9  * the terms of the GNU Lesser General Public License as published by the
10  * Free Software Foundation, either version 3 of the License, or (at your
11  * option) any later version.
12  *
13  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
14  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
15  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16  * License for more details.
17  *
18  * You should have received a copy of the GNU Lesser General Public
19  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
20 
21 namespace MoFEM {
22 
24  Interface &m_field)
25  : VolumeElementForcesAndSourcesCore(m_field, MBPRISM),
26  dataH1TrianglesOnly(MBPRISM), dataH1TroughThickness(MBPRISM),
27  opHOCoordsAndNormals(hoCoordsAtGaussPtsF3, nOrmals_at_GaussPtF3,
28  tAngent1_at_GaussPtF3, tAngent2_at_GaussPtF3,
29  hoCoordsAtGaussPtsF4, nOrmals_at_GaussPtF4,
30  tAngent1_at_GaussPtF4, tAngent2_at_GaussPtF4) {}
31 
34 
35  if (numeredEntFiniteElementPtr->getEntType() != MBPRISM)
38 
39  auto get_fe_coordinates = [&]() {
42  int num_nodes;
43  const EntityHandle *conn;
44  CHKERR mField.get_moab().get_connectivity(ent, conn, num_nodes, true);
45  coords.resize(num_nodes * 3, false);
46  CHKERR mField.get_moab().get_coords(conn, num_nodes,
47  &*coords.data().begin());
49  };
50 
51  auto calculate_area_of_triangles = [&] {
53  normal.resize(6, false);
56  aRea[0] = cblas_dnrm2(3, &normal[0], 1) * 0.5;
57  aRea[1] = cblas_dnrm2(3, &normal[3], 1) * 0.5;
59  };
60 
61  CHKERR get_fe_coordinates();
62  CHKERR calculate_area_of_triangles();
63 
67 
68  auto get_h1_base_data = [&](auto &dataH1) {
70  CHKERR getEntitySense<MBEDGE>(dataH1);
71  CHKERR getEntitySense<MBTRI>(dataH1);
72  CHKERR getEntitySense<MBQUAD>(dataH1);
73  CHKERR getEntityDataOrder<MBEDGE>(dataH1, H1);
74  CHKERR getEntityDataOrder<MBTRI>(dataH1, H1);
75  CHKERR getEntityDataOrder<MBQUAD>(dataH1, H1);
76  CHKERR getEntityDataOrder<MBPRISM>(dataH1, H1);
77  // Triangles only
78  CHKERR getEntitySense<MBEDGE>(dataH1TrianglesOnly);
79  CHKERR getEntitySense<MBTRI>(dataH1TrianglesOnly);
80  CHKERR getEntityDataOrder<MBEDGE>(dataH1TrianglesOnly, H1);
81  CHKERR getEntityDataOrder<MBTRI>(dataH1TrianglesOnly, H1);
82  // Through thickness
83  CHKERR getEntitySense<MBEDGE>(dataH1TroughThickness);
84  CHKERR getEntityDataOrder<MBEDGE>(dataH1TroughThickness, H1);
86  };
87 
88  // H1
89  if ((dataH1.spacesOnEntities[MBEDGE]).test(H1))
90  CHKERR get_h1_base_data(dataH1);
91 
92  // Hdiv
93  if ((dataH1.spacesOnEntities[MBTRI]).test(HDIV))
94  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented yet");
95 
96  // Hcurl
97  if ((dataH1.spacesOnEntities[MBEDGE]).test(HCURL))
98  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented yet");
99 
100  // L2
101  if ((dataH1.spacesOnEntities[MBPRISM]).test(L2))
102  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented yet");
103 
104  // get approx. on triangles, i.e. faces 3 and 4
105  auto set_gauss_points_on_triangle = [&](int &nb_gauss_pts_on_faces) {
107  int order_triangles_only = 1;
108  int valid_edges[] = {1, 1, 1, 0, 0, 0, 1, 1, 1};
109  for (unsigned int ee = 0; ee < 9; ee++) {
110  if (!valid_edges[ee])
111  continue;
112  order_triangles_only = std::max(
113  order_triangles_only,
114  dataH1TrianglesOnly.dataOnEntities[MBEDGE][ee].getDataOrder());
115  }
116  for (unsigned int ff = 3; ff <= 4; ff++) {
117  order_triangles_only = std::max(
118  order_triangles_only,
119  dataH1TrianglesOnly.dataOnEntities[MBTRI][ff].getDataOrder());
120  }
121  for (unsigned int qq = 0; qq < 3; qq++) {
122  order_triangles_only = std::max(
123  order_triangles_only,
124  dataH1TroughThickness.dataOnEntities[MBQUAD][qq].getDataOrder());
125  }
126  order_triangles_only = std::max(
127  order_triangles_only,
128  dataH1TroughThickness.dataOnEntities[MBPRISM][0].getDataOrder());
129 
130  // integration pts on the triangles surfaces
131  nb_gauss_pts_on_faces = 0;
132  int rule = getRuleTrianglesOnly(order_triangles_only);
133  if (rule >= 0) {
134  if (rule < QUAD_2D_TABLE_SIZE) {
135  if (QUAD_2D_TABLE[rule]->dim != 2) {
136  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
137  }
138  if (QUAD_2D_TABLE[rule]->order < rule) {
139  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
140  "wrong order %d != %d", QUAD_2D_TABLE[rule]->order, rule);
141  }
142  nb_gauss_pts_on_faces = QUAD_2D_TABLE[rule]->npoints;
143  gaussPtsTrianglesOnly.resize(3, nb_gauss_pts_on_faces, false);
144  cblas_dcopy(nb_gauss_pts_on_faces, &QUAD_2D_TABLE[rule]->points[1], 3,
145  &gaussPtsTrianglesOnly(0, 0), 1);
146  cblas_dcopy(nb_gauss_pts_on_faces, &QUAD_2D_TABLE[rule]->points[2], 3,
147  &gaussPtsTrianglesOnly(1, 0), 1);
148  cblas_dcopy(nb_gauss_pts_on_faces, QUAD_2D_TABLE[rule]->weights, 1,
149  &gaussPtsTrianglesOnly(2, 0), 1);
150  dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(
151  nb_gauss_pts_on_faces, 3, false);
152  double *shape_ptr = &*dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0]
153  .getN(NOBASE)
154  .data()
155  .begin();
156  cblas_dcopy(3 * nb_gauss_pts_on_faces, QUAD_2D_TABLE[rule]->points, 1,
157  shape_ptr, 1);
158  dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE).resize(
159  1, 6, false);
160  std::copy(Tools::diffShapeFunMBTRI.begin(),
162  &*dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0]
163  .getDiffN(NOBASE)
164  .data()
165  .begin());
166  } else
167  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
168  "rule > quadrature order %d < %d", rule, QUAD_2D_TABLE_SIZE);
169 
170  } else {
171  CHKERR setGaussPtsTrianglesOnly(order_triangles_only);
172  nb_gauss_pts_on_faces = gaussPtsTrianglesOnly.size2();
173  if (nb_gauss_pts_on_faces == 0)
175  dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(
176  nb_gauss_pts_on_faces, 3, false);
177  if (nb_gauss_pts_on_faces) {
179  .getN(NOBASE)
180  .data()
181  .begin(),
182  &gaussPtsTrianglesOnly(0, 0),
183  &gaussPtsTrianglesOnly(1, 0), nb_gauss_pts_on_faces);
184  dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE).resize(
185  1, 6, false);
186  std::copy(Tools::diffShapeFunMBTRI.begin(),
188  &*dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0]
189  .getDiffN(NOBASE)
190  .data()
191  .begin());
192  }
193  }
195  };
196 
197  // approx. trough prism thickness
198  auto set_gauss_points_through_thickness =
199  [&](int &nb_gauss_pts_through_thickness) {
201  nb_gauss_pts_through_thickness = 0;
202  int order_thickness = 1;
203  for (unsigned int ee = 3; ee <= 5; ee++) {
204  order_thickness = std::max(
205  order_thickness,
206  dataH1TroughThickness.dataOnEntities[MBEDGE][ee].getDataOrder());
207  }
208  for (unsigned int qq = 0; qq < 3; qq++) {
209  order_thickness = std::max(
210  order_thickness,
211  dataH1TroughThickness.dataOnEntities[MBQUAD][qq].getDataOrder());
212  }
213  order_thickness = std::max(
214  order_thickness,
215  dataH1TroughThickness.dataOnEntities[MBPRISM][0].getDataOrder());
216  // integration points
217  int rule = getRuleThroughThickness(order_thickness);
218  if (rule >= 0) {
219  if (rule < QUAD_1D_TABLE_SIZE) {
220  if (QUAD_1D_TABLE[rule]->dim != 1) {
221  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
222  "wrong dimension");
223  }
224  if (QUAD_1D_TABLE[rule]->order < rule) {
225  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
226  "wrong order %d != %d", QUAD_1D_TABLE[rule]->order,
227  rule);
228  }
229  nb_gauss_pts_through_thickness = QUAD_1D_TABLE[rule]->npoints;
230  gaussPtsThroughThickness.resize(2, nb_gauss_pts_through_thickness,
231  false);
232  cblas_dcopy(nb_gauss_pts_through_thickness,
233  &QUAD_1D_TABLE[rule]->points[1], 2,
234  &gaussPtsThroughThickness(0, 0), 1);
235  cblas_dcopy(nb_gauss_pts_through_thickness,
236  QUAD_1D_TABLE[rule]->weights, 1,
237  &gaussPtsThroughThickness(1, 0), 1);
238  } else {
239  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
240  "rule > quadrature order %d < %d", rule,
242  nb_gauss_pts_through_thickness = 0;
243  }
244  } else {
245  CHKERR setGaussPtsThroughThickness(order_thickness);
246  nb_gauss_pts_through_thickness = gaussPtsThroughThickness.size2();
247  }
249  };
250 
251  // Generate integration pts.
252  auto set_gauss_points_in_volume = [&](int nb_gauss_pts_on_faces,
253  int nb_gauss_pts_through_thickness,
254  int &nb_gauss_pts) {
256  nb_gauss_pts = nb_gauss_pts_on_faces * nb_gauss_pts_through_thickness;
257  gaussPts.resize(4, nb_gauss_pts, false);
258  int gg = 0;
259  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
260  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
261  gaussPts(0, gg) = gaussPtsTrianglesOnly(0, ggf);
262  gaussPts(1, gg) = gaussPtsTrianglesOnly(1, ggf);
263  gaussPts(2, gg) = gaussPtsThroughThickness(0, ggt);
264  gaussPts(3, gg) =
266  }
267  }
269  };
270 
271  int nb_gauss_pts, nb_gauss_pts_through_thickness, nb_gauss_pts_on_faces;
272  CHKERR set_gauss_points_on_triangle(nb_gauss_pts_on_faces);
273  if (!nb_gauss_pts_on_faces)
275  CHKERR set_gauss_points_through_thickness(nb_gauss_pts_through_thickness);
276  if (!nb_gauss_pts_through_thickness)
278  CHKERR set_gauss_points_in_volume(
279  nb_gauss_pts_on_faces, nb_gauss_pts_through_thickness, nb_gauss_pts);
280 
281  auto calc_coordinates_at_triangles = [&]() {
283  coordsAtGaussPtsTrianglesOnly.resize(nb_gauss_pts_on_faces, 6, false);
284  for (int gg = 0; gg < nb_gauss_pts_on_faces; gg++) {
285  for (int dd = 0; dd < 3; dd++) {
287  cblas_ddot(3,
288  &dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0].getN(
289  NOBASE)(gg, 0),
290  1, &coords[dd], 3);
292  cblas_ddot(3,
293  &dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0].getN(
294  NOBASE)(gg, 0),
295  1, &coords[9 + dd], 3);
296  }
297  }
299  };
300 
301  auto calc_vertex_base_on_prism = [&]() {
303  // Calculate "nobase" base functions on prism, this is cartesian product
304  // of base functions on triangles with base functions through thickness
305  // FIXME: This could be effectively implemented with tensors
306  dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nb_gauss_pts, 6,
307  false);
308  dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE).resize(nb_gauss_pts, 18,
309  false);
310  for (int dd = 0; dd != 6; dd++) {
311  int gg = 0;
312  for (int ggf = 0; ggf < nb_gauss_pts_on_faces; ggf++) {
313  int ddd = dd > 2 ? dd - 3 : dd;
314  double tri_n = dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0].getN(
315  NOBASE)(ggf, ddd);
316  double dksi_tri_n =
317  dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE)(
318  0, 2 * ddd + 0);
319  double deta_tri_n =
320  dataH1TrianglesOnly.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE)(
321  0, 2 * ddd + 1);
322  for (int ggt = 0; ggt < nb_gauss_pts_through_thickness; ggt++, gg++) {
323  double zeta = gaussPtsThroughThickness(0, ggt);
324  double dzeta, edge_shape;
325  if (dd < 3) {
326  dzeta = diffN_MBEDGE0;
327  edge_shape = N_MBEDGE0(zeta);
328  } else {
329  dzeta = diffN_MBEDGE1;
330  edge_shape = N_MBEDGE1(zeta);
331  }
332  dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE)(gg, dd) =
333  tri_n * edge_shape;
334  dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE)(gg, 3 * dd + 0) =
335  dksi_tri_n * edge_shape;
336  dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE)(gg, 3 * dd + 1) =
337  deta_tri_n * edge_shape;
338  dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE)(gg, 3 * dd + 2) =
339  tri_n * dzeta;
340  }
341  }
342  }
344  };
345 
346  auto calc_base_on_prism = [&]() {
348  // Calculate base functions on prism
349  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
350  if (dataH1.bAse.test(b)) {
351  switch (static_cast<FieldApproximationBase>(b)) {
354  if (dataH1.spacesOnEntities[MBVERTEX].test(H1)) {
356  gaussPts,
357  boost::shared_ptr<BaseFunctionCtx>(
362  static_cast<FieldApproximationBase>(b), NOBASE)));
363  }
364  if (dataH1.spacesOnEntities[MBTRI].test(HDIV)) {
365  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
366  "Not yet implemented");
367  }
368  if (dataH1.spacesOnEntities[MBEDGE].test(HCURL)) {
369  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
370  "Not yet implemented");
371  }
372  if (dataH1.spacesOnEntities[MBTET].test(L2)) {
373  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
374  "Not yet implemented");
375  }
376  break;
377  default:
378  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
379  "Not yet implemented");
380  }
381  }
382  }
384  };
385 
386  auto calc_coordinate_on_prism = [&]() {
388  /// Calculate coordinates at integration points
389  coordsAtGaussPts.resize(nb_gauss_pts, 3, false);
390  for (int gg = 0; gg < nb_gauss_pts; gg++) {
391  for (int dd = 0; dd < 3; dd++) {
393  6, &dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE)(gg, 0), 1,
394  &coords[dd], 3);
395  }
396  }
398  };
399 
400  auto calc_ho_triangle_face_normals = [&]() {
402  // Calculate ho-geometry tangents and normals
403  if (dataPtr->get<FieldName_mi_tag>().find(meshPositionsFieldName) !=
404  dataPtr->get<FieldName_mi_tag>().end()) {
405  hoCoordsAtGaussPtsF3.resize(nb_gauss_pts_on_faces, 3, false);
406  nOrmals_at_GaussPtF3.resize(nb_gauss_pts_on_faces, 3, false);
407  tAngent1_at_GaussPtF3.resize(nb_gauss_pts_on_faces, 3, false);
408  tAngent2_at_GaussPtF3.resize(nb_gauss_pts_on_faces, 3, false);
409  hoCoordsAtGaussPtsF4.resize(nb_gauss_pts_on_faces, 3, false);
410  nOrmals_at_GaussPtF4.resize(nb_gauss_pts_on_faces, 3, false);
411  tAngent1_at_GaussPtF4.resize(nb_gauss_pts_on_faces, 3, false);
412  tAngent2_at_GaussPtF4.resize(nb_gauss_pts_on_faces, 3, false);
415  MBEDGE);
417  MBEDGE);
420  } else {
421  hoCoordsAtGaussPtsF3.resize(0, 0, false);
422  nOrmals_at_GaussPtF3.resize(0, 0, false);
423  tAngent1_at_GaussPtF3.resize(0, 0, false);
424  tAngent2_at_GaussPtF3.resize(0, 0, false);
425  hoCoordsAtGaussPtsF4.resize(0, 0, false);
426  nOrmals_at_GaussPtF4.resize(0, 0, false);
427  tAngent1_at_GaussPtF4.resize(0, 0, false);
428  tAngent2_at_GaussPtF4.resize(0, 0, false);
429  }
431  };
432 
433  CHKERR calc_coordinates_at_triangles();
434  CHKERR calc_vertex_base_on_prism();
435  CHKERR calc_base_on_prism();
436  CHKERR calc_coordinate_on_prism();
437  CHKERR calc_ho_triangle_face_normals();
438 
439  // Iterate over operators
441 
443 }
444 
445 } // namespace MoFEM
Class used to pass element data to calculate base functions on fat prism.
MatrixDouble gaussPts
Matrix of integration points.
field with continuous normal traction
Definition: definitions.h:173
static QUAD *const QUAD_1D_TABLE[]
Definition: quad.h:164
virtual moab::Interface & get_moab()=0
MoFEMErrorCode createDataOnElement()
Create a entity data on element object.
DataForcesAndSourcesCore & dataH1
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
int npoints
Definition: quad.h:29
void cblas_dcopy(const int N, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_dcopy.c:11
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:477
#define N_MBEDGE1(x)
edge shape function
Definition: fem_tools.h:80
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
MoFEMErrorCode getValue(MatrixDouble &pts, boost::shared_ptr< BaseFunctionCtx > ctx_ptr)
boost::shared_ptr< const FEDofEntity_multiIndex > dataPtr
Pointer to finite element data dofs.
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
double cblas_dnrm2(const int N, const double *X, const int incX)
Definition: cblas_dnrm2.c:12
#define diffN_MBEDGE1
derivative of edge shape function
Definition: fem_tools.h:82
double cblas_ddot(const int N, const double *X, const int incX, const double *Y, const int incY)
Definition: cblas_ddot.c:12
#define QUAD_2D_TABLE_SIZE
Definition: quad.h:174
static QUAD *const QUAD_2D_TABLE[]
Definition: quad.h:175
MoFEMErrorCode loopOverOperators()
Iterate user data operators.
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
MoFEMErrorCode getNodesFieldData(DataForcesAndSourcesCore &data, const std::string &field_name) const
Get data on nodes.
#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
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:146
static MoFEMErrorCode getTriNormal(const double *coords, double *normal)
Get the Tri Normal objectGet triangle normal.
Definition: Tools.cpp:258
#define diffN_MBEDGE0
derivative of edge shape function
Definition: fem_tools.h:81
field with continuous tangents
Definition: definitions.h:172
Calculate base functions on tetrahedralFIXME: Need moab and mofem finite element structure to work (t...
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition: Tools.hpp:87
virtual MoFEMErrorCode setGaussPtsThroughThickness(int order_thickness)
#define CHKERR
Inline error check.
Definition: definitions.h:596
virtual MoFEMErrorCode opRhs(DataForcesAndSourcesCore &data, const bool do_vertices, const bool do_edges, const bool do_quads, const bool do_tris, const bool do_tets, const bool do_prisms, const bool error_if_no_base=true)
MoFEMErrorCode getEntityFieldData(DataForcesAndSourcesCore &data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition: fem_tools.c:194
MultiIndex Tag for field name.
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
virtual MoFEMErrorCode setGaussPtsTrianglesOnly(int order_triangles_only)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:407
continuous field
Definition: definitions.h:171
MoFEMErrorCode getSpacesAndBaseOnEntities(DataForcesAndSourcesCore &data) const
Get field approximation space and base on entities.
std::bitset< LASTBASE > bAse
bases on element
MoFEMErrorCode operator()()
function is run for every finite element
field with C-1 continuity
Definition: definitions.h:174
#define QUAD_1D_TABLE_SIZE
Definition: quad.h:163