v0.13.2
Loading...
Searching...
No Matches
NeoHookean.hpp
Go to the documentation of this file.
1/** \file NeoHookean.hpp
2 * \ingroup nonlinear_elastic_elem
3 * \brief Implementation of Neo-Hookean elastic material
4 * \example NeoHookean.hpp
5 */
6
7
8
9#ifndef __NEOHOOKEAN_HPP__
10#define __NEOHOOKEAN_HPP__
11
12/** \brief NeoHookan equation
13 * \ingroup nonlinear_elastic_elem
14 */
15template <typename TYPE>
18 TYPE> {
19
22
26
27 TYPE detC;
28 TYPE logJ;
29 MatrixBoundedArray<TYPE, 9> invC;
31
35
36 /** \brief calculate second Piola Kirchhoff
37 *
38 * \f$\mathbf{S} =
39 \mu(\mathbf{I}-\mathbf{C}^{-1})+\lambda(\ln{J})\mathbf{C}^{-1}\f$
40
41 For details look to: <br>
42 NONLINEAR CONTINUUM MECHANICS FOR FINITE ELEMENT ANALYSIS, Javier Bonet,
43 Richard D. Wood
44
45 */
46 MoFEMErrorCode NeoHooke_PiolaKirchhoffII() {
48
49 detC = determinantTensor3by3(this->C);
50 CHKERR invertTensor3by3(this->C, detC, invC);
51 this->J = determinantTensor3by3(this->F);
52
53 logJ = log(sqrt(this->J * this->J));
54 constexpr auto t_kd = FTensor::Kronecker_Delta<double>();
55
56 this->t_S(i, j) = this->mu * (t_kd(i, j) - t_invC(i, j)) +
57 (this->lambda * logJ) * t_invC(i, j);
58
60 }
61
62 virtual MoFEMErrorCode calculateP_PiolaKirchhoffI(
64 boost::shared_ptr<const NumeredEntFiniteElement> fe_ptr) {
66
67 this->lambda = LAMBDA(block_data.E, block_data.PoissonRatio);
68 this->mu = MU(block_data.E, block_data.PoissonRatio);
71
72 this->t_P(i, j) = this->t_F(i, k) * this->t_S(k, j);
73
75 }
76
77 /** \brief calculate elastic energy density
78 *
79
80 For details look to: <br>
81 NONLINEAR CONTINUUM MECHANICS FOR FINITE ELEMENT ANALYSIS, Javier Bonet,
82 Richard D. Wood
83
84 */
85 MoFEMErrorCode NeoHookean_ElasticEnergy() {
87 this->eNergy = this->t_C(i, i);
88 this->eNergy = 0.5 * this->mu * (this->eNergy - 3);
89 logJ = log(sqrt(this->J * this->J));
90 this->eNergy += -this->mu * logJ + 0.5 * this->lambda * pow(logJ, 2);
92 }
93
94 MoFEMErrorCode calculateElasticEnergy(
96 boost::shared_ptr<const NumeredEntFiniteElement> fe_ptr) {
98
99 this->lambda = LAMBDA(block_data.E, block_data.PoissonRatio);
100 this->mu = MU(block_data.E, block_data.PoissonRatio);
102 this->J = determinantTensor3by3(this->F);
104
106 }
107};
108
109#endif
Kronecker Delta class.
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
#define MU(E, NU)
Definition: fem_tools.h:23
#define LAMBDA(E, NU)
Definition: fem_tools.h:22
constexpr auto t_kd
NeoHookan equation.
Definition: NeoHookean.hpp:18
MoFEMErrorCode NeoHooke_PiolaKirchhoffII()
calculate second Piola Kirchhoff
Definition: NeoHookean.hpp:46
MoFEMErrorCode NeoHookean_ElasticEnergy()
calculate elastic energy density
Definition: NeoHookean.hpp:85
MoFEMErrorCode calculateElasticEnergy(const NonlinearElasticElement::BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
Calculate elastic energy density.
Definition: NeoHookean.hpp:94
FTensor::Index< 'i', 3 > i
Definition: NeoHookean.hpp:32
MatrixBoundedArray< TYPE, 9 > invC
Definition: NeoHookean.hpp:29
virtual MoFEMErrorCode calculateP_PiolaKirchhoffI(const NonlinearElasticElement::BlockData block_data, boost::shared_ptr< const NumeredEntFiniteElement > fe_ptr)
Function overload to implement user material.
Definition: NeoHookean.hpp:62
FTensor::Index< 'j', 3 > j
Definition: NeoHookean.hpp:33
FTensor::Index< 'k', 3 > k
Definition: NeoHookean.hpp:34
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_invC
Definition: NeoHookean.hpp:30
data for calculation heat conductivity and heat capacity elements
Implementation of elastic (non-linear) St. Kirchhoff equation.
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_C
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_F
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_P
FTensor::Tensor2< FTensor::PackPtr< TYPE *, 0 >, 3, 3 > t_S
static auto resizeAndSet(MatrixBoundedArray< TYPE, 9 > &m)
structure grouping operators and data used for calculation of nonlinear elastic element