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.
axisymmetric_updated_lagrangian_U_P_element.hpp
Go to the documentation of this file.
1 //
2 // Project Name: KratosSolidMechanicsApplication $
3 // Created by: $Author: JMCarbonell $
4 // Last modified by: $Co-Author: $
5 // Date: $Date: July 2013 $
6 // Revision: $Revision: 0.0 $
7 //
8 //
9 
10 #if !defined(KRATOS_AXISYMMETRIC_UPDATED_LAGRANGIAN_U_P_ELEMENT_H_INCLUDED)
11 #define KRATOS_AXISYMMETRIC_UPDATED_LAGRANGIAN_U_P_ELEMENT_H_INCLUDED
12 
13 // System includes
14 
15 // External includes
16 
17 // Project includes
19 
20 
21 namespace Kratos
22 {
37 
39 
45 class KRATOS_API(SOLID_MECHANICS_APPLICATION) AxisymmetricUpdatedLagrangianUPElement
47 {
48 public:
49 
55  typedef ConstitutiveLawType::Pointer ConstitutiveLawPointerType;
64 
67 
71 
73  AxisymmetricUpdatedLagrangianUPElement(IndexType NewId, GeometryType::Pointer pGeometry);
74 
75  AxisymmetricUpdatedLagrangianUPElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties);
76 
79 
82 
86 
89 
93 
104  Element::Pointer Create(IndexType NewId, NodesArrayType const& ThisNodes, PropertiesType::Pointer pProperties) const override;
105 
106 
114  Element::Pointer Clone(IndexType NewId, NodesArrayType const& ThisNodes) const override;
115 
116  //************* GETTING METHODS
117 
118  //SET
119 
123  void SetValuesOnIntegrationPoints(const Variable<double>& rVariable, const std::vector<double>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
124 
125 
126  //GET:
127 
131  void CalculateOnIntegrationPoints(const Variable<double>& rVariable, std::vector<double>& rValues, const ProcessInfo& rCurrentProcessInfo) override;
132 
133 
134  //************* STARTING - ENDING METHODS
135 
140  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
141 
142  //************* COMPUTING METHODS
143 
150  void CalculateMassMatrix(MatrixType& rMassMatrix,
151  const ProcessInfo& rCurrentProcessInfo) override;
152  //************************************************************************************
153  //************************************************************************************
161  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
162 
166 
177 
178 protected:
184 
188  std::vector< Matrix > mDeformationGradientF0;
189 
194 
199  {
200  }
201 
205 
206 
211  void CalculateAndAddLHS(LocalSystemComponents& rLocalSystem,
212  ElementDataType& rVariables,
213  double& rIntegrationWeight) override;
214 
219  void CalculateAndAddRHS(LocalSystemComponents& rLocalSystem,
220  ElementDataType& rVariables,
221  Vector& rVolumeForce,
222  double& rIntegrationWeight) override;
223 
227  double& CalculateTotalMass( double& rTotalMass, const ProcessInfo& rCurrentProcessInfo ) override;
228 
229 
233  void CalculateAndAddKuug(MatrixType& rK,
234  ElementDataType & rVariables,
235  double& rIntegrationWeight
236  ) override;
237 
238 
242  void CalculateAndAddKup (MatrixType& rK,
243  ElementDataType & rVariables,
244  double& rIntegrationWeight
245  ) override;
246 
250  void CalculateAndAddKpu(MatrixType& rK,
251  ElementDataType & rVariables,
252  double& rIntegrationWeight
253  ) override;
254 
258  void CalculateAndAddKpp(MatrixType& rK,
259  ElementDataType & rVariables,
260  double& rIntegrationWeight
261  ) override;
262 
266  void CalculateAndAddKppStab(MatrixType& rK,
267  ElementDataType & rVariables,
268  double& rIntegrationWeight
269  ) override;
270 
271 
275  void CalculateAndAddPressureForces(VectorType& rRightHandSideVector,
276  ElementDataType & rVariables,
277  double& rIntegrationWeight
278  ) override;
279 
280 
284  void CalculateAndAddStabilizedPressure(VectorType& rRightHandSideVector,
285  ElementDataType & rVariables,
286  double& rIntegrationWeight
287  ) override;
288 
292  void InitializeElementData(ElementDataType & rVariables,
293  const ProcessInfo& rCurrentProcessInfo) override;
294 
295 
296 
300  void FinalizeStepVariables(ElementDataType & rVariables,
301  const double& rPointNumber ) override;
302 
306  void CalculateKinematics(ElementDataType& rVariables,
307  const double& rPointNumber) override;
308 
309 
313  void CalculateRadius(double & rCurrentRadius,
314  double & rReferenceRadius,
315  const Vector& rN);
316 
320  void CalculateDeformationGradient(Matrix& rF,
321  const Matrix& rDN_DX,
322  const Matrix& rDeltaPosition,
323  const double & rCurrentRadius,
324  const double & rReferenceRadius);
325 
329  void CalculateDeformationMatrix(Matrix& rB,
330  const Matrix& rDN_DX,
331  const Vector& rN,
332  const double & rCurrentRadius);
333 
337  void GetHistoricalVariables( ElementDataType& rVariables,
338  const double& rPointNumber ) override;
339 
340 
344  void CalculateGreenLagrangeStrain(const Matrix& rF,
345  Vector& rStrainVector) override;
346 
350  void CalculateAlmansiStrain(const Matrix& rF,
351  Vector& rStrainVector) override;
352 
353 
357  double& CalculateVolumeChange(double& rVolumeChange, ElementDataType& rVariables) override;
358 
369 
370 private:
371 
377 
378 
382 
383 
387 
388 
393 
397  friend class Serializer;
398 
399  // A private default constructor necessary for serialization
400 
401  void save(Serializer& rSerializer) const override;
402 
403  void load(Serializer& rSerializer) override;
404 
405 
412 
413 }; // Class AxisymmetricUpdatedLagrangianUPElement
414 
422 
423 } // namespace Kratos.
424 #endif // KRATOS_AXISYMMETRIC_UPDATED_LAGRANGIAN_U_P_ELEMENT_H_INCLUDED defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Axisymmetric Updated Lagrangian U-P Element for 2D geometries. For Linear Triangles.
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:47
LargeDisplacementUPElement::ElementDataType ElementDataType
Type for element variables.
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:63
Vector mDeterminantF0
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:193
GeometryData::IntegrationMethod IntegrationMethod
Type definition for integration methods.
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:59
std::vector< Matrix > mDeformationGradientF0
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:188
ConstitutiveLawType::StressMeasure StressMeasureType
StressMeasure from constitutive laws.
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:57
AxisymmetricUpdatedLagrangianUPElement()
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:198
GeometryData::SizeType SizeType
Type for size.
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:61
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(AxisymmetricUpdatedLagrangianUPElement)
Counted pointer of AxisymmetricUpdatedLagrangianUPElement.
ConstitutiveLawType::Pointer ConstitutiveLawPointerType
Pointer type for constitutive laws.
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:55
ConstitutiveLaw ConstitutiveLawType
Definition: axisymmetric_updated_lagrangian_U_P_element.hpp:53
Definition: constitutive_law.h:47
StressMeasure
Definition: constitutive_law.h:69
std::size_t IndexType
Definition: flags.h:74
IntegrationMethod
Definition: geometry_data.h:76
std::size_t SizeType
Definition: geometry_data.h:173
Large Displacement Lagrangian U-P Element for 3D and 2D geometries. Linear Triangles and Tetrahedra (...
Definition: large_displacement_U_P_element.hpp:47
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
double CalculateRadius(const PairedCondition *pCondition, const Vector &rNSlave)
Calculates the radius of axisymmetry.
Definition: mortar_explicit_contribution_utilities.cpp:675
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
def load(f)
Definition: ode_solve.py:307
Definition: solid_element.hpp:83