v0.15.0
Loading...
Searching...
No Matches
VolumeElementForcesAndSourcesCore.cpp
Go to the documentation of this file.
1/** \file VolumeElementForcesAndSourcesCore.cpp
2
3\brief Implementation of volume element
4
5*/
6
7namespace MoFEM {
8
10 Interface &m_field)
11 : ForcesAndSourcesCore(m_field), vOlume(elementMeasure),
12 meshPositionsFieldName("MESH_NODE_POSITIONS"), coords(24), jAc(3, 3),
13 invJac(3, 3), opSetInvJacH1(invJac),
14 opContravariantPiolaTransform(elementMeasure, jAc),
15 opCovariantPiolaTransform(invJac), opSetInvJacHdivAndHcurl(invJac),
16 tJac(&jAc(0, 0), &jAc(0, 1), &jAc(0, 2), &jAc(1, 0), &jAc(1, 1),
17 &jAc(1, 2), &jAc(2, 0), &jAc(2, 1), &jAc(2, 2)),
18 tInvJac(&invJac(0, 0), &invJac(0, 1), &invJac(0, 2), &invJac(1, 0),
19 &invJac(1, 1), &invJac(1, 2), &invJac(2, 0), &invJac(2, 1),
20 &invJac(2, 2)) {}
21
25
29 int order_data = getMaxDataOrder();
30 int order_row = getMaxRowOrder();
31 int order_col = getMaxColOrder();
32 const auto type = numeredEntFiniteElementPtr->getEntType();
33
34 auto get_rule_by_type = [&]() {
35 switch (type) {
36 case MBHEX:
37 return getRule(order_row + 1, order_col + 1, order_data + 1);
38 default:
39 return getRule(order_row, order_col, order_data);
40 }
41 };
42
43 const int rule = get_rule_by_type();
44
45 auto calc_base_for_tet = [&]() {
47 const size_t nb_gauss_pts = gaussPts.size2();
48 auto &base = dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE);
49 auto &diff_base = dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
50 base.resize(nb_gauss_pts, 4, false);
51 diff_base.resize(nb_gauss_pts, 12, false);
52 CHKERR Tools::shapeFunMBTET(&*base.data().begin(), &gaussPts(0, 0),
53 &gaussPts(1, 0), &gaussPts(2, 0), nb_gauss_pts);
54 double *diff_shape_ptr = &*diff_base.data().begin();
55 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
56 for (int nn = 0; nn != 4; ++nn) {
57 for (int dd = 0; dd != 3; ++dd, ++diff_shape_ptr) {
58 *diff_shape_ptr = Tools::diffShapeFunMBTET[3 * nn + dd];
59 }
60 }
61 }
63 };
64
65 auto calc_base_for_hex = [&]() {
67 const size_t nb_gauss_pts = gaussPts.size2();
68 auto &base = dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE);
69 auto &diff_base = dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
70 base.resize(nb_gauss_pts, 8, false);
71 diff_base.resize(nb_gauss_pts, 24, false);
72 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
73 const double ksi = gaussPts(0, gg);
74 const double zeta = gaussPts(1, gg);
75 const double eta = gaussPts(2, gg);
76 base(gg, 0) = N_MBHEX0(ksi, zeta, eta);
77 base(gg, 1) = N_MBHEX1(ksi, zeta, eta);
78 base(gg, 2) = N_MBHEX2(ksi, zeta, eta);
79 base(gg, 3) = N_MBHEX3(ksi, zeta, eta);
80 base(gg, 4) = N_MBHEX4(ksi, zeta, eta);
81 base(gg, 5) = N_MBHEX5(ksi, zeta, eta);
82 base(gg, 6) = N_MBHEX6(ksi, zeta, eta);
83 base(gg, 7) = N_MBHEX7(ksi, zeta, eta);
84 diff_base(gg, 0 * 3 + 0) = diffN_MBHEX0x(zeta, eta);
85 diff_base(gg, 1 * 3 + 0) = diffN_MBHEX1x(zeta, eta);
86 diff_base(gg, 2 * 3 + 0) = diffN_MBHEX2x(zeta, eta);
87 diff_base(gg, 3 * 3 + 0) = diffN_MBHEX3x(zeta, eta);
88 diff_base(gg, 4 * 3 + 0) = diffN_MBHEX4x(zeta, eta);
89 diff_base(gg, 5 * 3 + 0) = diffN_MBHEX5x(zeta, eta);
90 diff_base(gg, 6 * 3 + 0) = diffN_MBHEX6x(zeta, eta);
91 diff_base(gg, 7 * 3 + 0) = diffN_MBHEX7x(zeta, eta);
92 diff_base(gg, 0 * 3 + 1) = diffN_MBHEX0y(ksi, eta);
93 diff_base(gg, 1 * 3 + 1) = diffN_MBHEX1y(ksi, eta);
94 diff_base(gg, 2 * 3 + 1) = diffN_MBHEX2y(ksi, eta);
95 diff_base(gg, 3 * 3 + 1) = diffN_MBHEX3y(ksi, eta);
96 diff_base(gg, 4 * 3 + 1) = diffN_MBHEX4y(ksi, eta);
97 diff_base(gg, 5 * 3 + 1) = diffN_MBHEX5y(ksi, eta);
98 diff_base(gg, 6 * 3 + 1) = diffN_MBHEX6y(ksi, eta);
99 diff_base(gg, 7 * 3 + 1) = diffN_MBHEX7y(ksi, eta);
100 diff_base(gg, 0 * 3 + 2) = diffN_MBHEX0z(ksi, zeta);
101 diff_base(gg, 1 * 3 + 2) = diffN_MBHEX1z(ksi, zeta);
102 diff_base(gg, 2 * 3 + 2) = diffN_MBHEX2z(ksi, zeta);
103 diff_base(gg, 3 * 3 + 2) = diffN_MBHEX3z(ksi, zeta);
104 diff_base(gg, 4 * 3 + 2) = diffN_MBHEX4z(ksi, zeta);
105 diff_base(gg, 5 * 3 + 2) = diffN_MBHEX5z(ksi, zeta);
106 diff_base(gg, 6 * 3 + 2) = diffN_MBHEX6z(ksi, zeta);
107 diff_base(gg, 7 * 3 + 2) = diffN_MBHEX7z(ksi, zeta);
108 }
110 };
111
112 auto set_integration_pts_for_hex = [&]() {
115 rule);
117 };
118
119 auto set_integration_pts_for_tet = [&]() {
121 if (rule < QUAD_3D_TABLE_SIZE) {
122 if (QUAD_3D_TABLE[rule]->dim != 3) {
123 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
124 }
125 if (QUAD_3D_TABLE[rule]->order < rule) {
127 "wrong order %d != %d", QUAD_3D_TABLE[rule]->order, rule);
128 }
129 size_t nb_gauss_pts = QUAD_3D_TABLE[rule]->npoints;
130 gaussPts.resize(4, nb_gauss_pts, false);
131 cblas_dcopy(nb_gauss_pts, &QUAD_3D_TABLE[rule]->points[1], 4,
132 &gaussPts(0, 0), 1);
133 cblas_dcopy(nb_gauss_pts, &QUAD_3D_TABLE[rule]->points[2], 4,
134 &gaussPts(1, 0), 1);
135 cblas_dcopy(nb_gauss_pts, &QUAD_3D_TABLE[rule]->points[3], 4,
136 &gaussPts(2, 0), 1);
137 cblas_dcopy(nb_gauss_pts, QUAD_3D_TABLE[rule]->weights, 1,
138 &gaussPts(3, 0), 1);
139
140 // CHKERR calc_base_for_tet();
141
142 auto &base = dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE);
143 auto &diff_base = dataH1.dataOnEntities[MBVERTEX][0].getDiffN(NOBASE);
144 base.resize(nb_gauss_pts, 4, false);
145 diff_base.resize(nb_gauss_pts, 12, false);
146 double *shape_ptr = &*base.data().begin();
147 cblas_dcopy(4 * nb_gauss_pts, QUAD_3D_TABLE[rule]->points, 1, shape_ptr,
148 1);
149 double *diff_shape_ptr = &*diff_base.data().begin();
150 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
151 for (int nn = 0; nn != 4; ++nn) {
152 for (int dd = 0; dd != 3; ++dd, ++diff_shape_ptr) {
153 *diff_shape_ptr = Tools::diffShapeFunMBTET[3 * nn + dd];
154 }
155 }
156 }
157
158 } else {
159 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
160 "rule > quadrature order %d < %d", rule, QUAD_3D_TABLE_SIZE);
161 }
163 };
164
165 if (rule >= 0) {
166 switch (type) {
167 case MBTET:
168 CHKERR set_integration_pts_for_tet();
169 break;
170 case MBHEX:
171 CHKERR set_integration_pts_for_hex();
172 CHKERR calc_base_for_hex();
173 break;
174 default:
175 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
176 "Element type not implemented: %d", type);
177 }
178 } else {
179 // If rule is negative, set user defined integration points
180 CHKERR setGaussPts(order_row, order_col, order_data);
181 const size_t nb_gauss_pts = gaussPts.size2();
182 if (nb_gauss_pts) {
183 switch (type) {
184 case MBTET: {
185 CHKERR calc_base_for_tet();
186 } break;
187 case MBHEX:
188 CHKERR calc_base_for_hex();
189 break;
190 default:
191 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
192 "Element type not implemented: %d", type);
193 }
194 }
195 }
196
198}
199
202 const auto ent = numeredEntFiniteElementPtr->getEnt();
203 const auto type = numeredEntFiniteElementPtr->getEntType();
204 CHKERR mField.get_moab().get_connectivity(ent, conn, num_nodes, true);
205 coords.resize(3 * num_nodes, false);
206 CHKERR mField.get_moab().get_coords(conn, num_nodes, &*coords.data().begin());
207
208 auto get_tet_t_diff_n = [&]() {
212 };
213
214 auto get_hex_t_diff_n = [&]() {
219 };
220
221 auto get_t_diff_n = [&]() {
222 if (type == MBTET)
223 return get_tet_t_diff_n();
224 return get_hex_t_diff_n();
225 };
226
227 auto t_diff_n = get_t_diff_n();
228
230 &coords[0], &coords[1], &coords[2]);
231 FTensor::Index<'i', 3> i;
232 FTensor::Index<'j', 3> j;
233 jAc.clear();
234
235 for (size_t n = 0; n != num_nodes; ++n) {
236 tJac(i, j) += t_coords(i) * t_diff_n(j);
237 ++t_coords;
238 ++t_diff_n;
239 }
242
243 if (type == MBTET)
244 vOlume *= G_TET_W1[0] / 6.;
245
246#ifndef NDEBUG
247 if (vOlume <= std::numeric_limits<double>::epsilon() || vOlume != vOlume) {
248 MOFEM_LOG("SELF", Sev::warning) << "Volume is zero " << vOlume;
249 MOFEM_LOG("SELF", Sev::warning) << "Coords " << coords << endl;
250 }
251#endif
252
254}
255
259 // Get coords at Gauss points
260 FTensor::Index<'i', 3> i;
261
262 auto &shape_functions = dataH1.dataOnEntities[MBVERTEX][0].getN(NOBASE);
263 double *shape_functions_ptr = &*shape_functions.data().begin();
264 const size_t nb_base_functions = shape_functions.size2();
265 const size_t nb_gauss_pts = gaussPts.size2();
266 coordsAtGaussPts.resize(nb_gauss_pts, 3, false);
267 coordsAtGaussPts.clear();
268 FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3> t_coords_at_gauss_ptr(
269 &coordsAtGaussPts(0, 0), &coordsAtGaussPts(0, 1),
270 &coordsAtGaussPts(0, 2));
272 shape_functions_ptr);
273 for (unsigned int gg = 0; gg != nb_gauss_pts; ++gg) {
275 &coords[0], &coords[1], &coords[2]);
276 for (int bb = 0; bb != nb_base_functions; ++bb) {
277 t_coords_at_gauss_ptr(i) += t_coords(i) * t_shape_functions;
278 ++t_coords;
279 ++t_shape_functions;
280 };
281 ++t_coords_at_gauss_ptr;
282 }
284}
285
289
292
293 auto type = numeredEntFiniteElementPtr->getEntType();
294 auto dim_type = CN::Dimension(type);
295
296 auto get_data_on_ents = [&](auto lower_dim, auto space) {
298 auto data = dataOnElement[space];
299 data->facesNodes = dataH1.facesNodes;
300 data->facesNodesOrder = dataH1.facesNodesOrder;
301 for (auto dd = dim_type; dd >= lower_dim; --dd) {
302 int nb_ents = moab::CN::NumSubEntities(type, dd);
303 for (int ii = 0; ii != nb_ents; ++ii) {
304 auto sub_ent_type = moab::CN::SubEntityType(type, dd, ii);
305 if ((dataH1.spacesOnEntities[sub_ent_type]).test(space)) {
306 auto &data_on_ent = data->dataOnEntities[sub_ent_type];
307 CHKERR getEntitySense(sub_ent_type, data_on_ent);
308 CHKERR getEntityDataOrder(sub_ent_type, space, data_on_ent);
309 data->spacesOnEntities[sub_ent_type].set(space);
310 }
311 }
312 }
314 };
315
316 CHKERR get_data_on_ents(1, H1); // H1
317 CHKERR get_data_on_ents(1, HCURL); // Hcurl
318 CHKERR get_data_on_ents(2, HDIV); // Hdiv
319 CHKERR get_data_on_ents(3, L2); // L2
320
322}
323
326
327 if (numeredEntFiniteElementPtr->getEntType() == MBTET) {
328 // OK that code is not nice.
329 MatrixDouble new_diff_n;
330 for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
331 FTensor::Index<'i', 3> i;
332 FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
334 if ((data.getDiffN(base).size1() == 4) &&
335 (data.getDiffN(base).size2() == 3)) {
336 const size_t nb_gauss_pts = gaussPts.size2();
337 const size_t nb_base_functions = 4;
338 new_diff_n.resize(nb_gauss_pts, 3 * nb_base_functions, false);
339 double *new_diff_n_ptr = &*new_diff_n.data().begin();
341 new_diff_n_ptr, &new_diff_n_ptr[1], &new_diff_n_ptr[2]);
342 double *t_diff_n_ptr = &*data.getDiffN(base).data().begin();
343 for (unsigned int gg = 0; gg != nb_gauss_pts; gg++) {
345 t_diff_n_ptr, &t_diff_n_ptr[1], &t_diff_n_ptr[2]);
346 for (unsigned int bb = 0; bb != nb_base_functions; bb++) {
347 t_new_diff_n(i) = t_diff_n(i);
348 ++t_new_diff_n;
349 ++t_diff_n;
350 }
351 }
352 data.getDiffN(base).resize(new_diff_n.size1(), new_diff_n.size2(),
353 false);
354 data.getDiffN(base).swap(new_diff_n);
355 }
356 }
357 }
358
359 if (dataH1.spacesOnEntities[MBVERTEX].test(H1))
361
362 std::array<std::bitset<LASTSPACE>, 4> spaces_by_dim{
363
364 std::bitset<LASTSPACE>(0), //
365 std::bitset<LASTSPACE>(0), //
366 std::bitset<LASTSPACE>(0), //
367 std::bitset<LASTSPACE>(0)
368
369 };
370 for (auto type = MBEDGE; type != MBENTITYSET; ++type) {
371 spaces_by_dim[CN::Dimension(type)] |= dataH1.spacesOnEntities[type];
372 }
373
374 if (
375
376 spaces_by_dim[1].test(HCURL) || //
377 spaces_by_dim[2].test(HCURL) || //
378 spaces_by_dim[3].test(HCURL)
379
380 ) {
383 }
384
385 if (spaces_by_dim[2].test(HDIV) || spaces_by_dim[3].test(HDIV)) {
387 // Fix for tetrahedrons
388 if (numeredEntFiniteElementPtr->getEntType() == MBTET) {
389 for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
390 FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
391 for (auto t : {MBTRI, MBTET}) {
392 for (auto &d : dataHdiv.dataOnEntities[t]) {
393 d.getN(base) /= 6;
394 d.getDiffN(base) /= 6;
395 }
396 }
397 }
398 }
400 }
401
402 if (spaces_by_dim[3].test(L2)) {
404 }
405
407}
408
409MoFEMErrorCode VolumeElementForcesAndSourcesCore::UserDataOperator::setPtrFE(
412 if (!(ptrFE = dynamic_cast<VolumeElementForcesAndSourcesCore *>(ptr)))
413 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
414 "User operator and finite element do not work together");
416}
417
420
421 const auto type = numeredEntFiniteElementPtr->getEntType();
422 if (type != lastEvaluatedElementEntityType) {
423 switch (type) {
424 case MBTET:
426 boost::shared_ptr<BaseFunction>(new TetPolynomialBase(this));
427 break;
428 case MBHEX:
430 boost::shared_ptr<BaseFunction>(new HexPolynomialBase());
431 break;
432 default:
434 }
437 }
438
442 if (gaussPts.size2() == 0)
447
449
450 // Iterate over operators
452
454}
455
456} // namespace MoFEM
MatrixDouble invJac
FieldApproximationBase
approximation base
Definition definitions.h:58
@ LASTBASE
Definition definitions.h:69
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base nme:nme847.
Definition definitions.h:60
@ 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_MBHEX3z(x, y)
Definition fem_tools.h:98
#define diffN_MBHEX7z(x, y)
Definition fem_tools.h:102
#define diffN_MBHEX7y(x, z)
Definition fem_tools.h:94
#define diffN_MBHEX6x(y, z)
Definition fem_tools.h:85
#define N_MBHEX7(x, y, z)
Definition fem_tools.h:78
#define N_MBHEX3(x, y, z)
Definition fem_tools.h:74
#define N_MBHEX5(x, y, z)
Definition fem_tools.h:76
#define diffN_MBHEX1y(x, z)
Definition fem_tools.h:88
#define diffN_MBHEX3y(x, z)
Definition fem_tools.h:90
static const double G_TET_W1[]
Definition fem_tools.h:1115
#define diffN_MBHEX4y(x, z)
Definition fem_tools.h:91
#define diffN_MBHEX5x(y, z)
Definition fem_tools.h:84
#define diffN_MBHEX1z(x, y)
Definition fem_tools.h:96
#define diffN_MBHEX4z(x, y)
Definition fem_tools.h:99
#define diffN_MBHEX6y(x, z)
Definition fem_tools.h:93
#define diffN_MBHEX6z(x, y)
Definition fem_tools.h:101
#define N_MBHEX4(x, y, z)
Definition fem_tools.h:75
#define N_MBHEX0(x, y, z)
Definition fem_tools.h:71
#define N_MBHEX6(x, y, z)
Definition fem_tools.h:77
#define diffN_MBHEX5z(x, y)
Definition fem_tools.h:100
#define N_MBHEX2(x, y, z)
Definition fem_tools.h:73
#define diffN_MBHEX0x(y, z)
Definition fem_tools.h:79
#define diffN_MBHEX1x(y, z)
Definition fem_tools.h:80
#define diffN_MBHEX2z(x, y)
Definition fem_tools.h:97
#define diffN_MBHEX5y(x, z)
Definition fem_tools.h:92
#define diffN_MBHEX4x(y, z)
Definition fem_tools.h:83
#define N_MBHEX1(x, y, z)
Definition fem_tools.h:72
#define diffN_MBHEX3x(y, z)
Definition fem_tools.h:82
#define diffN_MBHEX2y(x, z)
Definition fem_tools.h:89
#define diffN_MBHEX0y(x, z)
Definition fem_tools.h:87
#define diffN_MBHEX0z(x, y)
Definition fem_tools.h:95
#define diffN_MBHEX2x(y, z)
Definition fem_tools.h:81
#define diffN_MBHEX7x(y, z)
Definition fem_tools.h:86
double eta
#define MOFEM_LOG(channel, severity)
Log.
FTensor::Index< 'i', SPACE_DIM > i
const double n
refractive index of diffusive medium
FTensor::Index< 'j', 3 > j
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.
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
MoFEMErrorCode invertTensor3by3(ublas::matrix< T, L, A > &jac_data, ublas::vector< T, A > &det_data, ublas::matrix< T, L, A > &inv_jac_data)
Calculate inverse of tensor rank 2 at integration points.
static auto determinantTensor3by3(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
double zeta
Viscous hardening.
Definition plastic.cpp:130
constexpr double t
plate stiffness
Definition plate.cpp:58
#define QUAD_3D_TABLE_SIZE
Definition quad.h:186
static QUAD *const QUAD_3D_TABLE[]
Definition quad.h:187
virtual moab::Interface & get_moab()=0
virtual MPI_Comm & get_comm() const =0
virtual MoFEMErrorCode opRhs(EntitiesFieldData &data, const bool error_if_no_base=false)
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
MatrixDouble & getDiffN(const FieldApproximationBase base)
get derivatives of base functions
MatrixInt facesNodes
nodes on finite element faces
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
MatrixInt facesNodesOrder
order of face nodes on element
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
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.
MoFEMErrorCode getFaceNodes(EntitiesFieldData &data) const
Get nodes on faces.
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.
Calculate base functions on tetrahedral.
Calculate base functions on tetrahedral.
static MoFEMErrorCode outerProductOfEdgeIntegrationPtsForHex(MatrixDouble &pts, const int edge0, const int edge1, const int edge2)
Definition Tools.cpp:662
static constexpr std::array< double, 24 > diffShapeFunMBHEXAtCenter
Definition Tools.hpp:218
static MoFEMErrorCode shapeFunMBTET(double *shape, const double *ksi, const double *eta, const double *zeta, const double nb)
Calculate shape functions on tetrahedron.
Definition Tools.hpp:747
static constexpr std::array< double, 12 > diffShapeFunMBTET
Definition Tools.hpp:271
virtual MoFEMErrorCode transformBaseFunctions()
Transform base functions based on geometric element Jacobian.
virtual MoFEMErrorCode setIntegrationPts()
Set integration points.
MoFEMErrorCode operator()()
function is run for every finite element
virtual MoFEMErrorCode getSpaceBaseAndOrderOnElement()
Determine approximation space and order of base functions.
virtual MoFEMErrorCode calculateVolumeAndJacobian()
Calculate element volume and Jacobian.
virtual MoFEMErrorCode calculateCoordinatesAtGaussPts()
Calculate coordinate at integration points.
int order
Definition quad.h:28
int npoints
Definition quad.h:29