v0.8.23
Templates.hpp
Go to the documentation of this file.
1 /** \file Templates.hpp
2  * \brief Templates declarations
3  *
4  * MoFEM is free software: you can redistribute it and/or modify it under
5  * the terms of the GNU Lesser General Public License as published by the
6  * Free Software Foundation, either version 3 of the License, or (at your
7  * option) any later version.
8  *
9  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
11  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
12  * License for more details.
13  *
14  * You should have received a copy of the GNU Lesser General Public
15  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>
16  */
17 
18 #ifndef __TEMPLATES_HPP__
19 #define __TEMPLATES_HPP__
20 
21 namespace MoFEM {
22 
23 /**
24  * @brief Get Vector adaptor
25  *
26  * \code
27  *
28  * double *a;
29  * CHKERR VecGetArray(v,&a);
30  *
31  * for(int n = 0; n != nodes; ++n) {
32  *
33  * auto a = getVectorAdaptor(&a[3*n], 3);
34  * double dot = inner_prod(a, a);
35  *
36  * }
37  *
38  * CHKERR VecRetsoreArray(v,&a);
39  * \endcode
40  *
41  */
42 template <typename T1>
43 inline auto getVectorAdaptor(T1 ptr, const size_t n) {
44  typedef typename std::remove_pointer<T1>::type T;
46  ublas::shallow_array_adaptor<T>(n, ptr));
47 };
48 
49 /**
50  * @brief Get Matrix adaptor
51  *
52  * \code
53  *
54  * double *a;
55  * CHKERR VecGetArray(v,&a);
56  *
57  * for(int n = 0; n != nodes; ++n) {
58  *
59  * auto F = getMatrixAdaptor(&a[3*3*n], 3, 3);
60  * MatrixDouble C = prod(F, trans(F));
61  *
62  * }
63  *
64  * CHKERR VecRetsoreArray(v,&a);
65  * \endcode
66  *
67  */
68 template <typename T1>
69 inline auto getMatrixAdaptor(T1 ptr, const size_t n, const size_t m) {
70  typedef typename std::remove_pointer<T1>::type T;
72  n, m, ublas::shallow_array_adaptor<T>(n * m, ptr));
73 };
74 
75 /**
76  * This small utility that cascades two key extractors will be
77  * used throughout the boost example
78  * <a
79  * href=http://www.boost.org/doc/libs/1_53_0/libs/multi_index/example/complex_structs.cpp>
80  * http://www.boost.org/doc/libs/1_53_0/libs/multi_index/example/complex_structs.cpp
81  * </a>
82  */
83 template <class KeyExtractor1, class KeyExtractor2> struct KeyFromKey {
84 public:
85  typedef typename KeyExtractor1::result_type result_type;
86 
87  KeyFromKey(const KeyExtractor1 &key1_ = KeyExtractor1(),
88  const KeyExtractor2 &key2_ = KeyExtractor2())
89  : key1(key1_), key2(key2_) {}
90 
91  template <typename Arg> result_type operator()(Arg &arg) const {
92  return key1(key2(arg));
93  }
94 
95 private:
96  KeyExtractor1 key1;
97  KeyExtractor2 key2;
98 };
99 
100 template <typename id_type> struct LtBit {
101  inline bool operator()(const id_type &valueA, const id_type &valueB) const {
102  return valueA.to_ulong() < valueB.to_ulong();
103  }
104 };
105 
106 template <typename id_type> struct EqBit {
107  inline bool operator()(const id_type &valueA, const id_type &valueB) const {
108  return valueA.to_ulong() == valueB.to_ulong();
109  }
110 };
111 
112 template <typename id_type> struct HashBit {
113  inline unsigned int operator()(const id_type &value) const {
114  return value.to_ulong();
115  }
116 };
117 
118 template <class X> inline std::string toString(X x) {
119  std::ostringstream buffer;
120  buffer << x;
121  return buffer.str();
122 }
123 
124 /**
125 * \brief Get tensor rank 0 (scalar) form data vector
126 
127 Example how to use it.
128 \code
129 VectorDouble vec;
130 vec.resize(nb_gauss_pts,false);
131 vec.clear();
132 auto t0 = getFTensor0FromData(data);
133 for(int gg = 0;gg!=nb_gauss_pts;gg++) {
134 
135  ++t0;
136 }
137 \endcode
138 
139 */
140 template <class T, class A>
142 getFTensor0FromVec(ublas::vector<T, A> &data) {
143  static_assert(!std::is_same<T, T>::value, "not implemented");
145 }
146 
147 template <>
150  ublas::vector<double, DoubleAllocator> &data) {
151  return FTensor::Tensor0<FTensor::PackPtr<double *, 1>>(&*data.data().begin());
152 }
153 
154 /**
155  * \brief Get tensor rank 1 (vector) form data matrix
156  */
157 template <int Tensor_Dim, class T, class L, class A>
158 static inline FTensor::Tensor1<FTensor::PackPtr<T *, 1>, Tensor_Dim>
159 getFTensor1FromMat(ublas::matrix<T, L, A> &data) {
160  static_assert(!std::is_same<T, T>::value, "not implemented");
161 }
162 
163 /**
164  * \brief Get tensor rank 1 (vector) form data matrix (specialization)
165  */
166 template <int Tensor_Dim>
167 static inline FTensor::Tensor1<FTensor::PackPtr<double *, 1>, Tensor_Dim>
169  return getFTensor1FromMat<Tensor_Dim, double, ublas::row_major,
170  DoubleAllocator>(data);
171 }
172 
173 template <>
175 getFTensor1FromMat<3, double, ublas::row_major, DoubleAllocator>(
176  MatrixDouble &data) {
177  if (data.size1() != 3)
178  THROW_MESSAGE("getFTensor1FromMat<3>: wrong size of data matrix, number of "
179  "rows should be 3 but is %d" +
180  boost::lexical_cast<std::string>(data.size1()));
181 
183  &data(0, 0), &data(1, 0), &data(2, 0));
184 }
185 
186 template <>
188 getFTensor1FromMat<2, double, ublas::row_major, DoubleAllocator>(
189  MatrixDouble &data) {
190  if (data.size1() != 2)
191  THROW_MESSAGE("getFTensor1FromMat<2>: wrong size of data matrix, number of "
192  "rows should be 2 but is %d" +
193  boost::lexical_cast<std::string>(data.size1()));
194 
195  return FTensor::Tensor1<FTensor::PackPtr<double *, 1>, 2>(&data(0, 0),
196  &data(1, 0));
197 }
198 
199 /**
200  * \brief Get tensor rank 2 (matrix) form data matrix
201  */
202 template <int Tensor_Dim0, int Tensor_Dim1, class T, class L, class A>
203 static inline FTensor::Tensor2<FTensor::PackPtr<T *, 1>, Tensor_Dim0,
204  Tensor_Dim1>
205 getFTensor2FromMat(ublas::matrix<T, L, A> &data) {
206  static_assert(!std::is_same<T, T>::value,
207  "Such getFTensor2FromMat specialisation is not implemented");
208 }
209 
210 /**
211  * Template specialization for getFTensor2FromMat
212  *
213  */
214 template <>
217  if (data.size1() != 9)
218  THROW_MESSAGE("getFTensor2FromMat<3,3>: wrong size of data matrix; numer "
219  "of rows should be 9 but is " +
220  boost::lexical_cast<std::string>(data.size1()));
221 
223  &data(0, 0), &data(1, 0), &data(2, 0), &data(3, 0), &data(4, 0),
224  &data(5, 0), &data(6, 0), &data(7, 0), &data(8, 0));
225 }
226 
227 /**
228  * Template specialization for getFTensor2FromMat
229  */
230 template <>
233  if (data.size1() != 6)
234  THROW_MESSAGE("getFTensor2FromMat<3,3>: wrong size of data matrix, numer "
235  "of rows should be 6 but is " +
236  boost::lexical_cast<std::string>(data.size1()));
237 
239  &data(0, 0), &data(1, 0), &data(2, 0), &data(3, 0), &data(4, 0),
240  &data(5, 0));
241 }
242 
243 /**
244  * Template specialization for getFTensor2FromMat
245  */
246 template <>
249  if (data.size1() != 4)
250  THROW_MESSAGE("getFTensor2FromMat<2,2>: wrong size of data matrix, numer "
251  "of rows should be 4 but is " +
252  boost::lexical_cast<std::string>(data.size1()));
253 
255  &data(0, 0), &data(1, 0), &data(2, 0), &data(3, 0));
256 }
257 
258 /**
259  * \brief Get tensor rank 2 (matrix) form data matrix (specialization)
260  */
261 template <int Tensor_Dim0, int Tensor_Dim1>
262 static inline FTensor::Tensor2<FTensor::PackPtr<double *, 1>, Tensor_Dim0,
263  Tensor_Dim1>
265  return getFTensor2FromMat<Tensor_Dim0, Tensor_Dim1, double, ublas::row_major,
266  DoubleAllocator>(data);
267 }
268 
269 /**
270  * \brief Get symmetric tensor rank 2 (matrix) form data matrix
271  */
272 template <int Tensor_Dim, class T, class L, class A>
274 getFTensor2SymmetricFromMat(ublas::matrix<T, L, A> &data) {
275  static_assert(
276  !std::is_same<T, T>::value,
277  "Such getFTensor2SymmetricFromMat specialisation is not implemented");
278 }
279 
280 /**
281  * @brief Get symmetric tensor rank 2 form matrix of for dimension 3
282  *
283  * Specialisation for symmetric tensor 2
284  *
285  * @tparam
286  * @param data
287  * @return FTensor::Tensor2_symmetric<FTensor::PackPtr<double *, 1>, 3>
288  */
289 template <>
292  if (data.size1() != 6)
294  "getFTensor2SymmetricFromMat<3>: wrong size of data matrix, numer "
295  "of rows should be 6 but is " +
296  boost::lexical_cast<std::string>(data.size1()));
297 
299  &data(0, 0), &data(1, 0), &data(2, 0), &data(3, 0), &data(4, 0),
300  &data(5, 0));
301 }
302 
303 /**
304  * @brief Get symmetric tensor rank 2 form matrix
305  *
306  * Specialisation for symmetric tensor 2
307  *
308  * @tparam Tensor_Dim
309  * @param data
310  * @return FTensor::Tensor2_symmetric<FTensor::PackPtr<double *, 1>, Tensor_Dim>
311  */
312 template <int Tensor_Dim>
314  Tensor_Dim>
317  DoubleAllocator>(data);
318 }
319 
320 /**
321  * @brief Make Tensor1 from pointer
322  *
323  * @tparam DIM
324  * @param ptr
325  * @return FTensor::Tensor2<FTensor::PackPtr<double *, 3 * DIM>, 3, DIM>
326  */
327 template <int DIM>
329 getFTensor1FromPtr(double *ptr) {
330  static_assert(DIM != 3,
331  "Such getFTensor1FromPtr specialization is not implemented");
332 };
333 
334 template <>
336 getFTensor1FromPtr<3>(double *ptr) {
338  &ptr[HVEC0], &ptr[HVEC1], &ptr[HVEC2]);
339 };
340 
341 /**
342  * @brief Make Tensor2 from pointer
343  *
344  * @tparam DIM
345  * @param ptr
346  * @return FTensor::Tensor2<FTensor::PackPtr<double *, DIM1 * DIM2>, DIM1, DIM2>
347  */
348 template <int DIM1, int DIM2>
350 getFTensor2FromPtr(double *ptr) {
351  static_assert(DIM1 != 3, "Such getFTensor2FromPtr is not implemented");
352  static_assert(DIM2 >= 2 && DIM2 <= 3,
353  "Such getFTensor2FromPtr is not implemented");
354 };
355 
356 template <>
358 inline getFTensor2FromPtr<3, 2>(double *ptr) {
360  &ptr[HVEC0_0], &ptr[HVEC0_1],
361 
362  &ptr[HVEC1_0], &ptr[HVEC1_1],
363 
364  &ptr[HVEC2_0], &ptr[HVEC2_1]);
365 };
366 
367 template <>
369 inline getFTensor2FromPtr<3, 3>(double *ptr) {
371  &ptr[HVEC0_0], &ptr[HVEC0_1], &ptr[HVEC0_2],
372 
373  &ptr[HVEC1_0], &ptr[HVEC1_1], &ptr[HVEC1_2],
374 
375  &ptr[HVEC2_0], &ptr[HVEC2_1], &ptr[HVEC2_2]);
376 };
377 
378 /**
379  * @brief Calculate the determinant of a 3x3 matrix or a tensor of rank 2
380  *
381  * @tparam T
382  * @param t
383  * @return double
384  */
385 template <class T> static inline double dEterminant(T &t) {
386  return t(0, 0) * t(1, 1) * t(2, 2) + t(1, 0) * t(2, 1) * t(0, 2) +
387  t(2, 0) * t(0, 1) * t(1, 2) - t(0, 0) * t(2, 1) * t(1, 2) -
388  t(2, 0) * t(1, 1) * t(0, 2) - t(1, 0) * t(0, 1) * t(2, 2);
389 }
390 
391 /**
392  * \brief Calculate inverse of tensor rank 2 at integration points
393 
394  */
395 template <int Tensor_Dim, class T, class L, class A>
396 inline MoFEMErrorCode invertTensor3by3(ublas::matrix<T, L, A> &jac_data,
397  ublas::vector<T, A> &det_data,
398  ublas::matrix<T, L, A> &inv_jac_data) {
400  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
401  "Specialization for this template not yet implemented");
403 }
404 
405 template <>
406 inline MoFEMErrorCode
407 invertTensor3by3<3, double, ublas::row_major, DoubleAllocator>(
408  MatrixDouble &jac_data, VectorDouble &det_data, MatrixDouble &inv_jac_data);
409 
410 /**
411  * \brief Calculate determinant 3 by 3
412 
413  */
414 template <class T1, class T2>
415 inline MoFEMErrorCode determinantTensor3by3(T1 &t, T2 &det) {
417  det = +t(0, 0) * t(1, 1) * t(2, 2) + t(1, 0) * t(2, 1) * t(0, 2) +
418  t(2, 0) * t(0, 1) * t(1, 2) - t(0, 0) * t(2, 1) * t(1, 2) -
419  t(2, 0) * t(1, 1) * t(0, 2) - t(1, 0) * t(0, 1) * t(2, 2);
421 }
422 
423 /**
424  * \brief Calculate determinant 2 by 2
425 
426  */
427 template <class T1, class T2>
428 inline MoFEMErrorCode determinantTensor2by2(T1 &t, T2 &det) {
430  det = t(0, 0) * t(1, 1) - t(0, 1) * t(1, 0);
432 }
433 
434 /**
435  * \brief Calculate matrix inverse 3 by 3
436 
437  */
438 template <class T1, class T2, class T3>
439 inline MoFEMErrorCode invertTensor3by3(T1 &t, T2 &det, T3 &inv_t) {
441  const auto inv_det = 1. / det;
442  inv_t(0, 0) = (t(1, 1) * t(2, 2) - t(1, 2) * t(2, 1)) * inv_det;
443  inv_t(0, 1) = (t(0, 2) * t(2, 1) - t(0, 1) * t(2, 2)) * inv_det;
444  inv_t(0, 2) = (t(0, 1) * t(1, 2) - t(0, 2) * t(1, 1)) * inv_det;
445  inv_t(1, 0) = (t(1, 2) * t(2, 0) - t(1, 0) * t(2, 2)) * inv_det;
446  inv_t(1, 1) = (t(0, 0) * t(2, 2) - t(0, 2) * t(2, 0)) * inv_det;
447  inv_t(1, 2) = (t(0, 2) * t(1, 0) - t(0, 0) * t(1, 2)) * inv_det;
448  inv_t(2, 0) = (t(1, 0) * t(2, 1) - t(1, 1) * t(2, 0)) * inv_det;
449  inv_t(2, 1) = (t(0, 1) * t(2, 0) - t(0, 0) * t(2, 1)) * inv_det;
450  inv_t(2, 2) = (t(0, 0) * t(1, 1) - t(0, 1) * t(1, 0)) * inv_det;
452 }
453 
454 /**
455  * \brief Calculate matrix inverse 2 by 2
456 
457  */
458 template <class T1, class T2, class T3>
459 inline MoFEMErrorCode invertTensor2by2(T1 &t, T2 &det, T3 &inv_t) {
461  const auto inv_det = 1. / det;
462  inv_t(0, 0) = t(1, 1) * inv_det;
463  inv_t(0, 1) = -t(0, 1) * inv_det;
464  inv_t(1, 0) = -t(1, 0) * inv_det;
465  inv_t(1, 1) = t(0, 0) * inv_det;
467 }
468 
469 #ifdef WITH_ADOL_C
470 
471 /**
472  * \brief Calculate matrix inverse, specialization for adouble tensor
473 
474  */
475 template <>
476 inline MoFEMErrorCode invertTensor3by3<FTensor::Tensor2<adouble, 3, 3>, adouble,
481  inv_t(0, 0) = (t(1, 1) * t(2, 2) - t(1, 2) * t(2, 1)) / det;
482  inv_t(0, 1) = (t(0, 2) * t(2, 1) - t(0, 1) * t(2, 2)) / det;
483  inv_t(0, 2) = (t(0, 1) * t(1, 2) - t(0, 2) * t(1, 1)) / det;
484  inv_t(1, 0) = (t(1, 2) * t(2, 0) - t(1, 0) * t(2, 2)) / det;
485  inv_t(1, 1) = (t(0, 0) * t(2, 2) - t(0, 2) * t(2, 0)) / det;
486  inv_t(1, 2) = (t(0, 2) * t(1, 0) - t(0, 0) * t(1, 2)) / det;
487  inv_t(2, 0) = (t(1, 0) * t(2, 1) - t(1, 1) * t(2, 0)) / det;
488  inv_t(2, 1) = (t(0, 1) * t(2, 0) - t(0, 0) * t(2, 1)) / det;
489  inv_t(2, 2) = (t(0, 0) * t(1, 1) - t(0, 1) * t(1, 0)) / det;
491 }
492 
493 #endif
494 
495 /**
496  * \brief Calculate matrix inverse, specialization for symmetric tensor
497 
498  */
499 template <>
500 inline MoFEMErrorCode
501 invertTensor3by3<FTensor::Tensor2_symmetric<double, 3>, double,
506  const auto inv_det = 1. / det;
507  inv_t(0, 0) = (t(1, 1) * t(2, 2) - t(1, 2) * t(2, 1)) * inv_det;
508  inv_t(0, 1) = (t(0, 2) * t(2, 1) - t(0, 1) * t(2, 2)) * inv_det;
509  inv_t(0, 2) = (t(0, 1) * t(1, 2) - t(0, 2) * t(1, 1)) * inv_det;
510  inv_t(1, 1) = (t(0, 0) * t(2, 2) - t(0, 2) * t(2, 0)) * inv_det;
511  inv_t(1, 2) = (t(0, 2) * t(1, 0) - t(0, 0) * t(1, 2)) * inv_det;
512  inv_t(2, 2) = (t(0, 0) * t(1, 1) - t(0, 1) * t(1, 0)) * inv_det;
514 }
515 
516 #ifdef WITH_ADOL_C
517 
518 /**
519  * \brief Calculate matrix inverse, specialization for adouble symmetric tensor
520 
521  */
522 template <>
523 inline MoFEMErrorCode
524 invertTensor3by3<FTensor::Tensor2_symmetric<adouble, 3>, adouble,
529  inv_t(0, 0) = (t(1, 1) * t(2, 2) - t(1, 2) * t(2, 1)) / det;
530  inv_t(0, 1) = (t(0, 2) * t(2, 1) - t(0, 1) * t(2, 2)) / det;
531  inv_t(0, 2) = (t(0, 1) * t(1, 2) - t(0, 2) * t(1, 1)) / det;
532  inv_t(1, 1) = (t(0, 0) * t(2, 2) - t(0, 2) * t(2, 0)) / det;
533  inv_t(1, 2) = (t(0, 2) * t(1, 0) - t(0, 0) * t(1, 2)) / det;
534  inv_t(2, 2) = (t(0, 0) * t(1, 1) - t(0, 1) * t(1, 0)) / det;
536 }
537 
538 #endif
539 
540 /**
541  * \brief Calculate matrix inverse, specialization for symmetric (pointer)
542  tensor
543 
544  */
545 template <>
546 inline MoFEMErrorCode
547 invertTensor3by3<FTensor::Tensor2_symmetric<double, 3>, double,
552  const auto inv_det = 1. / det;
553  inv_t(0, 0) = (t(1, 1) * t(2, 2) - t(1, 2) * t(2, 1)) * inv_det;
554  inv_t(0, 1) = (t(0, 2) * t(2, 1) - t(0, 1) * t(2, 2)) * inv_det;
555  inv_t(0, 2) = (t(0, 1) * t(1, 2) - t(0, 2) * t(1, 1)) * inv_det;
556  inv_t(1, 1) = (t(0, 0) * t(2, 2) - t(0, 2) * t(2, 0)) * inv_det;
557  inv_t(1, 2) = (t(0, 2) * t(1, 0) - t(0, 0) * t(1, 2)) * inv_det;
558  inv_t(2, 2) = (t(0, 0) * t(1, 1) - t(0, 1) * t(1, 0)) * inv_det;
560 }
561 
562 /**
563  * @brief Extract entity handle form multi-index container
564  *
565  */
567  template <typename Iterator>
568  static inline EntityHandle extract(const Iterator &it) {
569  return (*it)->getRefEnt();
570  }
571 };
572 
573 /**
574  * @brief Insert ordered mofem multi-index into range
575  *
576  * \code
577  * auto hi_rit = refEntsPtr->upper_bound(start);
578  * auto hi_rit = refEntsPtr->upper_bound(end);
579  * Range to_erase;
580  * insertOrdered(to_erase, RefEntExtractor(), rit, hi_rit);
581  * \endcode
582  *
583  * @tparam Iterator
584  * @param r
585  * @param begin_iter
586  * @param end_iter
587  * @return moab::Range::iterator
588  */
589 template <typename Extractor, typename Iterator>
590 moab::Range::iterator insertOrdered(Range &r, Extractor, Iterator begin_iter,
591  Iterator end_iter) {
592  moab::Range::iterator hint = r.begin();
593  while (begin_iter != end_iter) {
594  size_t j = 0;
595  auto bi = Extractor::extract(begin_iter);
596  Iterator pj = begin_iter;
597  while (pj != end_iter && (bi + j) == Extractor::extract(pj)) {
598  ++pj;
599  ++j;
600  }
601  hint = r.insert(hint, bi, bi + (j - 1));
602  begin_iter = pj;
603  }
604  return hint;
605 };
606 
607 } // namespace MoFEM
608 
609 #endif //__TEMPLATES_HPP__
static FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:142
FTensor::Tensor2< FTensor::PackPtr< double *, DIM1 *DIM2 >, DIM1, DIM2 > getFTensor2FromPtr(double *ptr)
Make Tensor2 from pointer.
Definition: Templates.hpp:350
MoFEMErrorCode determinantTensor3by3(T1 &t, T2 &det)
Calculate determinant 3 by 3.
Definition: Templates.hpp:415
static EntityHandle extract(const Iterator &it)
Definition: Templates.hpp:568
FTensor::Tensor1< FTensor::PackPtr< double *, DIM >, DIM > getFTensor1FromPtr(double *ptr)
Make Tensor1 from pointer.
Definition: Templates.hpp:329
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:501
ublas::matrix< double, ublas::row_major, ublas::shallow_array_adaptor< double > > MatrixShallowArrayAdaptor
Definition: Types.hpp:107
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > getFTensor0FromVec< double, DoubleAllocator >(ublas::vector< double, DoubleAllocator > &data)
Definition: Templates.hpp:149
KeyFromKey(const KeyExtractor1 &key1_=KeyExtractor1(), const KeyExtractor2 &key2_=KeyExtractor2())
Definition: Templates.hpp:87
#define THROW_MESSAGE(a)
Throw MoFEM exception.
Definition: definitions.h:620
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:508
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
KeyExtractor1 key1
Definition: Templates.hpp:96
bool operator()(const id_type &valueA, const id_type &valueB) const
Definition: Templates.hpp:107
moab::Range::iterator insertOrdered(Range &r, Extractor, Iterator begin_iter, Iterator end_iter)
Insert ordered mofem multi-index into range.
Definition: Templates.hpp:590
MoFEMErrorCode determinantTensor2by2(T1 &t, T2 &det)
Calculate determinant 2 by 2.
Definition: Templates.hpp:428
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.
Definition: Templates.hpp:396
Extract entity handle form multi-index container.
Definition: Templates.hpp:566
KeyExtractor2 key2
Definition: Templates.hpp:97
FTensor::Tensor1< FTensor::PackPtr< double *, 3 >, 3 > getFTensor1FromPtr< 3 >(double *ptr)
Definition: Templates.hpp:336
unsigned int operator()(const id_type &value) const
Definition: Templates.hpp:113
FTensor::Tensor2< FTensor::PackPtr< double *, 9 >, 3, 3 > getFTensor2FromPtr< 3, 3 >(double *ptr)
Definition: Templates.hpp:369
MoFEMErrorCode invertTensor2by2(T1 &t, T2 &det, T3 &inv_t)
Calculate matrix inverse 2 by 2.
Definition: Templates.hpp:459
static double dEterminant(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
Definition: Templates.hpp:385
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
static FTensor::Tensor1< FTensor::PackPtr< T *, 1 >, Tensor_Dim > getFTensor1FromMat(ublas::matrix< T, L, A > &data)
Get tensor rank 1 (vector) form data matrix.
Definition: Templates.hpp:159
auto getVectorAdaptor(T1 ptr, const size_t n)
Get Vector adaptor.
Definition: Templates.hpp:43
std::string toString(X x)
Definition: Templates.hpp:118
bool operator()(const id_type &valueA, const id_type &valueB) const
Definition: Templates.hpp:101
FTensor::Tensor2< FTensor::PackPtr< double *, 6 >, 3, 2 > getFTensor2FromPtr< 3, 2 >(double *ptr)
Definition: Templates.hpp:358
ublas::vector< T, ublas::shallow_array_adaptor< T > > VectorShallowArrayAdaptor
Definition: Types.hpp:100
auto getMatrixAdaptor(T1 ptr, const size_t n, const size_t m)
Get Matrix adaptor.
Definition: Templates.hpp:69
result_type operator()(Arg &arg) const
Definition: Templates.hpp:91
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72
static FTensor::Tensor2< FTensor::PackPtr< T *, 1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2FromMat(ublas::matrix< T, L, A > &data)
Get tensor rank 2 (matrix) form data matrix.
Definition: Templates.hpp:205
KeyExtractor1::result_type result_type
Definition: Templates.hpp:85
std::vector< double, std::allocator< double > > DoubleAllocator
Definition: Types.hpp:70
static FTensor::Tensor2_symmetric< FTensor::PackPtr< T *, 1 >, Tensor_Dim > getFTensor2SymmetricFromMat(ublas::matrix< T, L, A > &data)
Get symmetric tensor rank 2 (matrix) form data matrix.
Definition: Templates.hpp:274