v0.13.1
HODataOperators.hpp
Go to the documentation of this file.
1/** \file HODataOperators.hpp
2 * \brief Operators managing HO geometry
3
4*/
5
6#ifndef __HO_DATA_OPERATORS_HPP__
7#define __HO_DATA_OPERATORS_HPP__
8
9namespace MoFEM {
10
11/**
12 * @brief Calculate jacobian on Hex or other volume which is not simplex.
13 *
14 * Using not a field data, but geometrical coordinates of element.
15 *
16 * \note Use OpCalculateVectorFieldGradient to calculate Jacobian from field
17 * data.
18 *
19 */
22
23 OpCalculateHOJacForVolume(boost::shared_ptr<MatrixDouble> jac_ptr);
24
27
28private:
29 boost::shared_ptr<MatrixDouble> jacPtr;
30};
31
32/** \deprecated use OpCalculateHOJacForVolume
33 */
35
36/**
37 * @brief Calculate HO coordinates at gauss points
38 *
39 */
40template <int FIELD_DIM = 3>
42
45
48};
49
50/**
51 * @deprecated This class should be DIM = 3 specialization when default
52 * parameter is removed
53 *
54 */
56
58
61};
62
63/**
64 * \brief Set inverse jacobian to base functions
65 *
66 * \deprecated Version with default variant DIM = 3 is deprecated, keep for
67 * back compatibility. Should be removed in future. Use of default variant make
68 * code implicit, what can be source of some hidden error. Explict interface is
69 * better, when user see and control outcome, and is aware of existing variants.
70 *
71 */
72template <int DIM = 3>
74 using OpSetHOInvJacToScalarBasesImpl::OpSetHOInvJacToScalarBasesImpl;
75};
76
77template <>
79 : public OpSetInvJacSpaceForFaceImpl<2, 1> {
80
81 using OpSetInvJacSpaceForFaceImpl<2, 1>::OpSetInvJacSpaceForFaceImpl;
82};
83
84/**
85 * \brief transform local reference derivatives of shape function to global
86 derivatives if higher order geometry is given
87 *
88
89 * \ingroup mofem_forces_and_sources
90*/
92
94 boost::shared_ptr<MatrixDouble> inv_jac_ptr)
95 : ForcesAndSourcesCore::UserDataOperator(space), invJacPtr(inv_jac_ptr) {
96
97 if (!invJacPtr)
99 "Pointer for invJacPtr not allocated");
100
101 doVertices = false;
102 if (space == HDIV)
103 doEdges = false;
104 }
105
108
109private:
110 boost::shared_ptr<MatrixDouble> invJacPtr;
113};
114
115/**
116 * @brief Modify integration weights on face to take in account higher-order
117 * geometry
118 * @ingroup mofem_forces_and_sources_tri_element
119 *
120 */
127};
128
129/**
130 * \brief Set inverse jacobian to base functions
131 *
132 */
134
135 OpSetHOWeights(boost::shared_ptr<VectorDouble> det_ptr)
137 detPtr(det_ptr) {
138 if (!detPtr)
140 "Pointer for detPtr not allocated");
141 }
142
145
146private:
147 boost::shared_ptr<VectorDouble> detPtr;
148};
149
150/** \brief Apply contravariant (Piola) transfer to Hdiv space for HO geometry
151
152* \ingroup mofem_forces_and_sources
153*/
156
158 boost::shared_ptr<VectorDouble> det_ptr,
159 boost::shared_ptr<MatrixDouble> jac_ptr)
161 jacPtr(jac_ptr) {
162 doVertices = false;
163 if (space == HDIV)
164 doEdges = false;
165 }
166
169
170private:
171 boost::shared_ptr<VectorDouble> detPtr;
172 boost::shared_ptr<MatrixDouble> jacPtr;
173
176};
177
178/** \brief Apply covariant (Piola) transfer to Hcurl space for HO geometry
179 */
182
184 boost::shared_ptr<MatrixDouble> jac_inv_ptr)
186 jacInvPtr(jac_inv_ptr) {
187 doVertices = false;
188 if (space == HDIV)
189 doEdges = false;
190 if (!jacInvPtr)
192 "Pointer for jacPtr not allocated");
193 }
194
197
198private:
199 boost::shared_ptr<MatrixDouble> jacInvPtr;
200
203};
204
205/** \brief Calculate jacobian for face element
206
207 It is assumed that face element is XY plane. Applied
208 only for 2d problems and 2d problems embedded in 3d space
209
210 \note If you push operators for HO normal befor this operator, HO geometry is
211 taken into account when you calculate jacobian.
212
213*/
214template <int DIM> struct OpCalculateHOJacForFaceImpl;
215
216template <>
219
220 OpCalculateHOJacForFaceImpl(boost::shared_ptr<MatrixDouble> jac_ptr);
221
222 MoFEMErrorCode doWork(int side, EntityType type,
224
225protected:
226 boost::shared_ptr<MatrixDouble> jacPtr;
227};
228
229template <>
231
232 using OpCalculateHOJacForFaceImpl<2>::OpCalculateHOJacForFaceImpl;
233
234 MoFEMErrorCode doWork(int side, EntityType type,
236};
237
240
241template <int DIM> struct OpCalculateHOJac;
242
243template <> struct OpCalculateHOJac<3> : public OpCalculateHOJacForVolume {
245};
246
247template <> struct OpCalculateHOJac<2> : public OpCalculateHOJacForFaceImpl<2> {
248 using OpCalculateHOJacForFaceImpl<2>::OpCalculateHOJacForFaceImpl;
249};
250
251/** \brief Calculate normals at Gauss points of triangle element
252 * \ingroup mofem_forces_and_source
253 */
254template <int FIELD_DIM = 3>
257
260
263};
264
265/** \brief transform Hdiv base fluxes from reference element to physical
266 * triangle \ingroup mofem_forces_and_sources
267 *
268 * \note Matrix which keeps normal is assumed to have three columns, and number
269 * of rows should be equal to number of integration points.
270 *
271 */
274
276 const FieldSpace space,
277 boost::shared_ptr<MatrixDouble> normals_at_gauss_pts = nullptr)
279 normalsAtGaussPts(normals_at_gauss_pts) {}
280
283
284private:
285 boost::shared_ptr<MatrixDouble> normalsAtGaussPts;
286};
287
288/** \brief transform Hcurl base fluxes from reference element to physical edge
289 * \ingroup mofem_forces_and_sources
290 */
293
295 const FieldSpace space = HCURL,
296 boost::shared_ptr<MatrixDouble> tangent_at_pts = nullptr)
298 tangentAtGaussPts(tangent_at_pts) {}
299
302
303private:
304 boost::shared_ptr<MatrixDouble> tangentAtGaussPts;
305};
306
307/** \brief transform Hcurl base fluxes from reference element to physical
308 * triangle \ingroup mofem_forces_and_sources
309 */
312
314 const FieldSpace space,
315 boost::shared_ptr<MatrixDouble> normals_at_pts = nullptr,
316 boost::shared_ptr<MatrixDouble> tangent1_at_pts = nullptr,
317 boost::shared_ptr<MatrixDouble> tangent2_at_pts = nullptr)
319 normalsAtPts(normals_at_pts), tangent1AtPts(tangent1_at_pts),
320 tangent2AtPts(tangent2_at_pts) {
321 if (normals_at_pts || tangent1_at_pts || tangent2_at_pts)
322 if (normals_at_pts && tangent1_at_pts && tangent2_at_pts)
324 "All elements in constructor have to set pointer");
325 }
326
329
330private:
331 boost::shared_ptr<MatrixDouble> normalsAtPts;
332 boost::shared_ptr<MatrixDouble> tangent1AtPts;
333 boost::shared_ptr<MatrixDouble> tangent2AtPts;
334
337};
338
339/** \brief Calculate tangent vector on edge form HO geometry approximation
340 * \ingroup mofem_forces_and_sources
341 */
342template <int FIELD_DIM = 3>
345
347 std::string field_name,
348 boost::shared_ptr<MatrixDouble> tangents_at_pts = nullptr)
350 tangentsAtPts(tangents_at_pts) {}
351
354
355private:
356 boost::shared_ptr<MatrixDouble> tangentsAtPts;
357};
358
359/**
360 * @brief Scale base functions by inverses of measure of element
361 *
362 * @tparam OP
363 */
364template <typename OP> struct OpScaleBaseBySpaceInverseOfMeasure : public OP {
365
367 boost::shared_ptr<VectorDouble> det_jac_ptr, const FieldSpace space = L2)
368 : OP(space), fieldSpace(space), detJacPtr(det_jac_ptr) {}
369
373
374 if (!detJacPtr) {
375 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "detJacPtr not set");
376 }
377
378 auto scale = [&]() {
379 for (int b = AINSWORTH_LEGENDRE_BASE; b != USER_BASE; b++) {
381 auto &base_fun = data.getN(base);
382 auto &diff_base_fun = data.getDiffN(base);
383 if (detJacPtr) {
384
385 auto &det_vec = *detJacPtr;
386 const auto nb_base_fun = base_fun.size2();
387 const auto nb_int_pts = base_fun.size1();
388
389 if (nb_int_pts != det_vec.size())
390 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
391 "Number of integration pts in detJacPtr does not mush "
392 "number of integration pts in base function");
393
394 auto get_row = [](auto &m, auto gg) {
395 return ublas::matrix_row<decltype(m)>(m, gg);
396 };
397
398 for (auto gg = 0; gg != nb_int_pts; ++gg)
399 get_row(base_fun) /= det_vec[gg];
400
401 if (diff_base_fun.size1() == nb_int_pts) {
402 for (auto gg = 0; gg != nb_int_pts; ++gg)
403 get_row(diff_base_fun) /= det_vec[gg];
404 }
405 }
406 }
407 };
408
409 if (this->getFEDim() == 3) {
410 switch (fieldSpace) {
411 case H1:
412 scale();
413 break;
414 case HCURL:
415 if (type >= MBEDGE)
416 scale();
417 break;
418 case HDIV:
419 if (type >= MBTRI)
420 scale();
421 break;
422 case L2:
423 if (type >= MBTET)
424 scale();
425 break;
426 default:
427 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "impossible case");
428 }
429 } else if (this->getFEDim() == 2) {
430 switch (fieldSpace) {
431 case H1:
432 scale();
433 break;
434 case HCURL:
435 if (type >= MBEDGE)
436 scale();
437 break;
438 case HDIV:
439 case L2:
440 if (type >= MBTRI)
441 scale();
442 break;
443 default:
444 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "impossible case");
445 }
446 } else if (this->getFEDim() == 1) {
447 switch (fieldSpace) {
448 case H1:
449 scale();
450 break;
451 case HCURL:
452 case L2:
453 if (type >= MBEDGE)
454 scale();
455 break;
456 default:
457 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "impossible case");
458 }
459 } else {
460 SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSIBLE_CASE, "impossible case");
461 }
462
464 }
465
466private:
468 boost::shared_ptr<VectorDouble> detJacPtr;
469};
470
471/**
472 * @brief Add operators pushing bases from local to physical configuration
473 *
474 * @tparam FE_DIM dimension of element
475 * @tparam PROBLEM_DIM problem dimension
476 * @tparam SPACE_DIM space dimension
477 */
478template <int FE_DIM, int PROBLEM_DIM, int SPACE_DIM> struct AddHOOps;
479
480template <> struct AddHOOps<2, 2, 2> {
481 AddHOOps() = delete;
482 static MoFEMErrorCode
483 add(boost::ptr_vector<ForcesAndSourcesCore::UserDataOperator> &pipeline,
484 std::vector<FieldSpace> spaces, std::string geom_field_name = "");
485};
486
487template <> struct AddHOOps<1, 2, 2> {
488 AddHOOps() = delete;
489 static MoFEMErrorCode
490 add(boost::ptr_vector<ForcesAndSourcesCore::UserDataOperator> &pipeline,
491 std::vector<FieldSpace> spaces, std::string geom_field_name = "");
492};
493template <> struct AddHOOps<1, 3, 3> {
494 AddHOOps() = delete;
495 static MoFEMErrorCode
496 add(boost::ptr_vector<ForcesAndSourcesCore::UserDataOperator> &pipeline,
497 std::vector<FieldSpace> space, std::string geom_field_name = "");
498};
499
500template <> struct AddHOOps<2, 3, 3> {
501 AddHOOps() = delete;
502 static MoFEMErrorCode
503 add(boost::ptr_vector<ForcesAndSourcesCore::UserDataOperator> &pipeline,
504 std::vector<FieldSpace> space, std::string geom_field_name = "");
505};
506
507template <> struct AddHOOps<3, 3, 3> {
508 AddHOOps() = delete;
509 static MoFEMErrorCode
510 add(boost::ptr_vector<ForcesAndSourcesCore::UserDataOperator> &pipeline,
511 std::vector<FieldSpace> space, std::string geom_field_name = "");
512};
513
514template <int FIELD_DIM>
520 const auto nb_dofs = data.getFieldData().size() / FIELD_DIM;
521 if (nb_dofs) {
522 if (type == MBVERTEX)
523 getCoordsAtGaussPts().clear();
524 auto t_base = data.getFTensor0N();
525 auto t_coords = getFTensor1CoordsAtGaussPts();
526 const auto nb_integration_pts = data.getN().size1();
527 const auto nb_base_functions = data.getN().size2();
528 for (auto gg = 0; gg != nb_integration_pts; ++gg) {
529 auto t_dof = data.getFTensor1FieldData<FIELD_DIM>();
530 size_t bb = 0;
531 for (; bb != nb_dofs; ++bb) {
532 t_coords(i) += t_base * t_dof(i);
533 ++t_dof;
534 ++t_base;
535 }
536 for (; bb < nb_base_functions; ++bb)
537 ++t_base;
538 ++t_coords;
539 }
540 }
542}
543
544template <int FIELD_DIM>
549
552
553 auto get_ftensor1 = [](MatrixDouble &m) {
555 &m(0, 0), &m(0, 1), &m(0, 2));
556 };
557
558 unsigned int nb_dofs = data.getFieldData().size();
559 if (nb_dofs != 0) {
560
561 int nb_gauss_pts = data.getN().size1();
562 auto &tangent1_at_gauss_pts = getTangent1AtGaussPts();
563 auto &tangent2_at_gauss_pts = getTangent2AtGaussPts();
564 tangent1_at_gauss_pts.resize(nb_gauss_pts, 3, false);
565 tangent2_at_gauss_pts.resize(nb_gauss_pts, 3, false);
566
567 switch (type) {
568 case MBVERTEX: {
569 tangent1_at_gauss_pts.clear();
570 tangent2_at_gauss_pts.clear();
571 }
572 case MBEDGE:
573 case MBTRI:
574 case MBQUAD: {
575
576#ifndef NDEBUG
577 if (2 * data.getN().size2() != data.getDiffN().size2()) {
578 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
579 "expected two direcatives in local element coordinates");
580 }
581 if (nb_dofs % FIELD_DIM != 0) {
582 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
583 "expected that number of dofs is multiplicative of field "
584 "dimension");
585 }
586#endif
587
588 if (nb_dofs > FIELD_DIM * data.getN().size2()) {
589 unsigned int nn = 0;
590 for (; nn != nb_dofs; nn++) {
591 if (!data.getFieldDofs()[nn]->getActive())
592 break;
593 }
594 if (nn > FIELD_DIM * data.getN().size2()) {
595 SETERRQ1(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
596 "Data inconsistency for base %s",
598 } else {
599 nb_dofs = nn;
600 if (!nb_dofs)
602 }
603 }
604 const int nb_base_functions = data.getN().size2();
605 auto t_base = data.getFTensor0N();
606 auto t_diff_base = data.getFTensor1DiffN<2>();
607 auto t_t1 = get_ftensor1(tangent1_at_gauss_pts);
608 auto t_t2 = get_ftensor1(tangent2_at_gauss_pts);
609 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
610 auto t_data = data.getFTensor1FieldData<FIELD_DIM>();
611 int bb = 0;
612 for (; bb != nb_dofs / FIELD_DIM; ++bb) {
614 t_t1(i) += t_data(i) * t_diff_base(N0);
615 t_t2(i) += t_data(i) * t_diff_base(N1);
616 ++t_data;
617 ++t_base;
618 ++t_diff_base;
619 }
620 for (; bb != nb_base_functions; ++bb) {
621 ++t_base;
622 ++t_diff_base;
623 }
624 ++t_t1;
625 ++t_t2;
626 }
627 } break;
628 default:
629 SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
630 }
631 }
632
633 if (moab::CN::Dimension(type) == 2) {
634
635 auto &normal_at_gauss_pts = getNormalsAtGaussPts();
636 auto &tangent1_at_gauss_pts = getTangent1AtGaussPts();
637 auto &tangent2_at_gauss_pts = getTangent2AtGaussPts();
638
639 const auto nb_gauss_pts = tangent1_at_gauss_pts.size1();
640 normal_at_gauss_pts.resize(nb_gauss_pts, 3, false);
641
642 auto t_normal = get_ftensor1(normal_at_gauss_pts);
643 auto t_t1 = get_ftensor1(tangent1_at_gauss_pts);
644 auto t_t2 = get_ftensor1(tangent2_at_gauss_pts);
645 for (unsigned int gg = 0; gg != nb_gauss_pts; ++gg) {
649 t_normal(j) = FTensor::levi_civita(i, j, k) * t_t1(k) * t_t2(i);
650 ++t_normal;
651 ++t_t1;
652 ++t_t2;
653 }
654 }
655
657}
658
659template <int FIELD_DIM>
664
665 auto get_tangent = [&]() -> MatrixDouble & {
666 if (tangentsAtPts)
667 return *tangentsAtPts;
668 else
669 return getTangentAtGaussPts();
670 };
671
672 auto &tangent = get_tangent();
673 int nb_gauss_pts = getGaussPts().size2();
674
675 if (type == MBVERTEX) {
676 tangent.resize(nb_gauss_pts, 3, false);
677 tangent.clear();
678 }
679
680 int nb_dofs = data.getFieldData().size();
681 if (nb_dofs != 0) {
682 const int nb_base_functions = data.getN().size2();
683 double *diff_base_function = &*data.getDiffN().data().begin();
684 auto tangent_at_gauss_pts =
685 getFTensor1FromPtr<FIELD_DIM, 3>(&*tangent.data().begin());
686
688 int size = nb_dofs / FIELD_DIM;
689 if (nb_dofs % FIELD_DIM) {
690 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
691 "Inconsistent number of dofs and filed dimension");
692 }
693 for (int gg = 0; gg != nb_gauss_pts; ++gg) {
694 auto field_data = data.getFTensor1FieldData<FIELD_DIM>();
695 int bb = 0;
696 tangent_at_gauss_pts(i) = 0;
697 for (; bb < size; ++bb) {
698 tangent_at_gauss_pts(i) += field_data(i) * (*diff_base_function);
699 ++field_data;
700 ++diff_base_function;
701 }
702 for (; bb != nb_base_functions; ++bb) {
703 ++diff_base_function;
704 }
705 ++tangent_at_gauss_pts;
706 }
707 }
709}
710
711/**
712 * @deprecated do not use this function, instead use AddHOOps.
713 *
714 * @tparam E
715 * @param field
716 * @param e
717 * @param h1
718 * @param hcurl
719 * @param hdiv
720 * @param l2
721 * @return MoFEMErrorCode
722 */
723template <typename E>
724DEPRECATED MoFEMErrorCode addHOOpsVol(const std::string field, E &e, bool h1,
725 bool hcurl, bool hdiv, bool l2) {
726 std::vector<FieldSpace> spaces;
727 if (h1)
728 spaces.push_back(H1);
729 if (hcurl)
730 spaces.push_back(HCURL);
731 if (hdiv)
732 spaces.push_back(HDIV);
733 if (l2)
734 spaces.push_back(L2);
735 return AddHOOps<3, 3, 3>::add(e.getOpPtrVector(), spaces, field);
736}
737
738/**
739 * @deprecated do not use this function, instead use AddHOOps.
740 *
741 * @tparam E
742 * @param field
743 * @param e
744 * @param hcurl
745 * @param hdiv
746 * @return MoFEMErrorCode
747 */
748template <typename E>
749DEPRECATED MoFEMErrorCode addHOOpsFace3D(const std::string field, E &e,
750 bool hcurl, bool hdiv) {
751 std::vector<FieldSpace> spaces;
752 if (hcurl)
753 spaces.push_back(HCURL);
754 if (hdiv)
755 spaces.push_back(HDIV);
756 return AddHOOps<2, 3, 3>::add(e.getOpPtrVector(), spaces, field);
757}
758
759}; // namespace MoFEM
760
761#endif
static Number< 1 > N1
static Number< 0 > N0
ForcesAndSourcesCore::UserDataOperator UserDataOperator
constexpr int FIELD_DIM
FieldApproximationBase
approximation base
Definition: definitions.h:58
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
@ USER_BASE
user implemented approximation base
Definition: definitions.h:68
#define CHK_THROW_MESSAGE(err, msg)
Check and throw MoFEM exception.
Definition: definitions.h:595
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
FieldSpace
approximation spaces
Definition: definitions.h:82
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ H1
continuous field
Definition: definitions.h:85
@ NOSPACE
Definition: definitions.h:83
@ 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_IMPOSIBLE_CASE
Definition: definitions.h:35
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
static const char *const ApproximationBaseNames[]
Definition: definitions.h:72
#define DEPRECATED
Definition: definitions.h:17
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
FTensor::Index< 'm', SPACE_DIM > m
FTensor::Index< 'i', SPACE_DIM > i
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: MoFEM.hpp:24
DEPRECATED MoFEMErrorCode addHOOpsVol(const std::string field, E &e, bool h1, bool hcurl, bool hdiv, bool l2)
DEPRECATED MoFEMErrorCode addHOOpsFace3D(const std::string field, E &e, bool hcurl, bool hdiv)
double scale
Definition: plastic.cpp:94
constexpr auto field_name
Add operators pushing bases from local to physical configuration.
bool & doVertices
\deprectaed If false skip vertices
bool & doEdges
\deprectaed If false skip edges
Data on single entity (This is passed as argument to DataOperator::doWork)
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0N(const FieldApproximationBase base)
Get base function as Tensor0.
FieldApproximationBase & getBase()
Get approximation base.
MatrixDouble & getDiffN(const FieldApproximationBase base)
get derivatives of base functions
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
const VectorDouble & getFieldData() const
get dofs values
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1FieldData()
Return FTensor of rank 1, i.e. vector from filed data coeffinects.
const VectorDofs & getFieldDofs() const
get dofs data stature FEDofEntity
@ OPROW
operator doWork function is executed on FE rows
@ OPSPACE
operator do Work is execute on space data
structure to get information form mofem into EntitiesFieldData
Calculate HO coordinates at gauss points.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculateHOCoords(const std::string field_name)
boost::shared_ptr< MatrixDouble > jacPtr
Calculate jacobian for face element.
Calculate jacobian on Hex or other volume which is not simplex.
boost::shared_ptr< MatrixDouble > jacPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpCalculateHOJacForVolume(boost::shared_ptr< MatrixDouble > jac_ptr)
Calculate normals at Gauss points of triangle element.
OpGetHONormalsOnFace(std::string field_name)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Calculate tangent vector on edge form HO geometry approximation.
boost::shared_ptr< MatrixDouble > tangentsAtPts
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpGetHOTangentsOnEdge(std::string field_name, boost::shared_ptr< MatrixDouble > tangents_at_pts=nullptr)
transform Hcurl base fluxes from reference element to physical edge
boost::shared_ptr< MatrixDouble > tangentAtGaussPts
OpHOSetContravariantPiolaTransformOnEdge3D(const FieldSpace space=HCURL, boost::shared_ptr< MatrixDouble > tangent_at_pts=nullptr)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
transform Hdiv base fluxes from reference element to physical triangle
boost::shared_ptr< MatrixDouble > normalsAtGaussPts
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpHOSetContravariantPiolaTransformOnFace3D(const FieldSpace space, boost::shared_ptr< MatrixDouble > normals_at_gauss_pts=nullptr)
transform Hcurl base fluxes from reference element to physical triangle
OpHOSetCovariantPiolaTransformOnFace3D(const FieldSpace space, boost::shared_ptr< MatrixDouble > normals_at_pts=nullptr, boost::shared_ptr< MatrixDouble > tangent1_at_pts=nullptr, boost::shared_ptr< MatrixDouble > tangent2_at_pts=nullptr)
boost::shared_ptr< MatrixDouble > normalsAtPts
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< MatrixDouble > tangent2AtPts
boost::shared_ptr< MatrixDouble > tangent1AtPts
Scale base functions by inverses of measure of element.
OpScaleBaseBySpaceInverseOfMeasure(boost::shared_ptr< VectorDouble > det_jac_ptr, const FieldSpace space=L2)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
boost::shared_ptr< VectorDouble > detJacPtr
Apply contravariant (Piola) transfer to Hdiv space for HO geometry.
boost::shared_ptr< VectorDouble > detPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< MatrixDouble > jacPtr
OpSetHOContravariantPiolaTransform(const FieldSpace space, boost::shared_ptr< VectorDouble > det_ptr, boost::shared_ptr< MatrixDouble > jac_ptr)
Apply covariant (Piola) transfer to Hcurl space for HO geometry.
OpSetHOCovariantPiolaTransform(const FieldSpace space, boost::shared_ptr< MatrixDouble > jac_inv_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< MatrixDouble > jacInvPtr
Set inverse jacobian to base functions.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
transform local reference derivatives of shape function to global derivatives if higher order geometr...
boost::shared_ptr< MatrixDouble > invJacPtr
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpSetHOInvJacVectorBase(const FieldSpace space, boost::shared_ptr< MatrixDouble > inv_jac_ptr)
Set inverse jacobian to base functions.
OpSetHOWeights(boost::shared_ptr< VectorDouble > det_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
boost::shared_ptr< VectorDouble > detPtr
Modify integration weights on face to take in account higher-order geometry.
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Transform local reference derivatives of shape functions to global derivatives.
OpSetInvJacToScalarBasesBasic(FieldSpace space, boost::shared_ptr< MatrixDouble > inv_jac_ptr)