KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
geo_structural_base_element.hpp
Go to the documentation of this file.
1 // KRATOS___
2 // // ) )
3 // // ___ ___
4 // // ____ //___) ) // ) )
5 // // / / // // / /
6 // ((____/ / ((____ ((___/ / MECHANICS
7 //
8 // License: geo_mechanics_application/license.txt
9 //
10 // Main authors: Vahid Galavi
11 //
12 
13 #if !defined(KRATOS_GEO_STRUCTURAL_BASE_ELEMENT_H_INCLUDED)
14 #define KRATOS_GEO_STRUCTURAL_BASE_ELEMENT_H_INCLUDED
15 
16 // Project includes
17 #include "includes/define.h"
18 #include "includes/element.h"
19 #include "includes/variables.h"
20 
21 // Application includes
22 
23 namespace Kratos
24 {
25 
26 template <unsigned int TDim, unsigned int TNumNodes>
27 class KRATOS_API(GEO_MECHANICS_APPLICATION) GeoStructuralBaseElement : public Element
28 {
29 public:
31  typedef std::size_t SizeType;
32 
34 
36 
39 
42  : Element(NewId, ThisNodes)
43  {
44  }
45 
47  GeoStructuralBaseElement(IndexType NewId, GeometryType::Pointer pGeometry)
48  : Element(NewId, pGeometry)
49  {
50  }
51 
53  GeoStructuralBaseElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
54  : Element(NewId, pGeometry, pProperties)
55  {
56  mThisIntegrationMethod = this->GetIntegrationMethod();
57  }
58 
61 
63 
64  Element::Pointer Create(IndexType NewId,
65  NodesArrayType const& ThisNodes,
66  PropertiesType::Pointer pProperties) const override;
67 
68  Element::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, PropertiesType::Pointer pProperties) const override;
69 
70  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
71 
72  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
73 
74  void GetDofList(DofsVectorType& rElementalDofList, const ProcessInfo& rCurrentProcessInfo) const override;
75 
76  GeometryData::IntegrationMethod GetIntegrationMethod() const override;
77 
79 
80  void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix,
81  VectorType& rRightHandSideVector,
82  const ProcessInfo& rCurrentProcessInfo) override;
83 
84  void CalculateLeftHandSide(MatrixType& rLeftHandSideMatrix, const ProcessInfo& rCurrentProcessInfo) override;
85 
86  void CalculateRightHandSide(VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
87 
88  void EquationIdVector(EquationIdVectorType& rResult, const ProcessInfo& rCurrentProcessInfo) const override;
89 
90  void CalculateMassMatrix(MatrixType& rMassMatrix, const ProcessInfo& rCurrentProcessInfo) override;
91 
92  void CalculateDampingMatrix(MatrixType& rDampingMatrix, const ProcessInfo& rCurrentProcessInfo) override;
93 
94  void GetValuesVector(Vector& rValues, int Step = 0) const override;
95 
96  void GetFirstDerivativesVector(Vector& rValues, int Step = 0) const override;
97 
98  void GetSecondDerivativesVector(Vector& rValues, int Step = 0) const override;
99 
101 
102  void SetValuesOnIntegrationPoints(const Variable<double>& rVariable,
103  const std::vector<double>& rValues,
104  const ProcessInfo& rCurrentProcessInfo) override;
105 
107 
108 protected:
109  static constexpr SizeType N_DOF_NODE = (TDim == 2 ? 3 : 6);
110  static constexpr SizeType N_DOF_ELEMENT = N_DOF_NODE * TNumNodes;
111  static constexpr SizeType VoigtSize = (TDim == 3 ? VOIGT_SIZE_3D : VOIGT_SIZE_2D_PLANE_STRESS);
112 
116 
118 
124 
126 
129 
133 
144  double detF;
146 
149  };
150 
152 
153  std::vector<ConstitutiveLaw::Pointer> mConstitutiveLawVector;
154  std::vector<Vector> mStressVector;
155 
156  virtual SizeType GetTotalNumberIntegrationPoints() const;
157  virtual SizeType GetCrossNumberIntegrationPoints() const;
158  virtual SizeType GetAlongNumberIntegrationPoints() const;
159 
160  virtual void InitializeElementVariables(ElementVariables& rVariables,
161  ConstitutiveLaw::Parameters& rConstitutiveParameters,
162  const GeometryType& Geom,
163  const PropertiesType& Prop,
164  const ProcessInfo& rCurrentProcessInfo) const;
165 
166  virtual void GetNodalDofValuesVector(Vector& rNodalVariableVector,
167  const GeometryType& Geom,
168  IndexType SolutionStepIndex = 0) const;
169 
171 
172  virtual void CalculateStiffnessMatrix(MatrixType& rStiffnessMatrix, const ProcessInfo& rCurrentProcessInfo);
173 
174  virtual void CalculateAll(MatrixType& rLeftHandSideMatrix,
175  VectorType& rRightHandSideVector,
176  const ProcessInfo& rCurrentProcessInfo,
177  const bool CalculateStiffnessMatrixFlag,
178  const bool CalculateResidualVectorFlag);
179 
180  virtual void CalculateRHS(VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo);
181 
182  virtual void CalculateNodalCrossDirection(Matrix& NodalCrossDirection) const;
183 
185 
186 private:
189 
192 
194 
195  friend class Serializer;
196 
197  void save(Serializer& rSerializer) const override
198  {
200  }
201 
202  void load(Serializer& rSerializer) override
203  {
205  }
206 
207 }; // Class GeoStructuralBaseElement
208 
209 } // namespace Kratos
210 
211 #endif // KRATOS_GEO_STRUCTURAL_BASE_ELEMENT_H_INCLUDED defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Base class for all Elements.
Definition: element.h:60
std::size_t SizeType
Definition: element.h:94
std::size_t IndexType
Definition: flags.h:74
Definition: geo_structural_base_element.hpp:28
GeoStructuralBaseElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: geo_structural_base_element.hpp:41
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(GeoStructuralBaseElement)
virtual ~GeoStructuralBaseElement()
Destructor.
Definition: geo_structural_base_element.hpp:60
std::vector< Vector > mStressVector
Definition: geo_structural_base_element.hpp:154
GeoStructuralBaseElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using Geometry.
Definition: geo_structural_base_element.hpp:47
std::size_t SizeType
The definition of the sizetype.
Definition: geo_structural_base_element.hpp:31
GeoStructuralBaseElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Constructor using Properties.
Definition: geo_structural_base_element.hpp:53
std::vector< ConstitutiveLaw::Pointer > mConstitutiveLawVector
Definition: geo_structural_base_element.hpp:153
GeometryData::IntegrationMethod mThisIntegrationMethod
Definition: geo_structural_base_element.hpp:151
GeoStructuralBaseElement(IndexType NewId=0)
Default Constructor.
Definition: geo_structural_base_element.hpp:38
IntegrationMethod
Definition: geometry_data.h:76
Geometry base class.
Definition: geometry.h:71
PointerVector is a container like stl vector but using a vector to store pointers to its data.
Definition: pointer_vector.h:72
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
Properties encapsulates data shared by different Elements or Conditions. It can store any type of dat...
Definition: properties.h:69
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
std::size_t IndexType
The definition of the index type.
Definition: key_hash.h:35
Matrix MatrixType
Definition: geometrical_transformation_utilities.h:55
Modeler::Pointer Create(const std::string &ModelerName, Model &rModel, const Parameters ModelParameters)
Checks if the modeler is registered.
Definition: modeler_factory.cpp:30
void SetValuesOnIntegrationPoints(TObject &dummy, const Variable< TDataType > &rVariable, const std::vector< TDataType > &values, const ProcessInfo &rCurrentProcessInfo)
Definition: add_mesh_to_python.cpp:185
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
ModelPart::NodesContainerType NodesArrayType
Definition: gid_gauss_point_container.h:42
constexpr SizeType VOIGT_SIZE_3D
Definition: geo_mechanics_application_constants.h:41
constexpr SizeType VOIGT_SIZE_2D_PLANE_STRESS
Definition: geo_mechanics_application_constants.h:42
def load(f)
Definition: ode_solve.py:307
Definition: constitutive_law.h:189
Member Variables.
Definition: geo_structural_base_element.hpp:114
double detF
Definition: geo_structural_base_element.hpp:144
double IntegrationCoefficient
Definition: geo_structural_base_element.hpp:136
Matrix B
Variables computed at each GP.
Definition: geo_structural_base_element.hpp:131
array_1d< double, TNumNodes *TDim > DisplacementVector
Properties variables.
Definition: geo_structural_base_element.hpp:120
Vector StrainVector
Constitutive Law parameters.
Definition: geo_structural_base_element.hpp:138
Matrix F
Definition: geo_structural_base_element.hpp:143
Matrix UVoigtMatrix
Auxiliary Variables.
Definition: geo_structural_base_element.hpp:148
Matrix TransformationMatrix
Definition: geo_structural_base_element.hpp:134
double HalfThickness
Definition: geo_structural_base_element.hpp:145
Matrix GradNe
Definition: geo_structural_base_element.hpp:142
Vector StressVector
Definition: geo_structural_base_element.hpp:139
Matrix NodalCrossDirection
General elemental variables.
Definition: geo_structural_base_element.hpp:128
BoundedMatrix< double, TDim, TNumNodes *TDim > NuTot
Definition: geo_structural_base_element.hpp:132
Matrix ConstitutiveMatrix
Definition: geo_structural_base_element.hpp:140
Vector DofValuesVector
Definition: geo_structural_base_element.hpp:125
array_1d< double, TNumNodes *TDim > UVector
Definition: geo_structural_base_element.hpp:123
array_1d< double, TNumNodes *TDim > VelocityVector
Definition: geo_structural_base_element.hpp:121
array_1d< double, TDim > GaussVolumeAcceleration
Definition: geo_structural_base_element.hpp:135
array_1d< double, TNumNodes *TDim > NodalVolumeAcceleration
Definition: geo_structural_base_element.hpp:122
Vector Nu
Definition: geo_structural_base_element.hpp:141