v0.8.20
EshelbianPlasticity.hpp
Go to the documentation of this file.
1 /**
2  * \file EshelbianPlasticity.hpp
3  * \brief Eshelbian plasticity interface
4  */
5 
6 #ifndef __ESHELBIAN_PLASTICITY_HPP__
7 #define __ESHELBIAN_PLASTICITY_HPP__
8 
9 namespace EshelbianPlasticity {
10 
12 
13 static const MOFEMuuid IDD_MOFEMEshelbianCrackInterface =
15 
16 typedef boost::shared_ptr<MatrixDouble> MatrixPtr;
17 typedef boost::shared_ptr<VectorDouble> VectorPtr;
18 
22  &m(0, 0), &m(1, 0), &m(2, 0), &m(3, 0), &m(4, 0), &m(5, 0), &m(6, 0),
23  &m(7, 0), &m(8, 0), &m(9, 0), &m(10, 0), &m(11, 0), &m(12, 0), &m(13, 0),
24  &m(14, 0), &m(15, 0), &m(16, 0), &m(17, 0), &m(18, 0), &m(19, 0),
25  &m(20, 0), &m(21, 0), &m(22, 0), &m(23, 0), &m(24, 0), &m(25, 0),
26  &m(26, 0));
27 }
28 
29 template <typename T>
36  t_R(i, j) = 0;
37  t_R(0, 0) = 1;
38  t_R(1, 1) = 1;
39  t_R(2, 2) = 1;
40 
41  const double angle = sqrt(t_omega(i) * t_omega(i));
42 
43  const double tol = 1e-12;
44  if (fabs(angle) < tol) {
45  return t_R;
46  }
47 
49  t_Omega(i, j) = levi_civita(i, j, k) * t_omega(k);
50  const double a = sin(angle) / angle;
51  const double ss_2 = sin(angle / 2.);
52  const double b = 2 * ss_2 * ss_2 / (angle * angle);
53  t_R(i, j) += a * t_Omega(i, j);
54  t_R(i, j) += b * t_Omega(i, k) * t_Omega(k, j);
55 
56  return t_R;
57 }
58 
59 template <typename T>
62 
67 
68  double angle = sqrt(t_omega(i) * t_omega(i));
69  const double tol = 1e-12;
70  if (fabs(angle) < tol) {
72  auto t_R = getRotationFormVector(t_omega);
73  t_diff_R(i, j, k) = levi_civita(i, l, k) * t_R(l, j);
74  return t_diff_R;
75  }
76 
79  t_Omega(i, j) = levi_civita(i, j, k) * t_omega(k);
80 
81  const double angle2 = angle * angle;
82  const double ss = sin(angle);
83  const double a = ss / angle;
84  const double ss_2 = sin(angle / 2.);
85  const double b = 2 * ss_2 * ss_2 / angle2;
86  t_diff_R(i, j, k) = 0;
87  t_diff_R(i, j, k) = a * levi_civita(i, j, k);
88  t_diff_R(i, j, k) += b * (t_Omega(i, l) * levi_civita(l, j, k) +
89  levi_civita(i, l, k) * t_Omega(l, j));
90 
91  const double cc = cos(angle);
92  const double cc_2 = cos(angle / 2.);
93  const double diff_a = (angle * cc - ss) / angle2;
94 
95  const double diff_b =
96  (2. * angle * ss_2 * cc_2 - 4. * ss_2 * ss_2) / (angle2 * angle);
97  t_diff_R(i, j, k) += diff_a * t_Omega(i, j) * (t_omega(k) / angle);
98  t_diff_R(i, j, k) +=
99  diff_b * t_Omega(i, l) * t_Omega(l, j) * (t_omega(k) / angle);
100 
101  return t_diff_R;
102 }
103 
108 
110  Mat Suu;
111  AO aoSuu;
112  Mat SuuSuu;
115  : Suu(PETSC_NULL), aoSuu(PETSC_NULL), SuuSuu(PETSC_NULL),
116  aoSuuSuu(PETSC_NULL) {}
118  if (Suu)
119  MatDestroy(&Suu);
120  if (aoSuu)
121  AODestroy(&aoSuu);
122  if (SuuSuu)
123  MatDestroy(&SuuSuu);
124  if (aoSuuSuu)
125  AODestroy(&aoSuuSuu);
126  }
127 
130  if (this->Suu)
131  MatDestroy(&this->Suu);
132  if (this->aoSuu)
133  AODestroy(&this->aoSuu);
134  this->Suu = Suu;
135  this->aoSuu = aoSuu;
136  PetscObjectReference((PetscObject)this->Suu);
137  PetscObjectReference((PetscObject)this->aoSuu);
139  }
140 
143  if (this->SuuSuu)
144  MatDestroy(&this->SuuSuu);
145  if (this->aoSuuSuu)
146  AODestroy(&this->aoSuuSuu);
147  this->SuuSuu = SuuSuu;
148  this->aoSuuSuu = aoSuuSuu;
149  PetscObjectReference((PetscObject)this->SuuSuu);
150  PetscObjectReference((PetscObject)this->aoSuuSuu);
152  }
153 };
154 
155 template <typename E> struct EpElement : public E, EpElementBase {
156  EpElement(MoFEM::Interface &m_field) : E(m_field), EpElementBase() {}
157 };
158 
159 template <> struct EpElement<BasicMethod> : public FEMethod, EpElementBase {
161 };
162 
163 template <> struct EpElement<FEMethod> : public FEMethod, EpElementBase {
165 };
166 
167 struct EpFEMethod : EpElement<FEMethod> {
171  CHKERR MatZeroEntries(Suu);
172  CHKERR MatZeroEntries(SuuSuu);
174  }
175 
178  CHKERR MatAssemblyBegin(Suu, MAT_FINAL_ASSEMBLY);
179  CHKERR MatAssemblyEnd(Suu, MAT_FINAL_ASSEMBLY);
180  CHKERR MatAssemblyBegin(SuuSuu, MAT_FINAL_ASSEMBLY);
181  CHKERR MatAssemblyEnd(SuuSuu, MAT_FINAL_ASSEMBLY);
182  // CHKERR MatView(SuuSuu, PETSC_VIEWER_DRAW_WORLD);
183  // std::string wait;
184  // std::cin >> wait;
185  // CHKERR MatView(Suu, PETSC_VIEWER_DRAW_WORLD);
186  // std::cin >> wait;
187  // CHKERR MatView(ts_B, PETSC_VIEWER_DRAW_WORLD);
188  // std::cin >> wait;
190  }
191 };
192 
194 
213 
222 
243 
244  // Not really data at integration points, used to calculate Schur complement
245 
250  std::map<EntityType, std::vector<MatrixDouble>> UP;
251  std::map<EntityType, std::vector<MatrixDouble>> invUP;
252  std::map<EntityType, std::vector<MatrixDouble>> invUx;
253 
259  std::map<EntityType, std::vector<MatrixDouble>> invPB;
260  std::map<EntityType, std::vector<MatrixDouble>> PB;
261  std::map<EntityType, std::vector<MatrixDouble>> invBx;
262 };
263 
264 struct OpJacobian;
265 
267 
272 
275 
276  PhysicalEquations() = delete;
277  PhysicalEquations(const int size_active, const int size_dependent)
278  : activeVariables(size_active, 0), dependentVariables(size_dependent, 0),
279  dependentVariablesDirevatives(size_dependent * size_active, 0) {}
280  virtual ~PhysicalEquations() {}
281 
282  virtual MoFEMErrorCode recordTape(const int tag) = 0;
283 
284  virtual OpJacobian *
285  returnOpJacobian(const std::string &field_name, const int tag,
286  const bool eval_rhs, const bool eval_lhs,
287  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
288  boost::shared_ptr<PhysicalEquations> &physics_ptr) = 0;
289 
290  std::vector<double> activeVariables;
291  std::vector<double> dependentVariables;
292  std::vector<double> dependentVariablesDirevatives;
293 
294  /** \name Active variables */
295 
296  /**@{*/
297 
298  template <int S>
299  inline static DTensor2Ptr get_VecTensor2(std::vector<double> &v) {
300  return DTensor2Ptr(&v[S + 0], &v[S + 1], &v[S + 2], &v[S + 3], &v[S + 4],
301  &v[S + 5], &v[S + 6], &v[S + 7], &v[S + 8]);
302  }
303 
304  template <int S>
305  inline static DTensor0Ptr get_VecTensor0(std::vector<double> &v) {
306  return DTensor0Ptr(&v[S + 0]);
307  }
308 
309  template <int S0, int S1, int S2>
310  inline static DTensor3Ptr get_vecTensor3(std::vector<double> &v) {
311  constexpr int A00 = 18 * S0 + 18 * 0 + 9 * S1 + 3 * S2;
312  constexpr int A01 = 18 * S0 + 18 * 1 + 9 * S1 + 3 * S2;
313  constexpr int A02 = 18 * S0 + 18 * 2 + 9 * S1 + 3 * S2;
314  constexpr int A10 = 18 * S0 + 18 * 3 + 9 * S1 + 3 * S2;
315  constexpr int A11 = 18 * S0 + 18 * 4 + 9 * S1 + 3 * S2;
316  constexpr int A12 = 18 * S0 + 18 * 5 + 9 * S1 + 3 * S2;
317  constexpr int A20 = 18 * S0 + 18 * 6 + 9 * S1 + 3 * S2;
318  constexpr int A21 = 18 * S0 + 18 * 7 + 9 * S1 + 3 * S2;
319  constexpr int A22 = 18 * S0 + 18 * 8 + 9 * S1 + 3 * S2;
320  return DTensor3Ptr(
321 
322  &v[A00 + 0], &v[A00 + 1], &v[A00 + 2], &v[A01 + 0], &v[A01 + 1],
323  &v[A01 + 2], &v[A02 + 0], &v[A02 + 1], &v[A02 + 2],
324 
325  &v[A10 + 0], &v[A10 + 1], &v[A10 + 2], &v[A11 + 0], &v[A11 + 1],
326  &v[A11 + 2], &v[A12 + 0], &v[A12 + 1], &v[A12 + 2],
327 
328  &v[A20 + 0], &v[A20 + 1], &v[A20 + 2], &v[A21 + 0], &v[A21 + 1],
329  &v[A21 + 2], &v[A22 + 0], &v[A22 + 1], &v[A22 + 2]
330 
331  );
332  }
333 
334  /**@}*/
335 
336  /** \name Active variables */
337 
338  /**@{*/
339 
340  inline DTensor2Ptr get_h() { return get_VecTensor2<0>(activeVariables); }
341  inline DTensor2Ptr get_H() { return get_VecTensor2<9>(activeVariables); }
342 
343  /**@}*/
344 
345  /** \name Dependent variables */
346 
347  /**@{*/
348 
349  inline DTensor2Ptr get_P() { return get_VecTensor2<0>(dependentVariables); }
351  return get_VecTensor2<9>(dependentVariables);
352  }
353  inline DTensor0Ptr get_Phi() {
354  return get_VecTensor0<18>(dependentVariables);
355  }
356  inline double &get_PhiRef() { return dependentVariables[18]; }
358  return get_VecTensor2<9 + 9 + 1>(dependentVariables);
359  }
360 
361  /**@}*/
362 
363  /** \name Derivatives of dependent variables */
364 
365  /**@{*/
366 
368  return get_vecTensor3<0, 0, 0>(dependentVariablesDirevatives);
369  }
371  return get_vecTensor3<0, 1, 0>(dependentVariablesDirevatives);
372  }
374  return get_vecTensor3<0, 0, 1>(dependentVariablesDirevatives);
375  }
377  return get_vecTensor3<0, 1, 1>(dependentVariablesDirevatives);
378  }
380  return get_vecTensor3<0, 0, 2>(dependentVariablesDirevatives);
381  }
383  return get_vecTensor3<0, 1, 2>(dependentVariablesDirevatives);
384  }
385 
387  return get_vecTensor3<9, 0, 0>(dependentVariablesDirevatives);
388  }
390  return get_vecTensor3<9, 1, 0>(dependentVariablesDirevatives);
391  }
393  return get_vecTensor3<9, 0, 1>(dependentVariablesDirevatives);
394  }
396  return get_vecTensor3<9, 1, 1>(dependentVariablesDirevatives);
397  }
399  return get_vecTensor3<9, 0, 2>(dependentVariablesDirevatives);
400  }
402  return get_vecTensor3<9, 1, 2>(dependentVariablesDirevatives);
403  }
404 
406  return get_VecTensor2<18 * (9 + 9) + 0>(dependentVariablesDirevatives);
407  }
409  return get_VecTensor2<18 * (9 + 9) + 9>(dependentVariablesDirevatives);
410  }
411 
413  return get_vecTensor3<9 + 9 + 1, 0, 0>(dependentVariablesDirevatives);
414  }
416  return get_vecTensor3<9 + 9 + 1, 1, 0>(dependentVariablesDirevatives);
417  }
419  return get_vecTensor3<9 + 9 + 1, 0, 1>(dependentVariablesDirevatives);
420  }
422  return get_vecTensor3<9 + 9 + 1, 1, 1>(dependentVariablesDirevatives);
423  }
425  return get_vecTensor3<9 + 9 + 1, 0, 2>(dependentVariablesDirevatives);
426  }
428  return get_vecTensor3<9 + 9 + 1, 1, 2>(dependentVariablesDirevatives);
429  }
430 
431  /**@}*/
432 };
433 
434 struct BcDisp {
435  BcDisp(std::string name, std::vector<double> &attr, Range &faces);
436  std::string blockName;
437  Range faces;
440 };
441 typedef std::vector<BcDisp> BcDispVec;
442 
443 struct BcRot {
444  BcRot(std::string name, std::vector<double> &attr, Range &faces);
445  std::string blockName;
446  Range faces;
448  double theta;
449 };
450 typedef std::vector<BcRot> BcRotVec;
451 
452 typedef std::vector<Range> TractionFreeBc;
453 
454 struct TractionBc {
455  TractionBc(std::string name, std::vector<double> &attr, Range &faces);
456  std::string blockName;
457  Range faces;
460 };
461 typedef std::vector<TractionBc> TractionBcVec;
462 
463 struct OpJacobian : public UserDataOperator {
464  const int tAg; ///< adol-c tape
465  const bool evalRhs;
466  const bool evalLhs;
467  boost::shared_ptr<DataAtIntegrationPts>
468  dataAtPts; ///< data at integration pts
469  boost::shared_ptr<PhysicalEquations>
470  physicsPtr; ///< material physical equations
471 
472  OpJacobian(const std::string &field_name, const int tag, const bool eval_rhs,
473  const bool eval_lhs,
474  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
475  boost::shared_ptr<PhysicalEquations> &physics_ptr)
476  : UserDataOperator(field_name, OPROW), tAg(tag), evalRhs(eval_rhs),
477  evalLhs(eval_lhs), dataAtPts(data_ptr), physicsPtr(physics_ptr) {}
478 
479  virtual MoFEMErrorCode evaluateRhs(EntData &data) = 0;
480  virtual MoFEMErrorCode evaluateLhs(EntData &data) = 0;
481 
482  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
483 };
484 
485 template <typename T> struct OpAssembleBasic : public T {
486 
487  boost::shared_ptr<DataAtIntegrationPts>
488  dataAtPts; ///< data at integration pts
489  const bool assembleSymmetry;
490  const bool assembleSuu;
491  const bool assembleSuuSuu;
492 
493  OpAssembleBasic(const std::string &field_name,
494  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
495  const char type)
496  : T(field_name, type), dataAtPts(data_ptr), assembleSymmetry(false),
497  assembleSuu(true), assembleSuuSuu(true) {}
498 
499  OpAssembleBasic(const std::string &row_field, const std::string &col_field,
500  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
501  const char type, const bool assemble_symmetry,
502  const bool assemble_suu, const bool assemble_suu_suu)
503  : T(row_field, col_field, type, false), dataAtPts(data_ptr),
504  assembleSymmetry(assemble_symmetry), assembleSuu(assemble_suu),
505  assembleSuuSuu(assemble_suu && assemble_suu_suu) {}
506 
507  VectorDouble nF; ///< local right hand side vector
508  MatrixDouble K; ///< local tangent matrix
510 
513  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
515  }
516  virtual MoFEMErrorCode assemble(EntData &data) {
518  double *vec_ptr = &*nF.begin();
519  int nb_dofs = data.getIndices().size();
520  int *ind_ptr = &*data.getIndices().begin();
521  CHKERR VecSetValues(T::getFEMethod()->ts_F, nb_dofs, ind_ptr, vec_ptr,
522  ADD_VALUES);
524  }
525 
526  virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) {
528  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
530  }
531 
532  virtual MoFEMErrorCode assemble(int row_side, int col_side,
533  EntityType row_type, EntityType col_type,
534  EntData &row_data, EntData &col_data) {
536 
537  double *mat_ptr = &*K.data().begin();
538  int *row_ind_ptr = &*row_data.getIndices().begin();
539  int *col_ind_ptr = &*col_data.getIndices().begin();
540  int row_nb_dofs = row_data.getIndices().size();
541  int col_nb_dofs = col_data.getIndices().size();
542  CHKERR MatSetValues(T::getFEMethod()->ts_B, row_nb_dofs, row_ind_ptr,
543  col_nb_dofs, col_ind_ptr, mat_ptr, ADD_VALUES);
544  if (assembleSymmetry) {
545  transposeK.resize(col_nb_dofs, row_nb_dofs);
546  noalias(transposeK) = trans(K);
547  double *mat_ptr = &*transposeK.data().begin();
548  CHKERR MatSetValues(T::getFEMethod()->ts_B, col_nb_dofs, col_ind_ptr,
549  row_nb_dofs, row_ind_ptr, mat_ptr, ADD_VALUES);
550  }
551 
552  Mat Suu = PETSC_NULL;
553  AO aoSuu = PETSC_NULL;
554  if (assembleSuu) {
555 
556  auto set_data = [&Suu, &aoSuu](auto fe) {
557  aoSuu = fe->aoSuu;
558  Suu = fe->Suu;
559  };
560 
561  if (auto ep_fe_vol_ptr = dynamic_cast<
563  T::getFEMethod()))
564  set_data(ep_fe_vol_ptr);
565  else if (auto ep_fe_face_ptr = dynamic_cast<
567  T::getFEMethod()))
568  set_data(ep_fe_face_ptr);
569  }
570 
571  if (Suu) {
572  auto row_ind = row_data.getIndices();
573  auto col_ind = col_data.getIndices();
574  CHKERR AOApplicationToPetsc(aoSuu, row_nb_dofs, &*row_ind.begin());
575  CHKERR AOApplicationToPetsc(aoSuu, col_nb_dofs, &*col_ind.begin());
576  double *mat_ptr = &*K.data().begin();
577  CHKERR MatSetValues(Suu, row_nb_dofs, &*row_ind.begin(), col_nb_dofs,
578  &*col_ind.begin(), mat_ptr, ADD_VALUES);
579  if (assembleSymmetry) {
580  double *mat_ptr = &*transposeK.data().begin();
581  CHKERR MatSetValues(Suu, col_nb_dofs, &*col_ind.begin(), row_nb_dofs,
582  &*row_ind.begin(), mat_ptr, ADD_VALUES);
583  }
584  }
585 
586  Mat SuuSuu = PETSC_NULL;
587  AO aoSuuSuu = PETSC_NULL;
588  if (assembleSuuSuu) {
589 
590  auto set_data = [&aoSuu, &SuuSuu, &aoSuuSuu](auto fe) {
591  aoSuu = fe->aoSuu;
592  SuuSuu = fe->SuuSuu;
593  aoSuuSuu = fe->aoSuuSuu;
594  };
595 
596  if (auto ep_fe_vol_ptr = dynamic_cast<
598  T::getFEMethod()))
599  set_data(ep_fe_vol_ptr);
600  else if (auto ep_fe_face_ptr = dynamic_cast<
602  T::getFEMethod()))
603  set_data(ep_fe_face_ptr);
604  }
605 
606  if (SuuSuu) {
607  auto row_ind = row_data.getIndices();
608  auto col_ind = col_data.getIndices();
609  CHKERR AOApplicationToPetsc(aoSuu, row_nb_dofs, &*row_ind.begin());
610  CHKERR AOApplicationToPetsc(aoSuu, col_nb_dofs, &*col_ind.begin());
611  CHKERR AOApplicationToPetsc(aoSuuSuu, row_nb_dofs, &*row_ind.begin());
612  CHKERR AOApplicationToPetsc(aoSuuSuu, col_nb_dofs, &*col_ind.begin());
613  double *mat_ptr = &*K.data().begin();
614  CHKERR MatSetValues(SuuSuu, row_nb_dofs, &*row_ind.begin(), col_nb_dofs,
615  &*col_ind.begin(), mat_ptr, ADD_VALUES);
616  if (assembleSymmetry) {
617  double *mat_ptr = &*transposeK.data().begin();
618  CHKERR MatSetValues(SuuSuu, col_nb_dofs, &*col_ind.begin(), row_nb_dofs,
619  &*row_ind.begin(), mat_ptr, ADD_VALUES);
620  }
621  }
622 
624  }
625 
626  MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
628  if (data.getIndices().empty())
630  nF.resize(data.getIndices().size(), false);
631  nF.clear();
632  CHKERR integrate(data);
633  CHKERR assemble(data);
635  }
636 
637  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
638  EntityType col_type, EntData &row_data,
639  EntData &col_data) {
641  if (row_data.getIndices().empty())
643  if (col_data.getIndices().empty())
645  K.resize(row_data.getIndices().size(), col_data.getIndices().size(), false);
646  K.clear();
647  CHKERR integrate(row_data, col_data);
648  CHKERR assemble(row_side, col_side, row_type, col_type, row_data, col_data);
650  }
651 };
652 
653 struct OpAssembleVolume : public OpAssembleBasic<VolUserDataOperator> {
654  OpAssembleVolume(const std::string &field,
655  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
656  const char type)
657  : OpAssembleBasic<VolUserDataOperator>(field, data_ptr, type) {}
658 
659  OpAssembleVolume(const std::string &row_field, const std::string &col_field,
660  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
661  const char type, const bool assemble_symmetry,
662  const bool assemble_suu, const bool assemble_suu_suu)
663  : OpAssembleBasic<VolUserDataOperator>(row_field, col_field, data_ptr,
664  type, assemble_symmetry,
665  assemble_suu, assemble_suu_suu) {}
666 };
667 
668 struct OpAssembleFace : public OpAssembleBasic<FaceUserDataOperator> {
669  OpAssembleFace(const std::string &field,
670  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
671  const char type)
672  : OpAssembleBasic<FaceUserDataOperator>(field, data_ptr, type) {}
673 
674  OpAssembleFace(const std::string &row_field, const std::string &col_field,
675  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
676  const char type, const bool assemble_symmetry,
677  const bool assemble_suu, const bool assemble_suu_suu)
678  : OpAssembleBasic<FaceUserDataOperator>(row_field, col_field, data_ptr,
679  type, assemble_symmetry,
680  assemble_suu, assemble_suu_suu) {}
681 };
682 
684  boost::shared_ptr<DataAtIntegrationPts>
685  dataAtPts; ///< data at integration pts
687  const std::string &field_name,
688  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
689  : VolUserDataOperator(field_name, OPROW), dataAtPts(data_ptr) {}
690  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
691 };
692 
694  OpSpatialEquilibrium(const std::string &field_name,
695  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
696  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
698 };
699 
701  OpSpatialRotation(const std::string &field_name,
702  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
703  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
705 };
706 
708  OpSpatialPhysical(const std::string &field_name,
709  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
710  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
712 };
713 
715  OpSpatialConsistencyP(const std::string &field_name,
716  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
717  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
719 };
720 
722  OpSpatialConsistencyBubble(const std::string &field_name,
723  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
724  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
726 };
727 
729  OpSpatialConsistencyDivTerm(const std::string &field_name,
730  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
731  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
733 };
734 
735 struct OpDispBc : public OpAssembleFace {
736  boost::shared_ptr<BcDispVec> bcDispPtr;
737  OpDispBc(const std::string &field_name,
738  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
739  boost::shared_ptr<BcDispVec> &bc_disp_ptr)
740  : OpAssembleFace(field_name, data_ptr, OPROW), bcDispPtr(bc_disp_ptr) {}
742 };
743 
744 struct OpDispBc_dx : public OpAssembleFace {
745  boost::shared_ptr<BcDispVec> bcDispPtr;
746  OpDispBc_dx(const std::string &row_field_name,
747  const std::string &col_field_name,
748  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
749  boost::shared_ptr<BcDispVec> &bc_disp_ptr)
750  : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL,
751  false, true, true),
752  bcDispPtr(bc_disp_ptr) {}
753  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
754 };
755 
756 struct OpRotationBc : public OpAssembleFace {
757  boost::shared_ptr<BcRotVec> bcRotPtr;
758  OpRotationBc(const std::string &field_name,
759  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
760  boost::shared_ptr<BcRotVec> &bc_rot_ptr)
761  : OpAssembleFace(field_name, data_ptr, OPROW), bcRotPtr(bc_rot_ptr) {}
763 };
764 
766  boost::shared_ptr<BcRotVec> bcRotPtr;
767  OpRotationBc_dx(const std::string &row_field_name,
768  const std::string &col_field_name,
769  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
770  boost::shared_ptr<BcRotVec> &bc_rot_ptr)
771  : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL,
772  false, true, true),
773  bcRotPtr(bc_rot_ptr) {}
774  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
775 };
776 
777 struct FeTractionBc;
778 
780  OpTractionBc(std::string row_field, FeTractionBc &bc_fe)
781  : FaceUserDataOperator(row_field, FaceUserDataOperator::OPROW),
782  bcFe(bc_fe) {}
783  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
784 
785 protected:
789 };
790 
791 struct FeTractionBc : public FEMethod {
793  FeTractionBc(MoFEM::Interface &m_field,const std::string field_name,
794  boost::shared_ptr<TractionBcVec> &bc)
795  : FEMethod(), mField(m_field), bcData(bc), fieldName(field_name) {}
797  MoFEMErrorCode postProcess() { return 0; }
798 
799  friend MoFEMErrorCode OpTractionBc::doWork(int side, EntityType type,
800  EntData &data);
801 
802 protected:
804  boost::shared_ptr<TractionBcVec> bcData;
805  std::string fieldName;
806 };
807 
808 template <>
810  EpElement(MoFEM::Interface &m_field, const std::string field_name,
811  boost::shared_ptr<TractionBcVec> &bc)
812  : FeTractionBc(m_field, field_name, bc), EpElementBase() {}
817  }
818 };
819 
821  OpSpatialEquilibrium_dw_dP(const std::string &row_field,
822  const std::string &col_field,
823  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
824  const bool assemble_off = false)
825  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
826  true, true) {
827  sYmm = false;
828  }
829  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
830 };
831 
833  OpSpatialPhysical_du_du(const std::string &row_field,
834  const std::string &col_field,
835  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
836  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false, false,
837  false) {
838  sYmm = false;
839  }
840  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
841  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
842  EntityType col_type, EntData &row_data,
843  EntData &col_data);
844  private:
847 };
848 
850  OpSpatialPhysical_du_dP(const std::string &row_field,
851  const std::string &col_field,
852  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
853  const bool assemble_off = false)
854  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
855  false, false) {
856  sYmm = false;
857  }
858 
859  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
860  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
861  EntityType col_type, EntData &row_data,
862  EntData &col_data);
863 };
864 
867  const std::string &row_field, const std::string &col_field,
868  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
869  const bool assemble_off = false)
870  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
871  false, false) {
872  sYmm = false;
873  }
874 
875  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
876 
877  /**
878  * @brief Assemble and saves UB, and invUB
879  *
880  * @param row_side
881  * @param col_side
882  * @param row_type
883  * @param col_type
884  * @param row_data
885  * @param col_data
886  * @return MoFEMErrorCode
887  */
888  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
889  EntityType col_type, EntData &row_data,
890  EntData &col_data);
891 };
892 
894  OpSpatialPhysical_du_domega(const std::string &row_field,
895  const std::string &col_field,
896  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
897  const bool assemble_off = false)
898  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
899  false, false) {
900  sYmm = false;
901  }
902  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
903  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
904  EntityType col_type, EntData &row_data,
905  EntData &col_data);
906 };
907 
909  OpSpatialPhysical_du_dx(const std::string &row_field,
910  const std::string &col_field,
911  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
912  const bool assemble_off = false)
913  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
914  false, false) {
915  sYmm = false;
916  }
917  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
918  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
919  EntityType col_type, EntData &row_data,
920  EntData &col_data);
921 };
922 
924  OpSpatialRotation_domega_dP(const std::string &row_field,
925  const std::string &col_field,
926  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
927  const bool assemble_off = false)
928  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
929  true, true) {
930  sYmm = false;
931  }
932  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
933 };
934 
937  const std::string &row_field, const std::string &col_field,
938  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
939  const bool assemble_off = false)
940  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
941  true, false) {
942  sYmm = false;
943  }
944  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
945 
946  /**
947  * @brief Assemble and calculate invOmegaB and OmegaB
948  *
949  * @param row_side
950  * @param col_side
951  * @param row_type
952  * @param col_type
953  * @param row_data
954  * @param col_data
955  * @return MoFEMErrorCode
956  */
957  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
958  EntityType col_type, EntData &row_data,
959  EntData &col_data);
960 };
961 
964  const std::string &row_field, const std::string &col_field,
965  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
966  const bool assemble_off = false)
967  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
968  true, true) {
969  sYmm = false;
970  }
971  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
972 };
973 
975  OpSpatialRotation_domega_dx(const std::string &row_field,
976  const std::string &col_field,
977  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
978  const bool assemble_off = false)
979  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
980  true, true) {
981  sYmm = false;
982  }
983  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
984 };
985 
988  const std::string &row_field, const std::string &col_field,
989  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
990  const bool assemble_off = false)
991  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
992  true, true) {
993  sYmm = false;
994  }
995  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
996 };
997 
1000  const std::string &row_field, const std::string &col_field,
1001  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1002  const bool assemble_off = false)
1003  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1004  true, false) {
1005  sYmm = false;
1006  }
1007  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1008 
1009  /**
1010  * @brief Assemble and calculate invBOmega and BOmega
1011  *
1012  * @param row_side
1013  * @param col_side
1014  * @param row_type
1015  * @param col_type
1016  * @param row_data
1017  * @param col_data
1018  * @return MoFEMErrorCode
1019  */
1020  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1021  EntityType col_type, EntData &row_data,
1022  EntData &col_data);
1023 };
1024 
1026  OpSpatialConsistency_dP_dx(const std::string &row_field,
1027  const std::string &col_field,
1028  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1029  const bool assemble_off = false)
1030  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1031  true, true) {
1032  sYmm = false;
1033  }
1034  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1035 };
1036 
1039  const std::string &row_field, const std::string &col_field,
1040  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1041  const bool assemble_off = false)
1042  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1043  true, false) {
1044  sYmm = false;
1045  }
1046  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1047 
1048  /**
1049  * @brief Assemble and calculate invBx and xB
1050  *
1051  * @param row_side
1052  * @param col_side
1053  * @param row_type
1054  * @param col_type
1055  * @param row_data
1056  * @param col_data
1057  * @return MoFEMErrorCode
1058  */
1059  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1060  EntityType col_type, EntData &row_data,
1061  EntData &col_data);
1062 };
1063 
1065 
1066  moab::Interface &postProcMesh;
1067  std::vector<EntityHandle> &mapGaussPts;
1068  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
1069 
1070  OpPostProcDataStructure(moab::Interface &post_proc_mesh,
1071  std::vector<EntityHandle> &map_gauss_pts,
1072  const std::string field_name,
1073  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1074  : VolUserDataOperator(field_name, UserDataOperator::OPROW),
1075  postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
1076  dataAtPts(data_ptr) {}
1077 
1078  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
1079 };
1080 
1082  OpSpatialPrj(const std::string &row_field,
1083  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1084  : OpAssembleVolume(row_field, data_ptr, OPROW) {}
1085  MoFEMErrorCode integrate(EntData &row_data);
1086 };
1087 
1089  OpSpatialPrj_dx_dx(const std::string &row_field, const std::string &col_field,
1090  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1091  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false, true,
1092  true) {
1093  // FIXME: That is symmetric
1094  sYmm = false;
1095  }
1096  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1097 };
1098 
1100  OpSpatialPrj_dx_dw(const std::string &row_field, const std::string &col_field,
1101  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1102  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false, true,
1103  true) {
1104  sYmm = false;
1105  }
1106  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1107 };
1108 
1111  const std::string &row_field, const std::string &col_field,
1112  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1113  const bool assemble_off = false)
1114  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1115  true, true) {
1116  sYmm = true;
1117  }
1118  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1119  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1120  EntityType col_type, EntData &row_data,
1121  EntData &col_data);
1122 };
1123 
1126  const std::string &row_field, const std::string &col_field,
1127  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1128  const bool assemble_off = false)
1129  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1130  true, false) {
1131  sYmm = false;
1132  }
1133  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1134  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1135  EntityType col_type, EntData &row_data,
1136  EntData &col_data);
1137 };
1138 
1141  const std::string &row_field, const std::string &col_field,
1142  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1143  const bool assemble_off = false)
1144  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1145  true, false) {
1146  sYmm = true;
1147  }
1148  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1149  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1150  EntityType col_type, EntData &row_data,
1151  EntData &col_data);
1152 
1153 private:
1156 
1157 };
1158 
1161  const std::string &row_field, const std::string &col_field,
1162  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1163  const bool assemble_off = false)
1164  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1165  true, false) {
1166  sYmm = false;
1167  }
1168  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1169  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1170  EntityType col_type, EntData &row_data,
1171  EntData &col_data);
1172 };
1173 
1176  const std::string &row_field, const std::string &col_field,
1177  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1178  const bool assemble_off = false)
1179  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1180  true, true) {
1181  sYmm = false;
1182  }
1183  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1184  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1185  EntityType col_type, EntData &row_data,
1186  EntData &col_data);
1187 };
1188 
1191  const std::string &row_field, const std::string &col_field,
1192  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1193  const bool assemble_off = false)
1194  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1195  true, true) {
1196  sYmm = false;
1197  }
1198  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1199  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1200  EntityType col_type, EntData &row_data,
1201  EntData &col_data);
1202 };
1203 
1206  const std::string &row_field, const std::string &col_field,
1207  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1208  const bool assemble_off = false)
1209  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1210  true, true) {
1211  sYmm = false;
1212  }
1213  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1214  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1215  EntityType col_type, EntData &row_data,
1216  EntData &col_data);
1217 };
1218 
1221  const std::string &row_field, const std::string &col_field,
1222  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1223  const bool assemble_off = false)
1224  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1225  true, true) {
1226  sYmm = false;
1227  }
1228  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1229  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1230  EntityType col_type, EntData &row_data,
1231  EntData &col_data);
1232 };
1233 
1236  const std::string &row_field, const std::string &col_field,
1237  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1238  const bool assemble_off = false)
1239  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1240  true, false) {
1241  sYmm = false;
1242  }
1243  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1244  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1245  EntityType col_type, EntData &row_data,
1246  EntData &col_data);
1247 };
1248 
1251  const std::string &row_field, const std::string &col_field,
1252  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1253  const bool assemble_off = false)
1254  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, assemble_off,
1255  true, true) {
1256  sYmm = false;
1257  }
1258  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) { return 0; }
1259  MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type,
1260  EntityType col_type, EntData &row_data,
1261  EntData &col_data);
1262 };
1263 
1265 
1266  OpCalculateStrainEnergy(const std::string field_name,
1267  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
1268  boost::shared_ptr<double> &e)
1269  : VolUserDataOperator(field_name, VolUserDataOperator::OPROW),
1270  dataAtPts(data_at_pts), energy(e) {}
1271 
1272  MoFEMErrorCode doWork(int side, EntityType type,
1274 
1275 private:
1276  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
1277  boost::shared_ptr<double> energy;
1278 };
1279 
1281 
1282  /**
1283  * \brief Getting interface of core database
1284  * @param uuid unique ID of interface
1285  * @param iface returned pointer to interface
1286  * @return error code
1287  */
1288  MoFEMErrorCode query_interface(const MOFEMuuid &uuid,
1289  UnknownInterface **iface) const;
1290 
1292 
1293  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
1294  boost::shared_ptr<PhysicalEquations> physicalEquations;
1295 
1296  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> elasticFeRhs;
1297  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> elasticFeLhs;
1298  boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> elasticBcLhs;
1299  boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> elasticBcRhs;
1300  boost::shared_ptr<EpFEMethod> schurAssembly;
1301 
1302  DM dM; ///< Coupled problem all fields
1303  DM dmElastic; ///< Elastic problem
1304  DM dmElasticSchurStreach; ///< Sub problem of dmElastic Schur
1305  DM dmElasticSchurBubble; ///< Sub problem of dmElastic Schur
1306  DM dmMaterial; ///< Material problem
1307  DM dmPrjSpatial; ///< Project spatial displacement
1308 
1309  const std::string piolaStress;
1310  const std::string eshelbyStress;
1311  const std::string spatialDisp;
1312  const std::string spatialPosition;
1313  const std::string materialDisp;
1314  const std::string streachTensor;
1315  const std::string rotAxis;
1316  const std::string materialGradient;
1317  const std::string tauField;
1318  const std::string lambdaField;
1319  const std::string bubbleField;
1320 
1321  const std::string elementVolumeName;
1322  const std::string naturalBcElement;
1323  const std::string essentialBcElement;
1324 
1325  EshelbianCore(MoFEM::Interface &m_field);
1326  virtual ~EshelbianCore();
1327 
1329  PetscBool ghostFrameOn;
1332  double alpha;
1333  double gamma;
1334 
1336 
1337  boost::shared_ptr<BcDispVec> bcSpatialDispVecPtr;
1338  boost::shared_ptr<BcRotVec> bcSpatialRotationVecPtr;
1339  boost::shared_ptr<TractionBcVec> bcSpatialTraction;
1340  boost::shared_ptr<TractionFreeBc> bcSpatialFreeTraction;
1341 
1342  template <typename BC>
1343  MoFEMErrorCode getBc(boost::shared_ptr<BC> &bc_vec_ptr,
1344  const std::string block_set_name,
1345  const int nb_attributes) {
1348  auto block_name = it->getName();
1349  if (block_name.compare(0, block_set_name.length(), block_set_name) == 0) {
1350  std::vector<double> block_attributes;
1351  CHKERR it->getAttributes(block_attributes);
1352  if (block_attributes.size() != nb_attributes) {
1353  SETERRQ3(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1354  "In block %s expected %d attributes, but given %d",
1355  it->getName().c_str(), nb_attributes,
1356  block_attributes.size());
1357  }
1358  Range faces;
1359  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
1360  true);
1361  bc_vec_ptr->emplace_back(block_name, block_attributes, faces);
1362  }
1363  }
1365  }
1366 
1368  bcSpatialDispVecPtr = boost::make_shared<BcDispVec>();
1369  return getBc(bcSpatialDispVecPtr, "SPATIAL_DISP_BC", 6);
1370  }
1371 
1373  bcSpatialRotationVecPtr = boost::make_shared<BcRotVec>();
1374  return getBc(bcSpatialRotationVecPtr, "SPATIAL_ROTATION_BC", 4);
1375  }
1376 
1378  bcSpatialTraction = boost::make_shared<TractionBcVec>();
1379  return getBc(bcSpatialTraction, "SPATIAL_TRACTION_BC", 6);
1380  }
1381 
1383  boost::shared_ptr<TractionFreeBc> &bc_ptr,
1384  const std::string disp_block_set_name,
1385  const std::string rot_block_set_name);
1386  inline MoFEMErrorCode
1389  boost::shared_ptr<TractionFreeBc>(new TractionFreeBc());
1390  return getTractionFreeBc(meshset, bcSpatialFreeTraction, "SPATIAL_DISP_BC",
1391  "SPATIAL_ROTATION_BC");
1392  }
1393 
1394  MoFEMErrorCode addFields(const EntityHandle meshset = 0);
1397  MoFEMErrorCode addDMs(const BitRefLevel bit = BitRefLevel().set(0));
1398 
1400  const double lambda,
1401  const double mu,
1402  const double sigma_y);
1403 
1404  MoFEMErrorCode addMaterial_HMHMooneyRivlin(const int tape, const double alpha,
1405  const double beta,
1406  const double lambda,
1407  const double sigma_y);
1408 
1410  const int tag, const bool do_rhs, const bool do_lhs,
1411  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> &fe);
1412 
1414  const int tag, const bool add_elastic, const bool add_material,
1415  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> &fe_rhs,
1416  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> &fe_lhs);
1417 
1419  const bool add_elastic, const bool add_material,
1420  boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> &fe_rhs,
1421  boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> &fe_lhs);
1422 
1423  MoFEMErrorCode setElasticElementOps(const int tag);
1425 
1426  MoFEMErrorCode setUpTSElastic(TS ts, Mat m, Vec f);
1427  MoFEMErrorCode solveElastic(TS ts, Vec x);
1428 
1429  MoFEMErrorCode postProcessResults(const int tag, const std::string file);
1430 };
1431 
1432 } // namespace EshelbianPlasticity
1433 
1434 #endif //__ESHELBIAN_PLASTICITY_HPP__
boost::shared_ptr< PhysicalEquations > physicsPtr
material physical equations
std::map< EntityType, std::vector< MatrixDouble > > invBx
MoFEMErrorCode setGenericFaceElementOps(const bool add_elastic, const bool add_material, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &fe_rhs, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &fe_lhs)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialPhysical_du_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpSpatialPhysical_du_dBubble(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpSpatialConsistencySchur_domega_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode setBaseVolumeElementOps(const int tag, const bool do_rhs, const bool do_lhs, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
VectorDouble nF
local right hand side vector
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
virtual moab::Interface & get_moab()=0
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpJacobian(const std::string &field_name, const int tag, const bool eval_rhs, const bool eval_lhs, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< PhysicalEquations > &physics_ptr)
MoFEMErrorCode setElasticElementOps(const int tag)
std::map< EntityType, std::vector< MatrixDouble > > invUx
OpSpatialPrj(const std::string &row_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpDispBc(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcDispVec > &bc_disp_ptr)
MoFEMErrorCode integrate(EntData &data)
boost::shared_ptr< EpFEMethod > schurAssembly
DM dmPrjSpatial
Project spatial displacement.
OpAssembleBasic(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type, const bool assemble_symmetry, const bool assemble_suu, const bool assemble_suu_suu)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
FTensor::Tensor3< adouble, 3, 3, 3 > ATensor3
OpRotationBc(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcRotVec > &bc_rot_ptr)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
boost::shared_ptr< BcRotVec > bcRotPtr
boost::shared_ptr< TractionBcVec > bcData
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialConsistencyDivTerm(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
boost::shared_ptr< BcRotVec > bcSpatialRotationVecPtr
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
PhysicalEquations(const int size_active, const int size_dependent)
OpSpatialRotation(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpDispBc_dx(const std::string &row_field_name, const std::string &col_field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcDispVec > &bc_disp_ptr)
OpAssembleVolume(const std::string &field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type)
std::vector< EntityHandle > & mapGaussPts
OpSpatialRotation_domega_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpAssembleFace(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type, const bool assemble_symmetry, const bool assemble_suu, const bool assemble_suu_suu)
MoFEMErrorCode addDMs(const BitRefLevel bit=BitRefLevel().set(0))
OpSpatialPhysical_du_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:77
base class for all interface classes
OpAssembleVolume(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type, const bool assemble_symmetry, const bool assemble_suu, const bool assemble_suu_suu)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
VectorBoundedArray< int, 3 > VectorInt3
Definition: Types.hpp:83
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialPhysical_du_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
std::vector< BcRot > BcRotVec
MoFEMErrorCode integrate(EntData &data)
virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Assemble and calculate invBOmega and BOmega.
EpElement(MoFEM::Interface &m_field)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Assemble and calculate invOmegaB and OmegaB.
FTensor::Tensor3< double, 3, 3, 3 > getDiffRotationFormVector(FTensor::Tensor1< T, 3 > &t_omega)
virtual MoFEMErrorCode assemble(EntData &data)
OpAssembleBasic(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpSpatialRotation_domega_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
EshelbianCore(MoFEM::Interface &m_field)
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
OpSpatialConsistencySchur_dP_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
std::map< EntityType, std::vector< MatrixDouble > > PB
MoFEMErrorCode integrate(EntData &data)
MoFEMErrorCode addFields(const EntityHandle meshset=0)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
std::map< EntityType, std::vector< MatrixDouble > > invPB
static const MOFEMuuid IDD_MOFEMEshelbianCrackInterface
OpSpatialConsistency_dBubble_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
std::map< EntityType, std::vector< MatrixDouble > > invUP
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
boost::shared_ptr< BcRotVec > bcRotPtr
std::vector< Range > TractionFreeBc
OpSpatialConsistencySchur_domega_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpCalculateStrainEnergy(const std::string field_name, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, boost::shared_ptr< double > &e)
MoFEMErrorCode integrate(EntData &row_data)
OpSpatialConsistency_dP_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpSpatialPhysical(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
boost::shared_ptr< VectorDouble > VectorPtr
std::vector< TractionBc > TractionBcVec
OpTractionBc(std::string row_field, FeTractionBc &bc_fe)
OpSpatialConsistency_dBubble_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
static DTensor2Ptr get_VecTensor2(std::vector< double > &v)
boost::shared_ptr< MatrixDouble > MatrixPtr
OpSpatialConsistency_dP_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpSpatialRotation_domega_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode addStreachSchurMatrix(Mat Suu, AO aoSuu)
virtual MoFEMErrorCode evaluateLhs(EntData &data)=0
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MatrixDouble K
local tangent matrix
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > getFTensor3FromMat(MatrixDouble &m)
OpSpatialConsistencyP(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
boost::shared_ptr< BcDispVec > bcSpatialDispVecPtr
double tol
virtual MoFEMErrorCode evaluateRhs(EntData &data)=0
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode postProcessResults(const int tag, const std::string file)
MoFEMErrorCode addBoundaryFiniteElement(const EntityHandle meshset=0)
MoFEMErrorCode addVolumeFiniteElement(const EntityHandle meshset=0)
OpCalculateRotationAndSpatialGradient(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode setGenericVolumeElementOps(const int tag, const bool add_elastic, const bool add_material, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe_rhs, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe_lhs)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
boost::shared_ptr< TractionFreeBc > bcSpatialFreeTraction
TractionBc(std::string name, std::vector< double > &attr, Range &faces)
MoFEMErrorCode setUpTSElastic(TS ts, Mat m, Vec f)
OpSpatialPrj_dx_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
BcRot(std::string name, std::vector< double > &attr, Range &faces)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
OpSpatialConsistencySchur_dBubble_dx(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode addMaterial_HMHHStVenantKirchhoff(const int tape, const double lambda, const double mu, const double sigma_y)
OpSpatialPhysical_du_du(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Assemble and calculate invBx and xB.
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Assemble and saves UB, and invUB.
MoFEMErrorCode getTractionFreeBc(const EntityHandle meshset, boost::shared_ptr< TractionFreeBc > &bc_ptr, const std::string disp_block_set_name, const std::string rot_block_set_name)
OpSpatialPrj_dx_dw(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
static DTensor0Ptr get_VecTensor0(std::vector< double > &v)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
OpSpatialConsistencySchur_dP_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
FeTractionBc(MoFEM::Interface &m_field, const std::string field_name, boost::shared_ptr< TractionBcVec > &bc)
OpSpatialConsistencySchur_domega_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
boost::shared_ptr< PhysicalEquations > physicalEquations
OpSpatialRotation_domega_dBubble(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode solveElastic(TS ts, Vec x)
FTensor::Tensor2< adouble, 3, 3 > ATensor2
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > DTensor0Ptr
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
boost::shared_ptr< BcDispVec > bcDispPtr
DataForcesAndSourcesCore::EntData EntData
OpSpatialConsistencySchur_dBubble_dBubble(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpPostProcDataStructure(moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, const std::string field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode getBc(boost::shared_ptr< BC > &bc_vec_ptr, const std::string block_set_name, const int nb_attributes)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field...
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcRhs
OpSpatialConsistencySchur_dP_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
#define CHKERR
Inline error check.
Definition: definitions.h:594
MoFEMErrorCode getSpatialTractionFreeBc(const EntityHandle meshset=0)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
DM dM
Coupled problem all fields.
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
FaceElementForcesAndSourcesCore::UserDataOperator FaceUserDataOperator
ForcesAndSourcesCore::UserDataOperator UserDataOperator
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
DM dmElasticSchurStreach
Sub problem of dmElastic Schur.
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:75
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:85
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
virtual MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
OpSpatialConsistencySchur_dBubble_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
EpElement(MoFEM::Interface &m_field, const std::string field_name, boost::shared_ptr< TractionBcVec > &bc)
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeLhs
OpSpatialConsistencySchur_dP_dBubble(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcLhs
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > DTensor3Ptr
const VectorInt & getIndices() const
Get global indices of dofs on entity.
boost::shared_ptr< TractionBcVec > bcSpatialTraction
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
virtual OpJacobian * returnOpJacobian(const std::string &field_name, const int tag, const bool eval_rhs, const bool eval_lhs, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< PhysicalEquations > &physics_ptr)=0
Data on single entity (This is passed as argument to DataOperator::doWork)
OpAssembleFace(const std::string &field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type)
OpSpatialConsistencyBubble(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialEquilibrium_dw_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
BcDisp(std::string name, std::vector< double > &attr, Range &faces)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
boost::shared_ptr< BcDispVec > bcDispPtr
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:405
OpRotationBc_dx(const std::string &row_field_name, const std::string &col_field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcRotVec > &bc_rot_ptr)
MoFEMErrorCode integrate(EntData &data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
static DTensor3Ptr get_vecTensor3(std::vector< double > &v)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
std::map< EntityType, std::vector< MatrixDouble > > UP
VolumeElementForcesAndSourcesCore::UserDataOperator VolUserDataOperator
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:76
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
DM dmElasticSchurBubble
Sub problem of dmElastic Schur.
virtual MoFEMErrorCode recordTape(const int tag)=0
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
moab::Interface & postProcMesh
std::vector< BcDisp > BcDispVec
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
FTensor::Tensor2< double, 3, 3 > getRotationFormVector(FTensor::Tensor1< T, 3 > &t_omega)
FTensor::Tensor1< adouble, 3 > ATensor1
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode addBubbleSchurMatrix(Mat SuuSuu, AO aoSuuSuu)
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeRhs
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Getting interface of core database.
OpSpatialEquilibrium(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
virtual MoFEMErrorCode integrate(EntData &data)
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, 3, 3 > DTensor2Ptr
MoFEMErrorCode addMaterial_HMHMooneyRivlin(const int tape, const double alpha, const double beta, const double lambda, const double sigma_y)
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use