v0.8.16
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 
124 
132 
153 
154 };
155 
156 struct OpJacobian;
157 
159 
163 
166 
167  PhysicalEquations() = delete;
168  PhysicalEquations(const int size_active, const int size_dependent)
169  : activeVariables(size_active, 0), dependentVariables(size_dependent, 0),
170  dependentVariablesDirevatives(size_dependent * size_active, 0) {}
171  virtual ~PhysicalEquations() {}
172 
173  virtual MoFEMErrorCode recordTape(const int tag) = 0;
174 
175  virtual OpJacobian *
176  returnOpJacobian(const std::string &field_name, const int tag,
177  const bool eval_rhs, const bool eval_lhs,
178  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
179  boost::shared_ptr<PhysicalEquations> &physics_ptr) = 0;
180 
181  std::vector<double> activeVariables;
182  std::vector<double> dependentVariables;
183  std::vector<double> dependentVariablesDirevatives;
184 
185  /** \name Active variables */
186 
187  /**@{*/
188 
189  template <int S>
190  inline static DTensor2Ptr get_VecTensor2(std::vector<double> &v) {
191  return DTensor2Ptr(&v[S + 0], &v[S + 1], &v[S + 2], &v[S + 3], &v[S + 4],
192  &v[S + 5], &v[S + 6], &v[S + 7], &v[S + 8]);
193  }
194 
195  template <int S>
196  inline static DTensor0Ptr get_VecTensor0(std::vector<double> &v) {
197  return DTensor0Ptr(&v[S + 0]);
198  }
199 
200  template <int S0, int S1, int S2>
201  inline static DTensor3Ptr get_vecTensor3(std::vector<double> &v) {
202  constexpr int A00 = 18 * S0 + 18 * 0 + 9 * S1 + 3 * S2;
203  constexpr int A01 = 18 * S0 + 18 * 1 + 9 * S1 + 3 * S2;
204  constexpr int A02 = 18 * S0 + 18 * 2 + 9 * S1 + 3 * S2;
205  constexpr int A10 = 18 * S0 + 18 * 3 + 9 * S1 + 3 * S2;
206  constexpr int A11 = 18 * S0 + 18 * 4 + 9 * S1 + 3 * S2;
207  constexpr int A12 = 18 * S0 + 18 * 5 + 9 * S1 + 3 * S2;
208  constexpr int A20 = 18 * S0 + 18 * 6 + 9 * S1 + 3 * S2;
209  constexpr int A21 = 18 * S0 + 18 * 7 + 9 * S1 + 3 * S2;
210  constexpr int A22 = 18 * S0 + 18 * 8 + 9 * S1 + 3 * S2;
211  return DTensor3Ptr(
212 
213  &v[A00 + 0], &v[A00 + 1], &v[A00 + 2], &v[A01 + 0], &v[A01 + 1],
214  &v[A01 + 2], &v[A02 + 0], &v[A02 + 1], &v[A02 + 2],
215 
216  &v[A10 + 0], &v[A10 + 1], &v[A10 + 2], &v[A11 + 0], &v[A11 + 1],
217  &v[A11 + 2], &v[A12 + 0], &v[A12 + 1], &v[A12 + 2],
218 
219  &v[A20 + 0], &v[A20 + 1], &v[A20 + 2], &v[A21 + 0], &v[A21 + 1],
220  &v[A21 + 2], &v[A22 + 0], &v[A22 + 1], &v[A22 + 2]
221 
222  );
223  }
224 
225  /**@}*/
226 
227  /** \name Active variables */
228 
229  /**@{*/
230 
231  inline DTensor2Ptr get_h() { return get_VecTensor2<0>(activeVariables); }
232  inline DTensor2Ptr get_H() { return get_VecTensor2<9>(activeVariables); }
233 
234  /**@}*/
235 
236  /** \name Dependent variables */
237 
238  /**@{*/
239 
240  inline DTensor2Ptr get_P() { return get_VecTensor2<0>(dependentVariables); }
242  return get_VecTensor2<9>(dependentVariables);
243  }
244  inline DTensor0Ptr get_Phi() {
245  return get_VecTensor0<18>(dependentVariables);
246  }
247  inline double &get_PhiRef() { return dependentVariables[18]; }
249  return get_VecTensor2<9 + 9 + 1>(dependentVariables);
250  }
251 
252  /**@}*/
253 
254  /** \name Derivatives of dependent variables */
255 
256  /**@{*/
257 
259  return get_vecTensor3<0, 0, 0>(dependentVariablesDirevatives);
260  }
262  return get_vecTensor3<0, 1, 0>(dependentVariablesDirevatives);
263  }
265  return get_vecTensor3<0, 0, 1>(dependentVariablesDirevatives);
266  }
268  return get_vecTensor3<0, 1, 1>(dependentVariablesDirevatives);
269  }
271  return get_vecTensor3<0, 0, 2>(dependentVariablesDirevatives);
272  }
274  return get_vecTensor3<0, 1, 2>(dependentVariablesDirevatives);
275  }
276 
278  return get_vecTensor3<9, 0, 0>(dependentVariablesDirevatives);
279  }
281  return get_vecTensor3<9, 1, 0>(dependentVariablesDirevatives);
282  }
284  return get_vecTensor3<9, 0, 1>(dependentVariablesDirevatives);
285  }
287  return get_vecTensor3<9, 1, 1>(dependentVariablesDirevatives);
288  }
290  return get_vecTensor3<9, 0, 2>(dependentVariablesDirevatives);
291  }
293  return get_vecTensor3<9, 1, 2>(dependentVariablesDirevatives);
294  }
295 
297  return get_VecTensor2<18 * (9 + 9) + 0>(dependentVariablesDirevatives);
298  }
300  return get_VecTensor2<18 * (9 + 9) + 9>(dependentVariablesDirevatives);
301  }
302 
304  return get_vecTensor3<9 + 9 + 1, 0, 0>(dependentVariablesDirevatives);
305  }
307  return get_vecTensor3<9 + 9 + 1, 1, 0>(dependentVariablesDirevatives);
308  }
310  return get_vecTensor3<9 + 9 + 1, 0, 1>(dependentVariablesDirevatives);
311  }
313  return get_vecTensor3<9 + 9 + 1, 1, 1>(dependentVariablesDirevatives);
314  }
316  return get_vecTensor3<9 + 9 + 1, 0, 2>(dependentVariablesDirevatives);
317  }
319  return get_vecTensor3<9 + 9 + 1, 1, 2>(dependentVariablesDirevatives);
320  }
321 
322  /**@}*/
323 };
324 
326  static double dispBcScale;
327  static double rotBcScale;
329 };
330 
331 struct BcDisp {
332  BcDisp(std::vector<double> &attr, Range &faces);
333  Range faces;
336 };
337 typedef std::vector<BcDisp> BcDispVec;
338 
339 struct BcRot {
340  BcRot(std::vector<double> &attr, Range &faces);
341  Range faces;
343  double theta;
344 };
345 typedef std::vector<BcRot> BcRotVec;
346 
347 typedef std::vector<Range> TractionFreeBc;
348 
350  const int tAg; ///< adol-c tape
351  const bool evalRhs;
352  const bool evalLhs;
353  boost::shared_ptr<DataAtIntegrationPts>
354  dataAtPts; ///< data at integration pts
355  boost::shared_ptr<PhysicalEquations>
356  physicsPtr; ///< material physical equations
357 
358  OpJacobian(const std::string &field_name, const int tag, const bool eval_rhs,
359  const bool eval_lhs,
360  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
361  boost::shared_ptr<PhysicalEquations> &physics_ptr)
362  : ForcesAndSourcesCore::UserDataOperator(field_name, OPROW), tAg(tag),
363  evalRhs(eval_rhs), evalLhs(eval_lhs), dataAtPts(data_ptr),
364  physicsPtr(physics_ptr) {}
365 
366  virtual MoFEMErrorCode evaluateRhs(EntData &data) = 0;
367  virtual MoFEMErrorCode evaluateLhs(EntData &data) = 0;
368 
369  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
370 };
371 
372 template <typename T> struct OpAssembleBasic : public T {
373 
374  boost::shared_ptr<DataAtIntegrationPts>
375  dataAtPts; ///< data at integration pts
376  const bool assembleSymmetry;
377 
378  OpAssembleBasic(const std::string &field_name,
379  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
380  const char type)
381  : T(field_name, type), dataAtPts(data_ptr), assembleSymmetry(false) {}
382 
383  OpAssembleBasic(const std::string &row_field, const std::string &col_field,
384  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
385  const char type, const bool assemble_symmetry = false)
386  : T(row_field, col_field, type, false), dataAtPts(data_ptr),
387  assembleSymmetry(assemble_symmetry) {}
388 
389  VectorDouble nF; ///< local right hand side vector
390  MatrixDouble K; ///< local tangent matrix
392 
395  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
397  }
398  virtual MoFEMErrorCode assemble(EntData &data) {
400  double *vec_ptr = &*nF.begin();
401  int nb_dofs = data.getIndices().size();
402  int *ind_ptr = &*data.getIndices().begin();
403  CHKERR VecSetValues(T::getFEMethod()->snes_f, nb_dofs, ind_ptr, vec_ptr,
404  ADD_VALUES);
406  }
407 
408  virtual MoFEMErrorCode integrate(EntData &row_data, EntData &col_data) {
410  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not yet implemented");
412  }
413 
414  virtual MoFEMErrorCode assemble(EntData &row_data, EntData &col_data) {
416  double *mat_ptr = &*K.data().begin();
417  int *row_ind_ptr = &*row_data.getIndices().begin();
418  int *col_ind_ptr = &*col_data.getIndices().begin();
419  int row_nb_dofs = row_data.getIndices().size();
420  int col_nb_dofs = col_data.getIndices().size();
421  CHKERR MatSetValues(T::getFEMethod()->snes_B, row_nb_dofs, row_ind_ptr,
422  col_nb_dofs, col_ind_ptr, mat_ptr, ADD_VALUES);
423  if(assembleSymmetry) {
424  transposeK.resize(col_nb_dofs, row_nb_dofs);
425  noalias(transposeK) = trans(K);
426  double *mat_ptr = &*transposeK.data().begin();
427  CHKERR MatSetValues(T::getFEMethod()->snes_B, col_nb_dofs, col_ind_ptr,
428  row_nb_dofs, row_ind_ptr, mat_ptr, ADD_VALUES);
429  }
431  }
432 
433  MoFEMErrorCode doWork(int side, EntityType type, EntData &data) {
435  if (data.getIndices().empty())
437  nF.resize(data.getIndices().size(), false);
438  nF.clear();
439  CHKERR integrate(data);
440  CHKERR assemble(data);
442  }
443 
444  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
445  EntityType col_type, EntData &row_data,
446  EntData &col_data) {
448  if (row_data.getIndices().empty())
450  if (col_data.getIndices().empty())
452  K.resize(row_data.getIndices().size(), col_data.getIndices().size(), false);
453  K.clear();
454  CHKERR integrate(row_data, col_data);
455  CHKERR assemble(row_data, col_data);
457  }
458 };
459 
460 struct OpAssembleVolume : public OpAssembleBasic<VolUserDataOperator> {
461  OpAssembleVolume(const std::string &field,
462  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
463  const char type)
464  : OpAssembleBasic<VolUserDataOperator>(field, data_ptr, type) {}
465 
466  OpAssembleVolume(const std::string &row_field, const std::string &col_field,
467  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
468  const char type, const bool assemble_symmetry = false)
469  : OpAssembleBasic<VolUserDataOperator>(row_field, col_field, data_ptr,
470  type, assemble_symmetry) {}
471 };
472 
473 struct OpAssembleFace : public OpAssembleBasic<FaceUserDataOperator> {
474  OpAssembleFace(const std::string &field,
475  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
476  const char type)
477  : OpAssembleBasic<FaceUserDataOperator>(field, data_ptr, type) {}
478 
479  OpAssembleFace(const std::string &row_field, const std::string &col_field,
480  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
481  const char type)
482  : OpAssembleBasic<FaceUserDataOperator>(row_field, col_field, data_ptr,
483  type) {}
484 };
485 
487  boost::shared_ptr<DataAtIntegrationPts>
488  dataAtPts; ///< data at integration pts
489  OpCalculateRotationAndSpatialGradient(const std::string &field_name,
490  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
491  : VolUserDataOperator(field_name, OPROW), dataAtPts(data_ptr) {}
492  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
493 };
494 
496  OpSpatialEquilibrium(const std::string &field_name,
497  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
498  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
500 };
501 
503  OpSpatialRotation(const std::string &field_name,
504  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
505  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
507 };
508 
510  OpSpatialPhysical(const std::string &field_name,
511  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
512  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
514 };
515 
517  OpSpatialConsistencyP(const std::string &field_name,
518  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
519  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
521 };
522 
524  OpSpatialConsistencyBubble(const std::string &field_name,
525  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
526  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
528 };
529 
531  OpSpatialConsistencyDivTerm(const std::string &field_name,
532  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
533  : OpAssembleVolume(field_name, data_ptr, OPROW) {}
535 };
536 
537 struct OpDispBc : public OpAssembleFace {
538  boost::shared_ptr<BcDispVec> bcDispPtr;
539  OpDispBc(const std::string &field_name,
540  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
541  boost::shared_ptr<BcDispVec> &bc_disp_ptr)
542  : OpAssembleFace(field_name, data_ptr, OPROW), bcDispPtr(bc_disp_ptr) {}
544 };
545 
546 struct OpRotationBc : public OpAssembleFace {
547  boost::shared_ptr<BcRotVec> bcRotPtr;
548  OpRotationBc(const std::string &field_name,
549  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
550  boost::shared_ptr<BcRotVec> &bc_rot_ptr)
551  : OpAssembleFace(field_name, data_ptr, OPROW), bcRotPtr(bc_rot_ptr) {}
553 };
554 
556  std::string fieldName;
557  boost::shared_ptr<TractionFreeBc> bcData;
559  FeTractionFreeBc(const std::string field_name,
560  boost::shared_ptr<TractionFreeBc> &bc_data)
561  : BasicMethod(), fieldName(field_name), bcData(bc_data) {}
564 };
565 
567  OpSpatialEquilibrium_dw_dP(const std::string &row_field,
568  const std::string &col_field,
569  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
570  const bool assemble_off = false)
571  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
572  assemble_off) {
573  sYmm = false;
574  }
575  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
576 };
577 
579  OpSpatialPhysical_du_du(const std::string &row_field,
580  const std::string &col_field,
581  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
582  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL) {
583  sYmm = false;
584  }
585  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
586 };
587 
589  OpSpatialPhysical_du_dP(const std::string &row_field,
590  const std::string &col_field,
591  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
592  const bool assemble_off = false)
593  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
594  assemble_off) {
595  sYmm = false;
596  }
597  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
598 };
599 
601  OpSpatialPhysical_du_dBubble(const std::string &row_field,
602  const std::string &col_field,
603  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
604  const bool assemble_off = false)
605  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
606  assemble_off) {
607  sYmm = false;
608  }
609  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
610 };
611 
613  OpSpatialPhysical_du_domega(const std::string &row_field,
614  const std::string &col_field,
615  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
616  const bool assemble_off = false)
617  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
618  assemble_off) {
619  sYmm = false;
620  }
621  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
622 };
623 
625  OpSpatialRotation_domega_dP(const std::string &row_field,
626  const std::string &col_field,
627  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
628  const bool assemble_off = false)
629  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
630  assemble_off) {
631  sYmm = false;
632  }
633  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
634 };
635 
637  OpSpatialRotation_domega_dBubble(const std::string &row_field,
638  const std::string &col_field,
639  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
640  const bool assemble_off = false)
641  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
642  assemble_off) {
643  sYmm = false;
644  }
645  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
646 };
647 
650  const std::string &row_field, const std::string &col_field,
651  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
652  const bool assemble_off = false)
653  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
654  assemble_off) {
655  sYmm = false;
656  }
657  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
658 };
659 
662  const std::string &row_field, const std::string &col_field,
663  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
664  const bool assemble_off = false)
665  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
666  assemble_off) {
667  sYmm = false;
668  }
669  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
670 };
671 
674  const std::string &row_field, const std::string &col_field,
675  boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
676  const bool assemble_off = false)
677  : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL,
678  assemble_off) {
679  sYmm = false;
680  }
681  MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
682 };
683 
684 // struct OpRhsMaterial : public OpAssembleVolume {
685 // OpRhsMaterial(const std::string &field_name,
686 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
687 // : OpAssembleVolume(field_name, data_ptr, OPROW) {}
688 // MoFEMErrorCode integrate(EntData &data);
689 // };
690 
691 // struct OpMaterial_dW_dx : public OpAssembleVolume {
692 // OpMaterial_dW_dx(const std::string &row_field, const std::string &col_field,
693 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
694 // : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL) {
695 // sYmm = false;
696 // }
697 // MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
698 // };
699 
700 // struct OpMaterial_dW_dW : public OpAssembleVolume {
701 // OpMaterial_dW_dW(const std::string &row_field,
702 // const std::string &col_field,
703 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
704 // : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL) {
705 // sYmm = false;
706 // }
707 // MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
708 // };
709 
710 // struct OpMaterial_dW_dTau : public OpAssembleVolume {
711 // OpMaterial_dW_dTau(const std::string &row_field,
712 // const std::string &col_field,
713 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
714 // : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL) {
715 // sYmm = false;
716 // }
717 // MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
718 // };
719 
720 // struct OpRhsFlow : public OpAssembleVolume {
721 // OpRhsFlow(const std::string &field_name,
722 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
723 // : OpAssembleVolume(field_name, data_ptr, OPROW) {}
724 // MoFEMErrorCode integrate(EntData &data);
725 // };
726 
727 // struct OpFlow_dTau_dx : public OpAssembleVolume {
728 // OpFlow_dTau_dx(const std::string &row_field, const std::string &col_field,
729 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
730 // : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL) {
731 // sYmm = false;
732 // }
733 // MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
734 // };
735 
736 // struct OpFlow_dTau_dW : public OpAssembleVolume {
737 // OpFlow_dTau_dW(const std::string &row_field, const std::string &col_field,
738 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
739 // : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL) {
740 // sYmm = false;
741 // }
742 // MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
743 // };
744 
745 // struct OpFlow_dTau_dLambda : public OpAssembleVolume {
746 // OpFlow_dTau_dLambda(const std::string &row_field, const std::string &col_field,
747 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
748 // : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL) {
749 // sYmm = false;
750 // }
751 // MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
752 // };
753 
754 // struct OpRhsPhi : public OpAssembleVolume {
755 // double mu;
756 // double gamma;
757 // OpRhsPhi(const std::string &field,
758 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr, const double mu,
759 // const double gamma)
760 // : OpAssembleVolume(field, data_ptr, OPROW), mu(mu), gamma(gamma) {}
761 // MoFEMErrorCode integrate(EntData &data);
762 // };
763 
764 // struct OpPhi_dLambda_dW : public OpAssembleVolume {
765 // double mu;
766 // double gamma;
767 // OpPhi_dLambda_dW(const std::string &row_field, const std::string &col_field,
768 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
769 // const double mu, double gamma)
770 // : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL), mu(mu),
771 // gamma(gamma) {
772 // sYmm = false;
773 // }
774 
775 // MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
776 // };
777 
778 // struct OpPhi_dLambda_dx : public OpAssembleVolume {
779 // double mu;
780 // double gamma;
781 // OpPhi_dLambda_dx(const std::string &row_field, const std::string &col_field,
782 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
783 // const double mu, double gamma)
784 // : OpAssembleVolume(row_field, col_field, data_ptr, OPROWCOL), mu(mu),
785 // gamma(gamma) {
786 // sYmm = false;
787 // }
788 // MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
789 // };
790 
791 // struct OpPhi_dLambda_dLambda : public OpAssembleVolume {
792 // double mu;
793 // double gamma;
794 // OpPhi_dLambda_dLambda(const std::string &field,
795 // boost::shared_ptr<DataAtIntegrationPts> &data_ptr,
796 // const double mu, double gamma)
797 // : OpAssembleVolume(field, data_ptr, OPROWCOL), mu(mu),
798 // gamma(gamma) {
799 // sYmm = false;
800 // }
801 
802 // MoFEMErrorCode integrate(EntData &row_data, EntData &col_data);
803 // };
804 
806 
807  moab::Interface &postProcMesh;
808  std::vector<EntityHandle> &mapGaussPts;
809  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
810 
811  OpPostProcDataStructure(moab::Interface &post_proc_mesh,
812  std::vector<EntityHandle> &map_gauss_pts,
813  const std::string field_name,
814  boost::shared_ptr<DataAtIntegrationPts> &data_ptr)
815  : VolUserDataOperator(field_name,
816  ForcesAndSourcesCore::UserDataOperator::OPROW),
817  postProcMesh(post_proc_mesh), mapGaussPts(map_gauss_pts),
818  dataAtPts(data_ptr) {}
819 
820  MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
821 };
822 
824 
825  /**
826  * \brief Getting interface of core database
827  * @param uuid unique ID of interface
828  * @param iface returned pointer to interface
829  * @return error code
830  */
831  MoFEMErrorCode query_interface(const MOFEMuuid &uuid,
832  UnknownInterface **iface) const;
833 
835 
836  boost::shared_ptr<DataAtIntegrationPts> dataAtPts;
837  boost::shared_ptr<PhysicalEquations> physicalEquations;
838 
839  boost::shared_ptr<ForcesAndSourcesCore> elasticFeRhs;
840  boost::shared_ptr<ForcesAndSourcesCore> elasticFeLhs;
841  boost::shared_ptr<ForcesAndSourcesCore> materialFeRhs;
842  boost::shared_ptr<ForcesAndSourcesCore> materialFeLhs;
843  boost::shared_ptr<ForcesAndSourcesCore> epFeRhs;
844  boost::shared_ptr<ForcesAndSourcesCore> epFeLhs;
845  boost::shared_ptr<ForcesAndSourcesCore> elasticBcRhs;
846  boost::shared_ptr<ForcesAndSourcesCore> materialBcRhs;
847  boost::shared_ptr<ForcesAndSourcesCore> epBcRhs;
848  boost::shared_ptr<FEMethod> tractionFreeBc;
849 
850  DM dM;
853 
854  const std::string piolaStress;
855  const std::string eshelbyStress;
856  const std::string spatialDisp;
857  const std::string materialDisp;
858  const std::string streachGradient;
859  const std::string rotAxis;
860  const std::string materialGradient;
861  const std::string materialGradient0;
862  const std::string tauField;
863  const std::string lambdaField;
864  const std::string bubbleField;
865 
866  const std::string elementVolumeName;
867  const std::string elementFaceName;
868 
870  virtual ~EshelbianCore();
871 
874  double E;
875  double nu;
876  double sigma_y;
877  double gamma;
878  double loadFactor;
880 
882 
883  boost::shared_ptr<MethodForForceScaling> scaleSpatialBc;
884  boost::shared_ptr<MethodForForceScaling> scaleMaterialBc;
885 
886  boost::shared_ptr<BcDispVec> bcSpatialDispVecPtr;
887  boost::shared_ptr<BcRotVec> bcSpatialRotationVecPtr;
888  boost::shared_ptr<BcDispVec> bcMaterialDispVecPtr;
889  boost::shared_ptr<TractionFreeBc> tractionSpatialFreeBc;
890  boost::shared_ptr<TractionFreeBc> tractionMaterialFreeBc;
891 
892  template <typename BC>
893  MoFEMErrorCode getDisplacementsBc(boost::shared_ptr<BC> &bc_vec_ptr,
894  const std::string block_set_name,
895  const int nb_attributes) {
898  if (it->getName().compare(0, block_set_name.length(), block_set_name) ==
899  0) {
900  std::vector<double> block_attributes;
901  CHKERR it->getAttributes(block_attributes);
902  if (block_attributes.size() != nb_attributes) {
903  SETERRQ3(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
904  "In block %s expected %d attributes, but given %d",
905  it->getName().c_str(), nb_attributes,
906  block_attributes.size());
907  }
908  Range faces;
909  CHKERR it->getMeshsetIdEntitiesByDimension(mField.get_moab(), 2, faces,
910  true);
911  bc_vec_ptr->emplace_back(block_attributes, faces);
912  }
913  }
915  }
916 
918  bcSpatialDispVecPtr = boost::shared_ptr<BcDispVec>(new BcDispVec());
919  return getDisplacementsBc(bcSpatialDispVecPtr, "SPATIAL_DISP_BC",6);
920  }
921 
923  bcSpatialRotationVecPtr = boost::shared_ptr<BcRotVec>(new BcRotVec());
924  return getDisplacementsBc(bcSpatialRotationVecPtr, "SPATIAL_ROTATION_BC",
925  4);
926  }
927 
929  bcMaterialDispVecPtr = boost::shared_ptr<BcDispVec>(new BcDispVec());
930  return getDisplacementsBc(bcMaterialDispVecPtr, "MATERIAL_DISP_BC",6);
931  }
932 
934  boost::shared_ptr<TractionFreeBc> &bc_ptr,
935  const std::string disp_block_set_name,
936  const std::string rot_block_set_name);
937  inline MoFEMErrorCode
940  boost::shared_ptr<TractionFreeBc>(new TractionFreeBc());
941  return getTractionFreeBc(meshset, tractionSpatialFreeBc, "SPATIAL_DISP_BC",
942  "SPATIAL_ROTATION_BC");
943  }
944  inline MoFEMErrorCode
947  boost::shared_ptr<TractionFreeBc>(new TractionFreeBc());
949  "MATERIAL_DISP_BC", "MATERIAL_ROTATION_BC");
950  }
951 
952  MoFEMErrorCode addFields(const EntityHandle meshset = 0);
955  MoFEMErrorCode addDMs(const BitRefLevel bit = BitRefLevel().set(0));
956 
958  const double lambda,
959  const double mu,
960  const double sigma_y);
961 
962  MoFEMErrorCode setBaseVolumeElement(const int tag, const bool do_rhs,
963  const bool do_lhs,
964  boost::shared_ptr<ForcesAndSourcesCore> &fe);
965 
967  setGenericVolumeElement(const int tag,
968  const bool add_elastic, const bool add_material,
969  boost::shared_ptr<ForcesAndSourcesCore> &fe_rhs,
970  boost::shared_ptr<ForcesAndSourcesCore> &fe_lhs);
971 
973  setGenericFaceElement(const bool add_elastic, const bool add_material,
974  boost::shared_ptr<ForcesAndSourcesCore> &fe);
975 
976  MoFEMErrorCode setElasticElement(const int tag);
978  MoFEMErrorCode setMaterialElement(const int tag);
980  MoFEMErrorCode setEpElement(const int tag);
982 
983  MoFEMErrorCode setUpSnesElastic(SNES snes, Mat m, Vec f);
984  MoFEMErrorCode solveElastic(SNES snes, Vec x);
985  MoFEMErrorCode setUpSnesMaterial(SNES snes, Mat m, Vec f);
986  MoFEMErrorCode solveMaterial(SNES snes, Vec x);
987  MoFEMErrorCode setUpSnesEp(SNES snes, Mat m, Vec f);
988  MoFEMErrorCode solveEp(SNES snes_elastic, SNES snes_ep, Vec x_elastic,
989  Vec x_ep, double delta_lambda, const int nb_steps);
990 
991  MoFEMErrorCode postProcessResults(const int tag, const std::string file);
992 };
993 
994 } // namespace EshelbianPlasticity
995 
996 #endif //__ESHELBIAN_PLASTICITY_HPP__
boost::shared_ptr< PhysicalEquations > physicsPtr
material physical equations
MoFEMErrorCode setUpSnesElastic(SNES snes, Mat m, Vec f)
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)
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Common.hpp:60
VectorDouble nF
local right hand side vector
MoFEMErrorCode setUpSnesMaterial(SNES snes, Mat m, Vec f)
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)
boost::shared_ptr< TractionFreeBc > bcData
OpDispBc(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcDispVec > &bc_disp_ptr)
MoFEMErrorCode integrate(EntData &data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
OpRotationBc(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, boost::shared_ptr< BcRotVec > &bc_rot_ptr)
boost::shared_ptr< MethodForForceScaling > scaleMaterialBc
Class used to scale loads, f.e. in arc-length control.
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)
OpAssembleVolume(const std::string &field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr, const char type)
std::vector< EntityHandle > & mapGaussPts
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:459
boost::shared_ptr< DataAtIntegrationPts > dataAtPts
data at integration pts
MoFEMErrorCode solveEp(SNES snes_elastic, SNES snes_ep, Vec x_elastic, Vec x_ep, double delta_lambda, const int nb_steps)
MoFEMErrorCode integrate(EntData &row_data, EntData &col_data)
std::vector< BcRot > BcRotVec
MoFEMErrorCode integrate(EntData &data)
MoFEMErrorCode solveElastic(SNES snes, Vec x)
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:490
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 scaleNf(const FEMethod *fe, VectorDouble &nf)
MoFEMErrorCode integrate(EntData &data)
MoFEMErrorCode addFields(const EntityHandle meshset=0)
BcRot(std::vector< double > &attr, Range &faces)
MoFEMErrorCode setElasticElement(const int tag)
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
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)
static DTensor2Ptr get_VecTensor2(std::vector< double > &v)
boost::shared_ptr< MatrixDouble > MatrixPtr
boost::shared_ptr< ForcesAndSourcesCore > epFeLhs
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
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 setMaterialElement(const int tag)
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Common.hpp:212
boost::shared_ptr< TractionFreeBc > tractionSpatialFreeBc
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)
OpSpatialPhysical_du_du(const std::string &row_field, const std::string &col_field, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode setEpElement(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)
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)
FTensor::Tensor2< adouble, 3, 3 > ATensor2
FTensor::Tensor0< FTensor::PackPtr< double *, 1 > > DTensor0Ptr
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:578
boost::shared_ptr< ForcesAndSourcesCore > epBcRhs
MoFEMErrorCode getSpatialTractionFreeBc(const EntityHandle meshset=0)
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.
MoFEMErrorCode solveMaterial(SNES snes, Vec x)
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
MoFEMErrorCode setGenericVolumeElement(const int tag, const bool add_elastic, const bool add_material, boost::shared_ptr< ForcesAndSourcesCore > &fe_rhs, boost::shared_ptr< ForcesAndSourcesCore > &fe_lhs)
OpSpatialConsistencyBubble(const std::string &field_name, boost::shared_ptr< DataAtIntegrationPts > &data_ptr)
MoFEMErrorCode setBaseVolumeElement(const int tag, const bool do_rhs, const bool do_lhs, boost::shared_ptr< ForcesAndSourcesCore > &fe)
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 setUpSnesEp(SNES snes, Mat m, Vec f)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:403
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 setMaterialElementToSnes(DM dm)
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
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)
MoFEMErrorCode setGenericFaceElement(const bool add_elastic, const bool add_material, boost::shared_ptr< ForcesAndSourcesCore > &fe)
virtual MoFEMErrorCode integrate(EntData &data)
FTensor::Tensor2< FTensor::PackPtr< double *, 1 >, 3, 3 > DTensor2Ptr
BcDisp(std::vector< double > &attr, Range &faces)
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)