v0.8.19
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 
126 
135 
156 
163 
164 };
165 
166 struct OpJacobian;
167 
169 
174 
177 
178  PhysicalEquations() = delete;
179  PhysicalEquations(const int size_active, const int size_dependent)
180  : activeVariables(size_active, 0), dependentVariables(size_dependent, 0),
181  dependentVariablesDirevatives(size_dependent * size_active, 0) {}
182  virtual ~PhysicalEquations() {}
183 
184  virtual MoFEMErrorCode recordTape(const int tag) = 0;
185 
186  virtual OpJacobian *
187  returnOpJacobian(const std::string &field_name, const int tag,
188  const bool eval_rhs, const bool eval_lhs,
189  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
190  boost::shared_ptr<PhysicalEquations> &physics_ptr) = 0;
191 
192  std::vector<double> activeVariables;
193  std::vector<double> dependentVariables;
194  std::vector<double> dependentVariablesDirevatives;
195 
196  /** \name Active variables */
197 
198  /**@{*/
199 
200  template <int S>
201  inline static DTensor2Ptr get_VecTensor2(std::vector<double> &v) {
202  return DTensor2Ptr(&v[S + 0], &v[S + 1], &v[S + 2], &v[S + 3], &v[S + 4],
203  &v[S + 5], &v[S + 6], &v[S + 7], &v[S + 8]);
204  }
205 
206  template <int S>
207  inline static DTensor0Ptr get_VecTensor0(std::vector<double> &v) {
208  return DTensor0Ptr(&v[S + 0]);
209  }
210 
211  template <int S0, int S1, int S2>
212  inline static DTensor3Ptr get_vecTensor3(std::vector<double> &v) {
213  constexpr int A00 = 18 * S0 + 18 * 0 + 9 * S1 + 3 * S2;
214  constexpr int A01 = 18 * S0 + 18 * 1 + 9 * S1 + 3 * S2;
215  constexpr int A02 = 18 * S0 + 18 * 2 + 9 * S1 + 3 * S2;
216  constexpr int A10 = 18 * S0 + 18 * 3 + 9 * S1 + 3 * S2;
217  constexpr int A11 = 18 * S0 + 18 * 4 + 9 * S1 + 3 * S2;
218  constexpr int A12 = 18 * S0 + 18 * 5 + 9 * S1 + 3 * S2;
219  constexpr int A20 = 18 * S0 + 18 * 6 + 9 * S1 + 3 * S2;
220  constexpr int A21 = 18 * S0 + 18 * 7 + 9 * S1 + 3 * S2;
221  constexpr int A22 = 18 * S0 + 18 * 8 + 9 * S1 + 3 * S2;
222  return DTensor3Ptr(
223 
224  &v[A00 + 0], &v[A00 + 1], &v[A00 + 2], &v[A01 + 0], &v[A01 + 1],
225  &v[A01 + 2], &v[A02 + 0], &v[A02 + 1], &v[A02 + 2],
226 
227  &v[A10 + 0], &v[A10 + 1], &v[A10 + 2], &v[A11 + 0], &v[A11 + 1],
228  &v[A11 + 2], &v[A12 + 0], &v[A12 + 1], &v[A12 + 2],
229 
230  &v[A20 + 0], &v[A20 + 1], &v[A20 + 2], &v[A21 + 0], &v[A21 + 1],
231  &v[A21 + 2], &v[A22 + 0], &v[A22 + 1], &v[A22 + 2]
232 
233  );
234  }
235 
236  /**@}*/
237 
238  /** \name Active variables */
239 
240  /**@{*/
241 
242  inline DTensor2Ptr get_h() { return get_VecTensor2<0>(activeVariables); }
243  inline DTensor2Ptr get_H() { return get_VecTensor2<9>(activeVariables); }
244 
245  /**@}*/
246 
247  /** \name Dependent variables */
248 
249  /**@{*/
250 
251  inline DTensor2Ptr get_P() { return get_VecTensor2<0>(dependentVariables); }
253  return get_VecTensor2<9>(dependentVariables);
254  }
255  inline DTensor0Ptr get_Phi() {
256  return get_VecTensor0<18>(dependentVariables);
257  }
258  inline double &get_PhiRef() { return dependentVariables[18]; }
260  return get_VecTensor2<9 + 9 + 1>(dependentVariables);
261  }
262 
263  /**@}*/
264 
265  /** \name Derivatives of dependent variables */
266 
267  /**@{*/
268 
270  return get_vecTensor3<0, 0, 0>(dependentVariablesDirevatives);
271  }
273  return get_vecTensor3<0, 1, 0>(dependentVariablesDirevatives);
274  }
276  return get_vecTensor3<0, 0, 1>(dependentVariablesDirevatives);
277  }
279  return get_vecTensor3<0, 1, 1>(dependentVariablesDirevatives);
280  }
282  return get_vecTensor3<0, 0, 2>(dependentVariablesDirevatives);
283  }
285  return get_vecTensor3<0, 1, 2>(dependentVariablesDirevatives);
286  }
287 
289  return get_vecTensor3<9, 0, 0>(dependentVariablesDirevatives);
290  }
292  return get_vecTensor3<9, 1, 0>(dependentVariablesDirevatives);
293  }
295  return get_vecTensor3<9, 0, 1>(dependentVariablesDirevatives);
296  }
298  return get_vecTensor3<9, 1, 1>(dependentVariablesDirevatives);
299  }
301  return get_vecTensor3<9, 0, 2>(dependentVariablesDirevatives);
302  }
304  return get_vecTensor3<9, 1, 2>(dependentVariablesDirevatives);
305  }
306 
308  return get_VecTensor2<18 * (9 + 9) + 0>(dependentVariablesDirevatives);
309  }
311  return get_VecTensor2<18 * (9 + 9) + 9>(dependentVariablesDirevatives);
312  }
313 
315  return get_vecTensor3<9 + 9 + 1, 0, 0>(dependentVariablesDirevatives);
316  }
318  return get_vecTensor3<9 + 9 + 1, 1, 0>(dependentVariablesDirevatives);
319  }
321  return get_vecTensor3<9 + 9 + 1, 0, 1>(dependentVariablesDirevatives);
322  }
324  return get_vecTensor3<9 + 9 + 1, 1, 1>(dependentVariablesDirevatives);
325  }
327  return get_vecTensor3<9 + 9 + 1, 0, 2>(dependentVariablesDirevatives);
328  }
330  return get_vecTensor3<9 + 9 + 1, 1, 2>(dependentVariablesDirevatives);
331  }
332 
333  /**@}*/
334 };
335 
336 struct BcDisp {
337  BcDisp(std::vector<double> &attr, Range &faces);
338  Range faces;
341 };
342 typedef std::vector<BcDisp> BcDispVec;
343 
344 struct BcRot {
345  BcRot(std::vector<double> &attr, Range &faces);
346  Range faces;
348  double theta;
349 };
350 typedef std::vector<BcRot> BcRotVec;
351 
352 typedef std::vector<Range> TractionFreeBc;
353 
355  const int tAg; ///< adol-c tape
356  const bool evalRhs;
357  const bool evalLhs;
358  boost::shared_ptr<DataAtIntegrationPts>
359  dataAtPts; ///< data at integration pts
360  boost::shared_ptr<PhysicalEquations>
361  physicsPtr; ///< material physical equations
362 
363  OpJacobian(const std::string &field_name, const int tag, const bool eval_rhs,
364  const bool eval_lhs,
365  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
366  boost::shared_ptr<PhysicalEquations> &physics_ptr)
367  : ForcesAndSourcesCore::UserDataOperator(field_name, OPROW), tAg(tag),
368  evalRhs(eval_rhs), evalLhs(eval_lhs), dataAtPts(data_ptr),
369  physicsPtr(physics_ptr) {}
370 
371  virtual MoFEMErrorCode evaluateRhs(EntData &data) = 0;
372  virtual MoFEMErrorCode evaluateLhs(EntData &data) = 0;
373 
374  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
375 };
376 
377 template <typename T> struct OpAssembleBasic : public T {
378 
379  boost::shared_ptr<DataAtIntegrationPts>
380  dataAtPts; ///< data at integration pts
381  const bool assembleSymmetry;
382 
383  OpAssembleBasic(const std::string &field_name,
384  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
385  const char type)
386  : T(field_name, type), dataAtPts(data_ptr), assembleSymmetry(false) {}
387 
388  OpAssembleBasic(const std::string &row_field, const std::string &col_field,
389  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
390  const char type, const bool assemble_symmetry = false)
391  : T(row_field, col_field, type, false), dataAtPts(data_ptr),
392  assembleSymmetry(assemble_symmetry) {}
393 
394  VectorDouble nF; ///< local right hand side vector
395  MatrixDouble K; ///< local tangent matrix
397 
400  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
402  }
403  virtual MoFEMErrorCode assemble(EntData &data) {
405  double *vec_ptr = &*nF.begin();
406  int nb_dofs = data.getIndices().size();
407  int *ind_ptr = &*data.getIndices().begin();
408  CHKERR VecSetValues(T::getFEMethod()->ts_F, nb_dofs, ind_ptr, vec_ptr,
409  ADD_VALUES);
411  }
412 
413  virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) {
415  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
417  }
418 
419  virtual MoFEMErrorCode assemble(EntData &row_data, EntData &col_data) {
421  double *mat_ptr = &*K.data().begin();
422  int *row_ind_ptr = &*row_data.getIndices().begin();
423  int *col_ind_ptr = &*col_data.getIndices().begin();
424  int row_nb_dofs = row_data.getIndices().size();
425  int col_nb_dofs = col_data.getIndices().size();
426  CHKERR MatSetValues(T::getFEMethod()->ts_B, row_nb_dofs, row_ind_ptr,
427  col_nb_dofs, col_ind_ptr, mat_ptr, ADD_VALUES);
428  if (assembleSymmetry) {
429  transposeK.resize(col_nb_dofs, row_nb_dofs);
430  noalias(transposeK) = trans(K);
431  double *mat_ptr = &*transposeK.data().begin();
432  CHKERR MatSetValues(T::getFEMethod()->ts_B, col_nb_dofs, col_ind_ptr,
433  row_nb_dofs, row_ind_ptr, mat_ptr, ADD_VALUES);
434  }
436  }
437 
438  MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
440  if (data.getIndices().empty())
442  nF.resize(data.getIndices().size(), false);
443  nF.clear();
444  CHKERR integrate(data);
445  CHKERR assemble(data);
447  }
448 
449  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
450  EntityType col_type, EntData &row_data,
451  EntData &col_data) {
453  if (row_data.getIndices().empty())
455  if (col_data.getIndices().empty())
457  K.resize(row_data.getIndices().size(), col_data.getIndices().size(), false);
458  K.clear();
459  CHKERR integrate(row_data, col_data);
460  CHKERR assemble(row_data, col_data);
462  }
463 };
464 
465 struct OpAssembleVolume : public OpAssembleBasic<VolUserDataOperator> {
466  OpAssembleVolume(const std::string &field,
467  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
468  const char type)
469  : OpAssembleBasic<VolUserDataOperator>(field, data_ptr, type) {}
470 
471  OpAssembleVolume(const std::string &row_field, const std::string &col_field,
472  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
473  const char type, const bool assemble_symmetry = false)
474  : OpAssembleBasic<VolUserDataOperator>(row_field, col_field, data_ptr,
475  type, assemble_symmetry) {}
476 };
477 
478 struct OpAssembleFace : public OpAssembleBasic<FaceUserDataOperator> {
479  OpAssembleFace(const std::string &field,
480  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
481  const char type)
482  : OpAssembleBasic<FaceUserDataOperator>(field, data_ptr, type) {}
483 
484  OpAssembleFace(const std::string &row_field, const std::string &col_field,
485  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
486  const char type)
487  : OpAssembleBasic<FaceUserDataOperator>(row_field, col_field, data_ptr,
488  type) {}
489 };
490 
492  boost::shared_ptr<DataAtIntegrationPts>
493  dataAtPts; ///< data at integration pts
495  const std::string &field_name,
496  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
497  : VolUserDataOperator(field_name, OPROW), dataAtPts(data_ptr) {}
498  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
499 };
500 
502  OpSpatialEquilibrium(const std::string &field_name,
503  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
504  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
506 };
507 
509  OpSpatialRotation(const std::string &field_name,
510  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
511  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
513 };
514 
516  OpSpatialPhysical(const std::string &field_name,
517  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
518  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
520 };
521 
523  OpSpatialConsistencyP(const std::string &field_name,
524  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
525  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
527 };
528 
530  OpSpatialConsistencyBubble(const std::string &field_name,
531  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
532  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
534 };
535 
537  OpSpatialConsistencyDivTerm(const std::string &field_name,
538  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
539  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
541 };
542 
543 struct OpDispBc : public OpAssembleFace {
544  boost::shared_ptr<BcDispVec> bcDispPtr;
545  OpDispBc(const std::string &field_name,
546  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
547  boost::shared_ptr<BcDispVec> &bc_disp_ptr)
548  : OpAssembleFace(field_name, data_ptr, OPROW), bcDispPtr(bc_disp_ptr) {}
550 };
551 
552 struct OpDispBc_dx : public OpAssembleFace {
553  boost::shared_ptr<BcDispVec> bcDispPtr;
554  OpDispBc_dx(const std::string &row_field_name,
555  const std::string &col_field_name,
556  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
557  boost::shared_ptr<BcDispVec> &bc_disp_ptr)
558  : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL),
559  bcDispPtr(bc_disp_ptr) {}
560  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
561 };
562 
563 struct OpRotationBc : public OpAssembleFace {
564  boost::shared_ptr<BcRotVec> bcRotPtr;
565  OpRotationBc(const std::string &field_name,
566  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
567  boost::shared_ptr<BcRotVec> &bc_rot_ptr)
568  : OpAssembleFace(field_name, data_ptr, OPROW), bcRotPtr(bc_rot_ptr) {}
570 };
571 
573  boost::shared_ptr<BcRotVec> bcRotPtr;
574  OpRotationBc_dx(const std::string &row_field_name,
575  const std::string &col_field_name,
576  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
577  boost::shared_ptr<BcRotVec> &bc_rot_ptr)
578  : OpAssembleFace(row_field_name, col_field_name, data_ptr, OPROWCOL),
579  bcRotPtr(bc_rot_ptr) {}
580  MoFEMErrorCode integrate(EntData &row_data,EntData &col_data);
581 };
582 
583 struct FeTractionFreeBc : public BasicMethod {
584  std::string fieldName;
585  boost::shared_ptr<TractionFreeBc> bcData;
587  FeTractionFreeBc(const std::string field_name,
588  boost::shared_ptr<TractionFreeBc> &bc_data)
589  : BasicMethod(), fieldName(field_name), bcData(bc_data) {}
592 };
593 
595  OpSpatialEquilibrium_dw_dP(const std::string &row_field,
596  const std::string &col_field,
597  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
598  const bool assemble_off = false)
599  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
600  assemble_off) {
601  sYmm = false;
602  }
603  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
604 };
605 
607  OpSpatialPhysical_du_du(const std::string &row_field,
608  const std::string &col_field,
609  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
610  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL) {
611  sYmm = false;
612  }
613  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
614 };
615 
617  OpSpatialPhysical_du_dP(const std::string &row_field,
618  const std::string &col_field,
619  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
620  const bool assemble_off = false)
621  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
622  assemble_off) {
623  sYmm = false;
624  }
625  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
626 };
627 
630  const std::string &row_field, const std::string &col_field,
631  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
632  const bool assemble_off = false)
633  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
634  assemble_off) {
635  sYmm = false;
636  }
637  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
638 };
639 
641  OpSpatialPhysical_du_domega(const std::string &row_field,
642  const std::string &col_field,
643  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
644  const bool assemble_off = false)
645  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
646  assemble_off) {
647  sYmm = false;
648  }
649  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
650 };
651 
653  OpSpatialPhysical_du_dx(const std::string &row_field,
654  const std::string &col_field,
655  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
656  const bool assemble_off = false)
657  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
658  assemble_off) {
659  sYmm = false;
660  }
661  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
662 };
663 
665  OpSpatialRotation_domega_dP(const std::string &row_field,
666  const std::string &col_field,
667  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
668  const bool assemble_off = false)
669  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
670  assemble_off) {
671  sYmm = false;
672  }
673  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
674 };
675 
678  const std::string &row_field, const std::string &col_field,
679  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
680  const bool assemble_off = false)
681  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
682  assemble_off) {
683  sYmm = false;
684  }
685  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
686 };
687 
690  const std::string &row_field, const std::string &col_field,
691  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
692  const bool assemble_off = false)
693  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
694  assemble_off) {
695  sYmm = false;
696  }
697  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
698 };
699 
702  const std::string &row_field, const std::string &col_field,
703  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
704  const bool assemble_off = false)
705  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
706  assemble_off) {
707  sYmm = false;
708  }
709  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
710 };
711 
714  const std::string &row_field, const std::string &col_field,
715  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
716  const bool assemble_off = false)
717  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
718  assemble_off) {
719  sYmm = false;
720  }
721  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
722 };
723 
726  const std::string &row_field, const std::string &col_field,
727  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
728  const bool assemble_off = false)
729  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
730  assemble_off) {
731  sYmm = false;
732  }
733  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
734 };
735 
738  const std::string &row_field, const std::string &col_field,
739  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
740  const bool assemble_off = false)
741  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
742  assemble_off) {
743  sYmm = false;
744  }
745  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
746 };
747 
750  const std::string &row_field, const std::string &col_field,
751  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
752  const bool assemble_off = false)
753  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
754  assemble_off) {
755  sYmm = false;
756  }
757  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
758 };
759 
761 
762  moab::Interface &postProcMesh;
763  std::vector<EntityHandle> &mapGaussPts;
764  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
765 
766  OpPostProcDataStructure(moab::Interface &post_proc_mesh,
767  std::vector<EntityHandle> &map_gauss_pts,
768  const std::string field_name,
769  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
770  : VolUserDataOperator(field_name,
771  ForcesAndSourcesCore::UserDataOperator::OPROW),
772  postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
773  dataAtPts(data_ptr) {}
774 
775  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
776 };
777 
780  const std::string &row_field,
781  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
782  : OpAssembleVolume(row_field, row_field, data_ptr, OPROW) {
783  }
784  MoFEMErrorCode integrate(EntData &row_data);
785 };
786 
789  const std::string &row_field, const std::string &col_field,
790  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
791  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
792  sYmm = false;
793  }
794  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
795 };
796 
799  const std::string &row_field, const std::string &col_field,
800  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
801  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
802  sYmm = false;
803  }
804  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
805 };
806 
809  const std::string &row_field, const std::string &col_field,
810  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
811  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
812  sYmm = false;
813  }
814  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
815 };
816 
819  const std::string &row_field, const std::string &col_field,
820  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
821  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL, false) {
822  sYmm = false;
823  }
824  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
825 };
826 
828 
829  /**
830  * \brief Getting interface of core database
831  * @param uuid unique ID of interface
832  * @param iface returned pointer to interface
833  * @return error code
834  */
835  MoFEMErrorCode query_interface(const MOFEMuuid &uuid,
836  UnknownInterface **iface) const;
837 
839 
840  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
841  boost::shared_ptr<PhysicalEquations> physicalEquations;
842 
843  boost::shared_ptr<ForcesAndSourcesCore> elasticFeRhs;
844  boost::shared_ptr<ForcesAndSourcesCore> elasticFeLhs;
845  boost::shared_ptr<ForcesAndSourcesCore> materialFeRhs;
846  boost::shared_ptr<ForcesAndSourcesCore> materialFeLhs;
847  boost::shared_ptr<ForcesAndSourcesCore> epFeRhs;
848  boost::shared_ptr<ForcesAndSourcesCore> epFeLhs;
849  boost::shared_ptr<ForcesAndSourcesCore> elasticBcLhs;
850  boost::shared_ptr<ForcesAndSourcesCore> elasticBcRhs;
851  boost::shared_ptr<ForcesAndSourcesCore> materialBcRhs;
852  boost::shared_ptr<ForcesAndSourcesCore> materialBcLhs;
853  boost::shared_ptr<ForcesAndSourcesCore> epBcRhs;
854  boost::shared_ptr<ForcesAndSourcesCore> epBcLhs;
855  boost::shared_ptr<FEMethod> tractionFreeBc;
856 
857  DM dM; ///< Coupled problem all fields
858  DM dmElastic; ///< Elastic problem
859  DM dmMaterial; ///< Material problem
860  DM dmPrjSpatial; ///< Project spatial displacement
861 
862  const std::string piolaStress;
863  const std::string eshelbyStress;
864  const std::string spatialDisp;
865  const std::string spatialPosition;
866  const std::string materialDisp;
867  const std::string streachGradient;
868  const std::string rotAxis;
869  const std::string materialGradient;
870  const std::string tauField;
871  const std::string lambdaField;
872  const std::string bubbleField;
873 
874  const std::string elementVolumeName;
875  const std::string elementFaceName;
876 
878  virtual ~EshelbianCore();
879 
882  double gamma;
883  double loadFactor;
885 
887 
888  boost::shared_ptr<MethodForForceScaling> scaleSpatialBc;
889  boost::shared_ptr<MethodForForceScaling> scaleMaterialBc;
890 
891  boost::shared_ptr<BcDispVec> bcSpatialDispVecPtr;
892  boost::shared_ptr<BcRotVec> bcSpatialRotationVecPtr;
893  boost::shared_ptr<BcDispVec> bcMaterialDispVecPtr;
894  boost::shared_ptr<TractionFreeBc> tractionSpatialFreeBc;
895  boost::shared_ptr<TractionFreeBc> tractionMaterialFreeBc;
896 
897  template <typename BC>
898  MoFEMErrorCode getDisplacementsBc(boost::shared_ptr<BC> &bc_vec_ptr,
899  const std::string block_set_name,
900  const int nb_attributes) {
903  if (it->getName().compare(0, block_set_name.length(), block_set_name) ==
904  0) {
905  std::vector<double> block_attributes;
906  CHKERR it->getAttributes(block_attributes);
907  if (block_attributes.size() != nb_attributes) {
908  SETERRQ3(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
909  "In block %s expected %d attributes, but given %d",
910  it->getName().c_str(), nb_attributes,
911  block_attributes.size());
912  }
913  Range faces;
914  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
915  true);
916  bc_vec_ptr->emplace_back(block_attributes, faces);
917  }
918  }
920  }
921 
923  bcSpatialDispVecPtr = boost::shared_ptr<BcDispVec>(new BcDispVec());
924  return getDisplacementsBc(bcSpatialDispVecPtr, "SPATIAL_DISP_BC", 6);
925  }
926 
928  bcSpatialRotationVecPtr = boost::shared_ptr<BcRotVec>(new BcRotVec());
929  return getDisplacementsBc(bcSpatialRotationVecPtr, "SPATIAL_ROTATION_BC",
930  4);
931  }
932 
934  bcMaterialDispVecPtr = boost::shared_ptr<BcDispVec>(new BcDispVec());
935  return getDisplacementsBc(bcMaterialDispVecPtr, "MATERIAL_DISP_BC", 6);
936  }
937 
939  boost::shared_ptr<TractionFreeBc> &bc_ptr,
940  const std::string disp_block_set_name,
941  const std::string rot_block_set_name);
942  inline MoFEMErrorCode
945  boost::shared_ptr<TractionFreeBc>(new TractionFreeBc());
946  return getTractionFreeBc(meshset, tractionSpatialFreeBc, "SPATIAL_DISP_BC",
947  "SPATIAL_ROTATION_BC");
948  }
949  inline MoFEMErrorCode
952  boost::shared_ptr<TractionFreeBc>(new TractionFreeBc());
954  "MATERIAL_DISP_BC", "MATERIAL_ROTATION_BC");
955  }
956 
957  MoFEMErrorCode addFields(const EntityHandle meshset = 0);
960  MoFEMErrorCode addDMs(const BitRefLevel bit = BitRefLevel().set(0));
961 
963  const double lambda,
964  const double mu,
965  const double sigma_y);
966 
967  MoFEMErrorCode addMaterial_HMHMooneyRivlin(const int tape, const double alpha,
968  const double beta,
969  const double lambda,
970  const double sigma_y);
971 
973  setBaseVolumeElementOps(const int tag, const bool do_rhs, const bool do_lhs,
974  boost::shared_ptr<ForcesAndSourcesCore> &fe);
975 
977  setGenericVolumeElementOps(const int tag, const bool add_elastic,
978  const bool add_material,
979  boost::shared_ptr<ForcesAndSourcesCore> &fe_rhs,
980  boost::shared_ptr<ForcesAndSourcesCore> &fe_lhs);
981 
983  setGenericFaceElementOps(const bool add_elastic, const bool add_material,
984  boost::shared_ptr<ForcesAndSourcesCore> &fe_rhs,
985  boost::shared_ptr<ForcesAndSourcesCore> &fe_lhs);
986 
987  MoFEMErrorCode setElasticElementOps(const int tag);
989  MoFEMErrorCode setMaterialElementOps(const int tag);
990 
991  MoFEMErrorCode setUpTSElastic(TS ts, Mat m, Vec f);
992  MoFEMErrorCode solveElastic(TS ts, Vec x);
993 
994  MoFEMErrorCode postProcessResults(const int tag, const std::string file);
995 };
996 
997 } // namespace EshelbianPlasticity
998 
999 #endif //__ESHELBIAN_PLASTICITY_HPP__
MoFEMErrorCode setMaterialElementOps(const int tag)
boost::shared_ptr< PhysicalEquations > physicsPtr
material physical equations
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)
OpAssembleBasic(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type, const bool assemble_symmetry=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)
boost::shared_ptr< ForcesAndSourcesCore > materialBcRhs
OpAssembleVolume(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type, const bool assemble_symmetry=false)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Common.hpp:60
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)
boost::shared_ptr< TractionFreeBc > bcData
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)
DM dmPrjSpatial
Project spatial displacement.
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)
boost::shared_ptr< BcRotVec > bcRotPtr
boost::shared_ptr< MethodForForceScaling > scaleMaterialBc
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)
boost::shared_ptr< ForcesAndSourcesCore > elasticFeRhs
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)
VectorBoundedArray< int, 3 > VectorInt3
Definition: Common.hpp:218
base class for all interface classes
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:475
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)
std::vector< BcRot > BcRotVec
boost::shared_ptr< ForcesAndSourcesCore > epBcLhs
MoFEMErrorCode integrate(EntData &data)
virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return() ...
Definition: definitions.h:506
boost::shared_ptr< MethodForForceScaling > scaleSpatialBc
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)
MoFEMErrorCode integrate(EntData &data)
MoFEMErrorCode addFields(const EntityHandle meshset=0)
BcRot(std::vector< double > &attr, Range &faces)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
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)
boost::shared_ptr< BcRotVec > bcRotPtr
std::vector< Range > TractionFreeBc
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)
std::bitset< BITREFLEVEL_SIZE > BitRefLevel
Bit structure attached to each entity identifying to what mesh entity is attached.
Definition: Common.hpp:147
boost::shared_ptr< VectorDouble > VectorPtr
OpAssembleFace(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type)
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
boost::shared_ptr< ForcesAndSourcesCore > epFeLhs
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)
virtual MoFEMErrorCode evaluateLhs(EntData &data)=0
MatrixDouble K
local tangent matrix
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
boost::shared_ptr< ForcesAndSourcesCore > elasticFeLhs
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
MoFEMErrorCode getMaterialTractionFreeBc(const EntityHandle meshset=0)
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 getDisplacementsBc(boost::shared_ptr< BC > &bc_vec_ptr, const std::string block_set_name, const int nb_attributes)
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< ForcesAndSourcesCore > &fe_rhs, boost::shared_ptr< ForcesAndSourcesCore > &fe_lhs)
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Common.hpp:212
boost::shared_ptr< TractionFreeBc > tractionSpatialFreeBc
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)
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Common.hpp:211
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
MoFEMErrorCode addMaterial_HMHHStVenantKirchhoff(const int tape, const double lambda, const double mu, const double sigma_y)
MoFEMErrorCode setGenericFaceElementOps(const bool add_elastic, const bool add_material, boost::shared_ptr< ForcesAndSourcesCore > &fe_rhs, boost::shared_ptr< ForcesAndSourcesCore > &fe_lhs)
boost::shared_ptr< ForcesAndSourcesCore > materialBcLhs
OpSpatialPhysical_du_du(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
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)
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)
static DTensor0Ptr get_VecTensor0(std::vector< double > &v)
virtual MoFEMErrorCode assemble(EntData &row_data, EntData &col_data)
boost::shared_ptr< FEMethod > tractionFreeBc
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
OpPostProcDataStructure(moab::Interface &post_proc_mesh, std::vector< EntityHandle > &map_gauss_pts, const std::string field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
#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...
#define CHKERR
Inline error check.
Definition: definitions.h:594
boost::shared_ptr< ForcesAndSourcesCore > epBcRhs
MoFEMErrorCode getSpatialTractionFreeBc(const EntityHandle meshset=0)
OpSpatialPrj_dx_du(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
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)
boost::shared_ptr< BcDispVec > bcMaterialDispVecPtr
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
FTensor::Tensor3< FTensor::PackPtr< double *, 1 >, 3, 3, 3 > DTensor3Ptr
const VectorInt & getIndices() const
Get global indices of dofs on entity.
boost::shared_ptr< ForcesAndSourcesCore > elasticBcLhs
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)
std::bitset< BITINTERFACEUID_SIZE > BitIntefaceId
Definition: Common.hpp:152
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)
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)
VolumeElementForcesAndSourcesCore::UserDataOperator VolUserDataOperator
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
boost::shared_ptr< ForcesAndSourcesCore > materialFeRhs
virtual MoFEMErrorCode recordTape(const int tag)=0
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
boost::shared_ptr< ForcesAndSourcesCore > epFeRhs
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
boost::shared_ptr< ForcesAndSourcesCore > elasticBcRhs
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
VectorBoundedArray< double, 3 > VectorDouble3
Definition: Common.hpp:220
MoFEMErrorCode setBaseVolumeElementOps(const int tag, const bool do_rhs, const bool do_lhs, boost::shared_ptr< ForcesAndSourcesCore > &fe)
boost::shared_ptr< TractionFreeBc > tractionMaterialFreeBc
boost::shared_ptr< ForcesAndSourcesCore > materialFeLhs
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)
virtual MoFEMErrorCode integrate(EntData &data)
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, 3, 3 > DTensor2Ptr
BcDisp(std::vector< double > &attr, Range &faces)
OpSpatialPrj_dx_domega(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
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
FeTractionFreeBc(const std::string field_name, boost::shared_ptr< TractionFreeBc > &bc_data)