v0.12.1
Public Member Functions | Public Attributes | List of all members
NeoHookean< TYPE > Struct Template Reference

NeoHookan equation. More...

#include <users_modules/basic_finite_elements/nonlinear_elastic_materials/src/NeoHookean.hpp>

Inheritance diagram for NeoHookean< TYPE >:
[legend]
Collaboration diagram for NeoHookean< TYPE >:
[legend]

Public Member Functions

 NeoHookean ()
 
virtual MoFEMErrorCode NeoHooke_PiolaKirchhoffII ()
 calculate second Piola Kirchhoff More...
 
virtual MoFEMErrorCode calculateP_PiolaKirchhoffI (const NonlinearElasticElement::BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Function overload to implement user material. More...
 
virtual MoFEMErrorCode NeoHookean_ElasticEnergy ()
 calculate elastic energy density More...
 
MoFEMErrorCode calculateElasticEnergy (const NonlinearElasticElement::BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
 Calculate elastic energy density. More...
 
- Public Member Functions inherited from NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >
virtual ~FunctionsToCalculatePiolaKirchhoffI ()
 
MoFEMErrorCode dEterminant (ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 >> &a, TYPE &det)
 Calculate determinant of 3x3 matrix. More...
 
MoFEMErrorCode iNvert (TYPE det, ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 >> &a, ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 >> &inv_a)
 Calculate inverse of 3x3 matrix. More...
 
MoFEMErrorCode calculateC_CauchyDeformationTensor ()
 
MoFEMErrorCode calculateE_GreenStrain ()
 
MoFEMErrorCode calculateS_PiolaKirchhoffII ()
 
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 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

TYPE detC
 
TYPE logJ
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > invC
 
- Public Attributes inherited from NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >
double lambda
 
double mu
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > F
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > C
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > E
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > S
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > invF
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > P
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > SiGma
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > h
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > H
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > invH
 
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > 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...
 

Detailed Description

template<typename TYPE>
struct NeoHookean< TYPE >

NeoHookan equation.

Examples
NeoHookean.hpp.

Definition at line 28 of file NeoHookean.hpp.

Constructor & Destructor Documentation

◆ NeoHookean()

template<typename TYPE >
NeoHookean< TYPE >::NeoHookean ( )
Examples
NeoHookean.hpp.

Definition at line 30 of file NeoHookean.hpp.

Member Function Documentation

◆ calculateElasticEnergy()

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

Calculate elastic energy density.

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

Reimplemented from NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >.

Examples
NeoHookean.hpp.

Definition at line 102 of file NeoHookean.hpp.

105  {
107 
108  this->lambda = LAMBDA(block_data.E,block_data.PoissonRatio);
109  this->mu = MU(block_data.E,block_data.PoissonRatio);
111  ierr = this->dEterminant(this->F,this->J); CHKERRG(ierr);
114  }
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:440
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:476
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:433
#define MU(E, NU)
Definition: fem_tools.h:33
#define LAMBDA(E, NU)
Definition: fem_tools.h:32
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
virtual MoFEMErrorCode NeoHookean_ElasticEnergy()
calculate elastic energy density
Definition: NeoHookean.hpp:90
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > F
MoFEMErrorCode dEterminant(ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 >> &a, TYPE &det)
Calculate determinant of 3x3 matrix.

◆ calculateP_PiolaKirchhoffI()

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

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 from NonlinearElasticElement::FunctionsToCalculatePiolaKirchhoffI< TYPE >.

Examples
NeoHookean.hpp.

Definition at line 66 of file NeoHookean.hpp.

69  {
71 
72  this->lambda = LAMBDA(block_data.E,block_data.PoissonRatio);
73  this->mu = MU(block_data.E,block_data.PoissonRatio);
76  this->P.resize(3,3);
77  noalias(this->P) = prod(this->F,this->S);
78  //std::cerr << "P: " << P << std::endl;
80  }
virtual MoFEMErrorCode NeoHooke_PiolaKirchhoffII()
calculate second Piola Kirchhoff
Definition: NeoHookean.hpp:45
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > S
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > P

◆ NeoHooke_PiolaKirchhoffII()

template<typename TYPE >
virtual MoFEMErrorCode NeoHookean< TYPE >::NeoHooke_PiolaKirchhoffII ( )
virtual

calculate second Piola Kirchhoff

\(\mathbf{S} = \mu(\mathbf{I}-\mathbf{C}^{-1})+\lambda(\ln{J})\mathbf{C}^{-1}\)

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

Examples
NeoHookean.hpp.

Definition at line 45 of file NeoHookean.hpp.

45  {
47 
48  invC.resize(3,3);
49  this->S.resize(3,3);
50  ierr = this->dEterminant(this->C,detC); CHKERRG(ierr);
51  ierr = this->iNvert(detC,this->C,invC); CHKERRG(ierr);
52  ierr = this->dEterminant(this->F,this->J); CHKERRG(ierr);
53  // if(this->J<=0) {
54  // cerr << this->J << endl;
55  // cerr << this->F << endl;
56  // }
57  logJ = log(sqrt(this->J*this->J));
58  for(int i = 0;i<3;i++) {
59  for(int j = 0;j<3;j++) {
60  this->S(i,j) = this->mu*( ((i==j) ? 1 : 0) - invC(i,j) ) + this->lambda*logJ*invC(i,j);
61  }
62  }
64  }
FTensor::Index< 'j', 3 > j
FTensor::Index< 'i', 3 > i
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > invC
Definition: NeoHookean.hpp:34
MoFEMErrorCode iNvert(TYPE det, ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 >> &a, ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 >> &inv_a)
Calculate inverse of 3x3 matrix.
ublas::matrix< TYPE, ublas::row_major, ublas::bounded_array< TYPE, 9 > > C

◆ NeoHookean_ElasticEnergy()

template<typename TYPE >
virtual MoFEMErrorCode NeoHookean< TYPE >::NeoHookean_ElasticEnergy ( )
virtual

calculate elastic energy density

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

Examples
NeoHookean.hpp.

Definition at line 90 of file NeoHookean.hpp.

90  {
92  this->eNergy = 0;
93  for(int ii = 0;ii<3;ii++) {
94  this->eNergy += this->C(ii,ii);
95  }
96  this->eNergy = 0.5*this->mu*(this->eNergy-3);
97  logJ = log(sqrt(this->J*this->J));
98  this->eNergy += -this->mu*logJ + 0.5*this->lambda*pow(logJ,2);
100  }

Member Data Documentation

◆ detC

template<typename TYPE >
TYPE NeoHookean< TYPE >::detC
Examples
NeoHookean.hpp.

Definition at line 32 of file NeoHookean.hpp.

◆ invC

template<typename TYPE >
ublas::matrix<TYPE,ublas::row_major,ublas::bounded_array<TYPE,9> > NeoHookean< TYPE >::invC
Examples
NeoHookean.hpp.

Definition at line 34 of file NeoHookean.hpp.

◆ logJ

template<typename TYPE >
TYPE NeoHookean< TYPE >::logJ
Examples
NeoHookean.hpp.

Definition at line 33 of file NeoHookean.hpp.


The documentation for this struct was generated from the following file: