v0.14.0
Public Member Functions | Public Attributes | Static Protected Member Functions | List of all members
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE > Struct Template Reference

Implementation of elastic (non-linear) St. Kirchhoff equation. More...

#include <users_modules/basic_finite_elements/src/NonLinearElasticElement.hpp>

Inheritance diagram for NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >:
[legend]
Collaboration diagram for NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >:
[legend]

Public Member Functions

 FunctionsToCalculatePiolaKirchhoffI ()
 
virtual ~FunctionsToCalculatePiolaKirchhoffI ()=default
 
MoFEMErrorCode calculateC_CauchyDeformationTensor ()
 
MoFEMErrorCode calculateE_GreenStrain ()
 
MoFEMErrorCode calculateS_PiolaKirchhoffII ()
 
virtual MoFEMErrorCode calculateP_PiolaKirchhoffI (const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Function overload to implement user material. More...
 
virtual MoFEMErrorCode calculateCauchyStress (const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Function overload to implement user material. More...
 
virtual MoFEMErrorCode setUserActiveVariables (int &nb_active_variables)
 add additional active variables More...
 
virtual MoFEMErrorCode setUserActiveVariables (VectorDouble &activeVariables)
 Add additional independent variables More complex physical models depend on gradient of defamation and some additional variables. For example can depend on temperature. This function adds additional independent variables to the model. More...
 
virtual MoFEMErrorCode calculateElasticEnergy (const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Calculate elastic energy density. More...
 
virtual MoFEMErrorCode calculatesIGma_EshelbyStress (const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Calculate Eshelby stress. More...
 
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. More...
 

Public Attributes

FTensor::Index< 'i', 3 > i
 
FTensor::Index< 'j', 3 > j
 
FTensor::Index< 'k', 3 > k
 
double lambda
 
double mu
 
MatrixBoundedArray< TYPE, 9 > F
 
MatrixBoundedArray< TYPE, 9 > C
 
MatrixBoundedArray< TYPE, 9 > E
 
MatrixBoundedArray< TYPE, 9 > S
 
MatrixBoundedArray< TYPE, 9 > invF
 
MatrixBoundedArray< TYPE, 9 > P
 
MatrixBoundedArray< TYPE, 9 > sIGma
 
MatrixBoundedArray< TYPE, 9 > h
 
MatrixBoundedArray< TYPE, 9 > H
 
MatrixBoundedArray< TYPE, 9 > invH
 
MatrixBoundedArray< TYPE, 9 > sigmaCauchy
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_F
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_C
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_E
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_S
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_invF
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_P
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_sIGma
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_h
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_H
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_invH
 
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_sigmaCauchy
 
TYPE J
 
TYPE eNergy
 
TYPE detH
 
TYPE detF
 
int gG
 Gauss point number. More...
 
CommonDatacommonDataPtr
 
MoFEM::VolumeElementForcesAndSourcesCore::UserDataOperatoropPtr
 pointer to finite element tetrahedral operator More...
 

Static Protected Member Functions

static auto resizeAndSet (MatrixBoundedArray< TYPE, 9 > &m)
 

Detailed Description

template<typename TYPE>
struct NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >

Implementation of elastic (non-linear) St. Kirchhoff equation.

Examples
NeoHookean.hpp.

Definition at line 79 of file NonLinearElasticElement.hpp.

Constructor & Destructor Documentation

◆ FunctionsToCalculatePiolaKirchhoffI()

◆ ~FunctionsToCalculatePiolaKirchhoffI()

template<typename TYPE >
virtual NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::~FunctionsToCalculatePiolaKirchhoffI ( )
virtualdefault

Member Function Documentation

◆ calculateC_CauchyDeformationTensor()

template<typename TYPE >
MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::calculateC_CauchyDeformationTensor ( )
inline
Examples
NeoHookean.hpp.

Definition at line 160 of file NonLinearElasticElement.hpp.

160  {
162  t_C(i, j) = t_F(k, i) * t_F(k, j);
164  }

◆ calculateCauchyStress()

template<typename TYPE >
virtual MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::calculateCauchyStress ( const BlockData  block_data,
boost::shared_ptr< const NumeredEntFiniteElement >  fe_ptr 
)
inlinevirtual

Function overload to implement user material.

Calculation of Piola Kirchhoff I is implemented by user. Tangent matrix user implemented physical equation is calculated using automatic differentiation.

\(\mathbf{S} = \lambda\textrm{tr}[\mathbf{E}]\mathbf{I}+2\mu\mathbf{E}\)

Notes:
Number of actual Gauss point is accessed from variable gG.
Access to operator data structures is available by variable opPtr.
Access to common data is by commonDataPtr.

Parameters
block_dataused to give access to material parameters
fe_ptrpointer to element data structures

For details look to:
NONLINEAR CONTINUUM MECHANICS FOR FINITE ELEMENT ANALYSIS, Javier Bonet, Richard D. Wood

Definition at line 240 of file NonLinearElasticElement.hpp.

242  {
244  sigmaCauchy.resize(3, 3);
245  t_sigmaCauchy(i, j) = t_P(i, k) * t_F(j, k);
248  }

◆ calculateE_GreenStrain()

template<typename TYPE >
MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::calculateE_GreenStrain ( )
inline

Definition at line 166 of file NonLinearElasticElement.hpp.

166  {
168  constexpr auto t_kd = FTensor::Kronecker_Delta<double>();
169  t_E(i, j) = 0.5 * (t_C(i, j) - t_kd(i, j));
171  }

◆ calculateElasticEnergy()

template<typename TYPE >
virtual MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::calculateElasticEnergy ( const BlockData  block_data,
boost::shared_ptr< const NumeredEntFiniteElement >  fe_ptr 
)
inlinevirtual

Calculate elastic energy density.

\[\Psi = \frac{1}{2}\lambda(\textrm{tr}[\mathbf{E}])^2+\mu\mathbf{E}:\mathbf{E}\]

Reimplemented in SmallStrainTranverslyIsotropic< TYPE >, SmallStrainTranverslyIsotropic< double >, SmallStrainTranverslyIsotropic< adouble >, NeoHookean< TYPE >, and Hooke< TYPE >.

Definition at line 292 of file NonLinearElasticElement.hpp.

294  {
296  lambda = LAMBDA(block_data.E, block_data.PoissonRatio);
297  mu = MU(block_data.E, block_data.PoissonRatio);
300  TYPE trace = 0;
301  eNergy = 0;
302  for (int ii = 0; ii < 3; ii++) {
303  trace += E(ii, ii);
304  for (int jj = 0; jj < 3; jj++) {
305  TYPE e = E(ii, jj);
306  eNergy += mu * e * e;
307  }
308  }
309  eNergy += 0.5 * lambda * trace * trace;
311  }

◆ calculateP_PiolaKirchhoffI()

template<typename TYPE >
virtual MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::calculateP_PiolaKirchhoffI ( const BlockData  block_data,
boost::shared_ptr< const NumeredEntFiniteElement >  fe_ptr 
)
inlinevirtual

Function overload to implement user material.

Calculation of Piola Kirchhoff I is implemented by user. Tangent matrix user implemented physical equation is calculated using automatic differentiation.

\(\mathbf{S} = \lambda\textrm{tr}[\mathbf{E}]\mathbf{I}+2\mu\mathbf{E}\)

Notes:
Number of actual Gauss point is accessed from variable gG.
Access to operator data structures is available by variable opPtr.
Access to common data is by commonDataPtr.

Parameters
block_dataused to give access to material parameters
fe_ptrpointer to element data structures

For details look to:
NONLINEAR CONTINUUM MECHANICS FOR FINITE ELEMENT ANALYSIS, Javier Bonet, Richard D. Wood

Reimplemented in SmallStrainTranverslyIsotropic< TYPE >, SmallStrainTranverslyIsotropic< double >, SmallStrainTranverslyIsotropic< adouble >, VolumeLengthQuality< TYPE >, NeoHookean< TYPE >, and Hooke< TYPE >.

Definition at line 204 of file NonLinearElasticElement.hpp.

206  {
208  lambda = LAMBDA(block_data.E, block_data.PoissonRatio);
209  mu = MU(block_data.E, block_data.PoissonRatio);
213  t_P(i, j) = t_F(i, k) * t_S(k, j);
215  }

◆ calculateS_PiolaKirchhoffII()

template<typename TYPE >
MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::calculateS_PiolaKirchhoffII ( )
inline

Definition at line 174 of file NonLinearElasticElement.hpp.

174  {
176  constexpr auto t_kd = FTensor::Kronecker_Delta<double>();
177  t_S(i, j) = (2 * mu) * t_E(i, j) + (lambda * t_E(k, k)) * t_kd(i, j);
179  }

◆ calculatesIGma_EshelbyStress()

template<typename TYPE >
virtual MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::calculatesIGma_EshelbyStress ( const BlockData  block_data,
boost::shared_ptr< const NumeredEntFiniteElement >  fe_ptr 
)
inlinevirtual

Calculate Eshelby stress.

Definition at line 315 of file NonLinearElasticElement.hpp.

317  {
319  CHKERR calculateP_PiolaKirchhoffI(block_data, fe_ptr);
320  CHKERR calculateElasticEnergy(block_data, fe_ptr);
321  constexpr auto t_kd = FTensor::Kronecker_Delta<double>();
322  t_sIGma(i, j) = t_kd(i, j) * eNergy - t_F(k, i) * t_P(k, j);
324  }

◆ getDataOnPostProcessor()

template<typename TYPE >
virtual MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::getDataOnPostProcessor ( std::map< std::string, std::vector< VectorDouble >> &  field_map,
std::map< std::string, std::vector< MatrixDouble >> &  grad_map 
)
inlinevirtual

Do operations when pre-process.

Definition at line 328 of file NonLinearElasticElement.hpp.

330  {
333  }

◆ resizeAndSet()

template<typename TYPE >
static auto NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::resizeAndSet ( MatrixBoundedArray< TYPE, 9 > &  m)
inlinestaticprotected
Examples
NeoHookean.hpp.

Definition at line 336 of file NonLinearElasticElement.hpp.

336  {
337  m.resize(3, 3, false);
339  };

◆ setUserActiveVariables() [1/2]

template<typename TYPE >
virtual MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::setUserActiveVariables ( int &  nb_active_variables)
inlinevirtual

add additional active variables

Note
This member function if used should be implement by template member function Specialization, different implementation needed for TYPE=double or TYPE=adouble

More complex physical models depend on gradient of defamation and some additional variables. For example can depend on temperature. This function adds additional independent variables to the model.

Parameters
nb_active_variablesnumber of active variables
Returns
error code

Reimplemented in SmallStrainTranverslyIsotropicADouble, and MyMat< TYPE >.

Definition at line 264 of file NonLinearElasticElement.hpp.

264  {
267  }

◆ setUserActiveVariables() [2/2]

template<typename TYPE >
virtual MoFEMErrorCode NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::setUserActiveVariables ( VectorDouble &  activeVariables)
inlinevirtual

Add additional independent variables More complex physical models depend on gradient of defamation and some additional variables. For example can depend on temperature. This function adds additional independent variables to the model.

/note First 9 elements are reserved for gradient of deformation.

Parameters
activeVariablesvector of deepened variables, values after index 9 should be add.
Returns
error code

Reimplemented in SmallStrainTranverslyIsotropicADouble, and MyMat< TYPE >.

Definition at line 282 of file NonLinearElasticElement.hpp.

282  {
285  }

Member Data Documentation

◆ C

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::C
Examples
NeoHookean.hpp.

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ commonDataPtr

template<typename TYPE >
CommonData* NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::commonDataPtr

common data shared between entities (f.e. field values at Gauss pts.)

Definition at line 155 of file NonLinearElasticElement.hpp.

◆ detF

template<typename TYPE >
TYPE NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::detF

Definition at line 152 of file NonLinearElasticElement.hpp.

◆ detH

template<typename TYPE >
TYPE NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::detH

Definition at line 152 of file NonLinearElasticElement.hpp.

◆ E

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::E

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ eNergy

template<typename TYPE >
TYPE NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::eNergy
Examples
NeoHookean.hpp.

Definition at line 152 of file NonLinearElasticElement.hpp.

◆ F

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::F
Examples
NeoHookean.hpp.

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ gG

template<typename TYPE >
int NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::gG

Gauss point number.

Definition at line 154 of file NonLinearElasticElement.hpp.

◆ h

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::h

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ H

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::H

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ i

Definition at line 132 of file NonLinearElasticElement.hpp.

◆ invF

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::invF

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ invH

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::invH

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ j

Definition at line 133 of file NonLinearElasticElement.hpp.

◆ J

template<typename TYPE >
TYPE NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::J
Examples
NeoHookean.hpp.

Definition at line 152 of file NonLinearElasticElement.hpp.

◆ k

Definition at line 134 of file NonLinearElasticElement.hpp.

◆ lambda

template<typename TYPE >
double NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::lambda
Examples
NeoHookean.hpp.

Definition at line 146 of file NonLinearElasticElement.hpp.

◆ mu

Examples
NeoHookean.hpp.

Definition at line 146 of file NonLinearElasticElement.hpp.

◆ opPtr

pointer to finite element tetrahedral operator

Definition at line 158 of file NonLinearElasticElement.hpp.

◆ P

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::P

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ S

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::S

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ sIGma

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::sIGma

Definition at line 147 of file NonLinearElasticElement.hpp.

◆ sigmaCauchy

template<typename TYPE >
MatrixBoundedArray<TYPE, 9> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::sigmaCauchy

Definition at line 148 of file NonLinearElasticElement.hpp.

◆ t_C

Examples
NeoHookean.hpp.

Definition at line 149 of file NonLinearElasticElement.hpp.

◆ t_E

Definition at line 149 of file NonLinearElasticElement.hpp.

◆ t_F

Examples
NeoHookean.hpp.

Definition at line 149 of file NonLinearElasticElement.hpp.

◆ t_h

Definition at line 150 of file NonLinearElasticElement.hpp.

◆ t_H

Definition at line 150 of file NonLinearElasticElement.hpp.

◆ t_invF

template<typename TYPE >
FTensor::Tensor2<FTensor::PackPtr<TYPE *, 0>, 3, 3> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::t_invF

Definition at line 150 of file NonLinearElasticElement.hpp.

◆ t_invH

template<typename TYPE >
FTensor::Tensor2<FTensor::PackPtr<TYPE *, 0>, 3, 3> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::t_invH

Definition at line 150 of file NonLinearElasticElement.hpp.

◆ t_P

Examples
NeoHookean.hpp.

Definition at line 150 of file NonLinearElasticElement.hpp.

◆ t_S

Examples
NeoHookean.hpp.

Definition at line 149 of file NonLinearElasticElement.hpp.

◆ t_sIGma

template<typename TYPE >
FTensor::Tensor2<FTensor::PackPtr<TYPE *, 0>, 3, 3> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::t_sIGma

Definition at line 150 of file NonLinearElasticElement.hpp.

◆ t_sigmaCauchy

template<typename TYPE >
FTensor::Tensor2<FTensor::PackPtr<TYPE *, 0>, 3, 3> NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >::t_sigmaCauchy

Definition at line 150 of file NonLinearElasticElement.hpp.


The documentation for this struct was generated from the following file:
MoFEMFunctionReturnHot
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_H
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_H
Definition: NonLinearElasticElement.hpp:150
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::eNergy
TYPE eNergy
Definition: NonLinearElasticElement.hpp:152
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::calculateP_PiolaKirchhoffI
virtual MoFEMErrorCode calculateP_PiolaKirchhoffI(const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
Function overload to implement user material.
Definition: NonLinearElasticElement.hpp:204
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_P
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_P
Definition: NonLinearElasticElement.hpp:150
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_E
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_E
Definition: NonLinearElasticElement.hpp:149
E
LAMBDA
#define LAMBDA(E, NU)
Definition: fem_tools.h:22
FTensor::Kronecker_Delta
Kronecker Delta class.
Definition: Kronecker_Delta.hpp:15
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::h
MatrixBoundedArray< TYPE, 9 > h
Definition: NonLinearElasticElement.hpp:147
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_sigmaCauchy
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_sigmaCauchy
Definition: NonLinearElasticElement.hpp:150
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::j
FTensor::Index< 'j', 3 > j
Definition: NonLinearElasticElement.hpp:133
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::sIGma
MatrixBoundedArray< TYPE, 9 > sIGma
Definition: NonLinearElasticElement.hpp:147
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::i
FTensor::Index< 'i', 3 > i
Definition: NonLinearElasticElement.hpp:132
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_h
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_h
Definition: NonLinearElasticElement.hpp:150
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_sIGma
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_sIGma
Definition: NonLinearElasticElement.hpp:150
FTensor::Number< 0 >
CHKERR
#define CHKERR
Inline error check.
Definition: definitions.h:535
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_invH
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_invH
Definition: NonLinearElasticElement.hpp:150
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::H
MatrixBoundedArray< TYPE, 9 > H
Definition: NonLinearElasticElement.hpp:147
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::F
MatrixBoundedArray< TYPE, 9 > F
Definition: NonLinearElasticElement.hpp:147
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_F
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_F
Definition: NonLinearElasticElement.hpp:149
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_S
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_S
Definition: NonLinearElasticElement.hpp:149
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::calculateElasticEnergy
virtual MoFEMErrorCode calculateElasticEnergy(const BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
Calculate elastic energy density.
Definition: NonLinearElasticElement.hpp:292
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::C
MatrixBoundedArray< TYPE, 9 > C
Definition: NonLinearElasticElement.hpp:147
t_kd
constexpr auto t_kd
Definition: free_surface.cpp:137
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_C
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_C
Definition: NonLinearElasticElement.hpp:149
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::invH
MatrixBoundedArray< TYPE, 9 > invH
Definition: NonLinearElasticElement.hpp:147
MoFEM::getFTensor2FromArray3by3
auto getFTensor2FromArray3by3(ublas::matrix< T, L, A > &data, const FTensor::Number< S > &, const size_t rr, const size_t cc=0)
Definition: Templates.hpp:1313
MoFEM::determinantTensor3by3
static auto determinantTensor3by3(T &t)
Calculate the determinant of a 3x3 matrix or a tensor of rank 2.
Definition: Templates.hpp:1511
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::S
MatrixBoundedArray< TYPE, 9 > S
Definition: NonLinearElasticElement.hpp:147
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::calculateE_GreenStrain
MoFEMErrorCode calculateE_GreenStrain()
Definition: NonLinearElasticElement.hpp:166
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::mu
double mu
Definition: NonLinearElasticElement.hpp:146
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::resizeAndSet
static auto resizeAndSet(MatrixBoundedArray< TYPE, 9 > &m)
Definition: NonLinearElasticElement.hpp:336
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::sigmaCauchy
MatrixBoundedArray< TYPE, 9 > sigmaCauchy
Definition: NonLinearElasticElement.hpp:148
m
FTensor::Index< 'm', 3 > m
Definition: shallow_wave.cpp:80
PlasticOps::trace
double trace(FTensor::Tensor2_symmetric< T, 2 > &t_stress)
Definition: PlasticOpsGeneric.hpp:80
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::calculateC_CauchyDeformationTensor
MoFEMErrorCode calculateC_CauchyDeformationTensor()
Definition: NonLinearElasticElement.hpp:160
MoFEMFunctionBeginHot
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::lambda
double lambda
Definition: NonLinearElasticElement.hpp:146
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::invF
MatrixBoundedArray< TYPE, 9 > invF
Definition: NonLinearElasticElement.hpp:147
MoFEMFunctionReturn
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::calculateS_PiolaKirchhoffII
MoFEMErrorCode calculateS_PiolaKirchhoffII()
Definition: NonLinearElasticElement.hpp:174
MoFEMFunctionBegin
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::P
MatrixBoundedArray< TYPE, 9 > P
Definition: NonLinearElasticElement.hpp:147
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::t_invF
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_invF
Definition: NonLinearElasticElement.hpp:150
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::E
MatrixBoundedArray< TYPE, 9 > E
Definition: NonLinearElasticElement.hpp:147
MU
#define MU(E, NU)
Definition: fem_tools.h:23
NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI::k
FTensor::Index< 'k', 3 > k
Definition: NonLinearElasticElement.hpp:134