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.
U_Pw_element.hpp
Go to the documentation of this file.
1 // | / |
2 // ' / __| _` | __| _ \ __|
3 // . \ | ( | | ( |\__ `
4 // _|\_\_| \__,_|\__|\___/ ____/
5 // Multi-Physics
6 //
7 // License: BSD License
8 // Kratos default license: kratos/license.txt
9 //
10 // Main authors: Ignasi de Pouplana
11 //
12 
13 
14 #if !defined(KRATOS_U_PW_ELEMENT_H_INCLUDED )
15 #define KRATOS_U_PW_ELEMENT_H_INCLUDED
16 
17 // Project includes
18 #include "containers/array_1d.h"
19 #include "includes/define.h"
20 #include "includes/element.h"
21 #include "includes/serializer.h"
22 #include "geometries/geometry.h"
23 #include "utilities/math_utils.h"
25 
26 // Application includes
29 
30 namespace Kratos
31 {
32 
33 template< unsigned int TDim, unsigned int TNumNodes >
34 class KRATOS_API(POROMECHANICS_APPLICATION) UPwElement : public Element
35 {
36 
37 public:
38 
40 
42 
44  UPwElement(IndexType NewId = 0) : Element( NewId ) {}
45 
47  UPwElement(IndexType NewId, const NodesArrayType& ThisNodes) : Element(NewId, ThisNodes) {}
48 
50  UPwElement(IndexType NewId, GeometryType::Pointer pGeometry) : Element( NewId, pGeometry ) {}
51 
53  UPwElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties) : Element( NewId, pGeometry, pProperties )
54  {
55  mThisIntegrationMethod = this->GetIntegrationMethod();
56  }
57 
59  virtual ~UPwElement() {}
60 
62 
63  Element::Pointer Create(IndexType NewId, NodesArrayType const& ThisNodes, PropertiesType::Pointer pProperties) const override;
64 
65  Element::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, PropertiesType::Pointer pProperties) const override;
66 
67  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
68 
69  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
70 
71  void GetDofList( DofsVectorType& rElementalDofList, const ProcessInfo& rCurrentProcessInfo ) const override;
72 
73  GeometryData::IntegrationMethod GetIntegrationMethod() const override;
74 
76 
77  void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix,VectorType& rRightHandSideVector,const ProcessInfo& rCurrentProcessInfo ) override;
78 
79  void CalculateLeftHandSide(MatrixType& rLeftHandSideMatrix,const ProcessInfo& rCurrentProcessInfo ) override;
80 
81  void CalculateRightHandSide(VectorType& rRightHandSideVector,const ProcessInfo& rCurrentProcessInfo ) override;
82 
83  void EquationIdVector(EquationIdVectorType& rResult, const ProcessInfo& rCurrentProcessInfo) const override;
84 
85  void CalculateMassMatrix(MatrixType& rMassMatrix, const ProcessInfo& rCurrentProcessInfo) override;
86 
87  void CalculateDampingMatrix(MatrixType& rDampingMatrix, const ProcessInfo& rCurrentProcessInfo) override;
88 
89  void GetValuesVector(Vector& rValues, int Step = 0) const override;
90 
91  void GetFirstDerivativesVector(Vector& rValues, int Step = 0) const override;
92 
93  void GetSecondDerivativesVector(Vector& rValues, int Step = 0) const override;
94 
96 
97  void SetValuesOnIntegrationPoints(const Variable<double>& rVariable, const std::vector<double>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
98 
99  void SetValuesOnIntegrationPoints(const Variable<Matrix>& rVariable, const std::vector<Matrix>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
100 
101 
102  void CalculateOnIntegrationPoints(const Variable<double>& rVariable, std::vector<double>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
103 
104  void CalculateOnIntegrationPoints(const Variable<array_1d<double,3>>& rVariable, std::vector<array_1d<double,3>>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
105 
106  void CalculateOnIntegrationPoints(const Variable<Matrix>& rVariable, std::vector<Matrix>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
107 
108  void CalculateOnIntegrationPoints( const Variable<ConstitutiveLaw::Pointer>& rVariable, std::vector<ConstitutiveLaw::Pointer>& rValues,const ProcessInfo& rCurrentProcessInfo ) override;
109 
111 
112  void AddExplicitContribution(
113  const VectorType& rRHSVector,
114  const Variable<VectorType>& rRHSVariable,
115  const Variable<double >& rDestinationVariable,
116  const ProcessInfo& rCurrentProcessInfo
117  ) override;
118 
119  void AddExplicitContribution(const VectorType& rRHSVector,
120  const Variable<VectorType>& rRHSVariable,
121  const Variable<array_1d<double, 3> >& rDestinationVariable,
122  const ProcessInfo& rCurrentProcessInfo
123  ) override;
124 
125 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
126 
127 
128 protected:
129 
131 
133 
134  std::vector<ConstitutiveLaw::Pointer> mConstitutiveLawVector;
135 
137 
138  std::vector<double> mImposedZStrainVector;
139 
141 
142  virtual void CalculateStiffnessMatrix( MatrixType& rStiffnessMatrix, const ProcessInfo& CurrentProcessInfo );
143 
144  virtual void CalculateAll( MatrixType& rLeftHandSideMatrix, VectorType& rRightHandSideVector, const ProcessInfo& CurrentProcessInfo );
145 
146  virtual void CalculateRHS( VectorType& rRightHandSideVector, const ProcessInfo& CurrentProcessInfo );
147 
148  void CalculateIntegrationCoefficient(double& rIntegrationCoefficient, const double& detJ, const double& weight);
149 
150 
151  virtual void CalculateFluxResidual (VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo);
152 
153  virtual void CalculateMixBodyForce (VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo);
154 
155  virtual void CalculateNegInternalForce (VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo);
156 
157  virtual void CalculateExplicitContributions (VectorType& rFluxResidual, VectorType& rBodyForce, VectorType& rNegInternalForces, const ProcessInfo& rCurrentProcessInfo);
158 
159  virtual void CalculateLumpedMassMatrix( MatrixType& rLeftHandSideMatrix, const ProcessInfo& rCurrentProcessInfo );
160 
161  virtual void CalculateDampingMatrixWithLumpedMass( MatrixType& rLeftHandSideMatrix, const ProcessInfo& rCurrentProcessInfo );
162 
163  virtual void CalculateInertialForce (VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo);
164 
165  virtual void CalculateDampingForce (VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo);
166 
168 
169 private:
170 
172 
173  friend class Serializer;
174 
175  void save(Serializer& rSerializer) const override
176  {
178  }
179 
180  void load(Serializer& rSerializer) override
181  {
183  }
184 
186  UPwElement & operator=(UPwElement const& rOther);
187 
189  UPwElement(UPwElement const& rOther);
190 
191 
192 }; // Class UPwElement
193 
194 } // namespace Kratos
195 
196 #endif // KRATOS_U_PW_ELEMENT_H_INCLUDED defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Base class for all Elements.
Definition: element.h:60
std::vector< std::size_t > EquationIdVectorType
Definition: element.h:98
std::size_t IndexType
Definition: flags.h:74
IntegrationMethod
Definition: geometry_data.h:76
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
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
Definition: U_Pw_element.hpp:35
std::vector< ConstitutiveLaw::Pointer > mConstitutiveLawVector
Definition: U_Pw_element.hpp:134
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &rCurrentProcessInfo) const override
UPwElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: U_Pw_element.hpp:47
virtual ~UPwElement()
Destructor.
Definition: U_Pw_element.hpp:59
Matrix mIntrinsicPermeability
Definition: U_Pw_element.hpp:136
UPwElement(IndexType NewId=0)
Default Constructor.
Definition: U_Pw_element.hpp:44
UPwElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using Geometry.
Definition: U_Pw_element.hpp:50
void CalculateIntegrationCoefficient(double &rIntegrationCoefficient, const double &detJ, const double &weight)
UPwElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Constructor using Properties.
Definition: U_Pw_element.hpp:53
std::vector< double > mImposedZStrainVector
Definition: U_Pw_element.hpp:138
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(UPwElement)
GeometryData::IntegrationMethod mThisIntegrationMethod
Member Variables.
Definition: U_Pw_element.hpp:132
#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
pybind11::list CalculateOnIntegrationPoints(TObject &dummy, const Variable< TDataType > &rVariable, const ProcessInfo &rProcessInfo)
Definition: add_mesh_to_python.cpp:142
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
ModelPart::NodesContainerType NodesArrayType
Definition: gid_gauss_point_container.h:42
def load(f)
Definition: ode_solve.py:307