v0.9.1
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);
39 
40  virtual ~DataOperator() = default;
41 
42  /** \brief Operator for bi-linear form, usually to calculate values on left
43  * hand side
44  */
45  virtual MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
46  EntityType col_type,
50  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
52  }
53 
54  virtual MoFEMErrorCode opLhs(DataForcesAndSourcesCore &row_data,
55  DataForcesAndSourcesCore &col_data);
56 
57  /** \brief Operator for linear form, usually to calculate values on right hand
58  * side
59  */
60  virtual MoFEMErrorCode doWork(int side, EntityType type,
63  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented");
65  }
66 
67  virtual MoFEMErrorCode opRhs(DataForcesAndSourcesCore &data,
68  const bool error_if_no_base = false);
69 
70  bool sYmm; ///< If true assume that matrix is symmetric structure
71 
72  std::array<bool, MBMAXTYPE>
73  doEntities; ///< If true operator is executed for entity.
74 
75  // Deprecated variables. Use doEntities instead. I keep them for back
76  // compatibility with some older modules. It will be removed in some future.
77 
78  bool &doVertices; ///< \deprectaed If false skip vertices
79  bool &doEdges; ///< \deprectaed If false skip edges
80  bool &doQuads; ///< \deprectaed
81  bool &doTris; ///< \deprectaed
82  bool &doTets; ///< \deprectaed
83  bool &doPrisms; ///< \deprectaed
84 
85  /**
86  * \brief Get if operator uses symmetry of DOFs or not
87  *
88  * If symmetry is used, only not repeating combinations of entities are
89  * looped. For an example pair of (Vertex, Edge_0) and (Edge_0, Vertex) will
90  * calculate the same matrices only transposed. Implementing that this can be
91  * exploited by integrating only one pair.
92  *
93  * @return true if symmetry
94  */
95  inline bool getSymm() const { return sYmm; }
96 
97  /// set if operator is executed taking in account symmetry
98  inline void setSymm() { sYmm = true; }
99 
100  /// unset if operator is executed for non symmetric problem
101  inline void unSetSymm() { sYmm = false; }
102 
103 private:
104 
105  template <bool Symm>
106  inline MoFEMErrorCode opLhs(DataForcesAndSourcesCore &row_data,
107  DataForcesAndSourcesCore &col_data);
108 
109  template <bool ErrorIfNoBase>
110  inline MoFEMErrorCode opRhs(DataForcesAndSourcesCore &data,
111  const std::array<bool, MBMAXTYPE> &do_entities);
112 };
113 
114 /**
115  * \brief Transform local reference derivatives of shape function to global
116  derivatives
117 
118  * \ingroup mofem_forces_and_sources
119  */
120 struct OpSetInvJacH1 : public DataOperator {
121 
125 
127  : tInvJac(&inv_jac(0, 0), &inv_jac(0, 1), &inv_jac(0, 2), &inv_jac(1, 0),
128  &inv_jac(1, 1), &inv_jac(1, 2), &inv_jac(2, 0), &inv_jac(2, 1),
129  &inv_jac(2, 2)) {}
130 
132  MoFEMErrorCode doWork(int side, EntityType type,
134 };
135 
136 /**
137  * \brief brief Transform local reference derivatives of shape function to
138  global derivatives
139 
140  * \ingroup mofem_forces_and_sources
141  */
143 
148 
150  : tInvJac(&inv_jac(0, 0), &inv_jac(0, 1), &inv_jac(0, 2), &inv_jac(1, 0),
151  &inv_jac(1, 1), &inv_jac(1, 2), &inv_jac(2, 0), &inv_jac(2, 1),
152  &inv_jac(2, 2)) {}
153 
155  MoFEMErrorCode doWork(int side, EntityType type,
157 };
158 
159 /**
160  * \brief transform local reference derivatives of shape function to global
161  derivatives if higher order geometry is given
162 
163  * \ingroup mofem_forces_and_sources
164 */
165 struct OpSetHoInvJacH1 : public DataOperator {
166 
170  OpSetHoInvJacH1(MatrixDouble &inv_ho_jac) : invHoJac(inv_ho_jac) {}
171 
173  MoFEMErrorCode doWork(int side, EntityType type,
175 };
176 
177 /**
178  * \brief transform local reference derivatives of shape function to global
179  derivatives if higher order geometry is given
180  *
181 
182  * \ingroup mofem_forces_and_sources
183 */
185 
190 
191  OpSetHoInvJacHdivAndHcurl(MatrixDouble &inv_ho_jac) : invHoJac(inv_ho_jac) {}
192 
194  MoFEMErrorCode doWork(int side, EntityType type,
196 };
197 
198 /** \brief apply contravariant (Piola) transfer to Hdiv space
199 
200 Contravariant Piola transformation
201 \f[
202 \psi_i|_t = \frac{1}{\textrm{det}(J)}J_{ij}\hat{\psi}_j\\
203 \left.\frac{\partial \psi_i}{\partial \xi_j}\right|_t
204 =
205 \frac{1}{\textrm{det}(J)}J_{ik}\frac{\partial \hat{\psi}_k}{\partial \xi_j}
206 \f]
207 
208 * \ingroup mofem_forces_and_sources
209 
210 */
212 
213  double &vOlume;
214  // MatrixDouble3by3 &jAc;
215 
220 
222  : vOlume(volume),
223  // jAc(jac),
224  tJac(&jac(0, 0), &jac(0, 1), &jac(0, 2), &jac(1, 0), &jac(1, 1),
225  &jac(1, 2), &jac(2, 0), &jac(2, 1), &jac(2, 2)) {}
226 
229 
230  MoFEMErrorCode doWork(int side, EntityType type,
232 };
233 
234 /** \brief Apply contravariant (Piola) transfer to Hdiv space for HO geometr
235 
236 * \ingroup mofem_forces_and_sources
237 */
239 
245 
247  : detHoJac(det_jac), hoJac(jac) {}
248 
251  MoFEMErrorCode doWork(int side, EntityType type,
253 };
254 
255 /** \brief Apply covariant (Piola) transfer to Hcurl space for HO geometry
256  */
258 
263 
264  OpSetHoCovariantPiolaTransform(MatrixDouble &inv_jac) : hoInvJac(inv_jac) {}
265 
268  MoFEMErrorCode doWork(int side, EntityType type,
270 };
271 
272 /** \brief apply covariant transfer to Hcurl space
273 
274 Contravariant Piola transformation
275 \f[
276 \psi_i|_t = \frac{1}{\textrm{det}(J)}J_{ij}\hat{\psi}_j\\
277 \left.\frac{\partial \psi_i}{\partial \xi_j}\right|_t
278 =
279 \frac{1}{\textrm{det}(J)}J_{ik}\frac{\partial \hat{\psi}_k}{\partial \xi_j}
280 \f]
281 
282 
283 * \ingroup mofem_forces_and_sources
284 */
286 
291 
293  : tInvJac(&inv_jac(0, 0), &inv_jac(0, 1), &inv_jac(0, 2), &inv_jac(1, 0),
294  &inv_jac(1, 1), &inv_jac(1, 2), &inv_jac(2, 0), &inv_jac(2, 1),
295  &inv_jac(2, 2)) {}
296 
299 
300  MoFEMErrorCode doWork(int side, EntityType type,
302 };
303 
304 /** \brief Get field values and gradients at Gauss points
305  * \ingroup mofem_forces_and_sources
306  */
307 template <int RANK, int DIM> struct OpGetDataAndGradient : public DataOperator {
308 
311 
313  MatrixDouble &data_grad_at_gauss_pt)
314  : dataAtGaussPts(data_at_gauss_pt),
315  dataGradAtGaussPts(data_grad_at_gauss_pt) {}
316 
317  /**
318  * Return tensor associated with matrix storing values
319  */
320  template <int R>
322  THROW_MESSAGE("Not implemented");
323  }
324 
325  /**
326  * Return tensor associated with matrix storing gradient values
327  */
328  template <int R, int D>
330  THROW_MESSAGE("Not implemented");
331  }
332 
333  /**
334  * \brief Calculate gradient and values at integration points
335  * @param side side of entity on element
336  * @param type type of entity
337  * @param data data stored on entity (dofs values, dofs indices, etc.)
338  * @return error code
339  */
343  const int nb_base_functions = data.getN().size2();
344  bool constant_diff = false;
345  if (type == MBVERTEX && data.getDiffN().size1() * data.getDiffN().size2() ==
346  DIM * nb_base_functions) {
347  constant_diff = true;
348  }
349  const int nb_dofs = data.getFieldData().size();
350  for (unsigned int gg = 0; gg < data.getN().size1(); gg++) {
351  double *data_ptr, *n_ptr, *diff_n_ptr;
352  n_ptr = &data.getN()(gg, 0);
353  if (constant_diff) {
354  diff_n_ptr = &data.getDiffN()(0, 0);
355  } else {
356  diff_n_ptr = &data.getDiffN()(gg, 0);
357  }
358  data_ptr = &*data.getFieldData().data().begin();
359  for (int rr = 0; rr < RANK; rr++, data_ptr++) {
360  dataAtGaussPts(gg, rr) +=
361  cblas_ddot(nb_dofs / RANK, n_ptr, 1, data_ptr, RANK);
362  double *diff_n_ptr2 = diff_n_ptr;
363  for (unsigned int dd = 0; dd < DIM; dd++, diff_n_ptr2++) {
364  dataGradAtGaussPts(gg, DIM * rr + dd) +=
365  cblas_ddot(nb_dofs / RANK, diff_n_ptr2, DIM, data_ptr, RANK);
366  }
367  }
368  }
370  }
371 
372  MoFEMErrorCode doWork(int side, EntityType type,
374 
376 
377  if (data.getFieldData().size() == 0) {
379  }
380 
381  unsigned int nb_dofs = data.getFieldData().size();
382  if (nb_dofs == 0)
384 
385  if (nb_dofs % RANK != 0) {
386  SETERRQ4(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
387  "data inconsistency, type %d, side %d, nb_dofs %d, rank %d",
388  type, side, nb_dofs, RANK);
389  }
390  if (nb_dofs / RANK > data.getN().size2()) {
391  std::cerr << side << " " << type << " "
392  << ApproximationBaseNames[data.getBase()] << std::endl;
393  std::cerr << data.getN() << std::endl;
394  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
395  "data inconsistency nb_dofs >= data.N.size2(), i.e. %u >= %u",
396  nb_dofs, data.getN().size2());
397  }
398 
399  if (type == MBVERTEX) {
400  dataAtGaussPts.resize(data.getN().size1(), RANK, false);
401  dataGradAtGaussPts.resize(data.getN().size1(), RANK * DIM, false);
402  dataAtGaussPts.clear();
403  dataGradAtGaussPts.clear();
404  }
405 
406  CHKERR calculateValAndGrad(side, type, data);
407 
409  }
410 };
411 
412 /**
413  * \brief Specialization for field with 3 coefficients in 3 dimension
414  */
415 template <>
416 template <>
418 OpGetDataAndGradient<3, 3>::getValAtGaussPtsTensor<3>(MatrixDouble &data);
419 
420 /**
421  * \brief Specialization for field with 3 coefficients in 3 dimension
422  */
423 template <>
424 template <>
426 OpGetDataAndGradient<3, 3>::getGradAtGaussPtsTensor<3, 3>(MatrixDouble &data);
427 
428 /**
429  * \brief Specialization for field with 3 coefficients in 3 dimension
430  */
431 template <>
433  int side, EntityType type, DataForcesAndSourcesCore::EntData &data);
434 
435 /**
436  * \brief Specialization for field with for scalar field in 3 dimension
437  */
438 template <>
439 MoFEMErrorCode OpGetDataAndGradient<1, 3>::calculateValAndGrad(
440  int side, EntityType type, DataForcesAndSourcesCore::EntData &data);
441 
442 /** \brief Calculate normals at Gauss points of triangle element
443  * \ingroup mofem_forces_and_source
444  */
446 
451 
453  MatrixDouble &normals_at_gausspt,
454  MatrixDouble &tangent1_at_gausspt,
455  MatrixDouble &tangent2_at_gausspt)
456  : cOords_at_GaussPt(coords_at_gausspt),
457  nOrmals_at_GaussPt(normals_at_gausspt),
458  tAngent1_at_GaussPt(tangent1_at_gausspt),
459  tAngent2_at_GaussPt(tangent2_at_gausspt) {}
460 
462  MoFEMErrorCode doWork(int side, EntityType type,
464 
465  MoFEMErrorCode calculateNormals();
466 };
467 
468 /** \brief calculate normals at Gauss points of triangle element
469  * \ingroup mofem_forces_and_sources
470  */
472 
481 
483  MatrixDouble &coords_at_gaussptf3, MatrixDouble &normals_at_gaussptf3,
484  MatrixDouble &tangent1_at_gaussptf3, MatrixDouble &tangent2_at_gaussptf3,
485  MatrixDouble &coords_at_gaussptf4, MatrixDouble &normals_at_gaussptf4,
486  MatrixDouble &tangent1_at_gaussptf4, MatrixDouble &tangent2_at_gaussptf4)
487  : cOords_at_GaussPtF3(coords_at_gaussptf3),
488  nOrmals_at_GaussPtF3(normals_at_gaussptf3),
489  tAngent1_at_GaussPtF3(tangent1_at_gaussptf3),
490  tAngent2_at_GaussPtF3(tangent2_at_gaussptf3),
491  cOords_at_GaussPtF4(coords_at_gaussptf4),
492  nOrmals_at_GaussPtF4(normals_at_gaussptf4),
493  tAngent1_at_GaussPtF4(tangent1_at_gaussptf4),
494  tAngent2_at_GaussPtF4(tangent2_at_gaussptf4) {}
495 
497  MoFEMErrorCode doWork(int side, EntityType type,
499 
500  MoFEMErrorCode calculateNormals();
501 };
502 
503 /** \brief transform Hdiv base fluxes from reference element to physical
504  * triangle \ingroup mofem_forces_and_sources
505  */
507 
510 
512  const MatrixDouble &normals_at_pts)
513  : nOrmal(normal), normalsAtGaussPts(normals_at_pts) {}
514 
515  MoFEMErrorCode doWork(int side, EntityType type,
517 };
518 
519 /** \brief transform Hcurl base fluxes from reference element to physical
520  * triangle \ingroup mofem_forces_and_sources
521  */
523 
530 
532  const MatrixDouble &normals_at_pts,
533  const VectorDouble &tangent0,
534  const MatrixDouble &tangent0_at_pts,
535  const VectorDouble &tangent1,
536  const MatrixDouble &tangent1_at_pts)
537  : nOrmal(normal), normalsAtGaussPts(normals_at_pts), tAngent0(tangent0),
538  tangent0AtGaussPt(tangent0_at_pts), tAngent1(tangent1),
539  tangent1AtGaussPt(tangent1_at_pts) {}
540 
541  MoFEMErrorCode doWork(int side, EntityType type,
543 };
544 
545 /** \brief Calculate tangent vector on edge form HO geometry approximation
546  * \ingroup mofem_forces_and_sources
547  */
549 
551 
552  OpGetHoTangentOnEdge(MatrixDouble &tangent) : tAngent(tangent) {}
553 
554  MoFEMErrorCode doWork(int side, EntityType type,
556 };
557 
558 /** \brief transform Hcurl base fluxes from reference element to physical edge
559  * \ingroup mofem_forces_and_sources
560  */
562 
565 
567  const MatrixDouble &tangent_at_pts)
568  : tAngent(tangent), tangentAtGaussPt(tangent_at_pts) {}
569 
570  MoFEMErrorCode doWork(int side, EntityType type,
572 };
573 
574 } // namespace MoFEM
575 
576 #endif //__DATAOPERATORS_HPP
577 
578 /**
579  * \defgroup mofem_forces_and_sources Forces and sources
580  **/
bool & doTris
\deprectaed
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)
FieldApproximationBase & getBase()
Get approximation base.
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,...
bool & doVertices
\deprectaed If false skip vertices
transform Hcurl base fluxes from reference element to physical triangle
FTensor::Index< 'i', 3 > i
FTensor::Index< 'i', 3 > i
void unSetSymm()
unset if operator is executed for non symmetric problem
bool & doPrisms
\deprectaed
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
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:74
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:482
OpSetContravariantPiolaTransform(double &volume, MatrixDouble3by3 &jac)
void setSymm()
set if operator is executed taking in account symmetry
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
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:625
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
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
FTensor::Tensor2< double *, R, D > getGradAtGaussPtsTensor(MatrixDouble &data)
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
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:163
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
DataForcesAndSourcesCore::EntData EntData
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
FTensor::Tensor1< double *, R > getValAtGaussPtsTensor(MatrixDouble &data)
OpSetContravariantPiolaTransformOnFace(const VectorDouble &normal, const MatrixDouble &normals_at_pts)
FTensor::Index< 'k', 3 > k
transform Hdiv base fluxes from reference element to physical triangle
bool & doEdges
\deprectaed If false skip edges
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:97
Apply contravariant (Piola) transfer to Hdiv space for HO geometr.
bool & doQuads
\deprectaed
transform local reference derivatives of shape function to global derivatives if higher order geometr...
#define CHKERR
Inline error check.
Definition: definitions.h:601
FTensor::Tensor2< double *, 3, 3 > tJac
Transform local reference derivatives of shape function to global derivatives.
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:412
FTensor::Index< 'j', 3 > j
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
Apply covariant (Piola) transfer to Hcurl space for HO geometry.
const VectorDouble & getFieldData() const
get dofs values
OpSetCovariantPiolaTransformOnFace(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.
bool & doTets
\deprectaed
MatrixDouble & getDiffN(const FieldApproximationBase base)
get derivatives of base functions