v0.10.0
EshelbianPlasticity.hpp
Go to the documentation of this file.
1 /**
2  * \file EshelbianPlasticity.hpp
3  * \brief Eshelbian plasticity interface
4  *
5  * \brief Problem implementation for mix element for large-strain elasticity
6  *
7  * \todo Implementation of plasticity
8  */
9 
10 #ifndef __ESHELBIAN_PLASTICITY_HPP__
11 #define __ESHELBIAN_PLASTICITY_HPP__
12 
13 namespace EshelbianPlasticity {
14 
16 
17 static const MOFEMuuid IDD_MOFEMEshelbianCrackInterface =
19 
20 typedef boost::shared_ptr<MatrixDouble> MatrixPtr;
21 typedef boost::shared_ptr<VectorDouble> VectorPtr;
22 
26  &m(0, 0), &m(1, 0), &m(2, 0), &m(3, 0), &m(4, 0), &m(5, 0), &m(6, 0),
27  &m(7, 0), &m(8, 0), &m(9, 0), &m(10, 0), &m(11, 0), &m(12, 0), &m(13, 0),
28  &m(14, 0), &m(15, 0), &m(16, 0), &m(17, 0), &m(18, 0), &m(19, 0),
29  &m(20, 0), &m(21, 0), &m(22, 0), &m(23, 0), &m(24, 0), &m(25, 0),
30  &m(26, 0));
31 }
32 
37 
38 struct EpElementBase {
39  SmartPetscObj<Mat> Suu;
40  SmartPetscObj<AO> aoSuu;
41  SmartPetscObj<Mat> SBubble;
42  SmartPetscObj<AO> aoSBubble;
43  SmartPetscObj<Mat> SOmega;
44  SmartPetscObj<AO> aoSOmega;
45  SmartPetscObj<Mat> Sw;
46  SmartPetscObj<AO> aoSw;
47  EpElementBase() = default;
48  virtual ~EpElementBase() = default;
49 
51  SmartPetscObj<AO> &aoSuu) {
53  this->Suu = Suu;
54  this->aoSuu = aoSuu;
56  }
57 
59  SmartPetscObj<AO> &aoSBubble) {
61  this->SBubble = SBubble;
62  this->aoSBubble = aoSBubble;
64  }
65 
67  SmartPetscObj<AO> &aoSw) {
69  this->Sw = Sw;
70  this->aoSw = aoSw;
72  }
73 
75  SmartPetscObj<AO> &aoSOmega) {
77  this->SOmega = SOmega;
78  this->aoSOmega = aoSOmega;
80  }
81 };
82 
83 template <typename E> struct EpElement : public E, EpElementBase {
84  EpElement(MoFEM::Interface &m_field) : E(m_field), EpElementBase() {}
85 };
86 
87 template <> struct EpElement<BasicMethod> : public FEMethod, EpElementBase {
89 };
90 
91 template <> struct EpElement<FEMethod> : public FEMethod, EpElementBase {
93 };
94 
95 struct EpFEMethod : EpElement<FEMethod> {
99  if (Suu)
100  CHKERR MatZeroEntries(Suu);
101  if (SBubble)
102  CHKERR MatZeroEntries(SBubble);
103  if (Sw)
104  CHKERR MatZeroEntries(Sw);
105  if (SOmega)
106  CHKERR MatZeroEntries(SOmega);
108  }
109 
112  auto assemble = [](Mat a) {
114  if (a) {
115  CHKERR MatAssemblyBegin(a, MAT_FINAL_ASSEMBLY);
116  CHKERR MatAssemblyEnd(a, MAT_FINAL_ASSEMBLY);
117  }
119  };
120  CHKERR assemble(Suu);
121  CHKERR assemble(SBubble);
122  CHKERR assemble(SOmega);
123  CHKERR assemble(Sw);
124  // std::string wait;
125  // CHKERR MatView(SOmega, PETSC_VIEWER_DRAW_WORLD);
126  // std::cin >> wait;
127  // CHKERR MatView(Sw, PETSC_VIEWER_DRAW_WORLD);
128  // std::cin >> wait;
129  // CHKERR MatView(SBubble, PETSC_VIEWER_DRAW_WORLD);
130  // std::cin >> wait;
131  // CHKERR MatView(Suu, PETSC_VIEWER_DRAW_WORLD);
132  // std::cin >> wait;
133  // CHKERR MatView(getTSB(), PETSC_VIEWER_DRAW_WORLD);
134  // std::cin >> wait;
136  }
137 };
138 
139 struct PhysicalEquations;
141  : public boost::enable_shared_from_this<DataAtIntegrationPts> {
142 
152 
163 
171 
192 
196 
198  return boost::shared_ptr<MatrixDouble>(shared_from_this(),
200  }
202  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &approxPAtPts);
203  }
204 
206  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divPAtPts);
207  }
208 
210  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &divSigmaAtPts);
211  }
212 
214  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wAtPts);
215  }
216 
218  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wDotAtPts);
219  }
220 
222  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &wDotDotAtPts);
223  }
224 
226  return boost::shared_ptr<MatrixDouble>(shared_from_this(),
228  }
229 
231  return boost::shared_ptr<MatrixDouble>(shared_from_this(),
233  }
234 
236  return boost::shared_ptr<MatrixDouble>(shared_from_this(),
238  }
239 
241  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &rotAxisAtPts);
242  }
243 
245  return boost::shared_ptr<MatrixDouble>(shared_from_this(),
246  &rotAxisDotAtPts);
247  }
248 
250  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
251  }
252 
254  return boost::shared_ptr<MatrixDouble>(shared_from_this(), &GAtPts);
255  }
256 
257  // Not really data at integration points, used to calculate Schur complement
258 
259  struct BlockMatData {
260 
261  std::string rowField;
262  std::string colField;
263  EntityType rowType;
264  EntityType colType;
265  int rowSide;
266  int colSide;
270 
272 
273  BlockMatData(const std::string row_field, const std::string col_field,
274  EntityType row_type, EntityType col_type, int row_side,
275  int col_side, const MatrixDouble &m, const VectorInt row_ind,
276  VectorInt col_ind)
277  : rowField(row_field), colField(col_field), rowType(row_type),
278  colType(col_type), rowSide(row_side), colSide(col_side),
279  setAtElement(true) {
280 
281  M.resize(m.size1(), m.size2(), false);
282  noalias(M) = m;
283  rowInd.resize(row_ind.size(), false);
284  noalias(rowInd) = row_ind;
285  colInd.resize(col_ind.size(), false);
286  noalias(colInd) = col_ind;
287  }
288 
289  void setInd(const VectorInt &row_ind, const VectorInt &col_ind) const {
290  auto &const_row_ind = const_cast<VectorInt &>(rowInd);
291  auto &const_col_ind = const_cast<VectorInt &>(colInd);
292  const_row_ind.resize(row_ind.size(), false);
293  noalias(const_row_ind) = row_ind;
294  const_col_ind.resize(col_ind.size(), false);
295  noalias(const_col_ind) = col_ind;
296  }
297 
298  void setMat(const MatrixDouble &m) const {
299  auto &const_m = const_cast<MatrixDouble &>(M);
300  const_m.resize(m.size1(), m.size2(), false);
301  noalias(const_m) = m;
302  }
303 
304  void addMat(const MatrixDouble &m) const {
305  auto &const_m = const_cast<MatrixDouble &>(M);
306  const_m += m;
307  }
308 
309  void clearMat() const {
310  auto &const_m = const_cast<MatrixDouble &>(M);
311  const_m.clear();
312  }
313 
314  void setSetAtElement() const {
315  bool &set = const_cast<bool &>(setAtElement);
316  set = true;
317  }
318 
319  void unSetAtElement() const {
320  bool &set = const_cast<bool &>(setAtElement);
321  set = false;
322  }
323  };
324 
325  typedef multi_index_container<
326  BlockMatData,
327  indexed_by<
328 
329  ordered_unique<
330 
331  composite_key<
332  BlockMatData,
333 
334  member<BlockMatData, std::string, &BlockMatData::rowField>,
335  member<BlockMatData, std::string, &BlockMatData::colField>,
336  member<BlockMatData, EntityType, &BlockMatData::rowType>,
337  member<BlockMatData, EntityType, &BlockMatData::colType>,
338  member<BlockMatData, int, &BlockMatData::rowSide>,
339  member<BlockMatData, int, &BlockMatData::colSide>
340 
341  >>,
342  ordered_non_unique<
343 
344  composite_key<
345  BlockMatData,
346 
347  member<BlockMatData, std::string, &BlockMatData::rowField>,
348  member<BlockMatData, std::string, &BlockMatData::colField>,
349  member<BlockMatData, EntityType, &BlockMatData::rowType>,
350  member<BlockMatData, EntityType, &BlockMatData::colType>
351 
352  >>,
353  ordered_non_unique<
354 
355  composite_key<
356  BlockMatData,
357 
358  member<BlockMatData, std::string, &BlockMatData::rowField>,
359  member<BlockMatData, std::string, &BlockMatData::colField>
360 
361  >>,
362  ordered_non_unique<
363  member<BlockMatData, std::string, &BlockMatData::rowField>>,
364  ordered_non_unique<
365  member<BlockMatData, std::string, &BlockMatData::colField>>
366 
367  >>
369 
373 
374  boost::shared_ptr<PhysicalEquations> physicsPtr;
375 };
376 
377 struct OpJacobian;
378 
380 
385 
388 
389  PhysicalEquations() = delete;
390  PhysicalEquations(const int size_active, const int size_dependent)
391  : activeVariables(size_active, 0), dependentVariables(size_dependent, 0),
392  dependentVariablesDirevatives(size_dependent * size_active, 0) {}
393  virtual ~PhysicalEquations() {}
394 
395  virtual MoFEMErrorCode recordTape(const int tag, DTensor2Ptr *t_h) = 0;
396 
397  virtual OpJacobian *
398  returnOpJacobian(const std::string &field_name, const int tag,
399  const bool eval_rhs, const bool eval_lhs,
400  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
401  boost::shared_ptr<PhysicalEquations> &physics_ptr) = 0;
402 
403  std::vector<double> activeVariables;
404  std::vector<double> dependentVariables;
405  std::vector<double> dependentVariablesDirevatives;
406 
407  /** \name Active variables */
408 
409  /**@{*/
410 
411  template <int S>
412  inline static DTensor2Ptr get_VecTensor2(std::vector<double> &v) {
413  return DTensor2Ptr(&v[S + 0], &v[S + 1], &v[S + 2], &v[S + 3], &v[S + 4],
414  &v[S + 5], &v[S + 6], &v[S + 7], &v[S + 8]);
415  }
416 
417  template <int S>
418  inline static DTensor0Ptr get_VecTensor0(std::vector<double> &v) {
419  return DTensor0Ptr(&v[S + 0]);
420  }
421 
422  template <int S0, int S1, int S2>
423  inline static DTensor3Ptr get_vecTensor3(std::vector<double> &v) {
424  constexpr int A00 = 18 * S0 + 18 * 0 + 9 * S1 + 3 * S2;
425  constexpr int A01 = 18 * S0 + 18 * 1 + 9 * S1 + 3 * S2;
426  constexpr int A02 = 18 * S0 + 18 * 2 + 9 * S1 + 3 * S2;
427  constexpr int A10 = 18 * S0 + 18 * 3 + 9 * S1 + 3 * S2;
428  constexpr int A11 = 18 * S0 + 18 * 4 + 9 * S1 + 3 * S2;
429  constexpr int A12 = 18 * S0 + 18 * 5 + 9 * S1 + 3 * S2;
430  constexpr int A20 = 18 * S0 + 18 * 6 + 9 * S1 + 3 * S2;
431  constexpr int A21 = 18 * S0 + 18 * 7 + 9 * S1 + 3 * S2;
432  constexpr int A22 = 18 * S0 + 18 * 8 + 9 * S1 + 3 * S2;
433  return DTensor3Ptr(
434 
435  &v[A00 + 0], &v[A00 + 1], &v[A00 + 2], &v[A01 + 0], &v[A01 + 1],
436  &v[A01 + 2], &v[A02 + 0], &v[A02 + 1], &v[A02 + 2],
437 
438  &v[A10 + 0], &v[A10 + 1], &v[A10 + 2], &v[A11 + 0], &v[A11 + 1],
439  &v[A11 + 2], &v[A12 + 0], &v[A12 + 1], &v[A12 + 2],
440 
441  &v[A20 + 0], &v[A20 + 1], &v[A20 + 2], &v[A21 + 0], &v[A21 + 1],
442  &v[A21 + 2], &v[A22 + 0], &v[A22 + 1], &v[A22 + 2]
443 
444  );
445  }
446 
447  /**@}*/
448 
449  /** \name Active variables */
450 
451  /**@{*/
452 
453  inline DTensor2Ptr get_h() { return get_VecTensor2<0>(activeVariables); }
454  inline DTensor2Ptr get_H() { return get_VecTensor2<9>(activeVariables); }
455 
456  /**@}*/
457 
458  /** \name Dependent variables */
459 
460  /**@{*/
461 
462  inline DTensor2Ptr get_P() { return get_VecTensor2<0>(dependentVariables); }
464  return get_VecTensor2<9>(dependentVariables);
465  }
466  inline DTensor0Ptr get_Phi() {
467  return get_VecTensor0<18>(dependentVariables);
468  }
469  inline double &get_PhiRef() { return dependentVariables[18]; }
471  return get_VecTensor2<9 + 9 + 1>(dependentVariables);
472  }
473 
474  /**@}*/
475 
476  /** \name Derivatives of dependent variables */
477 
478  /**@{*/
479 
481  return get_vecTensor3<0, 0, 0>(dependentVariablesDirevatives);
482  }
484  return get_vecTensor3<0, 1, 0>(dependentVariablesDirevatives);
485  }
487  return get_vecTensor3<0, 0, 1>(dependentVariablesDirevatives);
488  }
490  return get_vecTensor3<0, 1, 1>(dependentVariablesDirevatives);
491  }
493  return get_vecTensor3<0, 0, 2>(dependentVariablesDirevatives);
494  }
496  return get_vecTensor3<0, 1, 2>(dependentVariablesDirevatives);
497  }
498 
500  return get_vecTensor3<9, 0, 0>(dependentVariablesDirevatives);
501  }
503  return get_vecTensor3<9, 1, 0>(dependentVariablesDirevatives);
504  }
506  return get_vecTensor3<9, 0, 1>(dependentVariablesDirevatives);
507  }
509  return get_vecTensor3<9, 1, 1>(dependentVariablesDirevatives);
510  }
512  return get_vecTensor3<9, 0, 2>(dependentVariablesDirevatives);
513  }
515  return get_vecTensor3<9, 1, 2>(dependentVariablesDirevatives);
516  }
517 
519  return get_VecTensor2<18 * (9 + 9) + 0>(dependentVariablesDirevatives);
520  }
522  return get_VecTensor2<18 * (9 + 9) + 9>(dependentVariablesDirevatives);
523  }
524 
526  return get_vecTensor3<9 + 9 + 1, 0, 0>(dependentVariablesDirevatives);
527  }
529  return get_vecTensor3<9 + 9 + 1, 1, 0>(dependentVariablesDirevatives);
530  }
532  return get_vecTensor3<9 + 9 + 1, 0, 1>(dependentVariablesDirevatives);
533  }
535  return get_vecTensor3<9 + 9 + 1, 1, 1>(dependentVariablesDirevatives);
536  }
538  return get_vecTensor3<9 + 9 + 1, 0, 2>(dependentVariablesDirevatives);
539  }
541  return get_vecTensor3<9 + 9 + 1, 1, 2>(dependentVariablesDirevatives);
542  }
543 
544  /**@}*/
545 };
546 
547 struct BcDisp {
548  BcDisp(std::string name, std::vector<double> &attr, Range &faces);
549  std::string blockName;
550  Range faces;
553 };
554 typedef std::vector<BcDisp> BcDispVec;
555 
556 struct BcRot {
557  BcRot(std::string name, std::vector<double> &attr, Range &faces);
558  std::string blockName;
559  Range faces;
561  double theta;
562 };
563 typedef std::vector<BcRot> BcRotVec;
564 
565 typedef std::vector<Range> TractionFreeBc;
566 
567 struct TractionBc {
568  TractionBc(std::string name, std::vector<double> &attr, Range &faces);
569  std::string blockName;
570  Range faces;
573 };
574 typedef std::vector<TractionBc> TractionBcVec;
575 
576 struct OpJacobian : public UserDataOperator {
577  const int tAg; ///< adol-c tape
578  const bool evalRhs;
579  const bool evalLhs;
580  boost::shared_ptr<DataAtIntegrationPts>
581  dataAtPts; ///< data at integration pts
582  boost::shared_ptr<PhysicalEquations>
583  physicsPtr; ///< material physical equations
584 
585  OpJacobian(const std::string &field_name, const int tag, const bool eval_rhs,
586  const bool eval_lhs,
587  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
588  boost::shared_ptr<PhysicalEquations> &physics_ptr)
589  : UserDataOperator(field_name, OPROW), tAg(tag), evalRhs(eval_rhs),
590  evalLhs(eval_lhs), dataAtPts(data_ptr), physicsPtr(physics_ptr) {}
591 
592  virtual MoFEMErrorCode evaluateRhs(EntData &data) = 0;
593  virtual MoFEMErrorCode evaluateLhs(EntData &data) = 0;
594 
595  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
596 };
597 
598 template <typename T> struct OpAssembleBasic : public T {
599 
600  const bool assembleSymmetry;
601 
602  boost::shared_ptr<DataAtIntegrationPts>
603  dataAtPts; ///< data at integration pts
604 
605  OpAssembleBasic(const std::string &field_name,
606  boost::shared_ptr<DataAtIntegrationPts> data_ptr,
607  const char type)
608  : T(field_name, type), dataAtPts(data_ptr), assembleSymmetry(false) {}
609 
610  OpAssembleBasic(const std::string &row_field, const std::string &col_field,
611  boost::shared_ptr<DataAtIntegrationPts> data_ptr,
612  const char type, const bool assemble_symmetry)
613  : T(row_field, col_field, type, false), dataAtPts(data_ptr),
614  assembleSymmetry(assemble_symmetry) {}
615 
616  VectorDouble nF; ///< local right hand side vector
617  MatrixDouble K; ///< local tangent matrix
619 
622  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
624  }
625 
626  virtual MoFEMErrorCode integrate(int row_side, EntityType row_type,
627  EntData &data) {
629  CHKERR integrate(data);
631  }
632 
633  virtual MoFEMErrorCode assemble(EntData &data) {
635  double *vec_ptr = &*nF.begin();
636  int nb_dofs = data.getIndices().size();
637  int *ind_ptr = &*data.getIndices().begin();
638  CHKERR VecSetValues(this->getTSf(), nb_dofs, ind_ptr, vec_ptr, ADD_VALUES);
640  }
641 
642  virtual MoFEMErrorCode assemble(int row_side, EntityType row_type,
643  EntData &data) {
645  CHKERR assemble(data);
647  }
648 
649  virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) {
651  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
653  }
654 
655  virtual MoFEMErrorCode assemble(int row_side, int col_side,
656  EntityType row_type, EntityType col_type,
657  EntData &row_data, EntData &col_data) {
659 
660  auto &bmc = dataAtPts->blockMatContainor;
661  int *row_ind_ptr = &*row_data.getIndices().begin();
662  int *col_ind_ptr = &*col_data.getIndices().begin();
663  int row_nb_dofs = row_data.getIndices().size();
664  int col_nb_dofs = col_data.getIndices().size();
665 
666  auto add_block = [&](auto &row_name, auto &col_name, auto row_side,
667  auto col_side, auto row_type, auto col_type,
668  const auto &m, const auto &row_ind,
669  const auto &col_ind) {
670  auto it = bmc.get<0>().find(boost::make_tuple(
671  row_name, col_name, row_type, col_type, row_side, col_side));
672  if (it != bmc.get<0>().end()) {
673  it->setMat(m);
674  it->setInd(row_ind, col_ind);
675  it->setSetAtElement();
676  } else
678  row_name, col_name, row_type, col_type, row_side, col_side, m,
679  row_ind, col_ind));
680  };
681 
682  add_block(this->rowFieldName, this->colFieldName, row_side, col_side,
683  row_type, col_type, K, row_data.getIndices(),
684  col_data.getIndices());
685  if (assembleSymmetry) {
686  transposeK.resize(col_nb_dofs, row_nb_dofs, false);
687  noalias(transposeK) = trans(K);
688  add_block(this->colFieldName, this->rowFieldName, col_side, row_side,
689  col_type, row_type, transposeK, col_data.getIndices(),
690  row_data.getIndices());
691  }
692 
693  double *mat_ptr = &*K.data().begin();
694  CHKERR MatSetValues(this->getTSB(), row_nb_dofs, row_ind_ptr, col_nb_dofs,
695  col_ind_ptr, mat_ptr, ADD_VALUES);
696  if (assembleSymmetry) {
697  double *mat_ptr = &*transposeK.data().begin();
698  CHKERR MatSetValues(this->getTSB(), col_nb_dofs, col_ind_ptr, row_nb_dofs,
699  row_ind_ptr, mat_ptr, ADD_VALUES);
700  }
701 
703  }
704 
705  MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
707  if (data.getIndices().empty())
709  nF.resize(data.getIndices().size(), false);
710  nF.clear();
711  CHKERR integrate(side, type, data);
712  CHKERR assemble(side, type, data);
714  }
715 
716  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
717  EntityType col_type, EntData &row_data,
718  EntData &col_data) {
720  if (row_data.getIndices().empty())
722  if (col_data.getIndices().empty())
724  K.resize(row_data.getIndices().size(), col_data.getIndices().size(), false);
725  K.clear();
726  CHKERR integrate(row_data, col_data);
727  CHKERR assemble(row_side, col_side, row_type, col_type, row_data, col_data);
729  }
730 };
731 
732 struct OpAssembleVolume : public OpAssembleBasic<VolUserDataOperator> {
733  OpAssembleVolume(const std::string &field,
734  boost::shared_ptr<DataAtIntegrationPts> data_ptr,
735  const char type)
736  : OpAssembleBasic<VolUserDataOperator>(field, data_ptr, type) {}
737 
738  OpAssembleVolume(const std::string &row_field, const std::string &col_field,
739  boost::shared_ptr<DataAtIntegrationPts> data_ptr,
740  const char type, const bool assemble_symmetry)
741  : OpAssembleBasic<VolUserDataOperator>(row_field, col_field, data_ptr,
742  type, assemble_symmetry) {}
743 };
744 
745 struct OpAssembleFace : public OpAssembleBasic<FaceUserDataOperator> {
746  OpAssembleFace(const std::string &field,
747  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
748  const char type)
749  : OpAssembleBasic<FaceUserDataOperator>(field, data_ptr, type) {}
750 
751  OpAssembleFace(const std::string &row_field, const std::string &col_field,
752  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
753  const char type, const bool assemble_symmetry)
754  : OpAssembleBasic<FaceUserDataOperator>(row_field, col_field, data_ptr,
755  type, assemble_symmetry) {}
756 };
757 
759  boost::shared_ptr<DataAtIntegrationPts>
760  dataAtPts; ///< data at integration pts
762  const std::string &field_name,
763  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
764  : VolUserDataOperator(field_name, OPROW), dataAtPts(data_ptr) {}
765  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
766 };
767 
769  const double alphaW;
770  const double alphaRho;
771  OpSpatialEquilibrium(const std::string &field_name,
772  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
773  const double alpha, const double rho)
774  : OpAssembleVolume(field_name, data_ptr, OPROW), alphaW(alpha),
775  alphaRho(rho) {}
777 };
778 
780  OpSpatialRotation(const std::string &field_name,
781  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
782  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
784 };
785 
787  const double alphaU;
788  OpSpatialPhysical(const std::string &field_name,
789  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
790  const double alpha)
791  : OpAssembleVolume(field_name, data_ptr, OPROW), alphaU(alpha) {}
793 };
794 
796  OpSpatialConsistencyP(const std::string &field_name,
797  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
798  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
800 };
801 
803  OpSpatialConsistencyBubble(const std::string &field_name,
804  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
805  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
807 };
808 
810  OpSpatialConsistencyDivTerm(const std::string &field_name,
811  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
812  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
814 };
815 
816 struct OpDispBc : public OpAssembleFace {
817  boost::shared_ptr<BcDispVec> bcDispPtr;
818  OpDispBc(const std::string &field_name,
819  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
820  boost::shared_ptr<BcDispVec> &bc_disp_ptr)
821  : OpAssembleFace(field_name, data_ptr, OPROW), bcDispPtr(bc_disp_ptr) {}
823 };
824 
825 struct OpDispBc_dx : public OpAssembleFace {
826  boost::shared_ptr<BcDispVec> bcDispPtr;
827  OpDispBc_dx(const std::string &row_field_name,
828  const std::string &col_field_name,
829  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
830  boost::shared_ptr<BcDispVec> &bc_disp_ptr)
831  : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL,
832  false),
833  bcDispPtr(bc_disp_ptr) {}
834  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
835 };
836 
837 struct OpRotationBc : public OpAssembleFace {
838  boost::shared_ptr<BcRotVec> bcRotPtr;
839  OpRotationBc(const std::string &field_name,
840  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
841  boost::shared_ptr<BcRotVec> &bc_rot_ptr)
842  : OpAssembleFace(field_name, data_ptr, OPROW), bcRotPtr(bc_rot_ptr) {}
844 };
845 
847  boost::shared_ptr<BcRotVec> bcRotPtr;
848  OpRotationBc_dx(const std::string &row_field_name,
849  const std::string &col_field_name,
850  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
851  boost::shared_ptr<BcRotVec> &bc_rot_ptr)
852  : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL,
853  false),
854  bcRotPtr(bc_rot_ptr) {}
855  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
856 };
857 
858 struct FeTractionBc;
859 
861  OpTractionBc(std::string row_field, FeTractionBc &bc_fe)
862  : FaceUserDataOperator(row_field, FaceUserDataOperator::OPROW),
863  bcFe(bc_fe) {}
864  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
865 
866 protected:
870 };
871 
872 struct FeTractionBc : public FEMethod {
874  FeTractionBc(MoFEM::Interface &m_field, const std::string field_name,
875  boost::shared_ptr<TractionBcVec> &bc)
876  : FEMethod(), mField(m_field), bcData(bc), fieldName(field_name) {}
878  MoFEMErrorCode postProcess() { return 0; }
879 
880  friend MoFEMErrorCode OpTractionBc::doWork(int side, EntityType type,
881  EntData &data);
882 
883 protected:
885  boost::shared_ptr<TractionBcVec> bcData;
886  std::string fieldName;
887 };
888 
889 template <>
891  EpElement(MoFEM::Interface &m_field, const std::string field_name,
892  boost::shared_ptr<TractionBcVec> &bc)
893  : FeTractionBc(m_field, field_name, bc), EpElementBase() {}
898  }
899 };
900 
902  OpSpatialEquilibrium_dw_dP(const std::string &row_field,
903  const std::string &col_field,
904  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
905  const bool assemble_off = false)
906  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
907  assemble_off) {
908  sYmm = false;
909  }
910  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
911 };
912 
914  const double alphaW;
915  const double alphaRho;
916  OpSpatialEquilibrium_dw_dw(const std::string &row_field,
917  const std::string &col_field,
918  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
919  const double alpha, const double rho)
920  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false),
921  alphaW(alpha), alphaRho(rho) {
922  sYmm = true;
923  }
924  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
925 };
926 
928  const double alphaU;
929  OpSpatialPhysical_du_du(const std::string &row_field,
930  const std::string &col_field,
931  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
932  const double alpha)
933  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false),
934  alphaU(alpha) {
935  sYmm = false;
936  }
937  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
938 };
939 
941  OpSpatialPhysical_du_dP(const std::string &row_field,
942  const std::string &col_field,
943  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
944  const bool assemble_off = false)
945  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
946  assemble_off) {
947  sYmm = false;
948  }
949 
950  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
951 };
952 
955  const std::string &row_field, const std::string &col_field,
956  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
957  const bool assemble_off = false)
958  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
959  assemble_off) {
960  sYmm = false;
961  }
962 
963  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
964 };
965 
967  OpSpatialPhysical_du_domega(const std::string &row_field,
968  const std::string &col_field,
969  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
970  const bool assemble_off = false)
971  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
972  assemble_off) {
973  sYmm = false;
974  }
975  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
976 };
977 
979  OpSpatialPhysical_du_dx(const std::string &row_field,
980  const std::string &col_field,
981  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
982  const bool assemble_off = false)
983  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
984  assemble_off) {
985  sYmm = false;
986  }
987  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
988 };
989 
991  OpSpatialRotation_domega_dP(const std::string &row_field,
992  const std::string &col_field,
993  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
994  const bool assemble_off)
995  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
996  assemble_off) {
997  sYmm = false;
998  }
999  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1000 };
1001 
1004  const std::string &row_field, const std::string &col_field,
1005  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1006  const bool assemble_off)
1007  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
1008  assemble_off) {
1009  sYmm = false;
1010  }
1011  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1012 };
1013 
1016  const std::string &row_field, const std::string &col_field,
1017  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1018  const bool assemble_off = false)
1019  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
1020  assemble_off) {
1021  sYmm = false;
1022  }
1023  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1024 };
1025 
1028  const std::string &row_field, const std::string &col_field,
1029  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1030  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
1031  sYmm = false;
1032  }
1033  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1034 };
1035 
1038  const std::string &row_field, const std::string &col_field,
1039  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1040  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
1041  sYmm = false;
1042  }
1043  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1044 };
1045 
1047 
1049  std::vector<EntityHandle> &mapGaussPts;
1050  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
1051 
1053  std::vector<EntityHandle> &map_gauss_pts,
1054  const std::string field_name,
1055  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1056  : VolUserDataOperator(field_name, UserDataOperator::OPROW),
1057  postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
1058  dataAtPts(data_ptr) {}
1059 
1060  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
1061 };
1062 
1064  OpSpatialPrj(const std::string &row_field,
1065  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1066  : OpAssembleVolume(row_field, data_ptr, OPROW) {}
1068 };
1069 
1071  OpSpatialPrj_dx_dx(const std::string &row_field, const std::string &col_field,
1072  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1073  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
1074  // FIXME: That is symmetric
1075  sYmm = false;
1076  }
1077  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1078 };
1079 
1081  OpSpatialPrj_dx_dw(const std::string &row_field, const std::string &col_field,
1082  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1083  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
1084  sYmm = false;
1085  }
1086  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
1087 };
1088 
1090  OpSpatialSchurBegin(const std::string &row_field,
1091  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
1092  : OpAssembleVolume(row_field, data_ptr, OPROW) {
1093  sYmm = false;
1094  }
1095  MoFEMErrorCode integrate(EntData &row_data) { return 0; }
1096  MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data);
1097 };
1098 
1100  OpSpatialPreconditionMass(const std::string &row_field, MatrixPtr m_ptr)
1101  : OpAssembleVolume(row_field, nullptr, OPROW), mPtr(m_ptr) {
1102  sYmm = false;
1103  }
1104  MoFEMErrorCode integrate(EntData &row_data);
1105  MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data) {
1106  return 0;
1107  }
1108 
1109 private:
1111 };
1112 
1114  OpSpatialSchurEnd(const std::string &row_field,
1115  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
1116  const double eps)
1117  : OpAssembleVolume(row_field, data_ptr, OPROW), eps(eps) {
1118  sYmm = false;
1119  }
1120 
1121  MoFEMErrorCode integrate(EntData &row_data) { return 0; }
1122  MoFEMErrorCode assemble(int side, EntityType type, EntData &data);
1123  MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
1124  return assemble(side, type, data);
1125  }
1126 
1127 private:
1128  const double eps;
1129 
1133 
1134  map<std::string, MatrixDouble> invBlockMat;
1135  map<std::string, DataAtIntegrationPts::BlockMatContainor> blockMat;
1136 };
1137 
1139 
1140  OpCalculateStrainEnergy(const std::string field_name,
1141  boost::shared_ptr<DataAtIntegrationPts> data_at_pts,
1142  boost::shared_ptr<double> &e)
1143  : VolUserDataOperator(field_name, VolUserDataOperator::OPROW),
1144  dataAtPts(data_at_pts), energy(e) {}
1145 
1146  MoFEMErrorCode doWork(int side, EntityType type,
1148 
1149 private:
1150  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
1151  boost::shared_ptr<double> energy;
1152 };
1153 
1155 
1157 
1158  MoFEMErrorCode doWork(int side, EntityType type,
1160 };
1161 
1163 
1164  /**
1165  * \brief Getting interface of core database
1166  * @param uuid unique ID of interface
1167  * @param iface returned pointer to interface
1168  * @return error code
1169  */
1170  MoFEMErrorCode query_interface(const MOFEMuuid &uuid,
1171  UnknownInterface **iface) const;
1172 
1174 
1175  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
1176  boost::shared_ptr<PhysicalEquations> physicalEquations;
1177 
1178  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> elasticFeRhs;
1179  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> elasticFeLhs;
1180  boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> elasticBcLhs;
1181  boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> elasticBcRhs;
1182  boost::shared_ptr<EpFEMethod> schurAssembly;
1183 
1184  SmartPetscObj<DM> dM; ///< Coupled problem all fields
1185  SmartPetscObj<DM> dmElastic; ///< Elastic problem
1186  SmartPetscObj<DM> dmElasticSchurStreach; ///< Sub problem of dmElastic Schur
1187  SmartPetscObj<DM> dmElasticSchurBubble; ///< Sub problem of dmElastic Schur
1188  SmartPetscObj<DM>
1189  dmElasticSchurSpatialDisp; ///< Sub problem of dmElastic Schur
1190  SmartPetscObj<DM> dmElasticSchurOmega; ///< Sub problem of dmElastic Schur
1191  SmartPetscObj<DM> dmMaterial; ///< Material problem
1192 
1193  const std::string piolaStress;
1194  const std::string eshelbyStress;
1195  const std::string spatialDisp;
1196  const std::string materialDisp;
1197  const std::string streachTensor;
1198  const std::string rotAxis;
1199  const std::string materialGradient;
1200  const std::string tauField;
1201  const std::string lambdaField;
1202  const std::string bubbleField;
1203 
1204  const std::string elementVolumeName;
1205  const std::string naturalBcElement;
1206  const std::string essentialBcElement;
1207 
1208  EshelbianCore(MoFEM::Interface &m_field);
1209  virtual ~EshelbianCore();
1210 
1213  double alphaU;
1214  double alphaW;
1215  double alphaRho;
1216  double precEps;
1217 
1219 
1220  boost::shared_ptr<BcDispVec> bcSpatialDispVecPtr;
1221  boost::shared_ptr<BcRotVec> bcSpatialRotationVecPtr;
1222  boost::shared_ptr<TractionBcVec> bcSpatialTraction;
1223  boost::shared_ptr<TractionFreeBc> bcSpatialFreeTraction;
1224 
1225  template <typename BC>
1226  MoFEMErrorCode getBc(boost::shared_ptr<BC> &bc_vec_ptr,
1227  const std::string block_set_name,
1228  const int nb_attributes) {
1231  auto block_name = it->getName();
1232  if (block_name.compare(0, block_set_name.length(), block_set_name) == 0) {
1233  std::vector<double> block_attributes;
1234  CHKERR it->getAttributes(block_attributes);
1235  if (block_attributes.size() != nb_attributes) {
1236  SETERRQ3(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1237  "In block %s expected %d attributes, but given %d",
1238  it->getName().c_str(), nb_attributes,
1239  block_attributes.size());
1240  }
1241  Range faces;
1242  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
1243  true);
1244  bc_vec_ptr->emplace_back(block_name, block_attributes, faces);
1245  }
1246  }
1248  }
1249 
1251  bcSpatialDispVecPtr = boost::make_shared<BcDispVec>();
1252  return getBc(bcSpatialDispVecPtr, "SPATIAL_DISP_BC", 6);
1253  }
1254 
1256  bcSpatialRotationVecPtr = boost::make_shared<BcRotVec>();
1257  return getBc(bcSpatialRotationVecPtr, "SPATIAL_ROTATION_BC", 4);
1258  }
1259 
1261  bcSpatialTraction = boost::make_shared<TractionBcVec>();
1262  return getBc(bcSpatialTraction, "SPATIAL_TRACTION_BC", 6);
1263  }
1264 
1265  MoFEMErrorCode getTractionFreeBc(const EntityHandle meshset,
1266  boost::shared_ptr<TractionFreeBc> &bc_ptr,
1267  const std::string disp_block_set_name,
1268  const std::string rot_block_set_name);
1269  inline MoFEMErrorCode
1270  getSpatialTractionFreeBc(const EntityHandle meshset = 0) {
1272  boost::shared_ptr<TractionFreeBc>(new TractionFreeBc());
1273  return getTractionFreeBc(meshset, bcSpatialFreeTraction, "SPATIAL_DISP_BC",
1274  "SPATIAL_ROTATION_BC");
1275  }
1276 
1277  MoFEMErrorCode addFields(const EntityHandle meshset = 0);
1278  MoFEMErrorCode addVolumeFiniteElement(const EntityHandle meshset = 0);
1279  MoFEMErrorCode addBoundaryFiniteElement(const EntityHandle meshset = 0);
1280  MoFEMErrorCode addDMs(const BitRefLevel bit = BitRefLevel().set(0));
1281 
1283  const double lambda,
1284  const double mu,
1285  const double sigma_y);
1286 
1287  MoFEMErrorCode addMaterial_HMHMooneyRivlin(const int tape, const double alpha,
1288  const double beta,
1289  const double lambda,
1290  const double sigma_y);
1291 
1293  const int tag, const bool do_rhs, const bool do_lhs,
1294  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> &fe);
1295 
1297  const int tag, const bool add_elastic, const bool add_material,
1298  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> &fe_rhs,
1299  boost::shared_ptr<EpElement<VolumeElementForcesAndSourcesCore>> &fe_lhs);
1300 
1302  const bool add_elastic, const bool add_material,
1303  boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> &fe_rhs,
1304  boost::shared_ptr<EpElement<FaceElementForcesAndSourcesCore>> &fe_lhs);
1305 
1306  MoFEMErrorCode setElasticElementOps(const int tag);
1308 
1309  MoFEMErrorCode setUpTSElastic(TS ts, Mat m, Vec f);
1310  MoFEMErrorCode solveElastic(TS ts, Vec x);
1311 
1312  MoFEMErrorCode postProcessResults(const int tag, const std::string file);
1313 };
1314 
1315 } // namespace EshelbianPlasticity
1316 
1317 #endif //__ESHELBIAN_PLASTICITY_HPP__
static Index< 'm', 3 > m
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
@ L2
field with C-1 continuity
Definition: definitions.h:180
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
@ BLOCKSET
Definition: definitions.h:217
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:123
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:124
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#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.
const double T
const FTensor::Tensor2< T, Dim, Dim > Vec
std::vector< BcDisp > BcDispVec
DataForcesAndSourcesCore::EntData EntData
VolumeElementForcesAndSourcesCore::UserDataOperator VolUserDataOperator
boost::shared_ptr< VectorDouble > VectorPtr
std::vector< Range > TractionFreeBc
double f(const double v)
static const MOFEMuuid IDD_MOFEMEshelbianCrackInterface
std::vector< BcRot > BcRotVec
std::vector< TractionBc > TractionBcVec
boost::shared_ptr< MatrixDouble > MatrixPtr
FaceElementForcesAndSourcesCore::UserDataOperator FaceUserDataOperator
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > getFTensor3FromMat(MatrixDouble &m)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:73
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Types.hpp:56
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Types.hpp:51
VectorBoundedArray< int, 3 > VectorInt3
Definition: Types.hpp:87
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Types.hpp:92
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
DeprecatedCoreInterface Interface
Definition: Interface.hpp:1943
constexpr auto VecSetValues
constexpr auto MatSetValues
double rho
Definition: plastic.cpp:105
BcDisp(std::string name, std::vector< double > &attr, Range &faces)
BcRot(std::string name, std::vector< double > &attr, Range &faces)
BlockMatData(const std::string row_field, const std::string col_field, EntityType row_type, EntityType col_type, int row_side, int col_side, const MatrixDouble &m, const VectorInt row_ind, VectorInt col_ind)
void setInd(const VectorInt &row_ind, const VectorInt &col_ind) const
multi_index_container< BlockMatData, indexed_by< ordered_unique< composite_key< BlockMatData, member< BlockMatData, std::string, &BlockMatData::rowField >, member< BlockMatData, std::string, &BlockMatData::colField >, member< BlockMatData, EntityType, &BlockMatData::rowType >, member< BlockMatData, EntityType, &BlockMatData::colType >, member< BlockMatData, int, &BlockMatData::rowSide >, member< BlockMatData, int, &BlockMatData::colSide > > >, ordered_non_unique< composite_key< BlockMatData, member< BlockMatData, std::string, &BlockMatData::rowField >, member< BlockMatData, std::string, &BlockMatData::colField >, member< BlockMatData, EntityType, &BlockMatData::rowType >, member< BlockMatData, EntityType, &BlockMatData::colType > > >, ordered_non_unique< composite_key< BlockMatData, member< BlockMatData, std::string, &BlockMatData::rowField >, member< BlockMatData, std::string, &BlockMatData::colField > > >, ordered_non_unique< member< BlockMatData, std::string, &BlockMatData::rowField > >, ordered_non_unique< member< BlockMatData, std::string, &BlockMatData::colField > > > > BlockMatContainor
boost::shared_ptr< PhysicalEquations > physicsPtr
EpElement(MoFEM::Interface &m_field, const std::string field_name, boost::shared_ptr< TractionBcVec > &bc)
MoFEMErrorCode addBubbleSchurMatrix(SmartPetscObj< Mat > &SBubble, SmartPetscObj< AO > &aoSBubble)
MoFEMErrorCode addSpatialDispStressSchurMatrix(SmartPetscObj< Mat > &Sw, SmartPetscObj< AO > &aoSw)
MoFEMErrorCode addOmegaSchurMatrix(SmartPetscObj< Mat > &SOmega, SmartPetscObj< AO > &aoSOmega)
MoFEMErrorCode addStreachSchurMatrix(SmartPetscObj< Mat > &Suu, SmartPetscObj< AO > &aoSuu)
EpElement(MoFEM::Interface &m_field)
MoFEMErrorCode setUpTSElastic(TS ts, Mat m, Vec f)
SmartPetscObj< DM > dmMaterial
Material problem.
SmartPetscObj< DM > dmElasticSchurBubble
Sub problem of dmElastic Schur.
MoFEMErrorCode postProcessResults(const int tag, const std::string file)
MoFEMErrorCode addBoundaryFiniteElement(const EntityHandle meshset=0)
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeLhs
SmartPetscObj< DM > dmElasticSchurStreach
Sub problem of dmElastic Schur.
MoFEMErrorCode addDMs(const BitRefLevel bit=BitRefLevel().set(0))
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcRhs
MoFEMErrorCode addMaterial_HMHHStVenantKirchhoff(const int tape, const double lambda, const double mu, const double sigma_y)
MoFEMErrorCode setElasticElementOps(const int tag)
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)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
MoFEMErrorCode setBaseVolumeElementOps(const int tag, const bool do_rhs, const bool do_lhs, boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore >> &fe)
MoFEMErrorCode getSpatialTractionFreeBc(const EntityHandle meshset=0)
SmartPetscObj< DM > dmElasticSchurOmega
Sub problem of dmElastic Schur.
boost::shared_ptr< BcRotVec > bcSpatialRotationVecPtr
SmartPetscObj< DM > dmElasticSchurSpatialDisp
Sub problem of dmElastic Schur.
EshelbianCore(MoFEM::Interface &m_field)
MoFEMErrorCode addMaterial_HMHMooneyRivlin(const int tape, const double alpha, const double beta, const double lambda, const double sigma_y)
SmartPetscObj< DM > dM
Coupled problem all fields.
MoFEMErrorCode solveElastic(TS ts, Vec x)
MoFEMErrorCode addVolumeFiniteElement(const EntityHandle meshset=0)
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)
SmartPetscObj< DM > dmElastic
Elastic problem.
boost::shared_ptr< EpFEMethod > schurAssembly
boost::shared_ptr< BcDispVec > bcSpatialDispVecPtr
boost::shared_ptr< PhysicalEquations > physicalEquations
MoFEMErrorCode query_interface(const MOFEMuuid &uuid, UnknownInterface **iface) const
Getting interface of core database.
boost::shared_ptr< EpElement< VolumeElementForcesAndSourcesCore > > elasticFeRhs
boost::shared_ptr< TractionFreeBc > bcSpatialFreeTraction
boost::shared_ptr< TractionBcVec > bcSpatialTraction
MoFEMErrorCode addFields(const EntityHandle meshset=0)
MoFEMErrorCode setGenericFaceElementOps(const bool add_elastic, const bool add_material, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &fe_rhs, boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore >> &fe_lhs)
boost::shared_ptr< EpElement< FaceElementForcesAndSourcesCore > > elasticBcLhs
MoFEMErrorCode getBc(boost::shared_ptr< BC > &bc_vec_ptr, const std::string block_set_name, const int nb_attributes)
boost::shared_ptr< TractionBcVec > bcData
FeTractionBc(MoFEM::Interface &m_field, const std::string field_name, boost::shared_ptr< TractionBcVec > &bc)
virtual MoFEMErrorCode assemble(EntData &data)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
virtual MoFEMErrorCode integrate(int row_side, EntityType row_type, EntData &data)
virtual MoFEMErrorCode integrate(EntData &data)
virtual MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data)
MatrixDouble K
local tangent matrix
VectorDouble nF
local right hand side vector
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
virtual MoFEMErrorCode assemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpAssembleBasic(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type)
OpAssembleBasic(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type, const bool assemble_symmetry)
OpAssembleFace(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type, const bool assemble_symmetry)
OpAssembleFace(const std::string &field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type)
OpAssembleVolume(const std::string &field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type)
OpAssembleVolume(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > data_ptr, const char type, const bool assemble_symmetry)
OpCalculateRotationAndSpatialGradient(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
OpCalculateStrainEnergy(const std::string field_name, boost::shared_ptr< DataAtIntegrationPts > data_at_pts, boost::shared_ptr< double > &e)
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)
boost::shared_ptr< BcDispVec > bcDispPtr
boost::shared_ptr< BcDispVec > bcDispPtr
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< PhysicalEquations > physicsPtr
material physical equations
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)
virtual MoFEMErrorCode evaluateRhs(EntData &data)=0
virtual MoFEMErrorCode evaluateLhs(EntData &data)=0
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
MoFEMErrorCode doWork(int side, EntityType type, DataForcesAndSourcesCore::EntData &data)
moab::Interface & postProcMesh
OpPostProcDataStructure(moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, const std::string field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
std::vector< EntityHandle > & mapGaussPts
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
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 &row_data, EntData &col_data)
boost::shared_ptr< BcRotVec > bcRotPtr
MoFEMErrorCode integrate(EntData &data)
OpRotationBc(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcRotVec > &bc_rot_ptr)
boost::shared_ptr< BcRotVec > bcRotPtr
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialConsistency_dBubble_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialConsistency_dP_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpSpatialConsistencyBubble(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpSpatialConsistencyDivTerm(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpSpatialConsistencyP(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
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)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialEquilibrium_dw_dw(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha, const double rho)
OpSpatialEquilibrium(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha, const double rho)
OpSpatialPhysical_du_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)
OpSpatialPhysical_du_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off=false)
OpSpatialPhysical_du_domega(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)
OpSpatialPhysical_du_du(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
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)
OpSpatialPhysical(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double alpha)
MoFEMErrorCode integrate(EntData &data)
OpSpatialPreconditionMass(const std::string &row_field, MatrixPtr m_ptr)
MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialPrj_dx_dw(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
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)
MoFEMErrorCode integrate(EntData &row_data)
OpSpatialPrj(const std::string &row_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
OpSpatialRotation_domega_dBubble(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
OpSpatialRotation_domega_dP(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const bool assemble_off)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
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 integrate(EntData &row_data, EntData &col_data)
MoFEMErrorCode integrate(EntData &data)
OpSpatialRotation(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode assemble(int row_side, EntityType row_type, EntData &data)
MoFEMErrorCode integrate(EntData &row_data)
OpSpatialSchurBegin(const std::string &row_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode assemble(int side, EntityType type, EntData &data)
map< std::string, MatrixDouble > invBlockMat
MoFEMErrorCode integrate(EntData &row_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpSpatialSchurEnd(const std::string &row_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const double eps)
map< std::string, DataAtIntegrationPts::BlockMatContainor > blockMat
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpTractionBc(std::string row_field, FeTractionBc &bc_fe)
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > DTensor3Ptr
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, 3, 3 > DTensor2Ptr
FTensor::Tensor2< adouble, 3, 3 > ATensor2
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
static DTensor3Ptr get_vecTensor3(std::vector< double > &v)
PhysicalEquations(const int size_active, const int size_dependent)
static DTensor2Ptr get_VecTensor2(std::vector< double > &v)
FTensor::Tensor1< adouble, 3 > ATensor1
FTensor::Tensor3< adouble, 3, 3, 3 > ATensor3
virtual MoFEMErrorCode recordTape(const int tag, DTensor2Ptr *t_h)=0
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > DTensor0Ptr
static DTensor0Ptr get_VecTensor0(std::vector< double > &v)
TractionBc(std::string name, std::vector< double > &attr, Range &faces)
virtual moab::Interface & get_moab()=0
Data on single entity (This is passed as argument to DataOperator::doWork)
const VectorInt & getIndices() const
Get global indices of dofs on entity.
Deprecated interface functions.
base class for all interface classes