v0.8.23
DataOperators.hpp
Go to the documentation of this file.
1 /** \file DataOperators.hpp
2 
3  Number of structures preforming operations on integration points.
4 
5  For example:
6  - calculate Jacobian
7  - calculate Piola-Transform on shape functions
8  - calculate normals and tangent vectors at integration points on triangle
9 
10  */
11 
12 /* This file is part of MoFEM.
13  * MoFEM is free software: you can redistribute it and/or modify it under
14  * the terms of the GNU Lesser General Public License as published by the
15  * Free Software Foundation, either version 3 of the License, or (at your
16  * option) any later version.
17  *
18  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
19  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
21  * License for more details.
22  *
23  * You should have received a copy of the GNU Lesser General Public
24  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
25 
26 #ifndef __DATAOPERATORS_HPP
27 #define __DATAOPERATORS_HPP
28 
29 using namespace boost::numeric;
30 
31 namespace MoFEM {
32 
33 /** \brief base operator to do operations at Gauss Pt. level
34  * \ingroup mofem_forces_and_sources
35  */
36 struct DataOperator {
37 
38  DataOperator(const bool symm = true, const bool do_vertices = true,
39  const bool do_edges = true, const bool do_quads = true,
40  const bool do_tris = true, const bool do_tets = true,
41  const bool do_prisms = true)
42  : sYmm(symm), doVertices(do_vertices), doEdges(do_edges),
43  doQuads(do_quads), doTris(do_tris), doTets(do_tets),
44  doPrisms(do_prisms) {}
45 
46  virtual ~DataOperator() {}
47 
48  /** \brief Operator for bi-linear form, usually to calculate values on left
49  * hand side
50  */
51  virtual MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
52  EntityType col_type,
56  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
58  }
59 
60  virtual MoFEMErrorCode opLhs(DataForcesAndSourcesCore &row_data,
61  DataForcesAndSourcesCore &col_data,
62  bool symm = true);
63 
65  DataForcesAndSourcesCore &col_data) {
66  return opLhs(row_data, col_data, getSymm());
67  }
68 
69  /** \brief Operator for linear form, usually to calculate values on right hand
70  * side
71  */
72  virtual MoFEMErrorCode doWork(int side, EntityType type,
75  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
77  }
78 
79  virtual MoFEMErrorCode opRhs(DataForcesAndSourcesCore &data,
80  const bool do_vertices, const bool do_edges,
81  const bool do_quads, const bool do_tris,
82  const bool do_tets, const bool do_prisms,
83  const bool error_if_no_base = true);
84 
86  const bool error_if_no_base = true) {
87  return opRhs(data, doVertices, doEdges, doQuads, doTris, doTets, doPrisms,
88  error_if_no_base);
89  }
90 
91  bool sYmm; ///< If true assume that matrix is symmetric structure
92 
93  bool doVertices; ///< If false skip vertices
94  bool doEdges; ///< If false skip edges
95  bool doQuads;
96  bool doTris;
97  bool doTets;
98  bool doPrisms;
99 
100  /**
101  * \brief Get if operator uses symmetry of DOFs or not
102  *
103  * If symmetry is used, only not repeating combinations of entities are
104  * looped. For an example pair of (Vertex, Edge_0) and (Edge_0, Vertex) will
105  * calculate the same matrices only transposed. Implementing that this can be
106  * exploited by integrating only one pair.
107  *
108  * @return true if symmetry
109  */
110  inline bool getSymm() const { return sYmm; }
111 
112  /// set if operator is executed taking in account symmetry
113  inline void setSymm() { sYmm = true; }
114 
115  /// unset if operator is executed for non symmetric problem
116  inline void unSetSymm() { sYmm = false; }
117 };
118 
119 /**
120  * \brief Transform local reference derivatives of shape function to global
121  derivatives
122 
123  * \ingroup mofem_forces_and_sources
124  */
125 struct OpSetInvJacH1 : public DataOperator {
126 
130 
132  : tInvJac(&inv_jac(0, 0), &inv_jac(0, 1), &inv_jac(0, 2), &inv_jac(1, 0),
133  &inv_jac(1, 1), &inv_jac(1, 2), &inv_jac(2, 0), &inv_jac(2, 1),
134  &inv_jac(2, 2)) {}
135 
137  MoFEMErrorCode doWork(int side, EntityType type,
139 };
140 
141 /**
142  * \brief brief Transform local reference derivatives of shape function to
143  global derivatives
144 
145  * \ingroup mofem_forces_and_sources
146  */
148 
153 
155  : tInvJac(&inv_jac(0, 0), &inv_jac(0, 1), &inv_jac(0, 2), &inv_jac(1, 0),
156  &inv_jac(1, 1), &inv_jac(1, 2), &inv_jac(2, 0), &inv_jac(2, 1),
157  &inv_jac(2, 2)) {}
158 
160  MoFEMErrorCode doWork(int side, EntityType type,
162 };
163 
164 /**
165  * \brief transform local reference derivatives of shape function to global
166  derivatives if higher order geometry is given
167 
168  * \ingroup mofem_forces_and_sources
169 */
170 struct OpSetHoInvJacH1 : public DataOperator {
171 
175  OpSetHoInvJacH1(MatrixDouble &inv_ho_jac) : invHoJac(inv_ho_jac) {}
176 
178  MoFEMErrorCode doWork(int side, EntityType type,
180 };
181 
182 /**
183  * \brief transform local reference derivatives of shape function to global
184  derivatives if higher order geometry is given
185  *
186 
187  * \ingroup mofem_forces_and_sources
188 */
190 
195 
196  OpSetHoInvJacHdivAndHcurl(MatrixDouble &inv_ho_jac) : invHoJac(inv_ho_jac) {}
197 
199  MoFEMErrorCode doWork(int side, EntityType type,
201 };
202 
203 /** \brief apply contravariant (Piola) transfer to Hdiv space
204 
205 Contravariant Piola transformation
206 \f[
207 \psi_i|_t = \frac{1}{\textrm{det}(J)}J_{ij}\hat{\psi}_j\\
208 \left.\frac{\partial \psi_i}{\partial \xi_j}\right|_t
209 =
210 \frac{1}{\textrm{det}(J)}J_{ik}\frac{\partial \hat{\psi}_k}{\partial \xi_j}
211 \f]
212 
213 * \ingroup mofem_forces_and_sources
214 
215 */
217 
218  double &vOlume;
219  // MatrixDouble3by3 &jAc;
220 
225 
227  : vOlume(volume),
228  // jAc(jac),
229  tJac(&jac(0, 0), &jac(0, 1), &jac(0, 2), &jac(1, 0), &jac(1, 1),
230  &jac(1, 2), &jac(2, 0), &jac(2, 1), &jac(2, 2)) {}
231 
234 
235  MoFEMErrorCode doWork(int side, EntityType type,
237 };
238 
239 /** \brief Apply contravariant (Piola) transfer to Hdiv space for HO geometr
240 
241 * \ingroup mofem_forces_and_sources
242 */
244 
250 
252  : detHoJac(det_jac), hoJac(jac) {}
253 
256  MoFEMErrorCode doWork(int side, EntityType type,
258 };
259 
260 /** \brief Apply covariant (Piola) transfer to Hcurl space for HO geometry
261  */
263 
268 
269  OpSetHoCovariantPiolaTransform(MatrixDouble &inv_jac) : hoInvJac(inv_jac) {}
270 
273  MoFEMErrorCode doWork(int side, EntityType type,
275 };
276 
277 /** \brief apply covariant transfer to Hcurl space
278 
279 Contravariant Piola transformation
280 \f[
281 \psi_i|_t = \frac{1}{\textrm{det}(J)}J_{ij}\hat{\psi}_j\\
282 \left.\frac{\partial \psi_i}{\partial \xi_j}\right|_t
283 =
284 \frac{1}{\textrm{det}(J)}J_{ik}\frac{\partial \hat{\psi}_k}{\partial \xi_j}
285 \f]
286 
287 
288 * \ingroup mofem_forces_and_sources
289 */
291 
296 
298  : tInvJac(&inv_jac(0, 0), &inv_jac(0, 1), &inv_jac(0, 2), &inv_jac(1, 0),
299  &inv_jac(1, 1), &inv_jac(1, 2), &inv_jac(2, 0), &inv_jac(2, 1),
300  &inv_jac(2, 2)) {}
301 
304 
305  MoFEMErrorCode doWork(int side, EntityType type,
307 };
308 
309 /** \brief Get field values and gradients at Gauss points
310  * \ingroup mofem_forces_and_sources
311  */
312 template <int RANK, int DIM> struct OpGetDataAndGradient : public DataOperator {
313 
316 
318  MatrixDouble &data_grad_at_gauss_pt)
319  : dataAtGaussPts(data_at_gauss_pt),
320  dataGradAtGaussPts(data_grad_at_gauss_pt) {}
321 
322  /**
323  * Return tensor associated with matrix storing values
324  */
325  template <int R>
327  THROW_MESSAGE("Not implemented");
328  }
329 
330  /**
331  * Return tensor associated with matrix storing gradient values
332  */
333  template <int R, int D>
335  THROW_MESSAGE("Not implemented");
336  }
337 
338  /**
339  * \brief Calculate gradient and values at integration points
340  * @param side side of entity on element
341  * @param type type of entity
342  * @param data data stored on entity (dofs values, dofs indices, etc.)
343  * @return error code
344  */
345  MoFEMErrorCode calculateValAndGrad(int side, EntityType type,
348  const int nb_base_functions = data.getN().size2();
349  bool constant_diff = false;
350  if (type == MBVERTEX && data.getDiffN().size1() * data.getDiffN().size2() ==
351  DIM * nb_base_functions) {
352  constant_diff = true;
353  }
354  const int nb_dofs = data.getFieldData().size();
355  for (unsigned int gg = 0; gg < data.getN().size1(); gg++) {
356  double *data_ptr, *n_ptr, *diff_n_ptr;
357  n_ptr = &data.getN()(gg, 0);
358  if (constant_diff) {
359  diff_n_ptr = &data.getDiffN()(0, 0);
360  } else {
361  diff_n_ptr = &data.getDiffN()(gg, 0);
362  }
363  data_ptr = &*data.getFieldData().data().begin();
364  for (int rr = 0; rr < RANK; rr++, data_ptr++) {
365  dataAtGaussPts(gg, rr) +=
366  cblas_ddot(nb_dofs / RANK, n_ptr, 1, data_ptr, RANK);
367  double *diff_n_ptr2 = diff_n_ptr;
368  for (unsigned int dd = 0; dd < DIM; dd++, diff_n_ptr2++) {
369  dataGradAtGaussPts(gg, DIM * rr + dd) +=
370  cblas_ddot(nb_dofs / RANK, diff_n_ptr2, DIM, data_ptr, RANK);
371  }
372  }
373  }
375  }
376 
377  MoFEMErrorCode doWork(int side, EntityType type,
379 
381 
382  if (data.getFieldData().size() == 0) {
384  }
385 
386  unsigned int nb_dofs = data.getFieldData().size();
387  if (nb_dofs == 0)
389 
390  if (nb_dofs % RANK != 0) {
391  SETERRQ4(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
392  "data inconsistency, type %d, side %d, nb_dofs %d, rank %d",
393  type, side, nb_dofs, RANK);
394  }
395  if (nb_dofs / RANK > data.getN().size2()) {
396  std::cerr << side << " " << type << " "
397  << ApproximationBaseNames[data.getBase()] << std::endl;
398  std::cerr << data.getN() << std::endl;
399  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
400  "data inconsistency nb_dofs >= data.N.size2(), i.e. %u >= %u",
401  nb_dofs, data.getN().size2());
402  }
403 
404  if (type == MBVERTEX) {
405  dataAtGaussPts.resize(data.getN().size1(), RANK, false);
406  dataGradAtGaussPts.resize(data.getN().size1(), RANK * DIM, false);
407  dataAtGaussPts.clear();
408  dataGradAtGaussPts.clear();
409  }
410 
411  CHKERR calculateValAndGrad(side, type, data);
412 
414  }
415 };
416 
417 /**
418  * \brief Specialization for field with 3 coefficients in 3 dimension
419  */
420 template <>
421 template <>
423 OpGetDataAndGradient<3, 3>::getValAtGaussPtsTensor<3>(MatrixDouble &data);
424 
425 /**
426  * \brief Specialization for field with 3 coefficients in 3 dimension
427  */
428 template <>
429 template <>
431 OpGetDataAndGradient<3, 3>::getGradAtGaussPtsTensor<3, 3>(MatrixDouble &data);
432 
433 /**
434  * \brief Specialization for field with 3 coefficients in 3 dimension
435  */
436 template <>
438  int side, EntityType type, DataForcesAndSourcesCore::EntData &data);
439 
440 /**
441  * \brief Specialization for field with for scalar field in 3 dimension
442  */
443 template <>
444 MoFEMErrorCode OpGetDataAndGradient<1, 3>::calculateValAndGrad(
445  int side, EntityType type, DataForcesAndSourcesCore::EntData &data);
446 
447 /** \brief Calculate normals at Gauss points of triangle element
448  * \ingroup mofem_forces_and_source
449  */
451 
456 
458  MatrixDouble &normals_at_gausspt,
459  MatrixDouble &tangent1_at_gausspt,
460  MatrixDouble &tangent2_at_gausspt)
461  : cOords_at_GaussPt(coords_at_gausspt),
462  nOrmals_at_GaussPt(normals_at_gausspt),
463  tAngent1_at_GaussPt(tangent1_at_gausspt),
464  tAngent2_at_GaussPt(tangent2_at_gausspt) {}
465 
467  MoFEMErrorCode doWork(int side, EntityType type,
469 
470  MoFEMErrorCode calculateNormals();
471 };
472 
473 /** \brief calculate normals at Gauss points of triangle element
474  * \ingroup mofem_forces_and_sources
475  */
477 
486 
488  MatrixDouble &coords_at_gaussptf3, MatrixDouble &normals_at_gaussptf3,
489  MatrixDouble &tangent1_at_gaussptf3, MatrixDouble &tangent2_at_gaussptf3,
490  MatrixDouble &coords_at_gaussptf4, MatrixDouble &normals_at_gaussptf4,
491  MatrixDouble &tangent1_at_gaussptf4, MatrixDouble &tangent2_at_gaussptf4)
492  : cOords_at_GaussPtF3(coords_at_gaussptf3),
493  nOrmals_at_GaussPtF3(normals_at_gaussptf3),
494  tAngent1_at_GaussPtF3(tangent1_at_gaussptf3),
495  tAngent2_at_GaussPtF3(tangent2_at_gaussptf3),
496  cOords_at_GaussPtF4(coords_at_gaussptf4),
497  nOrmals_at_GaussPtF4(normals_at_gaussptf4),
498  tAngent1_at_GaussPtF4(tangent1_at_gaussptf4),
499  tAngent2_at_GaussPtF4(tangent2_at_gaussptf4) {}
500 
502  MoFEMErrorCode doWork(int side, EntityType type,
504 
505  MoFEMErrorCode calculateNormals();
506 };
507 
508 /** \brief transform Hdiv base fluxes from reference element to physical
509  * triangle \ingroup mofem_forces_and_sources
510  */
512 
515 
517  const MatrixDouble &normals_at_pts)
518  : nOrmal(normal), normalsAtGaussPts(normals_at_pts) {}
519 
520  MoFEMErrorCode doWork(int side, EntityType type,
522 };
523 
524 /** \brief transform Hcurl base fluxes from reference element to physical
525  * triangle \ingroup mofem_forces_and_sources
526  */
528 
535 
537  const MatrixDouble &normals_at_pts,
538  const VectorDouble &tangent0,
539  const MatrixDouble &tangent0_at_pts,
540  const VectorDouble &tangent1,
541  const MatrixDouble &tangent1_at_pts)
542  : nOrmal(normal), normalsAtGaussPts(normals_at_pts), tAngent0(tangent0),
543  tangent0AtGaussPt(tangent0_at_pts), tAngent1(tangent1),
544  tangent1AtGaussPt(tangent1_at_pts) {}
545 
546  MoFEMErrorCode doWork(int side, EntityType type,
548 };
549 
550 /** \brief Calculate tangent vector on edge form HO geometry approximation
551  * \ingroup mofem_forces_and_sources
552  */
554 
556 
557  OpGetHoTangentOnEdge(MatrixDouble &tangent) : tAngent(tangent) {}
558 
559  MoFEMErrorCode doWork(int side, EntityType type,
561 };
562 
563 /** \brief transform Hcurl base fluxes from reference element to physical edge
564  * \ingroup mofem_forces_and_sources
565  */
567 
570 
572  const MatrixDouble &tangent_at_pts)
573  : tAngent(tangent), tangentAtGaussPt(tangent_at_pts) {}
574 
575  MoFEMErrorCode doWork(int side, EntityType type,
577 };
578 
579 } // namespace MoFEM
580 
581 #endif //__DATAOPERATORS_HPP
582 
583 /***************************************************************************/ /**
584  * \defgroup mofem_forces_and_sources Forces and sources
585 ******************************************************************************/
Calculate tangent vector on edge form HO geometry approximation.
OpSetHoContravariantPiolaTransform(VectorDouble &det_jac, MatrixDouble &jac)
brief Transform local reference derivatives of shape function to global derivatives
MatrixDouble diffNinvJac
OpSetHoInvJacH1(MatrixDouble &inv_ho_jac)
transform Hdiv base fluxes from reference element to physical triangle
OpSetInvJacHdivAndHcurl(MatrixDouble3by3 &inv_jac)
base operator to do operations at Gauss Pt. level
FTensor::Index< 'j', 3 > j
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
FTensor::Index< 'i', 3 > i
FTensor::Index< 'i', 3 > i
void unSetSymm()
unset if operator is executed for non symmetric problem
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:500
OpGetDataAndGradient(MatrixDouble &data_at_gauss_pt, MatrixDouble &data_grad_at_gauss_pt)
OpSetHoInvJacHdivAndHcurl(MatrixDouble &inv_ho_jac)
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:77
transform Hcurl base fluxes from reference element to physical edge
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
OpSetContravariantPiolaTransform(double &volume, MatrixDouble3by3 &jac)
void setSymm()
set if operator is executed taking in account symmetry
FTensor::Tensor2< double *, 3, 3 > tInvJac
FTensor::Tensor2< double *, 3, 3 > tInvJac
virtual MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, DataForcesAndSourcesCore::EntData &row_data, DataForcesAndSourcesCore::EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
OpSetCovariantPiolaTransform(MatrixDouble3by3 &inv_jac)
FTensor::Tensor2< double *, 3, 3 > tInvJac
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:619
bool doVertices
If false skip vertices.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:507
virtual MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
OpGetCoordsAndNormalsOnFace(MatrixDouble &coords_at_gausspt, MatrixDouble &normals_at_gausspt, MatrixDouble &tangent1_at_gausspt, MatrixDouble &tangent2_at_gausspt)
FTensor::Index< 'i', 3 > i
virtual const MatrixDouble & getN(const FieldApproximationBase base) const
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
FTensor::Tensor2< double *, R, D > getGradAtGaussPtsTensor(MatrixDouble &data)
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
virtual MoFEMErrorCode opRhs(DataForcesAndSourcesCore &data, const bool error_if_no_base=true)
OpGetHoTangentOnEdge(MatrixDouble &tangent)
double cblas_ddot(const int N, const double *X, const int incX, const double *Y, const int incY)
Definition: cblas_ddot.c:12
FTensor::Index< 'i', 3 > i
static const char *const ApproximationBaseNames[]
Definition: definitions.h:157
virtual const MatrixDouble & getDiffN(const FieldApproximationBase base) const
get derivatives of base functions
Get field values and gradients at Gauss points.
FTensor::Index< 'j', 3 > j
apply contravariant (Piola) transfer to Hdiv space
FTensor::Index< 'k', 3 > k
bool sYmm
If true assume that matrix is symmetric structure.
bool getSymm() const
Get if operator uses symmetry of DOFs or not.
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
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
const VectorDouble & getFieldData() const
get dofs values
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
FTensor::Tensor1< double *, R > getValAtGaussPtsTensor(MatrixDouble &data)
DataForcesAndSourcesCore::EntData EntData
bool doEdges
If false skip edges.
FTensor::Index< 'k', 3 > k
DataOperator(const bool symm=true, const bool do_vertices=true, const bool do_edges=true, const bool do_quads=true, const bool do_tris=true, const bool do_tets=true, const bool do_prisms=true)
OpSetCovariantPiolaTransformOnEdge(const VectorDouble &tangent, const MatrixDouble &tangent_at_pts)
calculate normals at Gauss points of triangle element
OpSetInvJacH1(MatrixDouble3by3 &inv_jac)
OpGetCoordsAndNormalsOnPrism(MatrixDouble &coords_at_gaussptf3, MatrixDouble &normals_at_gaussptf3, MatrixDouble &tangent1_at_gaussptf3, MatrixDouble &tangent2_at_gaussptf3, MatrixDouble &coords_at_gaussptf4, MatrixDouble &normals_at_gaussptf4, MatrixDouble &tangent1_at_gaussptf4, MatrixDouble &tangent2_at_gaussptf4)
MatrixBoundedArray< double, 9 > MatrixDouble3by3
Definition: Types.hpp:98
Apply contravariant (Piola) transfer to Hdiv space for HO geometr.
transform local reference derivatives of shape function to global derivatives if higher order geometr...
#define CHKERR
Inline error check.
Definition: definitions.h:595
FTensor::Tensor2< double *, 3, 3 > tJac
Transform local reference derivatives of shape function to global derivatives.
FieldApproximationBase & getBase()
Get approximation base.
virtual MoFEMErrorCode opLhs(DataForcesAndSourcesCore &row_data, DataForcesAndSourcesCore &col_data)
Data on single entity (This is passed as argument to DataOperator::doWork)
OpSetHoCovariantPiolaTransform(MatrixDouble &inv_jac)
FTensor::Index< 'j', 3 > j
transform local reference derivatives of shape function to global derivatives if higher order geometr...
apply covariant transfer to Hcurl space
Calculate normals at Gauss points of triangle element.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406
transform Hcurl base fluxes from reference element to physical triangle
FTensor::Index< 'j', 3 > j
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:76
Apply covariant (Piola) transfer to Hcurl space for HO geometry.
OpSetContravariantPiolaTransformOnTriangle(const VectorDouble &normal, const MatrixDouble &normals_at_pts)
OpSetCovariantPiolaTransformOnTriangle(const VectorDouble &normal, const MatrixDouble &normals_at_pts, const VectorDouble &tangent0, const MatrixDouble &tangent0_at_pts, const VectorDouble &tangent1, const MatrixDouble &tangent1_at_pts)
MoFEMErrorCode calculateValAndGrad(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
Calculate gradient and values at integration points.