libMesh
Public Member Functions | Public Attributes | Private Member Functions | Static Private Member Functions | Private Attributes | List of all members
NonlinearNeoHookeCurrentConfig Class Reference

This class implements a constitutive formulation for an Neo-Hookean elastic solid in terms of the current configuration. More...

#include <nonlinear_neohooke_cc.h>

Public Member Functions

 NonlinearNeoHookeCurrentConfig (const std::vector< std::vector< RealGradient >> &dphi_in, GetPot &args, bool calculate_linearized_stiffness_in)
 
void init_for_qp (VectorValue< Gradient > &grad_u, unsigned int qp)
 Initialize the class for the given displacement gradient at the specified quadrature point. More...
 
void get_residual (DenseVector< Real > &residuum, unsigned int &i)
 Return the residual vector for the current state. More...
 
void get_linearized_stiffness (DenseMatrix< Real > &stiffness, unsigned int &i, unsigned int &j)
 Return the stiffness matrix for the current state. More...
 

Public Attributes

bool calculate_linearized_stiffness
 Flag to indicate if it is necessary to calculate values for stiffness matrix during initialization. More...
 

Private Member Functions

void build_b_0_mat (int i, DenseMatrix< Real > &b_l_mat)
 
void calculate_stress ()
 
void calculate_tangent ()
 

Static Private Member Functions

static void tensor_to_voigt (const RealTensor &tensor, DenseVector< Real > &vec)
 

Private Attributes

unsigned int current_qp
 
const std::vector< std::vector< RealGradient > > & dphi
 
DenseMatrix< Real > C_mat
 
Real E
 
Real nu
 
RealTensor F
 
RealTensor S
 
RealTensor tau
 
RealTensor sigma
 
DenseMatrix< Real > B_L
 
DenseMatrix< Real > B_K
 

Detailed Description

This class implements a constitutive formulation for an Neo-Hookean elastic solid in terms of the current configuration.

This implementation is suitable for computing large deformations. See e.g. "Nonlinear finite element methods" (P. Wriggers, 2008, Springer) for details.

Definition at line 40 of file nonlinear_neohooke_cc.h.

Constructor & Destructor Documentation

NonlinearNeoHookeCurrentConfig::NonlinearNeoHookeCurrentConfig ( const std::vector< std::vector< RealGradient >> &  dphi_in,
GetPot &  args,
bool  calculate_linearized_stiffness_in 
)

Definition at line 43 of file nonlinear_neohooke_cc.h.

45  :
46  calculate_linearized_stiffness(calculate_linearized_stiffness_in),
47  dphi(dphi_in)
48  {
49  E = args("material/neohooke/e_modulus", 10000.0);
50  nu = args("material/neohooke/nu", 0.3);
51  }
bool calculate_linearized_stiffness
Flag to indicate if it is necessary to calculate values for stiffness matrix during initialization...
const std::vector< std::vector< RealGradient > > & dphi

Member Function Documentation

void NonlinearNeoHookeCurrentConfig::build_b_0_mat ( int  i,
DenseMatrix< Real > &  b_l_mat 
)
private

Definition at line 141 of file nonlinear_neohooke_cc.C.

References current_qp, and dphi.

Referenced by get_linearized_stiffness(), and get_residual().

142 {
143  for (unsigned int ii = 0; ii < 3; ++ii)
144  b_0_mat(ii, ii) = dphi[i][current_qp](ii);
145 
146  b_0_mat(0, 3) = dphi[i][current_qp](1);
147  b_0_mat(1, 3) = dphi[i][current_qp](0);
148  b_0_mat(1, 4) = dphi[i][current_qp](2);
149  b_0_mat(2, 4) = dphi[i][current_qp](1);
150  b_0_mat(0, 5) = dphi[i][current_qp](2);
151  b_0_mat(2, 5) = dphi[i][current_qp](0);
152 }
const std::vector< std::vector< RealGradient > > & dphi
void NonlinearNeoHookeCurrentConfig::calculate_stress ( )
private

Definition at line 73 of file nonlinear_neohooke_cc.C.

References E, F, libMesh::TypeTensor< T >::inverse(), nu, libMesh::Real, S, sigma, and tau.

Referenced by init_for_qp().

74 {
75  double mu = E / (2.0 * (1.0 + nu));
76  double lambda = E * nu / ((1 + nu) * (1 - 2 * nu));
77 
78  Real detF = F.det();
79  RealTensor Ft = F.transpose();
80 
81  RealTensor C = Ft * F;
82  RealTensor b = F * Ft;
83  RealTensor identity;
84  identity(0, 0) = 1.0; identity(1, 1) = 1.0; identity(2, 2) = 1.0;
85  RealTensor invC = C.inverse();
86 
87  S = 0.5 * lambda * (detF * detF - 1) * invC + mu * (identity - invC);
88 
89  tau = (F * S) * Ft;
90  sigma = 1/detF * tau;
91 }
RealTensorValue RealTensor
TypeTensor< T > inverse() const
Definition: type_tensor.h:1022
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void NonlinearNeoHookeCurrentConfig::calculate_tangent ( )
private

Definition at line 52 of file nonlinear_neohooke_cc.C.

References C_mat, E, F, nu, libMesh::Real, and libMesh::DenseMatrix< T >::resize().

Referenced by init_for_qp().

53 {
54  Real mu = E / (2 * (1 + nu));
55  Real lambda = E * nu / ((1 + nu) * (1 - 2 * nu));
56 
57  Real detF = F.det();
58 
59  C_mat.resize(6, 6);
60  for (unsigned int i = 0; i < 3; ++i) {
61  for (unsigned int j = 0; j < 3; ++j) {
62  if (i == j) {
63  C_mat(i, j) = 2 * mu + lambda;
64  C_mat(i + 3, j + 3) = mu - 0.5 * lambda * (detF * detF - 1);
65  } else {
66  C_mat(i, j) = lambda * detF * detF;
67  }
68  }
69  }
70 }
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void resize(const unsigned int new_m, const unsigned int new_n)
Resize the matrix.
Definition: dense_matrix.h:776
void NonlinearNeoHookeCurrentConfig::get_linearized_stiffness ( DenseMatrix< Real > &  stiffness,
unsigned int i,
unsigned int j 
)

Return the stiffness matrix for the current state.

Definition at line 118 of file nonlinear_neohooke_cc.C.

References B_K, B_L, build_b_0_mat(), C_mat, current_qp, dphi, F, libMesh::DenseMatrix< T >::resize(), libMesh::DenseMatrix< T >::right_multiply(), libMesh::DenseMatrix< T >::right_multiply_transpose(), and sigma.

Referenced by SolidSystem::element_time_derivative().

121 {
122  stiffness.resize(3, 3);
123 
124  double G_IK = (sigma * dphi[i][current_qp]) * dphi[j][current_qp];
125  stiffness(0, 0) += G_IK;
126  stiffness(1, 1) += G_IK;
127  stiffness(2, 2) += G_IK;
128 
129  B_L.resize(3, 6);
130  this->build_b_0_mat(i, B_L);
131  B_K.resize(3, 6);
132  this->build_b_0_mat(j, B_K);
133 
136  B_L *= 1/F.det();
137 
138  stiffness += B_L;
139 }
void build_b_0_mat(int i, DenseMatrix< Real > &b_l_mat)
virtual void right_multiply(const DenseMatrixBase< T > &M2) libmesh_override
Performs the operation: (*this) <- (*this) * M3.
Definition: dense_matrix.C:210
void right_multiply_transpose(const DenseMatrix< T > &A)
Right multiplies by the transpose of the matrix A.
Definition: dense_matrix.C:257
void resize(const unsigned int new_m, const unsigned int new_n)
Resize the matrix.
Definition: dense_matrix.h:776
const std::vector< std::vector< RealGradient > > & dphi
void NonlinearNeoHookeCurrentConfig::get_residual ( DenseVector< Real > &  residuum,
unsigned int i 
)

Return the residual vector for the current state.

Definition at line 93 of file nonlinear_neohooke_cc.C.

References B_L, build_b_0_mat(), libMesh::DenseMatrix< T >::resize(), sigma, tensor_to_voigt(), and libMesh::DenseMatrix< T >::vector_mult().

Referenced by SolidSystem::element_time_derivative().

95 {
96  B_L.resize(3, 6);
97  DenseVector<Real> sigma_voigt(6);
98 
99  this->build_b_0_mat(i, B_L);
100 
101  tensor_to_voigt(sigma, sigma_voigt);
102 
103  B_L.vector_mult(residuum, sigma_voigt);
104 }
void build_b_0_mat(int i, DenseMatrix< Real > &b_l_mat)
static void tensor_to_voigt(const RealTensor &tensor, DenseVector< Real > &vec)
void resize(const unsigned int new_m, const unsigned int new_n)
Resize the matrix.
Definition: dense_matrix.h:776
void vector_mult(DenseVector< T > &dest, const DenseVector< T > &arg) const
Performs the matrix-vector multiplication, dest := (*this) * arg.
Definition: dense_matrix.C:382
void NonlinearNeoHookeCurrentConfig::init_for_qp ( VectorValue< Gradient > &  grad_u,
unsigned int  qp 
)

Initialize the class for the given displacement gradient at the specified quadrature point.

Definition at line 25 of file nonlinear_neohooke_cc.C.

References calculate_linearized_stiffness, calculate_stress(), calculate_tangent(), current_qp, F, libMesh::libmesh_real(), S, libMesh::TOLERANCE, and libMesh::TypeTensor< T >::zero().

Referenced by SolidSystem::element_time_derivative().

27 {
28  this->current_qp = qp;
29  F.zero();
30  S.zero();
31 
32  {
33  RealTensor invF;
34  invF.zero();
35  for (unsigned int i = 0; i < 3; ++i)
36  for (unsigned int j = 0; j < 3; ++j)
37  invF(i, j) += libmesh_real(grad_u(i)(j));
38 
39  F.add(invF.inverse());
40 
41  libmesh_assert_greater (F.det(), -TOLERANCE);
42  }
43 
45  this->calculate_tangent();
46 
47  this->calculate_stress();
48 }
T libmesh_real(T a)
bool calculate_linearized_stiffness
Flag to indicate if it is necessary to calculate values for stiffness matrix during initialization...
void zero()
Set all entries of the tensor to 0.
Definition: type_tensor.h:1248
RealTensorValue RealTensor
static const Real TOLERANCE
void NonlinearNeoHookeCurrentConfig::tensor_to_voigt ( const RealTensor &  tensor,
DenseVector< Real > &  vec 
)
staticprivate

Definition at line 106 of file nonlinear_neohooke_cc.C.

Referenced by get_residual().

108 {
109  vec(0) = tensor(0, 0);
110  vec(1) = tensor(1, 1);
111  vec(2) = tensor(2, 2);
112  vec(3) = tensor(0, 1);
113  vec(4) = tensor(1, 2);
114  vec(5) = tensor(0, 2);
115 
116 }

Member Data Documentation

DenseMatrix<Real> NonlinearNeoHookeCurrentConfig::B_K
private

Definition at line 93 of file nonlinear_neohooke_cc.h.

Referenced by get_linearized_stiffness().

DenseMatrix<Real> NonlinearNeoHookeCurrentConfig::B_L
private

Definition at line 92 of file nonlinear_neohooke_cc.h.

Referenced by get_linearized_stiffness(), and get_residual().

DenseMatrix<Real> NonlinearNeoHookeCurrentConfig::C_mat
private

Definition at line 88 of file nonlinear_neohooke_cc.h.

Referenced by calculate_tangent(), and get_linearized_stiffness().

bool NonlinearNeoHookeCurrentConfig::calculate_linearized_stiffness

Flag to indicate if it is necessary to calculate values for stiffness matrix during initialization.

Definition at line 77 of file nonlinear_neohooke_cc.h.

Referenced by init_for_qp().

unsigned int NonlinearNeoHookeCurrentConfig::current_qp
private

Definition at line 85 of file nonlinear_neohooke_cc.h.

Referenced by build_b_0_mat(), get_linearized_stiffness(), and init_for_qp().

const std::vector<std::vector<RealGradient> >& NonlinearNeoHookeCurrentConfig::dphi
private

Definition at line 86 of file nonlinear_neohooke_cc.h.

Referenced by build_b_0_mat(), and get_linearized_stiffness().

Real NonlinearNeoHookeCurrentConfig::E
private

Definition at line 89 of file nonlinear_neohooke_cc.h.

Referenced by calculate_stress(), and calculate_tangent().

RealTensor NonlinearNeoHookeCurrentConfig::F
private
Real NonlinearNeoHookeCurrentConfig::nu
private

Definition at line 90 of file nonlinear_neohooke_cc.h.

Referenced by calculate_stress(), and calculate_tangent().

RealTensor NonlinearNeoHookeCurrentConfig::S
private

Definition at line 91 of file nonlinear_neohooke_cc.h.

Referenced by calculate_stress(), and init_for_qp().

RealTensor NonlinearNeoHookeCurrentConfig::sigma
private
RealTensor NonlinearNeoHookeCurrentConfig::tau
private

Definition at line 91 of file nonlinear_neohooke_cc.h.

Referenced by calculate_stress().


The documentation for this class was generated from the following files: