v0.14.0
Loading...
Searching...
No Matches
FaceElementForcesAndSourcesCore.cpp
Go to the documentation of this file.
1/** \file FaceElementForcesAndSourcesCore.cpp
2
3\brief Implementation of face element
4
5*/
6
7namespace MoFEM {
8
10 Interface &m_field)
11 : ForcesAndSourcesCore(m_field),
12 meshPositionsFieldName("MESH_NODE_POSITIONS"), aRea(elementMeasure) {}
13
17
18 auto type = numeredEntFiniteElementPtr->getEntType();
19
23
24 auto get_ftensor_from_vec_3d = [](VectorDouble &v) {
26 &v[2]);
27 };
28
29 auto get_ftensor_n_diff = [&]() {
30 const auto &m = dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
32 &m(0, 1));
33 };
34
35 auto get_ftensor_from_mat_3d = [](MatrixDouble &m) {
37 &m(0, 0), &m(0, 1), &m(0, 2));
38 };
39
40 if (type == MBTRI) {
41
42 const size_t nb_gauss_pts = gaussPts.size2();
43 normalsAtGaussPts.resize(nb_gauss_pts, 3);
44 tangentOneAtGaussPts.resize(nb_gauss_pts, 3);
45 tangentTwoAtGaussPts.resize(nb_gauss_pts, 3);
46
47 auto t_tan1 = get_ftensor_from_mat_3d(tangentOneAtGaussPts);
48 auto t_tan2 = get_ftensor_from_mat_3d(tangentTwoAtGaussPts);
49 auto t_normal = get_ftensor_from_mat_3d(normalsAtGaussPts);
50
51 auto t_n =
54 &tangentOne[2]);
56 &tangentTwo[2]);
57
58 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
59 t_normal(i) = t_n(i);
60 t_tan1(i) = t_t1(i);
61 t_tan2(i) = t_t2(i);
62 ++t_tan1;
63 ++t_tan2;
64 ++t_normal;
65 }
66
67 } else if (type == MBQUAD) {
68
70 CHKERR mField.get_moab().get_connectivity(ent, conn, num_nodes, true);
71 coords.resize(num_nodes * 3, false);
72 CHKERR mField.get_moab().get_coords(conn, num_nodes,
73 &*coords.data().begin());
74
75 const size_t nb_gauss_pts = gaussPts.size2();
76 normalsAtGaussPts.resize(nb_gauss_pts, 3);
77 tangentOneAtGaussPts.resize(nb_gauss_pts, 3);
78 tangentTwoAtGaussPts.resize(nb_gauss_pts, 3);
79 normalsAtGaussPts.clear();
82
83 auto t_t1 = get_ftensor_from_mat_3d(tangentOneAtGaussPts);
84 auto t_t2 = get_ftensor_from_mat_3d(tangentTwoAtGaussPts);
85 auto t_normal = get_ftensor_from_mat_3d(normalsAtGaussPts);
86
89
90 auto t_diff = get_ftensor_n_diff();
91 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
92 auto t_coords = get_ftensor_from_vec_3d(coords);
93 for (int nn = 0; nn != num_nodes; ++nn) {
94 t_t1(i) += t_coords(i) * t_diff(N0);
95 t_t2(i) += t_coords(i) * t_diff(N1);
96 ++t_diff;
97 ++t_coords;
98 }
99 t_normal(j) = FTensor::levi_civita(i, j, k) * t_t1(k) * t_t2(i);
100
101 ++t_t1;
102 ++t_t2;
103 ++t_normal;
104 }
105 } else {
106 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
107 "Element type not implemented");
108 }
109
111}
112
115
117 CHKERR mField.get_moab().get_connectivity(ent, conn, num_nodes, true);
118 coords.resize(num_nodes * 3, false);
119 CHKERR mField.get_moab().get_coords(conn, num_nodes, &*coords.data().begin());
120 nOrmal.resize(3, false);
121 tangentOne.resize(3, false);
122 tangentTwo.resize(3, false);
123
124 auto calc_normal = [&](const double *diff_ptr) {
127 &coords[0], &coords[1], &coords[2]);
129 &nOrmal[0], &nOrmal[1], &nOrmal[2]);
131 &tangentOne[0], &tangentOne[1], &tangentOne[2]);
133 &tangentTwo[0], &tangentTwo[1], &tangentTwo[2]);
135 diff_ptr, &diff_ptr[1]);
136
140
143 t_t1(i) = 0;
144 t_t2(i) = 0;
145
146 for (int nn = 0; nn != num_nodes; ++nn) {
147 t_t1(i) += t_coords(i) * t_diff(N0);
148 t_t2(i) += t_coords(i) * t_diff(N1);
149 ++t_coords;
150 ++t_diff;
151 }
152 t_normal(j) = FTensor::levi_civita(i, j, k) * t_t1(k) * t_t2(i);
153 aRea = sqrt(t_normal(i) * t_normal(i));
155 };
156
157 const double *diff_ptr;
158 switch (numeredEntFiniteElementPtr->getEntType()) {
159 case MBTRI:
160 diff_ptr = Tools::diffShapeFunMBTRI.data();
161 CHKERR calc_normal(diff_ptr);
162 // FIXME: Normal should be divided not the area for triangle!!
163 aRea /= 2;
164 break;
165 case MBQUAD:
166 diff_ptr = Tools::diffShapeFunMBQUADAtCenter.data();
167 CHKERR calc_normal(diff_ptr);
168 break;
169 default:
170 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
171 "Element type not implemented");
172 }
173
175}
176
179 // Set integration points
180 int order_data = getMaxDataOrder();
181 int order_row = getMaxRowOrder();
182 int order_col = getMaxColOrder();
183
184 const auto type = numeredEntFiniteElementPtr->getEntType();
185
186 auto get_rule_by_type = [&]() {
187 switch (type) {
188 case MBQUAD:
189 return getRule(order_row + 1, order_col + 1, order_data + 1);
190 default:
191 return getRule(order_row, order_col, order_data);
192 }
193 };
194
195 const int rule = get_rule_by_type();
196
197 auto set_integration_pts_for_tri = [&]() {
199 if (rule < QUAD_2D_TABLE_SIZE) {
200 if (QUAD_2D_TABLE[rule]->dim != 2) {
201 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
202 }
203 if (QUAD_2D_TABLE[rule]->order < rule) {
204 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
205 "wrong order %d != %d", QUAD_2D_TABLE[rule]->order, rule);
206 }
207 const size_t nb_gauss_pts = QUAD_2D_TABLE[rule]->npoints;
208 gaussPts.resize(3, nb_gauss_pts, false);
209 cblas_dcopy(nb_gauss_pts, &QUAD_2D_TABLE[rule]->points[1], 3,
210 &gaussPts(0, 0), 1);
211 cblas_dcopy(nb_gauss_pts, &QUAD_2D_TABLE[rule]->points[2], 3,
212 &gaussPts(1, 0), 1);
213 cblas_dcopy(nb_gauss_pts, QUAD_2D_TABLE[rule]->weights, 1,
214 &gaussPts(2, 0), 1);
215 dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nb_gauss_pts, 3,
216 false);
217 double *shape_ptr =
218 &*dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin();
219 cblas_dcopy(3 * nb_gauss_pts, QUAD_2D_TABLE[rule]->points, 1, shape_ptr,
220 1);
221 dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE).resize(3, 2, false);
222 std::copy(
224 dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE).data().begin());
225 } else {
226 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
227 "rule > quadrature order %d < %d", rule, QUAD_2D_TABLE_SIZE);
228 }
230 };
231
232 auto calc_base_for_tri = [&]() {
234 const size_t nb_gauss_pts = gaussPts.size2();
235 auto &base = dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE);
236 auto &diff_base = dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
237 base.resize(nb_gauss_pts, 3, false);
238 diff_base.resize(3, 2, false);
239 CHKERR ShapeMBTRI(&*base.data().begin(), &gaussPts(0, 0), &gaussPts(1, 0),
240 nb_gauss_pts);
241 std::copy(
243 dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE).data().begin());
245 };
246
247 auto calc_base_for_quad = [&]() {
249 const size_t nb_gauss_pts = gaussPts.size2();
250 auto &base = dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE);
251 auto &diff_base = dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
252 base.resize(nb_gauss_pts, 4, false);
253 diff_base.resize(nb_gauss_pts, 8, false);
254 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
255 const double ksi = gaussPts(0, gg);
256 const double zeta = gaussPts(1, gg);
257 base(gg, 0) = N_MBQUAD0(ksi, zeta);
258 base(gg, 1) = N_MBQUAD1(ksi, zeta);
259 base(gg, 2) = N_MBQUAD2(ksi, zeta);
260 base(gg, 3) = N_MBQUAD3(ksi, zeta);
261 diff_base(gg, 0) = diffN_MBQUAD0x(zeta);
262 diff_base(gg, 1) = diffN_MBQUAD0y(ksi);
263 diff_base(gg, 2) = diffN_MBQUAD1x(zeta);
264 diff_base(gg, 3) = diffN_MBQUAD1y(ksi);
265 diff_base(gg, 4) = diffN_MBQUAD2x(zeta);
266 diff_base(gg, 5) = diffN_MBQUAD2y(ksi);
267 diff_base(gg, 6) = diffN_MBQUAD3x(zeta);
268 diff_base(gg, 7) = diffN_MBQUAD3y(ksi);
269 }
271 };
272
273 if (rule >= 0) {
274 switch (type) {
275 case MBTRI:
276 CHKERR set_integration_pts_for_tri();
277 break;
278 case MBQUAD:
280 rule);
281 CHKERR calc_base_for_quad();
282 break;
283 default:
284 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
285 "Element type not implemented: %d", type);
286 }
287
288 } else {
289 // If rule is negative, set user defined integration points
290 CHKERR setGaussPts(order_row, order_col, order_data);
291 const size_t nb_gauss_pts = gaussPts.size2();
292 if (nb_gauss_pts) {
293 switch (type) {
294 case MBTRI:
295 CHKERR calc_base_for_tri();
296 break;
297 case MBQUAD:
298 CHKERR calc_base_for_quad();
299 break;
300 default:
301 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
302 "Element type not implemented: %d", type);
303 }
304 }
305 }
307}
308
312 // Get spaces order/base and sense of entities.
313
315
316 // H1
317 if (dataH1.spacesOnEntities[MBEDGE].test(H1)) {
318 CHKERR getEntitySense<MBEDGE>(dataH1);
319 CHKERR getEntityDataOrder<MBEDGE>(dataH1, H1);
320 }
321 if (dataH1.spacesOnEntities[MBTRI].test(H1)) {
322 CHKERR getEntitySense<MBTRI>(dataH1);
323 CHKERR getEntityDataOrder<MBTRI>(dataH1, H1);
324 }
325 if (dataH1.spacesOnEntities[MBQUAD].test(H1)) {
326 CHKERR getEntitySense<MBQUAD>(dataH1);
327 CHKERR getEntityDataOrder<MBQUAD>(dataH1, H1);
328 }
329
330 // Hcurl
331 if (dataH1.spacesOnEntities[MBEDGE].test(HCURL)) {
332 CHKERR getEntitySense<MBEDGE>(dataHcurl);
333 CHKERR getEntityDataOrder<MBEDGE>(dataHcurl, HCURL);
334 dataHcurl.spacesOnEntities[MBEDGE].set(HCURL);
335 }
336 if (dataH1.spacesOnEntities[MBTRI].test(HCURL)) {
337 CHKERR getEntitySense<MBTRI>(dataHcurl);
338 CHKERR getEntityDataOrder<MBTRI>(dataHcurl, HCURL);
339 dataHcurl.spacesOnEntities[MBTRI].set(HCURL);
340 }
341 if (dataH1.spacesOnEntities[MBQUAD].test(HCURL)) {
342 CHKERR getEntitySense<MBQUAD>(dataHcurl);
343 CHKERR getEntityDataOrder<MBQUAD>(dataHcurl, HCURL);
344 dataHcurl.spacesOnEntities[MBQUAD].set(HCURL);
345 }
346
347 // Hdiv
348 if (dataH1.spacesOnEntities[MBTRI].test(HDIV)) {
349 CHKERR getEntitySense<MBTRI>(dataHdiv);
350 CHKERR getEntityDataOrder<MBTRI>(dataHdiv, HDIV);
351 dataHdiv.spacesOnEntities[MBTRI].set(HDIV);
352 }
353 if (dataH1.spacesOnEntities[MBQUAD].test(HDIV)) {
354 CHKERR getEntitySense<MBQUAD>(dataHdiv);
355 CHKERR getEntityDataOrder<MBQUAD>(dataHdiv, HDIV);
356 dataHdiv.spacesOnEntities[MBQUAD].set(HDIV);
357 }
358
359 // L2
360 if (dataH1.spacesOnEntities[MBTRI].test(L2)) {
361 CHKERR getEntitySense<MBTRI>(dataL2);
362 CHKERR getEntityDataOrder<MBTRI>(dataL2, L2);
363 dataL2.spacesOnEntities[MBTRI].set(L2);
364 }
365 if (dataH1.spacesOnEntities[MBQUAD].test(L2)) {
366 CHKERR getEntitySense<MBQUAD>(dataL2);
367 CHKERR getEntityDataOrder<MBQUAD>(dataL2, L2);
368 dataL2.spacesOnEntities[MBQUAD].set(L2);
369 }
370
372}
373
377
378 const size_t nb_nodes =
379 dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE).size2();
380 double *shape_functions =
381 &*dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin();
382 const size_t nb_gauss_pts = gaussPts.size2();
383 coordsAtGaussPts.resize(nb_gauss_pts, 3, false);
384 for (int gg = 0; gg != nb_gauss_pts; ++gg)
385 for (int dd = 0; dd != 3; ++dd)
386 coordsAtGaussPts(gg, dd) = cblas_ddot(
387 nb_nodes, &shape_functions[nb_nodes * gg], 1, &coords[dd], 3);
388
390}
391
392MoFEMErrorCode FaceElementForcesAndSourcesCore::UserDataOperator::setPtrFE(
395 if (!(ptrFE = dynamic_cast<FaceElementForcesAndSourcesCore *>(ptr)))
396 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
397 "User operator and finite element do not work together");
399}
400
403
404 const auto type = numeredEntFiniteElementPtr->getEntType();
405 if (type != lastEvaluatedElementEntityType) {
406 switch (type) {
407 case MBTRI:
409 boost::shared_ptr<BaseFunction>(new TriPolynomialBase());
410 break;
411 case MBQUAD:
413 boost::shared_ptr<BaseFunction>(new QuadPolynomialBase());
414 break;
415 default:
417 }
420 }
421
422 // Calculate normal and tangent vectors for face geometry
425
427 if (gaussPts.size2() == 0)
429
434
435 // Iterate over operators
437
439}
440
442FaceElementForcesAndSourcesCore::UserDataOperator::loopSideVolumes(
443 const string fe_name, VolumeElementForcesAndSourcesCoreOnSide &fe_method) {
444 return loopSide(fe_name, &fe_method, 3);
445}
446
450
451#ifndef NDEBUG
452 if (toElePtr->gaussPts.size1() != getGaussPts().size1()) {
453 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
454 "Inconsistent numer of weights %d != %d",
455 toElePtr->gaussPts.size1(), getGaussPts().size1());
456 }
457 if (toElePtr->gaussPts.size2() != getGaussPts().size2()) {
458 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
459 "Inconsistent numer of integtaion pts %d != %d",
460 toElePtr->gaussPts.size2(), getGaussPts().size2());
461 }
462#endif
463
464 // TODO: add support for quad element types
465 switch (getFEType()) {
466 case MBTRI:
467 break;
468 default:
469 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
470 "Element type not implemented");
471 }
472
473 auto get_ftensor_from_mat_3d = [](MatrixDouble &m) {
475 &m(0, 0), &m(0, 1), &m(0, 2));
476 };
477
478 // get local coordinates, i.e. local coordinates on child element using partent local
479 // cooridinates
480 auto get_local_coords_triangle = [&]() {
481 double ref_gauss_pts[2][3] = {{0, 1, 0}, {0, 0, 1}};
482 std::array<double, 3> ksi0 = {0, 1, 0};
483 std::array<double, 3> ksi1 = {0, 0, 1};
484 std::array<double, 9> ref_shapes;
485 CHKERR Tools::shapeFunMBTRI<1>(ref_shapes.data(), ksi0.data(), ksi1.data(),
486 3);
487 auto &node_coords = getCoords();
488 auto &glob_coords = toElePtr->coords;
489 std::array<double, 6> local_coords;
491 &*node_coords.begin(), &*glob_coords.begin(), 3, local_coords.data());
492 return local_coords;
493 };
494
495 // get derivative of shape functions
496 auto get_diff_triangle = [&]() {
497 auto diff_ptr = Tools::diffShapeFunMBTRI.data();
499 diff_ptr, &diff_ptr[1]);
500 };
501
502 // get jaconbian to map local coordinates of parent to child element
503 auto get_jac = [&](auto &&local_coords, auto &&t_diff) {
507 auto t_local_coords = getFTensor1FromPtr<2>(local_coords.data());
508 t_jac(I, J) = 0;
509 for (int nn = 0; nn != 3; ++nn) {
510 t_jac(I, J) += t_local_coords(I) * t_diff(J);
511 ++t_local_coords;
512 ++t_diff;
513 }
514 return t_jac;
515 };
516
517 // get tangent vectors tensor
518 auto t_mat_tangent = [&](auto &t1, auto &t2) {
520 &t1(0), &t1(1), &t1(2), &t2(0), &t2(1), &t2(2)};
521 };
522
523 // transform tangent vectors to child element tangents
524 auto transfrom = [&](auto &&t_mat_t, auto &&t_mat_out_t, auto &&t_inv_jac) {
530 for (auto gg = 0; gg != getGaussPts().size2(); ++gg) {
532 t_mat_out_t(J, i) = t_mat_t(I, i) * t_inv_jac(I, J);
533 ++t_mat_t;
534 ++t_mat_out_t;
535 }
536 };
537
538 // calculate normal vector on child element
539 auto calc_normal = [&](auto &n, auto &t1, auto &t2) {
543 auto t_t1 = get_ftensor_from_mat_3d(t1);
544 auto t_t2 = get_ftensor_from_mat_3d(t2);
545 auto t_n = get_ftensor_from_mat_3d(n);
546 for (auto gg = 0; gg != getGaussPts().size2(); ++gg) {
547 t_n(j) = FTensor::levi_civita(i, j, k) * t_t1(k) * t_t2(i);
548 ++t_t1;
549 ++t_t2;
550 ++t_n;
551 }
552 };
553
554 transfrom(
555
556 t_mat_tangent(getTangent1AtGaussPts(), getTangent2AtGaussPts()),
557 t_mat_tangent(toElePtr->tangentOneAtGaussPts,
558 toElePtr->tangentTwoAtGaussPts),
559 get_jac(get_local_coords_triangle(), get_diff_triangle())
560
561 );
562 calc_normal(toElePtr->normalsAtGaussPts, toElePtr->tangentOneAtGaussPts,
563 toElePtr->tangentTwoAtGaussPts);
564
566}
567
568} // namespace MoFEM
static Index< 'J', 3 > J
static Number< 1 > N1
static Number< 0 > N0
static MoFEMErrorCode get_jac(EntitiesFieldData::EntData &col_data, int gg, MatrixDouble &jac_stress, MatrixDouble &jac)
@ NOBASE
Definition: definitions.h:59
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ H1
continuous field
Definition: definitions.h:85
@ HCURL
field with continuous tangents
Definition: definitions.h:86
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
constexpr int order
const int dim
#define diffN_MBQUAD2y(x)
Definition: fem_tools.h:66
#define diffN_MBQUAD1x(y)
Definition: fem_tools.h:63
#define N_MBQUAD3(x, y)
quad shape function
Definition: fem_tools.h:60
#define diffN_MBQUAD0x(y)
Definition: fem_tools.h:61
#define diffN_MBQUAD1y(x)
Definition: fem_tools.h:64
#define diffN_MBQUAD3y(x)
Definition: fem_tools.h:68
#define diffN_MBQUAD0y(x)
Definition: fem_tools.h:62
#define N_MBQUAD0(x, y)
quad shape function
Definition: fem_tools.h:57
#define diffN_MBQUAD3x(y)
Definition: fem_tools.h:67
#define diffN_MBQUAD2x(y)
Definition: fem_tools.h:65
#define N_MBQUAD2(x, y)
quad shape function
Definition: fem_tools.h:59
#define N_MBQUAD1(x, y)
quad shape function
Definition: fem_tools.h:58
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition: fem_tools.c:182
FTensor::Index< 'n', SPACE_DIM > n
FTensor::Index< 'm', SPACE_DIM > m
FTensor::Index< 'i', SPACE_DIM > i
const double v
phase velocity of light in medium (cm/ns)
FTensor::Index< 'j', 3 > j
FTensor::Index< 'k', 3 > k
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
constexpr IntegrationType I
double zeta
Viscous hardening.
Definition: plastic.cpp:177
#define QUAD_2D_TABLE_SIZE
Definition: quad.h:174
static QUAD *const QUAD_2D_TABLE[]
Definition: quad.h:175
virtual moab::Interface & get_moab()=0
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
virtual MoFEMErrorCode calculateAreaAndNormal()
Calculate element area and normal of the face.
MoFEMErrorCode operator()()
function is run for every finite element
virtual MoFEMErrorCode getSpaceBaseAndOrderOnElement()
Determine approximation space and order of base functions.
virtual MoFEMErrorCode calculateAreaAndNormalAtIntegrationPts()
Calculate element area and normal of the face at integration points.
virtual MoFEMErrorCode setIntegrationPts()
Set integration points.
virtual MoFEMErrorCode calculateCoordinatesAtGaussPts()
Calculate coordinate at integration points.
structure to get information form mofem into EntitiesFieldData
MoFEMErrorCode loopOverOperators()
Iterate user data operators.
int getMaxRowOrder() const
Get max order of approximation for field in rows.
MoFEMErrorCode getSpacesAndBaseOnEntities(EntitiesFieldData &data) const
Get field approximation space and base on entities.
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
auto & getElementPolynomialBase()
Get the Entity Polynomial Base object.
MatrixDouble coordsAtGaussPts
coordinated at gauss points
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement()
Calculate base functions.
MatrixDouble gaussPts
Matrix of integration points.
MoFEMErrorCode calBernsteinBezierBaseFunctionsOnElement()
Calculate Bernstein-Bezier base.
int getMaxColOrder() const
Get max order of approximation for field in columns.
EntityType lastEvaluatedElementEntityType
Last evaluated type of element entity.
int getMaxDataOrder() const
Get max order of approximation for data fields.
MoFEMErrorCode createDataOnElement(EntityType type)
Create a entity data on element object.
Copy geometry-related data from one element to other.
Calculate base functions on triangle.
static MoFEMErrorCode getLocalCoordinatesOnReferenceTriNodeTri(const T *elem_coords, const T *glob_coords, const int nb_nodes, T *local_coords)
Get the local coordinates on reference three node tri object.
static MoFEMErrorCode outerProductOfEdgeIntegrationPtsForQuad(MatrixDouble &pts, const int edge0, const int edge1)
Definition: Tools.cpp:612
static constexpr std::array< double, 8 > diffShapeFunMBQUADAtCenter
Definition: Tools.hpp:212
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition: Tools.hpp:104
Calculate base functions on triangle.
Base volume element used to integrate on skeleton.
int npoints
Definition: quad.h:29