v0.14.0
HODataOperators.cpp
Go to the documentation of this file.
1 /** \file HODataOperators.cpp
2 
3 \brief Set of operators for high-order geometry approximation.
4 
5 */
6 
7 namespace MoFEM {
8 
10  boost::shared_ptr<MatrixDouble> jac_ptr)
12  jacPtr(jac_ptr) {
13 
14  for (auto t = MBEDGE; t != MBMAXTYPE; ++t)
15  doEntities[t] = false;
16 
17  if (!jacPtr)
18  THROW_MESSAGE("Jac pointer not allocated");
19 }
20 
25 
26  const auto nb_base_functions = data.getN(NOBASE).size2();
27  if (nb_base_functions) {
28 
29  const auto nb_gauss_pts = getGaussPts().size2();
30  auto t_diff_base = data.getFTensor1DiffN<3>(NOBASE);
31  auto &coords = getCoords();
32 
33 #ifndef NDEBUG
34  if (nb_gauss_pts != data.getDiffN().size1())
35  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
36  "Inconsistent number base functions and gauss points");
37  if (nb_base_functions != data.getDiffN().size2() / 3)
38  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
39  "Inconsistent number of base functions");
40  if (coords.size() != 3 * nb_base_functions)
41  SETERRQ2(
42  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
43  "Number of vertex coordinates and base functions is inconsistent "
44  "%d != %d",
45  coords.size() / 3, nb_base_functions);
46 #endif
47 
48  jacPtr->resize(9, nb_gauss_pts, false);
49  jacPtr->clear();
50  auto t_jac = getFTensor2FromMat<3, 3>(*jacPtr);
51  auto t_vol_inv_jac = getInvJac();
52 
56 
57  for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
59  &coords[0], &coords[1], &coords[2]);
60  for (size_t bb = 0; bb != nb_base_functions; ++bb) {
61  t_jac(i, j) += t_coords(i) * (t_vol_inv_jac(k, j) * t_diff_base(k));
62  ++t_diff_base;
63  ++t_coords;
64  }
65  ++t_jac;
66  }
67  }
68 
70 }
71 
76 
77  if (getFEDim() == 3) {
78 
79  auto transform_base = [&](MatrixDouble &diff_n) {
81  return applyTransform<3, 3, 3, 3>(diff_n);
83  };
84 
85  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
86  FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
87  CHKERR transform_base(data.getDiffN(base));
88  }
89 
90  switch (type) {
91  case MBVERTEX:
92  for (auto &m : data.getBBDiffNMap())
93  CHKERR transform_base(*(m.second));
94  break;
95  default:
96  for (auto &ptr : data.getBBDiffNByOrderArray())
97  if (ptr)
98  CHKERR transform_base(*ptr);
99  }
100 
101  } else {
102  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Use different operator");
103  }
104 
106 }
107 
112 
116 
117  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
118 
119  FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
120 
121  diffHdivInvJac.resize(data.getDiffN(base).size1(),
122  data.getDiffN(base).size2(), false);
123 
124  unsigned int nb_gauss_pts = data.getDiffN(base).size1();
125  unsigned int nb_base_functions = data.getDiffN(base).size2() / 9;
126  if (nb_base_functions == 0)
127  continue;
128 
129  auto t_diff_n = data.getFTensor2DiffN<3, 3>(base);
130  double *t_inv_diff_n_ptr = &*diffHdivInvJac.data().begin();
132  t_inv_diff_n_ptr, &t_inv_diff_n_ptr[HVEC0_1],
133  &t_inv_diff_n_ptr[HVEC0_2], &t_inv_diff_n_ptr[HVEC1_0],
134  &t_inv_diff_n_ptr[HVEC1_1], &t_inv_diff_n_ptr[HVEC1_2],
135  &t_inv_diff_n_ptr[HVEC2_0], &t_inv_diff_n_ptr[HVEC2_1],
136  &t_inv_diff_n_ptr[HVEC2_2]);
137  auto t_inv_jac = getFTensor2FromMat<3, 3>(*invJacPtr);
138  for (unsigned int gg = 0; gg != nb_gauss_pts; ++gg) {
139  for (unsigned int bb = 0; bb != nb_base_functions; ++bb) {
140  t_inv_diff_n(i, j) = t_inv_jac(k, j) * t_diff_n(i, k);
141  ++t_diff_n;
142  ++t_inv_diff_n;
143  }
144  ++t_inv_jac;
145  }
146 
147  data.getDiffN(base).swap(diffHdivInvJac);
148  }
149 
151 }
152 
156  const size_t nb_int_pts = getGaussPts().size2();
157  if (getNormalsAtGaussPts().size1()) {
158  if (getNormalsAtGaussPts().size1() == nb_int_pts) {
159  double a = getMeasure();
160  if (getNumeredEntFiniteElementPtr()->getEntType() == MBTRI)
161  a *= 2;
162  auto t_w = getFTensor0IntegrationWeight();
163  auto t_normal = getFTensor1NormalsAtGaussPts();
165  for (size_t gg = 0; gg != nb_int_pts; ++gg) {
166  t_w *= sqrt(t_normal(i) * t_normal(i)) / a;
167  ++t_w;
168  ++t_normal;
169  }
170  } else {
171  SETERRQ2(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE,
172  "Number of rows in getNormalsAtGaussPts should be equal to "
173  "number of integration points, but is not, i.e. %d != %d",
174  getNormalsAtGaussPts().size1(), nb_int_pts);
175  }
176  }
178 }
179 
183  const size_t nb_int_pts = getGaussPts().size2();
184  if (getTangentAtGaussPts().size1()) {
185  if (getTangentAtGaussPts().size1() == nb_int_pts) {
186  double a = getMeasure();
187  auto t_w = getFTensor0IntegrationWeight();
188  auto t_tangent = getFTensor1TangentAtGaussPts();
190  for (size_t gg = 0; gg != nb_int_pts; ++gg) {
191  t_w *= sqrt(t_tangent(i) * t_tangent(i)) / a;
192  ++t_w;
193  ++t_tangent;
194  }
195  } else {
196  SETERRQ2(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE,
197  "Number of rows in getTangentAtGaussPts should be equal to "
198  "number of integration points, but is not, i.e. %d != %d",
199  getTangentAtGaussPts().size1(), nb_int_pts);
200  }
201  }
203 }
204 
208 
209  const auto nb_integration_pts = detPtr->size();
210 
211 #ifndef NDEBUG
212  if (nb_integration_pts != getGaussPts().size2())
213  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
214  "Inconsistent number of data points");
215 #endif
216 
217  auto t_w = getFTensor0IntegrationWeight();
218  auto t_det = getFTensor0FromVec(*detPtr);
219  for (size_t gg = 0; gg != nb_integration_pts; ++gg) {
220  t_w *= t_det;
221  ++t_w;
222  ++t_det;
223  }
224 
226 }
227 
232 
236 
237 #ifndef NDEBUG
238  if (!detPtr)
239  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
240  "Pointer for detPtr not allocated");
241  if (!jacPtr)
242  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
243  "Pointer for jacPtr not allocated");
244 #endif
245 
246  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; ++b) {
247 
248  FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
249 
250  auto nb_gauss_pts = data.getNSharedPtr(base) ? data.getN(base).size1() : 0;
251  auto nb_base_functions =
252  data.getNSharedPtr(base) ? data.getN(base).size2() / 3 : 0;
253  auto nb_diff_base_functions =
254  data.getDiffNSharedPtr(base) ? data.getDiffN(base).size2() / 9 : 0;
255 
256 #ifndef NDEBUG
257  if (nb_diff_base_functions) {
258  if (nb_diff_base_functions != nb_base_functions)
259  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
260  "Wrong number base functions %d != %d", nb_diff_base_functions,
261  nb_base_functions);
262  if (data.getDiffN(base).size1() != nb_gauss_pts)
263  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
264  "Wrong number integration points");
265  }
266 #endif
267 
268  if (nb_gauss_pts && nb_base_functions) {
269 
270  piolaN.resize(nb_gauss_pts, data.getN(base).size2(), false);
271 
272  auto t_n = data.getFTensor1N<3>(base);
273  double *t_transformed_n_ptr = &*piolaN.data().begin();
275  t_transformed_n_ptr, // HVEC0
276  &t_transformed_n_ptr[HVEC1], &t_transformed_n_ptr[HVEC2]);
277 
278  auto t_det = getFTensor0FromVec(*detPtr);
279  auto t_jac = getFTensor2FromMat<3, 3>(*jacPtr);
280 
281  for (unsigned int gg = 0; gg != nb_gauss_pts; ++gg) {
282  for (unsigned int bb = 0; bb != nb_base_functions; ++bb) {
283  const double a = 1. / t_det;
284  t_transformed_n(i) = a * t_jac(i, k) * t_n(k);
285  ++t_n;
286  ++t_transformed_n;
287  }
288  ++t_det;
289  ++t_jac;
290  }
291 
292  data.getN(base).swap(piolaN);
293  }
294 
295  if (nb_gauss_pts && nb_diff_base_functions) {
296 
297  piolaDiffN.resize(nb_gauss_pts, data.getDiffN(base).size2(), false);
298 
299  auto t_diff_n = data.getFTensor2DiffN<3, 3>(base);
300  double *t_transformed_diff_n_ptr = &*piolaDiffN.data().begin();
302  t_transformed_diff_n(t_transformed_diff_n_ptr,
303  &t_transformed_diff_n_ptr[HVEC0_1],
304  &t_transformed_diff_n_ptr[HVEC0_2],
305  &t_transformed_diff_n_ptr[HVEC1_0],
306  &t_transformed_diff_n_ptr[HVEC1_1],
307  &t_transformed_diff_n_ptr[HVEC1_2],
308  &t_transformed_diff_n_ptr[HVEC2_0],
309  &t_transformed_diff_n_ptr[HVEC2_1],
310  &t_transformed_diff_n_ptr[HVEC2_2]);
311 
312  auto t_det = getFTensor0FromVec(*detPtr);
313  auto t_jac = getFTensor2FromMat<3, 3>(*jacPtr);
314 
315  for (unsigned int gg = 0; gg != nb_gauss_pts; ++gg) {
316  for (unsigned int bb = 0; bb != nb_diff_base_functions; ++bb) {
317  const double a = 1. / t_det;
318  t_transformed_diff_n(i, k) = a * t_jac(i, j) * t_diff_n(j, k);
319  ++t_diff_n;
320  ++t_transformed_diff_n;
321  }
322  ++t_det;
323  ++t_jac;
324  }
325 
326  data.getDiffN(base).swap(piolaDiffN);
327  }
328  }
329 
331 }
332 
337 
341 
342  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
343 
344  FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
345 
346  unsigned int nb_gauss_pts = data.getN(base).size1();
347  unsigned int nb_base_functions = data.getN(base).size2() / 3;
348  piolaN.resize(nb_gauss_pts, data.getN(base).size2(), false);
349  piolaDiffN.resize(nb_gauss_pts, data.getDiffN(base).size2(), false);
350 
351  auto t_n = data.getFTensor1N<3>(base);
352  double *t_transformed_n_ptr = &*piolaN.data().begin();
354  t_transformed_n_ptr, // HVEC0
355  &t_transformed_n_ptr[HVEC1], &t_transformed_n_ptr[HVEC2]);
356  auto t_diff_n = data.getFTensor2DiffN<3, 3>(base);
357  double *t_transformed_diff_n_ptr = &*piolaDiffN.data().begin();
358  FTensor::Tensor2<FTensor::PackPtr<double *, 9>, 3, 3> t_transformed_diff_n(
359  t_transformed_diff_n_ptr, &t_transformed_diff_n_ptr[HVEC0_1],
360  &t_transformed_diff_n_ptr[HVEC0_2], &t_transformed_diff_n_ptr[HVEC1_0],
361  &t_transformed_diff_n_ptr[HVEC1_1], &t_transformed_diff_n_ptr[HVEC1_2],
362  &t_transformed_diff_n_ptr[HVEC2_0], &t_transformed_diff_n_ptr[HVEC2_1],
363  &t_transformed_diff_n_ptr[HVEC2_2]);
364 
365  auto t_inv_jac = getFTensor2FromMat<3, 3>(*jacInvPtr);
366 
367  for (unsigned int gg = 0; gg != nb_gauss_pts; ++gg) {
368  for (unsigned int bb = 0; bb != nb_base_functions; ++bb) {
369  t_transformed_n(i) = t_inv_jac(k, i) * t_n(k);
370  t_transformed_diff_n(i, k) = t_inv_jac(j, i) * t_diff_n(j, k);
371  ++t_n;
372  ++t_transformed_n;
373  ++t_diff_n;
374  ++t_transformed_diff_n;
375  }
376  ++t_inv_jac;
377  }
378 
379  data.getN(base).swap(piolaN);
380  data.getDiffN(base).swap(piolaDiffN);
381  }
382 
384 }
385 
387  boost::shared_ptr<MatrixDouble> jac_ptr)
389  jacPtr(jac_ptr) {}
390 
394 
396 
397  const auto nb_gauss_pts = getGaussPts().size2();
398 
399  jacPtr->resize(4, nb_gauss_pts, false);
400  auto t_jac = getFTensor2FromMat<2, 2>(*jacPtr);
401 
404 
405  auto t_t1 = getFTensor1Tangent1AtGaussPts();
406  auto t_t2 = getFTensor1Tangent2AtGaussPts();
409 
410  for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
411  t_jac(i, N0) = t_t1(i);
412  t_jac(i, N1) = t_t2(i);
413  ++t_t1;
414  ++t_t2;
415  ++t_jac;
416  }
417 
419 };
420 
425 
429 
430  size_t nb_gauss_pts = getGaussPts().size2();
431  jacPtr->resize(9, nb_gauss_pts, false);
432 
433  auto t_t1 = getFTensor1Tangent1AtGaussPts();
434  auto t_t2 = getFTensor1Tangent2AtGaussPts();
435  auto t_normal = getFTensor1NormalsAtGaussPts();
436 
440 
441  auto t_jac = getFTensor2FromMat<3, 3>(*jacPtr);
442  for (size_t gg = 0; gg != nb_gauss_pts; ++gg, ++t_jac) {
443 
444  t_jac(i, N0) = t_t1(i);
445  t_jac(i, N1) = t_t2(i);
446 
447  const double l = sqrt(t_normal(j) * t_normal(j));
448 
449  t_jac(i, N2) = t_normal(i) / l;
450 
451  ++t_t1;
452  ++t_t2;
453  ++t_normal;
454  }
455 
457 }
458 
460  int side, EntityType type, EntitiesFieldData::EntData &data) {
463 
464  if (moab::CN::Dimension(type) != 2)
466 
467  auto get_normals_ptr = [&]() {
468  if (normalsAtGaussPts)
469  return &*normalsAtGaussPts->data().begin();
470  else
471  return &*getNormalsAtGaussPts().data().begin();
472  };
473 
474  auto apply_transform_nonlinear_geometry = [&](auto base, auto nb_gauss_pts,
475  auto nb_base_functions) {
477 
478  auto ptr = get_normals_ptr();
480  &ptr[0], &ptr[1], &ptr[2]);
481 
482  auto t_base = data.getFTensor1N<3>(base);
483  for (int gg = 0; gg != nb_gauss_pts; ++gg) {
484  const auto l2 = t_normal(i) * t_normal(i);
485  for (int bb = 0; bb != nb_base_functions; ++bb) {
486  const auto v = t_base(0);
487  t_base(i) = (v / l2) * t_normal(i);
488  ++t_base;
489  }
490  ++t_normal;
491  }
493  };
494 
495  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
496 
497  FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
498  const auto &base_functions = data.getN(base);
499  const auto nb_gauss_pts = base_functions.size1();
500 
501  if (nb_gauss_pts) {
502 
503  const auto nb_base_functions = base_functions.size2() / 3;
504  CHKERR apply_transform_nonlinear_geometry(base, nb_gauss_pts,
505  nb_base_functions);
506  }
507  }
508 
510 }
511 
513  int side, EntityType type, EntitiesFieldData::EntData &data) {
515 
516  const auto type_dim = moab::CN::Dimension(type);
517  if (type_dim != 1 && type_dim != 2)
519 
523 
524  auto get_jac = [&]() {
526  double *ptr_n = &*normalsAtPts->data().begin();
527  double *ptr_t1 = &*tangent1AtPts->data().begin();
528  double *ptr_t2 = &*tangent2AtPts->data().begin();
530  &ptr_t1[0], &ptr_t2[0], &ptr_n[0],
531 
532  &ptr_t1[1], &ptr_t2[1], &ptr_n[1],
533 
534  &ptr_t1[2], &ptr_t2[2], &ptr_n[2]);
535  } else {
536  double *ptr_n = &*getNormalsAtGaussPts().data().begin();
537  double *ptr_t1 = &*getTangent1AtGaussPts().data().begin();
538  double *ptr_t2 = &*getTangent2AtGaussPts().data().begin();
540  &ptr_t1[0], &ptr_t2[0], &ptr_n[0],
541 
542  &ptr_t1[1], &ptr_t2[1], &ptr_n[1],
543 
544  &ptr_t1[2], &ptr_t2[2], &ptr_n[2]);
545  }
546  };
547 
548  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; ++b) {
549 
550  FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
551 
552  auto &baseN = data.getN(base);
553  auto &diffBaseN = data.getDiffN(base);
554 
555  int nb_dofs = baseN.size2() / 3;
556  int nb_gauss_pts = baseN.size1();
557 
558  piolaN.resize(baseN.size1(), baseN.size2());
559  diffPiolaN.resize(diffBaseN.size1(), diffBaseN.size2());
560 
561  if (nb_dofs > 0 && nb_gauss_pts > 0) {
562 
563  auto t_h_curl = data.getFTensor1N<3>(base);
564  auto t_diff_h_curl = data.getFTensor2DiffN<3, 2>(base);
565 
566  FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3> t_transformed_h_curl(
567  &piolaN(0, HVEC0), &piolaN(0, HVEC1), &piolaN(0, HVEC2));
568 
570  t_transformed_diff_h_curl(
573  &diffPiolaN(0, HVEC2_0), &diffPiolaN(0, HVEC2_1));
574 
575  int cc = 0;
576  double det;
578 
579  // HO geometry is set, so jacobian is different at each gauss point
580  auto t_m_at_pts = get_jac();
581  for (int gg = 0; gg != nb_gauss_pts; ++gg) {
582  CHKERR determinantTensor3by3(t_m_at_pts, det);
583  CHKERR invertTensor3by3(t_m_at_pts, det, t_inv_m);
584  for (int ll = 0; ll != nb_dofs; ll++) {
585  t_transformed_h_curl(i) = t_inv_m(j, i) * t_h_curl(j);
586  t_transformed_diff_h_curl(i, k) = t_inv_m(j, i) * t_diff_h_curl(j, k);
587  ++t_h_curl;
588  ++t_transformed_h_curl;
589  ++t_diff_h_curl;
590  ++t_transformed_diff_h_curl;
591  ++cc;
592  }
593  ++t_m_at_pts;
594  }
595 
596 #ifndef NDEBUG
597  if (cc != nb_gauss_pts * nb_dofs)
598  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "Data inconsistency");
599 #endif
600 
601  baseN.swap(piolaN);
602  diffBaseN.swap(diffPiolaN);
603  }
604  }
605 
607 }
608 
610  int side, EntityType type, EntitiesFieldData::EntData &data) {
613 
614  if (type != MBEDGE)
616 
617  const auto nb_gauss_pts = getGaussPts().size2();
618 
619  if (tangentAtGaussPts)
620  if (tangentAtGaussPts->size1() != nb_gauss_pts)
621  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
622  "Wrong number of integration points %d!=%d",
623  tangentAtGaussPts->size1(), nb_gauss_pts);
624 
625  auto get_base_at_pts = [&](auto base) {
627  &data.getN(base)(0, HVEC0), &data.getN(base)(0, HVEC1),
628  &data.getN(base)(0, HVEC2));
629  return t_h_curl;
630  };
631 
632  auto get_tangent_ptr = [&]() {
633  if (tangentAtGaussPts) {
634  return &*tangentAtGaussPts->data().begin();
635  } else {
636  return &*getTangentAtGaussPts().data().begin();
637  }
638  };
639 
640  auto get_tangent_at_pts = [&]() {
641  auto ptr = get_tangent_ptr();
643  &ptr[0], &ptr[1], &ptr[2]);
644  return t_m_at_pts;
645  };
646 
647  auto calculate_squared_edge_length = [&]() {
648  std::vector<double> l1;
649  if (nb_gauss_pts) {
650  l1.resize(nb_gauss_pts);
651  auto t_m_at_pts = get_tangent_at_pts();
652  for (size_t gg = 0; gg != nb_gauss_pts; ++gg) {
653  l1[gg] = t_m_at_pts(i) * t_m_at_pts(i);
654  ++t_m_at_pts;
655  }
656  }
657  return l1;
658  };
659 
660  auto l1 = calculate_squared_edge_length();
661 
662  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
663 
664  FieldApproximationBase base = static_cast<FieldApproximationBase>(b);
665  const auto nb_dofs = data.getN(base).size2() / 3;
666 
667  if (nb_gauss_pts && nb_dofs) {
668  auto t_h_curl = get_base_at_pts(base);
669  int cc = 0;
670  auto t_m_at_pts = get_tangent_at_pts();
671  for (int gg = 0; gg != nb_gauss_pts; ++gg) {
672  const double l0 = l1[gg];
673  for (int ll = 0; ll != nb_dofs; ++ll) {
674  const double val = t_h_curl(0);
675  const double a = val / l0;
676  t_h_curl(i) = t_m_at_pts(i) * a;
677  ++t_h_curl;
678  ++cc;
679  }
680  ++t_m_at_pts;
681  }
682 
683 #ifndef NDEBUG
684  if (cc != nb_gauss_pts * nb_dofs)
685  SETERRQ(PETSC_COMM_SELF, MOFEM_IMPOSSIBLE_CASE, "Data inconsistency");
686 #endif
687  }
688  }
689 
691 }
692 
694  boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pipeline,
695  std::vector<FieldSpace> spaces, std::string geom_field_name,
696  boost::shared_ptr<MatrixDouble> jac_ptr,
697  boost::shared_ptr<VectorDouble> det_ptr,
698  boost::shared_ptr<MatrixDouble> inv_jac_ptr) {
700 
701  if (!jac_ptr)
702  jac_ptr = boost::make_shared<MatrixDouble>();
703  if (!det_ptr)
704  det_ptr = boost::make_shared<VectorDouble>();
705  if (!inv_jac_ptr)
706  inv_jac_ptr = boost::make_shared<MatrixDouble>();
707 
708  if (geom_field_name.empty()) {
709 
710  pipeline.push_back(new OpCalculateHOJac<2>(jac_ptr));
711 
712  } else {
713 
714  pipeline.push_back(new OpCalculateHOCoords<2>(geom_field_name));
715  pipeline.push_back(
716  new OpCalculateVectorFieldGradient<2, 2>(geom_field_name, jac_ptr));
717  pipeline.push_back(new OpGetHONormalsOnFace<2>(geom_field_name));
718  }
719 
720  pipeline.push_back(new OpInvertMatrix<2>(jac_ptr, det_ptr, inv_jac_ptr));
721  pipeline.push_back(new OpSetHOWeightsOnFace());
722 
723  for (auto s : spaces) {
724  switch (s) {
725  case NOSPACE:
726  break;
727  case H1:
728  case L2:
729  pipeline.push_back(new OpSetHOInvJacToScalarBases<2>(s, inv_jac_ptr));
730  break;
731  case HCURL:
732  pipeline.push_back(new OpSetCovariantPiolaTransformOnFace2D(inv_jac_ptr));
733  pipeline.push_back(new OpSetInvJacHcurlFace(inv_jac_ptr));
734  break;
735  case HDIV:
736  pipeline.push_back(new OpMakeHdivFromHcurl());
737  pipeline.push_back(new OpSetContravariantPiolaTransformOnFace2D(jac_ptr));
738  pipeline.push_back(new OpSetInvJacHcurlFace(inv_jac_ptr));
739  break;
740  default:
741  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
742  "Space %s not yet implemented", FieldSpaceNames[s]);
743  }
744  }
745 
747 }
748 
750  boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pipeline,
751  std::vector<FieldSpace> spaces, std::string geom_field_name,
752  boost::shared_ptr<MatrixDouble> jac_ptr,
753  boost::shared_ptr<VectorDouble> det_ptr,
754  boost::shared_ptr<MatrixDouble> inv_jac_ptr) {
756 
757  if (!jac_ptr)
758  jac_ptr = boost::make_shared<MatrixDouble>();
759  if (!det_ptr)
760  det_ptr = boost::make_shared<VectorDouble>();
761  if (!inv_jac_ptr)
762  inv_jac_ptr = boost::make_shared<MatrixDouble>();
763 
764  if (geom_field_name.empty()) {
765 
766  } else {
767 
768  pipeline.push_back(new OpCalculateHOCoords<3>(geom_field_name));
769  pipeline.push_back(new OpGetHONormalsOnFace<3>(geom_field_name));
770  }
771 
772  pipeline.push_back(new OpCalculateHOJacForFaceEmbeddedIn3DSpace(jac_ptr));
773  pipeline.push_back(new OpInvertMatrix<3>(jac_ptr, det_ptr, inv_jac_ptr));
774  pipeline.push_back(new OpSetHOWeightsOnFace());
775 
776  for (auto s : spaces) {
777  switch (s) {
778  case NOSPACE:
779  break;
780  case H1:
781  pipeline.push_back(
782  new OpSetInvJacH1ForFaceEmbeddedIn3DSpace(inv_jac_ptr));
783  break;
784  case HDIV:
785  pipeline.push_back(new OpMakeHdivFromHcurl());
786  pipeline.push_back(
788  jac_ptr));
789  pipeline.push_back(
790  new OpSetInvJacHcurlFaceEmbeddedIn3DSpace(inv_jac_ptr));
791  break;
792  case L2:
793  pipeline.push_back(
794  new OpSetInvJacL2ForFaceEmbeddedIn3DSpace(inv_jac_ptr));
795  break;
796  default:
797  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
798  "Space %s not yet implemented", FieldSpaceNames[s]);
799  }
800  }
801 
803 }
804 
806  boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pipeline,
807  std::vector<FieldSpace> spaces, std::string geom_field_name) {
809 
810  if (geom_field_name.empty()) {
811 
812  } else {
813 
814  pipeline.push_back(new OpCalculateHOCoords<2>(geom_field_name));
815  pipeline.push_back(new OpGetHOTangentsOnEdge<2>(geom_field_name));
816  }
817 
818  for (auto s : spaces) {
819  switch (s) {
820  case NOSPACE:
821  break;
822  case HCURL:
823  pipeline.push_back(new OpHOSetContravariantPiolaTransformOnEdge3D(HCURL));
824  break;
825  case HDIV:
826  pipeline.push_back(new OpSetContravariantPiolaTransformOnEdge2D());
827  break;
828  default:
829  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
830  "Space %s not yet implemented", FieldSpaceNames[s]);
831  }
832  }
833 
835 }
836 
838  boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pipeline,
839  std::vector<FieldSpace> spaces, std::string geom_field_name) {
841 
842  if (geom_field_name.empty()) {
843 
844  } else {
845 
846  pipeline.push_back(new OpCalculateHOCoords<3>(geom_field_name));
847  pipeline.push_back(new OpGetHOTangentsOnEdge<3>(geom_field_name));
848  }
849 
850  for (auto s : spaces) {
851  switch (s) {
852  case HCURL:
853  pipeline.push_back(new OpHOSetContravariantPiolaTransformOnEdge3D(HCURL));
854  break;
855  default:
856  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
857  "Space %s not yet implemented", FieldSpaceNames[s]);
858  }
859  }
860 
862 }
863 
865  boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pipeline,
866  std::vector<FieldSpace> spaces, std::string geom_field_name) {
868 
869  if (geom_field_name.empty()) {
870  } else {
871 
872  pipeline.push_back(new OpCalculateHOCoords<3>(geom_field_name));
873  pipeline.push_back(new OpGetHONormalsOnFace<3>(geom_field_name));
874  }
875 
876  for (auto s : spaces) {
877  switch (s) {
878  case NOSPACE:
879  break;
880  case HCURL:
881  pipeline.push_back(new OpHOSetCovariantPiolaTransformOnFace3D(HCURL));
882  break;
883  case HDIV:
884  pipeline.push_back(new OpHOSetContravariantPiolaTransformOnFace3D(HDIV));
885  break;
886  default:
887  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
888  "Space %s not yet implemented", FieldSpaceNames[s]);
889  break;
890  }
891  }
892 
894 }
895 
897  boost::ptr_deque<ForcesAndSourcesCore::UserDataOperator> &pipeline,
898  std::vector<FieldSpace> spaces, std::string geom_field_name,
899  boost::shared_ptr<MatrixDouble> jac, boost::shared_ptr<VectorDouble> det,
900  boost::shared_ptr<MatrixDouble> inv_jac) {
902 
903  if (!jac)
904  jac = boost::make_shared<MatrixDouble>();
905  if (!det)
906  det = boost::make_shared<VectorDouble>();
907  if (!inv_jac)
908  inv_jac = boost::make_shared<MatrixDouble>();
909 
910  if (geom_field_name.empty()) {
911 
912  pipeline.push_back(new OpCalculateHOJac<3>(jac));
913 
914  } else {
915 
916  pipeline.push_back(new OpCalculateHOCoords<3>(geom_field_name));
917  pipeline.push_back(
918  new OpCalculateVectorFieldGradient<3, 3>(geom_field_name, jac));
919  }
920 
921  pipeline.push_back(new OpInvertMatrix<3>(jac, det, inv_jac));
922  pipeline.push_back(new OpSetHOWeights(det));
923 
924  for (auto s : spaces) {
925  switch (s) {
926  case NOSPACE:
927  break;
928  case H1:
929  pipeline.push_back(new OpSetHOInvJacToScalarBases<3>(H1, inv_jac));
930  break;
931  case HCURL:
932  pipeline.push_back(new OpSetHOCovariantPiolaTransform(HCURL, inv_jac));
933  pipeline.push_back(new OpSetHOInvJacVectorBase(HCURL, inv_jac));
934  break;
935  case HDIV:
936  pipeline.push_back(
938  pipeline.push_back(new OpSetHOInvJacVectorBase(HDIV, inv_jac));
939  break;
940  case L2:
941  pipeline.push_back(new OpSetHOInvJacToScalarBases<3>(L2, inv_jac));
942  break;
943  default:
944  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
945  "Space %s not yet implemented", FieldSpaceNames[s]);
946  break;
947  }
948  }
949 
951 }
952 
953 } // namespace MoFEM
NOSPACE
@ NOSPACE
Definition: definitions.h:83
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
MoFEM::OpSetHOCovariantPiolaTransform::piolaN
MatrixDouble piolaN
Definition: HODataOperators.hpp:227
MoFEM::EntitiesFieldData::EntData
Data on single entity (This is passed as argument to DataOperator::doWork)
Definition: EntitiesFieldData.hpp:127
MoFEM::OpSetHOInvJacToScalarBasesImpl::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:73
MoFEM::OpSetHOContravariantPiolaTransform::piolaDiffN
MatrixDouble piolaDiffN
Definition: HODataOperators.hpp:201
MoFEM::OpSetHOWeights::detPtr
boost::shared_ptr< VectorDouble > detPtr
Definition: HODataOperators.hpp:173
MoFEM::OpGetHOTangentsOnEdge
Calculate tangent vector on edge form HO geometry approximation.
Definition: HODataOperators.hpp:369
MoFEM::EntitiesFieldData::EntData::getBBDiffNMap
virtual std::map< std::string, boost::shared_ptr< MatrixDouble > > & getBBDiffNMap()
get hash map of direvarives base function for BB base, key is a field name
Definition: EntitiesFieldData.cpp:869
MoFEM::OpSetInvJacHcurlFace
OpSetInvJacHcurlFaceImpl< 2 > OpSetInvJacHcurlFace
Definition: UserDataOperators.hpp:2982
H1
@ H1
continuous field
Definition: definitions.h:85
MoFEM::OpSetHOContravariantPiolaTransform::jacPtr
boost::shared_ptr< MatrixDouble > jacPtr
Definition: HODataOperators.hpp:198
LASTBASE
@ LASTBASE
Definition: definitions.h:69
FTensor::Tensor1
Definition: Tensor1_value.hpp:8
MoFEM::OpCalculateHOJacForFaceEmbeddedIn3DSpace
OpCalculateHOJacForFaceImpl< 3 > OpCalculateHOJacForFaceEmbeddedIn3DSpace
Definition: HODataOperators.hpp:265
HVEC0_2
@ HVEC0_2
Definition: definitions.h:198
MoFEM::OpSetContravariantPiolaTransformOnFace2D
OpSetContravariantPiolaTransformOnFace2DImpl< 2 > OpSetContravariantPiolaTransformOnFace2D
Definition: UserDataOperators.hpp:3080
MoFEM::OpHOSetCovariantPiolaTransformOnFace3D::normalsAtPts
boost::shared_ptr< MatrixDouble > normalsAtPts
Definition: HODataOperators.hpp:357
NOBASE
@ NOBASE
Definition: definitions.h:59
L2
@ L2
field with C-1 continuity
Definition: definitions.h:88
MoFEM::OpHOSetContravariantPiolaTransformOnEdge3D::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:609
MoFEM::OpSetHOContravariantPiolaTransform
Apply contravariant (Piola) transfer to Hdiv space for HO geometry.
Definition: HODataOperators.hpp:180
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::Types::MatrixDouble
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:77
MoFEM::OpSetHOCovariantPiolaTransform::piolaDiffN
MatrixDouble piolaDiffN
Definition: HODataOperators.hpp:228
HVEC0_1
@ HVEC0_1
Definition: definitions.h:195
MoFEM::OpHOSetContravariantPiolaTransformOnFace3D
transform Hdiv base fluxes from reference element to physical triangle
Definition: HODataOperators.hpp:298
MoFEM::OpSetHOInvJacToScalarBases< 2 >
Definition: HODataOperators.hpp:78
MoFEM::OpSetHOContravariantPiolaTransform::piolaN
MatrixDouble piolaN
Definition: HODataOperators.hpp:200
MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperator::getCoords
VectorDouble & getCoords()
nodal coordinates
Definition: VolumeElementForcesAndSourcesCore.hpp:180
MoFEM::OpCalculateHOJacForVolume::OpCalculateHOJacForVolume
OpCalculateHOJacForVolume(boost::shared_ptr< MatrixDouble > jac_ptr)
Definition: HODataOperators.cpp:9
MoFEM::EntitiesFieldData::EntData::getDiffN
MatrixDouble & getDiffN(const FieldApproximationBase base)
get derivatives of base functions
Definition: EntitiesFieldData.hpp:1316
MoFEM::OpSetCovariantPiolaTransformOnFace2D
OpSetCovariantPiolaTransformOnFace2DImpl< 2 > OpSetCovariantPiolaTransformOnFace2D
Definition: UserDataOperators.hpp:3032
MoFEM::EntitiesFieldData::EntData::getBBDiffNByOrderArray
virtual std::array< boost::shared_ptr< MatrixDouble >, MaxBernsteinBezierOrder > & getBBDiffNByOrderArray()
Definition: EntitiesFieldData.cpp:902
HVEC0_0
@ HVEC0_0
Definition: definitions.h:192
HVEC1_1
@ HVEC1_1
Definition: definitions.h:196
MOFEM_IMPOSSIBLE_CASE
@ MOFEM_IMPOSSIBLE_CASE
Definition: definitions.h:35
MoFEM::OpCalculateHOCoords
Calculate HO coordinates at gauss points.
Definition: HODataOperators.hpp:41
THROW_MESSAGE
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561
MoFEM::OpCalculateHOJacForFaceImpl
Calculate jacobian for face element.
Definition: HODataOperators.hpp:240
HVEC1
@ HVEC1
Definition: definitions.h:186
FTensor::Tensor2
Definition: Tensor2_value.hpp:16
MoFEM::OpCalculateHOJac< 2 >
Definition: HODataOperators.hpp:273
MoFEM::ForcesAndSourcesCore::UserDataOperator::getFTensor0IntegrationWeight
auto getFTensor0IntegrationWeight()
Get integration weights.
Definition: ForcesAndSourcesCore.hpp:1239
MoFEM::ForcesAndSourcesCore::UserDataOperator::getMeasure
double getMeasure() const
get measure of element
Definition: ForcesAndSourcesCore.hpp:1274
MoFEM::invertTensor3by3
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:1566
MoFEM::EdgeElementForcesAndSourcesCore::UserDataOperator::getTangentAtGaussPts
MatrixDouble & getTangentAtGaussPts()
get tangent vector to edge curve at integration points
Definition: EdgeElementForcesAndSourcesCore.hpp:204
MoFEM::ForcesAndSourcesCore::UserDataOperator::getGaussPts
MatrixDouble & getGaussPts()
matrix of integration (Gauss) points for Volume Element
Definition: ForcesAndSourcesCore.hpp:1235
MoFEM::OpCalculateHOJacForVolume::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:22
HVEC2_1
@ HVEC2_1
Definition: definitions.h:197
MoFEM::EntitiesFieldData::EntData::getFTensor2DiffN
FTensor::Tensor2< FTensor::PackPtr< double *, Tensor_Dim0 *Tensor_Dim1 >, Tensor_Dim0, Tensor_Dim1 > getFTensor2DiffN(FieldApproximationBase base)
Get derivatives of base functions for Hdiv space.
Definition: EntitiesFieldData.cpp:660
MoFEM::OpSetHOContravariantPiolaTransform::detPtr
boost::shared_ptr< VectorDouble > detPtr
Definition: HODataOperators.hpp:197
FTensor::Number< 0 >
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator::getNormalsAtGaussPts
MatrixDouble & getNormalsAtGaussPts()
if higher order geometry return normals at Gauss pts.
Definition: FaceElementForcesAndSourcesCore.hpp:290
MoFEM::OpSetInvJacL2ForFaceEmbeddedIn3DSpace
Definition: UserDataOperators.hpp:2944
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator::getTangent2AtGaussPts
MatrixDouble & getTangent2AtGaussPts()
if higher order geometry return tangent vector to triangle at Gauss pts.
Definition: FaceElementForcesAndSourcesCore.hpp:310
HVEC2_2
@ HVEC2_2
Definition: definitions.h:200
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator::getTangent1AtGaussPts
MatrixDouble & getTangent1AtGaussPts()
if higher order geometry return tangent vector to triangle at Gauss pts.
Definition: FaceElementForcesAndSourcesCore.hpp:304
HVEC1_2
@ HVEC1_2
Definition: definitions.h:199
MoFEM::OpSetContravariantPiolaTransformOnFace2DEmbeddedIn3DSpace
OpSetContravariantPiolaTransformOnFace2DImpl< 3 > OpSetContravariantPiolaTransformOnFace2DEmbeddedIn3DSpace
Definition: UserDataOperators.hpp:3082
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
a
constexpr double a
Definition: approx_sphere.cpp:30
MoFEM::OpSetHOCovariantPiolaTransform
Apply covariant (Piola) transfer to Hcurl space for HO geometry.
Definition: HODataOperators.hpp:206
MoFEM::ForcesAndSourcesCore::UserDataOperator
Definition: ForcesAndSourcesCore.hpp:549
MoFEM::OpSetHOInvJacToScalarBases
Set inverse jacobian to base functions.
Definition: HODataOperators.hpp:73
MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperator::getInvJac
FTensor::Tensor2< double *, 3, 3 > & getInvJac()
get element inverse Jacobian
Definition: VolumeElementForcesAndSourcesCore.hpp:175
convert.type
type
Definition: convert.py:64
MoFEM::OpSetHOWeights::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:205
MoFEM::OpHOSetCovariantPiolaTransformOnFace3D::piolaN
MatrixDouble piolaN
Definition: HODataOperators.hpp:361
MoFEM::getFTensor0FromVec
static auto getFTensor0FromVec(ublas::vector< T, A > &data)
Get tensor rank 0 (scalar) form data vector.
Definition: Templates.hpp:135
MoFEM::OpSetHOContravariantPiolaTransform::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:229
MoFEM::OpSetContravariantPiolaTransformOnEdge2D
Definition: UserDataOperators.hpp:3090
MoFEM::ForcesAndSourcesCore::UserDataOperator::getNumeredEntFiniteElementPtr
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
Return raw pointer to NumeredEntFiniteElement.
Definition: ForcesAndSourcesCore.hpp:999
MoFEM::OpSetHOCovariantPiolaTransform::jacInvPtr
boost::shared_ptr< MatrixDouble > jacInvPtr
Definition: HODataOperators.hpp:225
MoFEM::OpSetHOInvJacVectorBase::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:109
MoFEM::OpSetHOWeightsOnFace
Modify integration weights on face to take in account higher-order geometry.
Definition: HODataOperators.hpp:122
MoFEM::FaceElementForcesAndSourcesCore
Face finite element.
Definition: FaceElementForcesAndSourcesCore.hpp:23
MoFEM::EntitiesFieldData::EntData::getDiffNSharedPtr
virtual boost::shared_ptr< MatrixDouble > & getDiffNSharedPtr(const FieldApproximationBase base)
Definition: EntitiesFieldData.cpp:37
MoFEM::OpCalculateHOJac< 3 >
Definition: HODataOperators.hpp:269
get_jac
static MoFEMErrorCode get_jac(EntitiesFieldData::EntData &col_data, int gg, MatrixDouble &jac_stress, MatrixDouble &jac)
Definition: NonLinearElasticElement.cpp:731
MoFEM::OpGetHONormalsOnFace
Calculate normals at Gauss points of triangle element.
Definition: HODataOperators.hpp:281
MoFEM::OpHOSetCovariantPiolaTransformOnFace3D::diffPiolaN
MatrixDouble diffPiolaN
Definition: HODataOperators.hpp:362
MoFEM::OpHOSetCovariantPiolaTransformOnFace3D::tangent1AtPts
boost::shared_ptr< MatrixDouble > tangent1AtPts
Definition: HODataOperators.hpp:358
t
constexpr double t
plate stiffness
Definition: plate.cpp:59
MoFEM::VolumeElementForcesAndSourcesCore
Volume finite element base.
Definition: VolumeElementForcesAndSourcesCore.hpp:26
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
MoFEM::OpSetHOInvJacVectorBase::invJacPtr
boost::shared_ptr< MatrixDouble > invJacPtr
Definition: HODataOperators.hpp:111
FTensor::Index< 'i', 3 >
MoFEM::AddHOOps
Add operators pushing bases from local to physical configuration.
Definition: HODataOperators.hpp:503
MoFEM::OpCalculateVectorFieldGradient
Get field gradients at integration pts for scalar filed rank 0, i.e. vector field.
Definition: UserDataOperators.hpp:1536
MoFEM::determinantTensor3by3
static auto determinantTensor3by3(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
Definition: Templates.hpp:1518
v
const double v
phase velocity of light in medium (cm/ns)
Definition: initial_diffusion.cpp:40
MoFEM::OpCalculateHOJacForVolume::jacPtr
boost::shared_ptr< MatrixDouble > jacPtr
Definition: HODataOperators.hpp:29
MoFEM::OpSetHOWeightsOnEdge::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:180
MoFEM::DataOperator::doEntities
std::array< bool, MBMAXTYPE > doEntities
If true operator is executed for entity.
Definition: DataOperators.hpp:85
FieldSpaceNames
const static char *const FieldSpaceNames[]
Definition: definitions.h:92
MoFEM::OpHOSetCovariantPiolaTransformOnFace3D
transform Hcurl base fluxes from reference element to physical triangle
Definition: HODataOperators.hpp:336
MoFEM::EntitiesFieldData::EntData::getFTensor1DiffN
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1DiffN(const FieldApproximationBase base)
Get derivatives of base functions.
Definition: EntitiesFieldData.cpp:526
MoFEM::OpHOSetContravariantPiolaTransformOnEdge3D::tangentAtGaussPts
boost::shared_ptr< MatrixDouble > tangentAtGaussPts
Definition: HODataOperators.hpp:330
MoFEM::OpHOSetCovariantPiolaTransformOnFace3D::tangent2AtPts
boost::shared_ptr< MatrixDouble > tangent2AtPts
Definition: HODataOperators.hpp:359
MoFEM::EntitiesFieldData::EntData::getN
MatrixDouble & getN(const FieldApproximationBase base)
get base functions this return matrix (nb. of rows is equal to nb. of Gauss pts, nb....
Definition: EntitiesFieldData.hpp:1305
MoFEM::OpSetHOCovariantPiolaTransform::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:334
MoFEM::ForcesAndSourcesCore::UserDataOperator::getFEDim
int getFEDim() const
Get dimension of finite element.
Definition: ForcesAndSourcesCore.hpp:1007
MoFEM::OpSetHOInvJacVectorBase::diffHdivInvJac
MatrixDouble diffHdivInvJac
Definition: HODataOperators.hpp:112
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
MoFEM::EntitiesFieldData::EntData::getFTensor1N
FTensor::Tensor1< FTensor::PackPtr< double *, Tensor_Dim >, Tensor_Dim > getFTensor1N(FieldApproximationBase base)
Get base functions for Hdiv/Hcurl spaces.
Definition: EntitiesFieldData.cpp:640
MoFEM::OpSetHOWeights
Set inverse jacobian to base functions.
Definition: HODataOperators.hpp:159
AINSWORTH_LEGENDRE_BASE
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
HCURL
@ HCURL
field with continuous tangents
Definition: definitions.h:86
MOFEM_DATA_INCONSISTENCY
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
MoFEM::OpInvertMatrix
Definition: UserDataOperators.hpp:3254
FieldApproximationBase
FieldApproximationBase
approximation base
Definition: definitions.h:58
HVEC2_0
@ HVEC2_0
Definition: definitions.h:194
MoFEM::FaceElementForcesAndSourcesCore::UserDataOperator::getFTensor1NormalsAtGaussPts
auto getFTensor1NormalsAtGaussPts()
get normal at integration points
Definition: FaceElementForcesAndSourcesCore.hpp:316
MoFEM::OpHOSetContravariantPiolaTransformOnFace3D::normalsAtGaussPts
boost::shared_ptr< MatrixDouble > normalsAtGaussPts
Definition: HODataOperators.hpp:311
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
MoFEM::OpSetHOInvJacVectorBase
transform local reference derivatives of shape function to global derivatives if higher order geometr...
Definition: HODataOperators.hpp:91
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
MoFEM::OpMakeHdivFromHcurl
Make Hdiv space from Hcurl space in 2d.
Definition: UserDataOperators.hpp:2989
HVEC2
@ HVEC2
Definition: definitions.h:186
MoFEM::OpSetInvJacHcurlFaceEmbeddedIn3DSpace
OpSetInvJacHcurlFaceImpl< 3 > OpSetInvJacHcurlFaceEmbeddedIn3DSpace
Definition: UserDataOperators.hpp:2983
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
MoFEM::OpHOSetContravariantPiolaTransformOnFace3D::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:459
MoFEM::OpSetInvJacH1ForFaceEmbeddedIn3DSpace
Definition: UserDataOperators.hpp:2936
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
HDIV
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
MOFEM_NOT_IMPLEMENTED
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
MoFEM::OpSetHOWeightsOnFace::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:153
MoFEM::EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1TangentAtGaussPts
FTensor::Tensor1< FTensor::PackPtr< double *, 3 >, DIM > getFTensor1TangentAtGaussPts()
Get tangent vector at integration points.
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
l
FTensor::Index< 'l', 3 > l
Definition: matrix_function.cpp:21
MoFEM::OpHOSetCovariantPiolaTransformOnFace3D::doWork
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
Operator for linear form, usually to calculate values on right hand side.
Definition: HODataOperators.cpp:512
HVEC1_0
@ HVEC1_0
Definition: definitions.h:193
HVEC0
@ HVEC0
Definition: definitions.h:186
MoFEM::OpHOSetContravariantPiolaTransformOnEdge3D
transform Hcurl base fluxes from reference element to physical edge
Definition: HODataOperators.hpp:317
MoFEM::EntitiesFieldData::EntData::getNSharedPtr
virtual boost::shared_ptr< MatrixDouble > & getNSharedPtr(const FieldApproximationBase base, const BaseDerivatives direvatie)
Definition: EntitiesFieldData.cpp:27