v0.13.0
NonLinearElasticElement.hpp
Go to the documentation of this file.
1 /** \file NonLinearElasticElement.hpp
2  * \ingroup nonlinear_elastic_elem
3  * \brief Operators and data structures for non-linear elastic analysis
4  *
5  * Implementation of nonlinear elastic element.
6  */
7 
8 /*
9  * This file is part of MoFEM.
10  * MoFEM is free software: you can redistribute it and/or modify it under
11  * the terms of the GNU Lesser General Public License as published by the
12  * Free Software Foundation, either version 3 of the License, or (at your
13  * option) any later version.
14  *
15  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
16  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
17  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
18  * License for more details.
19  *
20  * You should have received a copy of the GNU Lesser General Public
21  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
22 
23 #ifndef __NONLINEAR_ELASTIC_HPP
24 #define __NONLINEAR_ELASTIC_HPP
25 
26 #ifndef WITH_ADOL_C
27 #error "MoFEM need to be compiled with ADOL-C"
28 #endif
29 
30 /** \brief structure grouping operators and data used for calculation of
31  * nonlinear elastic element \ingroup nonlinear_elastic_elem
32  *
33  * In order to assemble matrices and right hand vectors, the loops over
34  * elements, entities over that elements and finally loop over integration
35  * points are executed.
36  *
37  * Following implementation separate those three categories of loops and to each
38  * loop attach operator.
39  *
40  */
42 
43  /// \brief definition of volume element
45 
46  Mat A;
47  Vec F;
48 
49  int addToRule;
50 
51  MyVolumeFE(MoFEM::Interface &m_field);
52  virtual ~MyVolumeFE() = default;
53 
54  /** \brief it is used to calculate nb. of Gauss integration points
55  *
56  * for more details pleas look
57  * Reference:
58  *
59  * Albert Nijenhuis, Herbert Wilf,
60  * Combinatorial Algorithms for Computers and Calculators,
61  * Second Edition,
62  * Academic Press, 1978,
63  * ISBN: 0-12-519260-6,
64  * LC: QA164.N54.
65  *
66  * More details about algorithm
67  * http://people.sc.fsu.edu/~jburkardt/cpp_src/gm_rule/gm_rule.html
68  **/
69  int getRule(int order);
70 
71  SmartPetscObj<Vec> V;
72  double eNergy;
73 
76  };
77 
78  MyVolumeFE feRhs; ///< calculate right hand side for tetrahedral elements
79  MyVolumeFE &getLoopFeRhs() { return feRhs; } ///< get rhs volume element
80  MyVolumeFE feLhs; //< calculate left hand side for tetrahedral elements
81  MyVolumeFE &getLoopFeLhs() { return feLhs; } ///< get lhs volume element
82 
83  MyVolumeFE feEnergy; ///< calculate elastic energy
84  MyVolumeFE &getLoopFeEnergy() { return feEnergy; } ///< get energy fe
85 
87  short int tAg;
88 
89  NonlinearElasticElement(MoFEM::Interface &m_field, short int tag);
90  virtual ~NonlinearElasticElement() = default;
91 
92  template <typename TYPE> struct FunctionsToCalculatePiolaKirchhoffI;
93 
94  /** \brief data for calculation het conductivity and heat capacity elements
95  * \ingroup nonlinear_elastic_elem
96  */
97  struct BlockData {
98  int iD;
99  double E;
100  double PoissonRatio;
101  // Eberlein Fibres stiffness properties
102  double k1, k2;
103  Range tEts; ///< constrains elements in block set
104  boost::shared_ptr<FunctionsToCalculatePiolaKirchhoffI<adouble>>
106  boost::shared_ptr<FunctionsToCalculatePiolaKirchhoffI<double>>
110  };
111 
112  std::map<int, BlockData>
113  setOfBlocks; ///< maps block set id with appropriate BlockData
114 
115  /** \brief common data used by volume elements
116  * \ingroup nonlinear_elastic_elem
117  */
118  struct CommonData {
119 
120  std::map<std::string, std::vector<VectorDouble>> dataAtGaussPts;
121  std::map<std::string, std::vector<MatrixDouble>> gradAtGaussPts;
124  std::vector<MatrixDouble3by3> sTress;
125  std::vector<MatrixDouble>
126  jacStress; ///< this is simply material tangent operator
127 
128  // This part can be used to calculate stress directly from potential
129 
130  std::vector<double> eNergy;
131  std::vector<VectorDouble> jacEnergy;
132  std::vector<VectorDouble> hessianEnergy;
133  };
135 
139 
140  /** \brief Implementation of elastic (non-linear) St. Kirchhoff equation
141  * \ingroup nonlinear_elastic_elem
142  */
143  template <typename TYPE> struct FunctionsToCalculatePiolaKirchhoffI {
144 
148 
155  }
156 
158 
159  double lambda, mu;
160  MatrixBoundedArray<TYPE, 9> F, C, E, S, invF, P, sIGma, h, H, invH,
164 
166 
167  int gG; ///< Gauss point number
168  CommonData *commonDataPtr; ///< common data shared between entities (f.e.
169  ///< field values at Gauss pts.)
171  *opPtr; ///< pointer to finite element tetrahedral operator
172 
175  t_C(i, j) = t_F(k, i) * t_F(k, j);
177  }
178 
181  constexpr auto t_kd = FTensor::Kronecker_Delta<double>();
182  t_E(i, j) = 0.5 * (t_C(i, j) - t_kd(i, j));
184  }
185 
186  // St. Venant–Kirchhoff Material
189  constexpr auto t_kd = FTensor::Kronecker_Delta<double>();
190  t_S(i, j) = (2 * mu) * t_E(i, j) + (lambda * t_E(k, k)) * t_kd(i, j);
192  }
193 
194  /** \brief Function overload to implement user material
195  *
196 
197  * Calculation of Piola Kirchhoff I is implemented by user. Tangent matrix
198  * user implemented physical equation is calculated using automatic
199  * differentiation.
200 
201  * \f$\mathbf{S} =
202  \lambda\textrm{tr}[\mathbf{E}]\mathbf{I}+2\mu\mathbf{E}\f$
203 
204  * Notes: <br>
205  * Number of actual Gauss point is accessed from variable gG. <br>
206  * Access to operator data structures is available by variable opPtr. <br>
207  * Access to common data is by commonDataPtr. <br>
208 
209  * \param block_data used to give access to material parameters
210  * \param fe_ptr pointer to element data structures
211 
212  For details look to: <br>
213  NONLINEAR CONTINUUM MECHANICS FOR FINITE ELEMENT ANALYSIS, Javier Bonet,
214  Richard D. Wood
215 
216  */
218  const BlockData block_data,
219  boost::shared_ptr<const NumeredEntFiniteElement> fe_ptr) {
221  lambda = LAMBDA(block_data.E, block_data.PoissonRatio);
222  mu = MU(block_data.E, block_data.PoissonRatio);
226  t_P(i, j) = t_F(i, k) * t_S(k, j);
228  }
229 
230  /** \brief Function overload to implement user material
231  *
232 
233  * Calculation of Piola Kirchhoff I is implemented by user. Tangent matrix
234  * user implemented physical equation is calculated using automatic
235  * differentiation.
236 
237  * \f$\mathbf{S} =
238  \lambda\textrm{tr}[\mathbf{E}]\mathbf{I}+2\mu\mathbf{E}\f$
239 
240  * Notes: <br>
241  * Number of actual Gauss point is accessed from variable gG. <br>
242  * Access to operator data structures is available by variable opPtr. <br>
243  * Access to common data is by commonDataPtr. <br>
244 
245  * \param block_data used to give access to material parameters
246  * \param fe_ptr pointer to element data structures
247 
248  For details look to: <br>
249  NONLINEAR CONTINUUM MECHANICS FOR FINITE ELEMENT ANALYSIS, Javier Bonet,
250  Richard D. Wood
251 
252  */
254  const BlockData block_data,
255  boost::shared_ptr<const NumeredEntFiniteElement> fe_ptr) {
257  sigmaCauchy.resize(3, 3);
258  t_sigmaCauchy(i, j) = t_P(i, k) * t_F(j, k);
261  }
262 
263  /**
264  * \brief add additional active variables
265  *
266  * \note This member function if used should be implement by template member
267  * function Specialization, different implementation needed for TYPE=double
268  * or TYPE=adouble
269  *
270  * More complex physical models depend on gradient of defamation and some
271  * additional variables. For example can depend on temperature. This
272  * function adds additional independent variables to the model.
273  *
274  * @param nb_active_variables number of active variables
275  * @return error code
276  */
277  virtual MoFEMErrorCode setUserActiveVariables(int &nb_active_variables) {
280  }
281 
282  /**
283  * \brief Add additional independent variables
284  * More complex physical models depend on gradient of defamation and some
285  * additional variables. For example can depend on temperature. This
286  * function adds additional independent variables to the model.
287  *
288  * /note First 9 elements are reserved for gradient of deformation.
289  * @param activeVariables vector of deepened variables, values after index
290  * 9 should be add.
291  *
292  * @return error code
293  */
294  virtual MoFEMErrorCode
298  }
299 
300  /** \brief Calculate elastic energy density
301  *
302  * \f[\Psi =
303  * \frac{1}{2}\lambda(\textrm{tr}[\mathbf{E}])^2+\mu\mathbf{E}:\mathbf{E}\f]
304  */
306  const BlockData block_data,
307  boost::shared_ptr<const NumeredEntFiniteElement> fe_ptr) {
309  lambda = LAMBDA(block_data.E, block_data.PoissonRatio);
310  mu = MU(block_data.E, block_data.PoissonRatio);
313  TYPE trace = 0;
314  eNergy = 0;
315  for (int ii = 0; ii < 3; ii++) {
316  trace += E(ii, ii);
317  for (int jj = 0; jj < 3; jj++) {
318  TYPE e = E(ii, jj);
319  eNergy += mu * e * e;
320  }
321  }
322  eNergy += 0.5 * lambda * trace * trace;
324  }
325 
326  /** \brief Calculate Eshelby stress
327  */
329  const BlockData block_data,
330  boost::shared_ptr<const NumeredEntFiniteElement> fe_ptr) {
332  CHKERR calculateP_PiolaKirchhoffI(block_data, fe_ptr);
333  CHKERR calculateElasticEnergy(block_data, fe_ptr);
334  constexpr auto t_kd = FTensor::Kronecker_Delta<double>();
335  t_sIGma(i, j) = t_kd(i, j) * eNergy - t_F(k, i) * t_P(k, j);
337  }
338 
339  /** \brief Do operations when pre-process
340  */
342  std::map<std::string, std::vector<VectorDouble>> &field_map,
343  std::map<std::string, std::vector<MatrixDouble>> &grad_map) {
346  }
347 
348  protected:
349  inline static auto resizeAndSet(MatrixBoundedArray<TYPE, 9> &m) {
350  m.resize(3, 3, false);
352  };
353  };
354 
357 
358  std::vector<VectorDouble> &valuesAtGaussPts;
359  std::vector<MatrixDouble> &gradientAtGaussPts;
360  const EntityType zeroAtType;
361 
362  OpGetDataAtGaussPts(const std::string field_name,
363  std::vector<VectorDouble> &values_at_gauss_pts,
364  std::vector<MatrixDouble> &gradient_at_gauss_pts);
365 
366  /** \brief operator calculating deformation gradient
367  *
368  * temperature gradient is calculated multiplying derivatives of shape
369  * functions by degrees of freedom
370  */
371  MoFEMErrorCode doWork(int side, EntityType type,
373  };
374 
376  OpGetCommonDataAtGaussPts(const std::string field_name,
377  CommonData &common_data);
378  };
379 
380  /**
381  * \brief Operator performs automatic differentiation.
382  */
385 
386  BlockData &dAta; ///< Structure keeping data about problem, like material
387  ///< parameters
388  CommonData &commonData; ///< Structure keeping data abut this particular
389  ///< element, e.g. gradient of deformation at
390  ///< integration points
391  int tAg; //,lastId; ///< ADOL-C tag used for recording operations
392  int adlocReturnValue; ///< return value from ADOL-C, if non-zero that is
393  ///< error.
394  bool jAcobian; ///< if true Jacobian is calculated
395  bool fUnction; ///< if true stress i calculated
396  bool aLe; ///< true if arbitrary Lagrangian-Eulerian formulation
397  bool fieldDisp; ///< true if field of displacements is given, usually
398  ///< spatial positions are given.
399 
400  /**
401  \brief Construct operator to calculate Piola-Kirchhoff stress or its
402  derivatives over gradient deformation
403 
404  \param field_name approximation field name of spatial positions or
405  displacements \param data reference to block data (what is Young modulus,
406  Poisson ratio or what elements are part of the block) \param tag adol-c
407  unique tag of the tape \param jacobian if true derivative of Piola Stress
408  is calculated otherwise just stress is calculated \param field_disp if
409  true approximation field keeps displacements not spatial positions
410 
411  */
412  OpJacobianPiolaKirchhoffStress(const std::string field_name,
413  BlockData &data, CommonData &common_data,
414  int tag, bool jacobian, bool ale,
415  bool field_disp);
416 
419 
420  std::vector<MatrixDouble> *ptrh;
421  std::vector<MatrixDouble> *ptrH;
422 
426 
427  /**
428  * \brief Calculate Paola-Kirchhoff I stress
429  * @return error code
430  */
431  virtual MoFEMErrorCode calculateStress(const int gg);
432 
433  /**
434  * \brief Record ADOL-C tape
435  * @return error code
436  */
437  virtual MoFEMErrorCode recordTag(const int gg);
438 
439  /**
440  * \brief Play ADOL-C tape
441  * @return error code
442  */
443  virtual MoFEMErrorCode playTag(const int gg);
444 
445  /**
446  * \brief Cgeck if tape is recorded for given integration point
447  * @param gg integration point
448  * @return true if tag is recorded
449  */
450  virtual bool recordTagForIntegrationPoint(const int gg) {
451  // return true;
452  if (gg == 0)
453  return true;
454  return false;
455  }
456 
457  /**
458  * \brief Calculate stress or jacobian at gauss points
459  *
460  * @param row_side
461  * @param row_type
462  * @param row_data
463  * @return error code
464  */
465  MoFEMErrorCode doWork(int row_side, EntityType row_type,
466  EntitiesFieldData::EntData &row_data);
467  };
468 
469  /**
470  * \brief Calculate explicit derivative of free energy
471  */
474 
477 
478  int tAg; ///< tape tag
479  bool gRadient; ///< if set true gradient of energy is calculated
480  bool hEssian; ///< if set true hessian of energy is calculated
481  bool aLe; ///< true if arbitrary Lagrangian-Eulerian formulation
482  bool fieldDisp; ///< true if displacements instead spatial positions used
483 
484  OpJacobianEnergy(const std::string field_name, ///< field name for spatial
485  ///< positions or
486  ///< displacements
487  BlockData &data, CommonData &common_data, int tag,
488  bool gradient, bool hessian, bool ale, bool field_disp);
489 
492 
493  std::vector<MatrixDouble> *ptrh;
494  std::vector<MatrixDouble> *ptrH;
495 
499 
500  /**
501  * \brief Check if tape is recorded for given integration point
502  * @param gg integration point
503  * @return true if tag is recorded
504  */
505  virtual bool recordTagForIntegrationPoint(const int gg) {
506  if (gg == 0)
507  return true;
508  return false;
509  }
510 
511  /**
512  * \brief Calculate Paola-Kirchhoff I stress
513  * @return error code
514  */
515  virtual MoFEMErrorCode calculateEnergy(const int gg);
516 
517  /**
518  * \brief Record ADOL-C tape
519  * @return error code
520  */
521  virtual MoFEMErrorCode recordTag(const int gg);
522 
523  /**
524  * \brief Play ADOL-C tape
525  * @return error code
526  */
527  virtual MoFEMErrorCode playTag(const int gg);
528 
529  MoFEMErrorCode doWork(int row_side, EntityType row_type,
530  EntitiesFieldData::EntData &row_data);
531  };
532 
535 
538  bool fieldDisp;
539  bool aLe;
540 
541  ublas::vector<int> iNdices;
542  OpRhsPiolaKirchhoff(const std::string field_name, BlockData &data,
543  CommonData &common_data);
544 
546  MoFEMErrorCode doWork(int row_side, EntityType row_type,
547  EntitiesFieldData::EntData &row_data);
548 
549  virtual MoFEMErrorCode aSemble(int row_side, EntityType row_type,
550  EntitiesFieldData::EntData &row_data);
551  };
552 
553  struct OpEnergy
555 
558  SmartPetscObj<Vec> ghostVec;
559  bool fieldDisp;
560 
561  OpEnergy(const std::string field_name, BlockData &data,
562  CommonData &common_data, SmartPetscObj<Vec> ghost_vec,
563  bool field_disp);
564 
565  MoFEMErrorCode doWork(int row_side, EntityType row_type,
566  EntitiesFieldData::EntData &row_data);
567  };
568 
571 
574  int tAg;
575  bool aLe;
576 
577  ublas::vector<int> rowIndices;
578  ublas::vector<int> colIndices;
579 
580  OpLhsPiolaKirchhoff_dx(const std::string vel_field,
581  const std::string field_name, BlockData &data,
582  CommonData &common_data);
583 
585 
586  /**
587  \brief Directive of Piola Kirchhoff stress over spatial DOFs
588 
589  This project derivative \f$\frac{\partial P}{\partial F}\f$, that is
590  \f[
591  \frac{\partial P}{\partial x_\textrm{DOF}} = \frac{\partial P}{\partial
592  F}\frac{\partial F}{\partial x_\textrm{DOF}}, \f] where second therm
593  \f$\frac{\partial F}{\partial x_\textrm{DOF}}\f$ is derivative of shape
594  function
595 
596  */
598  int gg);
599 
600  virtual MoFEMErrorCode aSemble(int row_side, int col_side,
601  EntityType row_type, EntityType col_type,
602  EntitiesFieldData::EntData &row_data,
603  EntitiesFieldData::EntData &col_data);
604 
605  MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
606  EntityType col_type,
607  EntitiesFieldData::EntData &row_data,
608  EntitiesFieldData::EntData &col_data);
609  };
610 
612 
613  OpLhsPiolaKirchhoff_dX(const std::string vel_field,
614  const std::string field_name, BlockData &data,
615  CommonData &common_data);
616 
617  /// \brief Derivative of Piola Kirchhoff stress over material DOFs
619 
620  MoFEMErrorCode aSemble(int row_side, int col_side, EntityType row_type,
621  EntityType col_type,
622  EntitiesFieldData::EntData &row_data,
623  EntitiesFieldData::EntData &col_data);
624  };
625 
627 
628  OpJacobianEshelbyStress(const std::string field_name, BlockData &data,
629  CommonData &common_data, int tag, bool jacobian,
630  bool ale);
631 
632  MoFEMErrorCode calculateStress(const int gg);
633  };
634 
636 
637  OpRhsEshelbyStress(const std::string field_name, BlockData &data,
638  CommonData &common_data);
639  };
640 
641  /**
642  * \deprecated name with spelling mistake
643  */
645 
647 
648  OpLhsEshelby_dx(const std::string vel_field, const std::string field_name,
649  BlockData &data, CommonData &common_data);
650 
652  };
653 
655 
656  OpLhsEshelby_dX(const std::string vel_field, const std::string field_name,
657  BlockData &data, CommonData &common_data);
658 
660  };
661 
664  materialDoublePtr,
666  materialAdoublePtr);
667 
669  const std::string element_name,
670  const std::string spatial_position_field_name,
671  const std::string material_position_field_name = "MESH_NODE_POSITIONS",
672  const bool ale = false);
673 
674  /** \brief Set operators to calculate left hand tangent matrix and right hand
675  * residual
676  *
677  * \param fun class needed to calculate Piola Kirchhoff I Stress tensor
678  * \param spatial_position_field_name name of approximation field
679  * \param material_position_field_name name of field to define geometry
680  * \param ale true if arbitrary Lagrangian Eulerian formulation
681  * \param field_disp true if approximation field represents displacements
682  * otherwise it is field of spatial positions
683  */
685  const std::string spatial_position_field_name,
686  const std::string material_position_field_name = "MESH_NODE_POSITIONS",
687  const bool ale = false, const bool field_disp = false);
688 };
689 
690 #endif //__NONLINEAR_ELASTIC_HPP
691 
692 /**
693  * \defgroup nonlinear_elastic_elem NonLinear Elastic Element
694  * \ingroup user_modules
695  * \defgroup user_modules User modules
696  **/
EntitiesFieldData::EntData EntData
Kronecker Delta class.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
#define DEPRECATED
Definition: definitions.h:30
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
#define MU(E, NU)
Definition: fem_tools.h:33
#define LAMBDA(E, NU)
Definition: fem_tools.h:32
constexpr auto t_kd
const FTensor::Tensor2< T, Dim, Dim > Vec
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:88
UBlasVector< double > VectorDouble
Definition: Types.hpp:79
auto getFTensor2FromArray3by3(ublas::matrix< T, L, A > &data, const FTensor::Number< S > &, const size_t rr)
Definition: Templates.hpp:940
static auto determinantTensor3by3(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
Definition: Templates.hpp:1161
double trace(FTensor::Tensor2_symmetric< T, SPACE_DIM > &t_stress)
Definition: PlasticOps.hpp:351
FTensor::Index< 'm', 3 > m
Deprecated interface functions.
FTensor::Tensor2< double *, 3, 3 > & getJac()
get element Jacobian
data for calculation heat conductivity and heat capacity elements
boost::shared_ptr< FunctionsToCalculatePiolaKirchhoffI< adouble > > materialAdoublePtr
boost::shared_ptr< FunctionsToCalculatePiolaKirchhoffI< double > > materialDoublePtr
Range tEts
constrains elements in block set
common data used by volume elements
std::vector< MatrixDouble > jacStress
this is simply material tangent operator
std::map< std::string, std::vector< VectorDouble > > dataAtGaussPts
std::vector< MatrixDouble3by3 > sTress
std::map< std::string, std::vector< MatrixDouble > > gradAtGaussPts
Implementation of elastic (non-linear) St. Kirchhoff equation.
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_h
virtual MoFEMErrorCode calculatesIGma_EshelbyStress(const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
Calculate Eshelby stress.
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_invH
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_C
virtual MoFEMErrorCode calculateP_PiolaKirchhoffI(const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
Function overload to implement user material.
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_sIGma
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_F
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_sigmaCauchy
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_P
virtual MoFEMErrorCode calculateElasticEnergy(const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
Calculate elastic energy density.
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_H
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_E
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_invF
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_S
virtual MoFEMErrorCode setUserActiveVariables(int &nb_active_variables)
add additional active variables
MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperator * opPtr
pointer to finite element tetrahedral operator
static auto resizeAndSet(MatrixBoundedArray< TYPE, 9 > &m)
virtual MoFEMErrorCode getDataOnPostProcessor(std::map< std::string, std::vector< VectorDouble >> &field_map, std::map< std::string, std::vector< MatrixDouble >> &grad_map)
Do operations when pre-process.
virtual MoFEMErrorCode calculateCauchyStress(const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
Function overload to implement user material.
virtual MoFEMErrorCode setUserActiveVariables(VectorDouble &activeVariables)
Add additional independent variables More complex physical models depend on gradient of defamation an...
MoFEMErrorCode preProcess()
function is run at the beginning of loop
MoFEMErrorCode postProcess()
function is run at the end of loop
int getRule(int order)
it is used to calculate nb. of Gauss integration points
OpEnergy(const std::string field_name, BlockData &data, CommonData &common_data, SmartPetscObj< Vec > ghost_vec, bool field_disp)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
OpGetCommonDataAtGaussPts(const std::string field_name, CommonData &common_data)
MoFEMErrorCode doWork(int side, EntityType type, EntitiesFieldData::EntData &data)
operator calculating deformation gradient
OpGetDataAtGaussPts(const std::string field_name, std::vector< VectorDouble > &values_at_gauss_pts, std::vector< MatrixDouble > &gradient_at_gauss_pts)
Calculate explicit derivative of free energy.
bool gRadient
if set true gradient of energy is calculated
bool fieldDisp
true if displacements instead spatial positions used
virtual bool recordTagForIntegrationPoint(const int gg)
Check if tape is recorded for given integration point.
virtual MoFEMErrorCode playTag(const int gg)
Play ADOL-C tape.
bool hEssian
if set true hessian of energy is calculated
OpJacobianEnergy(const std::string field_name, BlockData &data, CommonData &common_data, int tag, bool gradient, bool hessian, bool ale, bool field_disp)
virtual MoFEMErrorCode recordTag(const int gg)
Record ADOL-C tape.
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
virtual MoFEMErrorCode calculateEnergy(const int gg)
Calculate Paola-Kirchhoff I stress.
bool aLe
true if arbitrary Lagrangian-Eulerian formulation
MoFEMErrorCode calculateStress(const int gg)
Calculate Paola-Kirchhoff I stress.
OpJacobianEshelbyStress(const std::string field_name, BlockData &data, CommonData &common_data, int tag, bool jacobian, bool ale)
virtual MoFEMErrorCode recordTag(const int gg)
Record ADOL-C tape.
bool aLe
true if arbitrary Lagrangian-Eulerian formulation
virtual MoFEMErrorCode playTag(const int gg)
Play ADOL-C tape.
OpJacobianPiolaKirchhoffStress(const std::string field_name, BlockData &data, CommonData &common_data, int tag, bool jacobian, bool ale, bool field_disp)
Construct operator to calculate Piola-Kirchhoff stress or its derivatives over gradient deformation.
virtual bool recordTagForIntegrationPoint(const int gg)
Cgeck if tape is recorded for given integration point.
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
Calculate stress or jacobian at gauss points.
virtual MoFEMErrorCode calculateStress(const int gg)
Calculate Paola-Kirchhoff I stress.
OpLhsEshelby_dX(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data)
OpLhsEshelby_dx(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data)
OpLhsPiolaKirchhoff_dX(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data)
MoFEMErrorCode aSemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
virtual MoFEMErrorCode aSemble(int row_side, int col_side, EntityType row_type, EntityType col_type, EntitiesFieldData::EntData &row_data, EntitiesFieldData::EntData &col_data)
OpLhsPiolaKirchhoff_dx(const std::string vel_field, const std::string field_name, BlockData &data, CommonData &common_data)
OpRhsEshelbyStress(const std::string field_name, BlockData &data, CommonData &common_data)
virtual MoFEMErrorCode aSemble(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
MoFEMErrorCode doWork(int row_side, EntityType row_type, EntitiesFieldData::EntData &row_data)
OpRhsPiolaKirchhoff(const std::string field_name, BlockData &data, CommonData &common_data)
structure grouping operators and data used for calculation of nonlinear elastic element
MyVolumeFE feRhs
calculate right hand side for tetrahedral elements
NonlinearElasticElement(MoFEM::Interface &m_field, short int tag)
MoFEMErrorCode addElement(const std::string element_name, const std::string spatial_position_field_name, const std::string material_position_field_name="MESH_NODE_POSITIONS", const bool ale=false)
std::map< int, BlockData > setOfBlocks
maps block set id with appropriate BlockData
FTensor::Index< 'k', 3 > k
DEPRECATED typedef OpRhsEshelbyStress OpRhsEshelbyStrees
MyVolumeFE & getLoopFeEnergy()
get energy fe
FTensor::Index< 'j', 3 > j
MoFEMErrorCode setBlocks(boost::shared_ptr< FunctionsToCalculatePiolaKirchhoffI< double >> materialDoublePtr, boost::shared_ptr< FunctionsToCalculatePiolaKirchhoffI< adouble >> materialAdoublePtr)
MyVolumeFE & getLoopFeLhs()
get lhs volume element
virtual ~NonlinearElasticElement()=default
MyVolumeFE feEnergy
calculate elastic energy
MyVolumeFE & getLoopFeRhs()
get rhs volume element
FTensor::Index< 'i', 3 > i
MoFEMErrorCode setOperators(const std::string spatial_position_field_name, const std::string material_position_field_name="MESH_NODE_POSITIONS", const bool ale=false, const bool field_disp=false)
Set operators to calculate left hand tangent matrix and right hand residual.
@ TYPE
Definition: inflate.h:32