libMesh
nonlinear_neohooke_cc.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 // \author Robert Weidlich
21 // \date Copyright 2012
22 
23 #ifndef NONLINEAR_NEOHOOKE_CC_H_
24 #define NONLINEAR_NEOHOOKE_CC_H_
25 
26 #include "libmesh/dense_vector.h"
27 #include "libmesh/dense_matrix.h"
28 #include "libmesh/vector_value.h"
29 #include "libmesh/tensor_value.h"
30 #include "libmesh/getpot.h"
31 
32 using namespace libMesh;
33 
41 {
42 public:
43  NonlinearNeoHookeCurrentConfig(const std::vector<std::vector<RealGradient>> & dphi_in,
44  GetPot & args,
45  bool calculate_linearized_stiffness_in) :
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  }
52 
57  void init_for_qp(VectorValue<Gradient> & grad_u,
58  unsigned int qp);
59 
63  void get_residual(DenseVector<Real> & residuum,
64  unsigned int & i);
65 
69  void get_linearized_stiffness(DenseMatrix<Real> & stiffness,
70  unsigned int & i,
71  unsigned int & j);
72 
78 private:
79  void build_b_0_mat(int i, DenseMatrix<Real> & b_l_mat);
80  void calculate_stress();
81  void calculate_tangent();
82  static void tensor_to_voigt(const RealTensor & tensor,
83  DenseVector<Real> & vec);
84 
85  unsigned int current_qp;
86  const std::vector<std::vector<RealGradient>> & dphi;
87 
91  RealTensor F, S, tau, sigma;
94 };
95 
96 #endif // NONLINEAR_NEOHOOKE_CC_H_
bool calculate_linearized_stiffness
Flag to indicate if it is necessary to calculate values for stiffness matrix during initialization...
This class implements a constitutive formulation for an Neo-Hookean elastic solid in terms of the cur...
This class defines a vector in LIBMESH_DIM dimensional Real or Complex space.
The libMesh namespace provides an interface to certain functionality in the library.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
NonlinearNeoHookeCurrentConfig(const std::vector< std::vector< RealGradient >> &dphi_in, GetPot &args, bool calculate_linearized_stiffness_in)
Defines a dense vector for use in Finite Element-type computations.
const std::vector< std::vector< RealGradient > > & dphi
This class defines a tensor in LIBMESH_DIM dimensional Real or Complex space.