v0.15.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
20 FTensor::Index<'i', 3> i;
21 FTensor::Index<'j', 3> j;
22 FTensor::Index<'k', 3> k;
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
137 FTensor::Index<'i', 3> i;
138 FTensor::Index<'j', 3> j;
139 FTensor::Index<'k', 3> k;
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 SETERRQ(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 SETERRQ(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 SETERRQ(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 SETERRQ(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 auto type = numeredEntFiniteElementPtr->getEntType();
317 auto dim_type = CN::Dimension(type);
318
319 auto get_data_on_ents = [&](auto lower_dim, auto space) {
321 auto data = dataOnElement[space];
322 for (auto dd = dim_type; dd >= lower_dim; --dd) {
323 int nb_ents = moab::CN::NumSubEntities(type, dd);
324 for (int ii = 0; ii != nb_ents; ++ii) {
325 auto sub_ent_type = moab::CN::SubEntityType(type, dd, ii);
326 if ((dataH1.spacesOnEntities[sub_ent_type]).test(space)) {
327 auto &data_on_ent = data->dataOnEntities[sub_ent_type];
328 CHKERR getEntitySense(sub_ent_type, data_on_ent);
329 CHKERR getEntityDataOrder(sub_ent_type, space, data_on_ent);
330 data->spacesOnEntities[sub_ent_type].set(space);
331 }
332 }
333 }
335 };
336
337 CHKERR get_data_on_ents(1, H1); // H1
338 CHKERR get_data_on_ents(1, HCURL); // Hcurl
339 CHKERR get_data_on_ents(2, HDIV); // Hdiv
340 CHKERR get_data_on_ents(2, L2); // L2
341
343}
344
348
349 const size_t nb_nodes =
350 dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE).size2();
351 double *shape_functions =
352 &*dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin();
353 const size_t nb_gauss_pts = gaussPts.size2();
354 coordsAtGaussPts.resize(nb_gauss_pts, 3, false);
355 for (int gg = 0; gg != nb_gauss_pts; ++gg)
356 for (int dd = 0; dd != 3; ++dd)
357 coordsAtGaussPts(gg, dd) = cblas_ddot(
358 nb_nodes, &shape_functions[nb_nodes * gg], 1, &coords[dd], 3);
359
361}
362
363MoFEMErrorCode FaceElementForcesAndSourcesCore::UserDataOperator::setPtrFE(
366 if (!(ptrFE = dynamic_cast<FaceElementForcesAndSourcesCore *>(ptr)))
367 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
368 "User operator and finite element do not work together");
370}
371
374
375 const auto type = numeredEntFiniteElementPtr->getEntType();
376 if (type != lastEvaluatedElementEntityType) {
377 switch (type) {
378 case MBTRI:
380 boost::shared_ptr<BaseFunction>(new TriPolynomialBase());
381 break;
382 case MBQUAD:
384 boost::shared_ptr<BaseFunction>(new QuadPolynomialBase());
385 break;
386 default:
388 }
391 }
392
393 // Calculate normal and tangent vectors for face geometry
396
398 if (gaussPts.size2() == 0)
400
405
406 // Iterate over operators
408
410}
411
413FaceElementForcesAndSourcesCore::UserDataOperator::loopSideVolumes(
414 const string fe_name, VolumeElementForcesAndSourcesCoreOnSide &fe_method) {
415 return loopSide(fe_name, &fe_method, 3);
416}
417
421
422#ifndef NDEBUG
423 if (toElePtr->gaussPts.size1() != getGaussPts().size1()) {
424 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
425 "Inconsistent numer of weights %ld != %ld",
426 toElePtr->gaussPts.size1(), getGaussPts().size1());
427 }
428 if (toElePtr->gaussPts.size2() != getGaussPts().size2()) {
429 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
430 "Inconsistent numer of integtaion pts %ld != %ld",
431 toElePtr->gaussPts.size2(), getGaussPts().size2());
432 }
433#endif
434
435 // TODO: add support for quad element types
436 switch (getFEType()) {
437 case MBTRI:
438 break;
439 default:
440 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
441 "Element type not implemented");
442 }
443
444 auto get_ftensor_from_mat_3d = [](MatrixDouble &m) {
446 &m(0, 0), &m(0, 1), &m(0, 2));
447 };
448
449 // get local coordinates, i.e. local coordinates on child element using parent
450 // local coordinates
451 auto get_local_coords_triangle = [&]() {
452 std::array<double, 3> ksi0 = {0, 1, 0};
453 std::array<double, 3> ksi1 = {0, 0, 1};
454 std::array<double, 9> ref_shapes;
455 CHKERR Tools::shapeFunMBTRI<1>(ref_shapes.data(), ksi0.data(), ksi1.data(),
456 3);
457 auto &node_coords = getCoords();
458 auto &glob_coords = toElePtr->coords;
459 std::array<double, 6> local_coords;
461 &*node_coords.begin(), &*glob_coords.begin(), 3, local_coords.data());
462 return local_coords;
463 };
464
465 // get derivative of shape functions
466 auto get_diff_triangle = [&]() {
467 auto diff_ptr = Tools::diffShapeFunMBTRI.data();
469 diff_ptr, &diff_ptr[1]);
470 };
471
472 // get jacobian to map local coordinates of parent to child element
473 auto get_jac = [&](auto &&local_coords, auto &&t_diff) {
474 FTensor::Index<'I', 2> I;
475 FTensor::Index<'J', 2> J;
477 auto t_local_coords = getFTensor1FromPtr<2>(local_coords.data());
478 t_jac(I, J) = 0;
479 for (int nn = 0; nn != 3; ++nn) {
480 t_jac(I, J) += t_local_coords(I) * t_diff(J);
481 ++t_local_coords;
482 ++t_diff;
483 }
484 return t_jac;
485 };
486
487 // get tangent vectors tensor
488 auto t_mat_tangent = [&](auto &t1, auto &t2) {
490 &t1(0), &t1(1), &t1(2), &t2(0), &t2(1), &t2(2)};
491 };
492
493 // transform tangent vectors to child element tangents
494 auto transform = [&](auto &&t_mat_t, auto &&t_mat_out_t, auto &&t_inv_jac) {
495 FTensor::Index<'I', 2> I;
496 FTensor::Index<'J', 2> J;
497 FTensor::Index<'i', 3> i;
500 for (auto gg = 0; gg != getGaussPts().size2(); ++gg) {
502 t_mat_out_t(J, i) = t_mat_t(I, i) * t_inv_jac(I, J);
503 ++t_mat_t;
504 ++t_mat_out_t;
505 }
506 };
507
508 // calculate normal vector on child element
509 auto calc_normal = [&](auto &n, auto &t1, auto &t2) {
510 FTensor::Index<'i', 3> i;
511 FTensor::Index<'j', 3> j;
512 FTensor::Index<'k', 3> k;
513 auto t_t1 = get_ftensor_from_mat_3d(t1);
514 auto t_t2 = get_ftensor_from_mat_3d(t2);
515 auto t_n = get_ftensor_from_mat_3d(n);
516 for (auto gg = 0; gg != getGaussPts().size2(); ++gg) {
517 t_n(j) = FTensor::levi_civita(i, j, k) * t_t1(k) * t_t2(i);
518 ++t_t1;
519 ++t_t2;
520 ++t_n;
521 }
522 };
523
524 transform(
525
526 t_mat_tangent(getTangent1AtGaussPts(), getTangent2AtGaussPts()),
527 t_mat_tangent(toElePtr->tangentOneAtGaussPts,
528 toElePtr->tangentTwoAtGaussPts),
529 get_jac(get_local_coords_triangle(), get_diff_triangle())
530
531 );
532 calc_normal(toElePtr->normalsAtGaussPts, toElePtr->tangentOneAtGaussPts,
533 toElePtr->tangentTwoAtGaussPts);
534
536}
537
538} // namespace MoFEM
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()
@ 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 ...
@ 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()
#define CHKERR
Inline error check.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
constexpr int order
#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< 'i', SPACE_DIM > i
const double v
phase velocity of light in medium (cm/ns)
const double n
refractive index of diffusive medium
FTensor::Index< 'J', DIM1 > J
Definition level_set.cpp:30
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.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
FTensor::Tensor1< FTensor::PackPtr< double *, S >, DIM > getFTensor1FromPtr(double *ptr)
Make Tensor1 from pointer.
constexpr IntegrationType I
double zeta
Viscous hardening.
Definition plastic.cpp:130
#define QUAD_2D_TABLE_SIZE
Definition quad.h:174
static QUAD *const QUAD_2D_TABLE[]
Definition quad.h:175
FTensor::Index< 'm', 3 > m
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 getEntitySense(const EntityType type, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
get sense (orientation) of entity
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
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.
MoFEMErrorCode getEntityDataOrder(const EntityType type, const FieldSpace space, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
Get the entity data order.
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 outerProductOfEdgeIntegrationPtsForQuad(MatrixDouble &pts, const int edge0, const int edge1)
Definition Tools.cpp:613
static constexpr std::array< double, 8 > diffShapeFunMBQUADAtCenter
Definition Tools.hpp:212
static constexpr std::array< double, 6 > diffShapeFunMBTRI
Definition Tools.hpp:104
static MoFEMErrorCode getLocalCoordinatesOnReferenceThreeNodeTri(const double *elem_coords, const double *glob_coords, const int nb_nodes, double *local_coords)
Get the local coordinates on reference three node tri object.
Definition Tools.cpp:188
static MoFEMErrorCode shapeFunMBTRI(double *shape, const double *ksi, const double *eta, const int nb)
Calculate shape functions on triangle.
Definition Tools.hpp:716
Calculate base functions on triangle.
int npoints
Definition quad.h:29