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.
transient_Pw_interface_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 #pragma once
14 
15 // Project includes
16 #include "includes/serializer.h"
17 
18 // Application includes
19 #include "custom_elements/U_Pw_small_strain_interface_element.hpp"
20 
21 namespace Kratos
22 {
23 
24 template <unsigned int TDim, unsigned int TNumNodes>
25 class KRATOS_API(GEO_MECHANICS_APPLICATION) TransientPwInterfaceElement
26  : public UPwSmallStrainInterfaceElement<TDim, TNumNodes>
27 {
28 public:
30 
32 
33  using IndexType = std::size_t;
35  using NodeType = Node;
38  using VectorType = Vector;
39  using MatrixType = Matrix;
40 
43 
45  using SizeType = std::size_t;
46 
47  using BaseType::CalculateRetentionResponse;
48  using BaseType::mRetentionLawVector;
49  using BaseType::mThisIntegrationMethod;
50 
51  using InterfaceElementVariables = typename BaseType::InterfaceElementVariables;
52  using SFGradAuxVariables = typename BaseType::SFGradAuxVariables;
53 
56  : UPwSmallStrainInterfaceElement<TDim, TNumNodes>(NewId)
57  {
58  }
59 
62  : UPwSmallStrainInterfaceElement<TDim, TNumNodes>(NewId, ThisNodes)
63  {
64  }
65 
67  TransientPwInterfaceElement(IndexType NewId, GeometryType::Pointer pGeometry)
68  : UPwSmallStrainInterfaceElement<TDim, TNumNodes>(NewId, pGeometry)
69  {
70  }
71 
73  TransientPwInterfaceElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
74  : UPwSmallStrainInterfaceElement<TDim, TNumNodes>(NewId, pGeometry, pProperties)
75  {
76  }
77 
78  ~TransientPwInterfaceElement() override = default;
83 
84  Element::Pointer Create(IndexType NewId,
85  NodesArrayType const& ThisNodes,
86  PropertiesType::Pointer pProperties) const override;
87 
88  Element::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, PropertiesType::Pointer pProperties) const override;
89 
90  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
91 
92  void GetDofList(DofsVectorType& rElementalDofList, const ProcessInfo& rCurrentProcessInfo) const override;
93 
94  void EquationIdVector(EquationIdVectorType& rResult, const ProcessInfo& rCurrentProcessInfo) const override;
95 
96  void GetValuesVector(Vector& rValues, int Step = 0) const override;
97 
98  void GetFirstDerivativesVector(Vector& rValues, int Step = 0) const override;
99 
100  void GetSecondDerivativesVector(Vector& rValues, int Step = 0) const override;
101 
102  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
103 
104  void CalculateMassMatrix(MatrixType& rMassMatrix, const ProcessInfo& rCurrentProcessInfo) override;
105 
106  void InitializeSolutionStep(const ProcessInfo& rCurrentProcessInfo) override;
107  void FinalizeSolutionStep(const ProcessInfo& rCurrentProcessInfo) override;
108 
109  void CalculateOnIntegrationPoints(const Variable<Matrix>& rVariable,
110  std::vector<Matrix>& rValues,
111  const ProcessInfo& rCurrentProcessInfo) override;
112 
113  void CalculateOnIntegrationPoints(const Variable<double>& rVariable,
114  std::vector<double>& rValues,
115  const ProcessInfo& rCurrentProcessInfo) override;
116 
118  std::vector<array_1d<double, 3>>& rValues,
119  const ProcessInfo& rCurrentProcessInfo) override;
120 
121 protected:
122  void CalculateOnLobattoIntegrationPoints(const Variable<array_1d<double, 3>>& rVariable,
123  std::vector<array_1d<double, 3>>& rOutput,
124  const ProcessInfo& rCurrentProcessInfo) override;
125 
126  void CalculateOnLobattoIntegrationPoints(const Variable<Matrix>& rVariable,
127  std::vector<Matrix>& rOutput,
128  const ProcessInfo& rCurrentProcessInfo) override;
129 
130  void CalculateAll(MatrixType& rLeftHandSideMatrix,
131  VectorType& rRightHandSideVector,
132  const ProcessInfo& CurrentProcessInfo,
133  const bool CalculateStiffnessMatrixFlag,
134  const bool CalculateResidualVectorFlag) override;
135 
136  void InitializeElementVariables(InterfaceElementVariables& rVariables,
137  const GeometryType& Geom,
138  const PropertiesType& Prop,
139  const ProcessInfo& CurrentProcessInfo) override;
140 
141  void CalculateAndAddLHS(MatrixType& rLeftHandSideMatrix, InterfaceElementVariables& rVariables) override;
142 
143  void CalculateAndAddCompressibilityMatrix(MatrixType& rLeftHandSideMatrix,
144  InterfaceElementVariables& rVariables) override;
145 
146  void CalculateAndAddPermeabilityMatrix(MatrixType& rLeftHandSideMatrix,
147  InterfaceElementVariables& rVariables) override;
148 
149  void CalculateAndAddRHS(VectorType& rRightHandSideVector,
150  InterfaceElementVariables& rVariables,
151  unsigned int GPoint) override;
152 
153  void CalculateAndAddCompressibilityFlow(VectorType& rRightHandSideVector,
154  InterfaceElementVariables& rVariables) override;
155 
156  void CalculateAndAddPermeabilityFlow(VectorType& rRightHandSideVector,
157  InterfaceElementVariables& rVariables) override;
158 
159  void CalculateAndAddFluidBodyFlow(VectorType& rRightHandSideVector,
160  InterfaceElementVariables& rVariables) override;
161 
162  unsigned int GetNumberOfDOF() const override;
163 
164 private:
166  friend class Serializer;
167 
168  void save(Serializer& rSerializer) const override
169  {
171  }
172 
173  void load(Serializer& rSerializer) override
174  {
176  }
177 
178 }; // Class TransientPwInterfaceElement
179 
180 } // namespace Kratos
Base class for all Elements.
Definition: element.h:60
std::size_t SizeType
Definition: element.h:94
std::vector< DofType::Pointer > DofsVectorType
Definition: element.h:100
std::vector< std::size_t > EquationIdVectorType
Definition: element.h:98
std::size_t IndexType
Definition: flags.h:74
This defines the geometrical object, base definition of the element and condition entities.
Definition: geometrical_object.h:58
Geometry base class.
Definition: geometry.h:71
PointerVector< TPointType > PointsArrayType
Definition: geometry.h:118
This class defines the node.
Definition: node.h:65
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
Definition: transient_Pw_interface_element.hpp:27
typename BaseType::SFGradAuxVariables SFGradAuxVariables
Definition: transient_Pw_interface_element.hpp:52
TransientPwInterfaceElement(const TransientPwInterfaceElement &)=delete
TransientPwInterfaceElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Constructor using Properties.
Definition: transient_Pw_interface_element.hpp:73
typename BaseType::InterfaceElementVariables InterfaceElementVariables
Definition: transient_Pw_interface_element.hpp:51
TransientPwInterfaceElement(IndexType NewId=0)
Default Constructor.
Definition: transient_Pw_interface_element.hpp:55
TransientPwInterfaceElement & operator=(const TransientPwInterfaceElement &)=delete
TransientPwInterfaceElement(TransientPwInterfaceElement &&)=delete
~TransientPwInterfaceElement() override=default
TransientPwInterfaceElement & operator=(TransientPwInterfaceElement &&)=delete
TransientPwInterfaceElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: transient_Pw_interface_element.hpp:61
TransientPwInterfaceElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using Geometry.
Definition: transient_Pw_interface_element.hpp:67
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(TransientPwInterfaceElement)
Definition: U_Pw_small_strain_interface_element.hpp:30
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
Modeler::Pointer Create(const std::string &ModelerName, Model &rModel, const Parameters ModelParameters)
Checks if the modeler is registered.
Definition: modeler_factory.cpp:30
pybind11::list CalculateOnIntegrationPoints(TObject &dummy, const Variable< TDataType > &rVariable, const ProcessInfo &rProcessInfo)
Definition: add_mesh_to_python.cpp:142
void InitializeSolutionStep(ConstructionUtility &rThisUtil, std::string ThermalSubModelPartName, std::string MechanicalSubModelPartName, std::string HeatFluxSubModelPartName, std::string HydraulicPressureSubModelPartName, bool thermal_conditions, bool mechanical_conditions, int phase)
Definition: add_custom_utilities_to_python.cpp:45
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
Internals::Matrix< double, AMatrix::dynamic, 1 > Vector
Definition: amatrix_interface.h:472
Internals::Matrix< double, AMatrix::dynamic, AMatrix::dynamic > Matrix
Definition: amatrix_interface.h:470
def load(f)
Definition: ode_solve.py:307