LifeV
NeoHookean.hpp
Go to the documentation of this file.
1 #ifndef NEOHOOKEAN_H
2 #define NEOHOOKEAN_H 1
3 
4 //@HEADER
5 /*
6  *******************************************************************************
7 
8  Copyright (C) 2004, 2005, 2007 EPFL, Politecnico di Milano, INRIA
9  Copyright (C) 2010 EPFL, Politecnico di Milano, Emory University
10 
11  This file is part of LifeV.
12 
13  LifeV is free software; you can redistribute it and/or modify
14  it under the terms of the GNU Lesser General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  LifeV is distributed in the hope that it will be useful,
19  but WITHOUT ANY WARRANTY; without even the implied warranty of
20  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
21  Lesser General Public License for more details.
22 
23  You should have received a copy of the GNU Lesser General Public License
24  along with LifeV. If not, see <http://www.gnu.org/licenses/>.
25 
26  *******************************************************************************
27  */
28 //@HEADER
29 /*!
30  @file
31  @brief Implementation of Neohookean structural model
32  @author Davide Forti <davide.forti@epfl.ch>
33  @maintainer Antonello Gerbi <antonello.gerbi@epfl.ch>
34  @date 25-05-2016
35 
36  This file contains an ETA implementation of a neohookean structural model. This class
37  temporary as it will disappear when the new structure model will be released.
38 
39 */
40 
41 #include <lifev/core/LifeV.hpp>
42 
43 #include <lifev/core/array/VectorEpetra.hpp>
44 
45 #include <lifev/core/array/MatrixEpetra.hpp>
46 
47 #include <lifev/core/fem/BCManage.hpp>
48 
49 #include <lifev/eta/fem/ETFESpace.hpp>
50 
51 #include <lifev/core/fem/FESpace.hpp>
52 
53 #include <lifev/eta/expression/Integrate.hpp>
54 
55 #include <lifev/fsi_blocks/solver/AssemblyElementalStructure.hpp>
56 
57 #include <lifev/fsi_blocks/solver/ExpressionDefinitions.hpp>
58 
59 namespace LifeV
60 {
61 
62 class NeoHookean
63 {
64 
65  typedef Epetra_Comm comm_Type;
66 
67  typedef std::shared_ptr< comm_Type > commPtr_Type;
68 
69  typedef VectorEpetra vector_Type;
70 
71  typedef std::shared_ptr<vector_Type> vectorPtr_Type;
72 
73  typedef MatrixEpetra<Real> matrix_Type;
74 
75  typedef std::shared_ptr<matrix_Type> matrixPtr_Type;
76 
77  typedef std::shared_ptr<BCHandler> bcPtr_Type;
78 
79  typedef RegionMesh<LinearTetra> mesh_Type;
80 
81  typedef std::shared_ptr<mesh_Type> meshPtr_Type;
82 
83  typedef MapEpetra map_Type;
84 
85  typedef ETFESpace< mesh_Type, map_Type, 3, 3 > ETFESpace_displacement;
86 
88 
89  typedef ExpressionDefinitions::determinantTensorF_Type determinantF_Type;
90 
92 
94 
95  typedef ExpressionDefinitions::traceTensor_Type traceTensor_Type;
96 
97 public:
98 
99  // Constructor
100  /*!
101  * @param communicator Epetra communicator
102  */
103  NeoHookean( const commPtr_Type& communicator );
104 
105  // empty destructor
106  ~NeoHookean();
107 
108  //! Set the physical parameters of the structure
109  /*!
110  * @param density Density of the structure
111  * @param young Young modulus of the structure
112  * @param poisson Poisson ratio of the structure
113  */
114  void setCoefficients ( const Real density, const Real young, const Real poisson);
115 
116  //! Setup: build the FE space and the ETA fespace
117  /*!
118  * @param mesh Computational mesh
119  * @param dOrder Degree of the finite element used
120  */
121  void setup( const meshPtr_Type& mesh, const std::string dOrder );
122 
123  //! Evaluate the residual of the problem
124  /*!
125  * @param solution Solution vector at the previous Newton iterate
126  * @param coefficient Coefficient which goes in front of the mass matrix depending on the time discretization scheme used
127  * @param csi Vector containing the solution at the previous time steps which depends on the time discretization scheme used
128  * @param residual Residual vector to be assembled
129  */
130  void evaluate_residual( const vectorPtr_Type& solution, const Real& coefficient, const vectorPtr_Type& csi, vectorPtr_Type& residual );
131 
132  //! Updates the Jacobian matrix
133  /*!
134  * @param solution Solution vector at the previous Newton iterate
135  * @param coefficient Coefficient which goes in front of the mass matrix depending on the time discretization scheme used
136  * @param jacobian Jacobian matrix to be updated
137  */
138  void update_jacobian(const vectorPtr_Type& solution, const Real& coefficient, matrixPtr_Type& jacobian );
139 
140  //! Getter of the standard FE space
141  /*!
142  * @return M_displacementFESpace FE space used for the displacement
143  *
144  */
145  const std::shared_ptr<FESpace<mesh_Type, map_Type> >& fespace ( ) const { return M_displacementFESpace; };
146 
147  //! Getter of the ETA FE space
148  /*!
149  * @return M_displacementFESpace_ETA ETA FE space used for the displacement
150  *
151  */
152  const std::shared_ptr<ETFESpace_displacement >& et_fespace ( ) const { return M_displacementFESpace_ETA; };
153 
154 private:
155 
156  // communicator
157  commPtr_Type M_comm;
158 
159  // Coefficients for the structure
160  Real M_density;
161 
162  // Lame coefficients computed on M_young and M_poisson
163  Real M_bulk;
164  Real M_mu;
165  Real M_offset;
166 
167  matrixSmall_Type M_identity;
168 
169  // order FE space displacement
170  std::string M_dOrder;
171 
172  // FE space
173  std::shared_ptr<FESpace<mesh_Type, map_Type> > M_displacementFESpace;
174 
175  // ET FE Space
176  std::shared_ptr<ETFESpace_displacement > M_displacementFESpace_ETA;
177 
178  // jacobian matrix
179  matrixPtr_Type M_jacobian;
180 
181 };
182 
183 } // end namespace LifeV
184 
185 #endif
ExpressionMinusTransposed< ExpressionAddition< ExpressionInterpolateGradient< MeshType, MapEpetra, 3, 3 >, ExpressionMatrix< 3, 3 > > > minusTransposedTensor_Type
void importFromHDF5(std::string const &fileName, std::string const &matrixName="matrix")
Read a matrix from a HDF5 (.h5) file.
MatrixSmall< 3, 3 > matrixSmall_Type
Epetra_Import const & importer()
Getter for the Epetra_Import.
Definition: MapEpetra.cpp:394
ExpressionProduct< ExpressionTranspose< ExpressionAddition< ExpressionInterpolateGradient< MeshType, MapEpetra, 3, 3 >, ExpressionMatrix< 3, 3 > > >, ExpressionAddition< ExpressionInterpolateGradient< MeshType, MapEpetra, 3, 3 >, ExpressionMatrix< 3, 3 > > > rightCauchyGreenTensor_Type
double Real
Generic real data.
Definition: LifeV.hpp:175
ExpressionAddition< ExpressionInterpolateGradient< MeshType, MapEpetra, 3, 3 >, ExpressionMatrix< 3, 3 > > deformationGradient_Type
ExpressionTrace< ExpressionProduct< ExpressionTranspose< ExpressionAddition< ExpressionInterpolateGradient< MeshType, MapEpetra, 3, 3 >, ExpressionMatrix< 3, 3 > > >, ExpressionAddition< ExpressionInterpolateGradient< MeshType, MapEpetra, 3, 3 >, ExpressionMatrix< 3, 3 > > > > traceTensor_Type
ExpressionDeterminant< ExpressionAddition< ExpressionInterpolateGradient< MeshType, MapEpetra, 3, 3 >, ExpressionMatrix< 3, 3 > > > determinantTensorF_Type
class ETFESpace A light, templated version of the FESpace
Definition: ETFESpace.hpp:93