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.
three_step_updated_lagrangian_element.h
Go to the documentation of this file.
1 //
2 // Project Name: KratosFluidDynamicsApplication $
3 // Last modified by: $Author: AFranci $
4 // Date: $Date: June 2021 $
5 // Revision: $Revision: 0.0 $
6 //
7 //
8 
9 #if !defined(KRATOS_THREE_STEP_UPDATED_LAGRANGIAN_ELEMENT_H_INCLUDED)
10 #define KRATOS_THREE_STEP_UPDATED_LAGRANGIAN_ELEMENT_H_INCLUDED
11 
12 // System includes
13 #include <string>
14 #include <iostream>
15 
16 // External includes
17 
18 // Project includes
19 #include "containers/array_1d.h"
20 #include "includes/define.h"
21 #include "includes/serializer.h"
22 #include "geometries/geometry.h"
23 #include "utilities/math_utils.h"
24 #include "utilities/geometry_utilities.h"
25 
27 
29 
30 #include "includes/model_part.h"
31 
32 namespace Kratos
33 {
34 
37 
40 
44 
48 
52 
56 
58 
61  template <unsigned int TDim>
63  {
64  public:
67 
70 
72 
74  typedef Node NodeType;
75 
78 
81 
83  typedef Vector VectorType;
84 
86  typedef Matrix MatrixType;
87 
88  typedef std::size_t IndexType;
89 
90  typedef std::size_t SizeType;
91 
92  typedef std::vector<std::size_t> EquationIdVectorType;
93 
94  typedef std::vector<Dof<double>::Pointer> DofsVectorType;
95 
97 
100 
103 
106 
109 
111  typedef ConstitutiveLawType::Pointer ConstitutiveLawPointerType;
112 
114 
115  typedef typename BaseType::PropertiesType::Pointer pPropertiesType;
116 
117  typedef typename BaseType::ElementalVariables ElementalVariables;
118 
122 
123  // Constructors.
124 
126 
130  {
131  }
132 
134 
138  ThreeStepUpdatedLagrangianElement(IndexType NewId, const NodesArrayType &ThisNodes) : BaseType(NewId, ThisNodes)
139  {
140  }
141 
143 
147  ThreeStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry) : BaseType(NewId, pGeometry)
148  {
149  }
150 
152 
157  ThreeStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties) : BaseType(NewId, pGeometry, pProperties)
158  {
159  }
160 
162 
164 
167  {
168  }
169 
173 
177 
179 
186  Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes,
187  pPropertiesType pProperties) const override
188  {
189  return Element::Pointer(new ThreeStepUpdatedLagrangianElement(NewId, this->GetGeometry().Create(ThisNodes), pProperties));
190  }
191 
192  Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override;
193 
194  void Initialize(const ProcessInfo &rCurrentProcessInfo) override{};
195 
197  void InitializeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override{};
198 
199  void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override{};
200 
202  virtual void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix,
203  VectorType &rRightHandSideVector,
204  const ProcessInfo &rCurrentProcessInfo) override{};
205 
206  void AddMomentumMassTerm(Matrix &rMassMatrix,
207  const ShapeFunctionsType &rN,
208  const double Weight);
209 
210  void AddMomentumRHSTerms(Vector &rRHSVector,
211  const double Density,
212  const array_1d<double, 3> &rBodyForce,
213  const double OldPressure,
214  const ShapeFunctionsType &rN,
215  const ShapeFunctionDerivativesType &rDN_DX,
216  const double Weight);
217 
218  void AddViscousTerm(MatrixType &rDampingMatrix,
219  const ShapeFunctionDerivativesType &rShapeDeriv,
220  const double Weight,
221  const double theta_velocity);
222 
223  void ComputeBoundLHSMatrix(MatrixType &BoundLHSMatrix,
224  const ShapeFunctionsType &rN,
225  const double Weight);
226 
227  void ComputeBoundRHSVector(VectorType &BoundRHSVector,
228  const ShapeFunctionsType &rN,
229  const double TimeStep,
230  const double BoundRHSCoeffAcc,
231  const double BoundRHSCoeffDev);
232 
234  const double TimeStep,
235  const double BoundRHSCoeffAcc,
236  const double BoundRHSCoeffDev,
237  const VectorType SpatialDefRate);
238 
239  void ComputeStabLaplacianMatrix(MatrixType &StabLaplacianMatrix,
240  const ShapeFunctionDerivativesType &rShapeDeriv,
241  const double Weight);
242 
243  void CalculateTauFIC(double &TauOne,
244  double ElemSize,
245  const double Density,
246  const double Viscosity,
247  const ProcessInfo &rCurrentProcessInfo);
248 
249  double ElementSize();
250 
251  void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix,
252  const ProcessInfo &rCurrentProcessInfo) override
253  {
254  KRATOS_TRY;
255  KRATOS_THROW_ERROR(std::logic_error, "ThreeStepUpdatedLagrangianElement::CalculateLeftHandSide not implemented", "");
256  KRATOS_CATCH("");
257  }
258 
264  virtual void Calculate(const Variable<array_1d<double, 3>> &rVariable,
265  array_1d<double, 3> &rOutput,
266  const ProcessInfo &rCurrentProcessInfo) override{};
267 
268  void CalculateRightHandSide(VectorType &rRightHandSideVector,
269  const ProcessInfo &rCurrentProcessInfo) override{};
270 
271  // The following methods have different implementations depending on TDim
273 
280  const ProcessInfo &rCurrentProcessInfo) const override;
281 
283 
287  void GetDofList(DofsVectorType &rElementalDofList,
288  const ProcessInfo &rCurrentProcessInfo) const override;
289 
290  virtual void UpdateCauchyStress(unsigned int g, const ProcessInfo &rCurrentProcessInfo) override{};
291 
292  void InitializeElementalVariables(ElementalVariables &rElementalVariables) override;
293 
297 
301 
303 
311  int Check(const ProcessInfo &rCurrentProcessInfo) const override { return 0; };
312 
316 
320 
322  std::string Info() const override
323  {
324  std::stringstream buffer;
325  buffer << "ThreeStepUpdatedLagrangianElement #" << this->Id();
326  return buffer.str();
327  }
328 
330  void PrintInfo(std::ostream &rOStream) const override {}
331 
332  // /// Print object's data.
333  // virtual void PrintData(std::ostream& rOStream) const;
334 
338 
340  protected:
343 
347 
348  ConstitutiveLaw::Pointer mpConstitutiveLaw = nullptr;
349 
353 
357 
358  virtual void CalculateLocalMomentumEquations(MatrixType &rLeftHandSideMatrix,
359  VectorType &rRightHandSideVector,
360  const ProcessInfo &rCurrentProcessInfo) override{};
361 
362  void CalculatePSPGPressureSystem(MatrixType &rLeftHandSideMatrix,
363  VectorType &rRightHandSideVector,
364  const ProcessInfo &rCurrentProcessInfo);
365 
366  // void CalculateFICPressureSystem(MatrixType &rLeftHandSideMatrix,
367  // VectorType &rRightHandSideVector,
368  // const ProcessInfo &rCurrentProcessInfo);
369 
371  VectorType &rRightHandSideVector,
372  const ShapeFunctionsType &rN,
373  const double lagMultiplier);
374 
375  void CalculateTauPSPG(double &TauOne,
376  double ElemSize,
377  const double Density,
378  const double Viscosity,
379  const ProcessInfo &rCurrentProcessInfo);
380 
381  virtual double GetThetaMomentum() override
382  {
383  std::cout << "I SHOULD NOT ENTER HERE!" << std::endl;
384  return 0.0;
385  };
386 
387  virtual double GetThetaContinuity() override
388  {
389  std::cout << "I SHOULD NOT ENTER HERE!" << std::endl;
390  return 0.0;
391  };
392 
393  void CalculateMassMatrix(Matrix &rMassMatrix,
394  const ProcessInfo &rCurrentProcessInfo) override{};
395 
396  void ComputeLumpedMassMatrix(Matrix &rMassMatrix,
397  const double Weight);
398 
399  void AddExternalForces(Vector &rRHSVector,
400  const double Density,
401  const ShapeFunctionsType &rN,
402  const double Weight) override;
403 
405  ElementalVariables &rElementalVariables,
406  const unsigned int g,
407  const Vector& rN,
408  const ProcessInfo &rCurrentProcessInfo,
409  double &Density,
410  double &DeviatoricCoeff,
411  double &VolumetricCoeff) override
412  {};
413 
417 
421 
425 
427 
428  private:
431 
435 
439 
440  friend class Serializer;
441 
442  void save(Serializer &rSerializer) const override
443  {
445  }
446 
447  void load(Serializer &rSerializer) override
448  {
450  }
451 
455 
459 
463 
467 
471 
474 
475  /* /// Copy constructor. */
476  /* ThreeStepUpdatedLagrangianElement(ThreeStepUpdatedLagrangianElement const& rOther); */
477 
479 
480  }; // Class ThreeStepUpdatedLagrangianElement
481 
483 
486 
490 
492  template <unsigned int TDim>
493  inline std::istream &operator>>(std::istream &rIStream,
495  {
496  return rIStream;
497  }
498 
500  template <unsigned int TDim>
501  inline std::ostream &operator<<(std::ostream &rOStream,
503  {
504  rThis.PrintInfo(rOStream);
505  rOStream << std::endl;
506  rThis.PrintData(rOStream);
507 
508  return rOStream;
509  }
510 
511 } // namespace Kratos.
512 
513 #endif // KRATOS_THREE_STEP_UPDATED_LAGRANGIAN_ELEMENT defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Definition: constitutive_law.h:47
Base class for all Elements.
Definition: element.h:60
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
GeometryType & GetGeometry()
Returns the reference of the geometry.
Definition: geometrical_object.h:158
Geometry base class.
Definition: geometry.h:71
This object defines an indexed object.
Definition: indexed_object.h:54
IndexType Id() const
Definition: indexed_object.h:107
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
A sorted associative container similar to an STL set, but uses a vector to store pointers to its data...
Definition: pointer_vector_set.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
A stabilized element for the incompressible Navier-Stokes equations.
Definition: three_step_updated_lagrangian_element.h:63
A stabilized element for the incompressible Navier-Stokes equations.
Definition: updated_lagrangian_element.h:64
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
#define KRATOS_THROW_ERROR(ExceptionType, ErrorMessage, MoreInfo)
Definition: define.h:77
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_TRY
Definition: define.h:109
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: three_step_updated_lagrangian_element.h:80
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, pPropertiesType pProperties) const override
Create a new element of this type.
Definition: three_step_updated_lagrangian_element.h:186
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:251
void GetDofList(DofsVectorType &rElementalDofList, const ProcessInfo &rCurrentProcessInfo) const override
Returns a list of the element's Dofs.
Definition: three_step_updated_lagrangian_element.cpp:63
UpdatedLagrangianElement< TDim > BaseType
Definition: three_step_updated_lagrangian_element.h:71
ConstitutiveLaw ConstitutiveLawType
Reference type definition for constitutive laws.
Definition: three_step_updated_lagrangian_element.h:108
virtual void Calculate(const Variable< array_1d< double, 3 >> &rVariable, array_1d< double, 3 > &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:264
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: three_step_updated_lagrangian_element.h:102
Vector VectorType
Vector type for local contributions to the linear system.
Definition: three_step_updated_lagrangian_element.h:83
double ElementSize()
Definition: three_step_updated_lagrangian_element.cpp:244
virtual void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Calculate the element's local contribution to the system for the current step.
Definition: three_step_updated_lagrangian_element.h:202
void CalculateTauPSPG(double &TauOne, double ElemSize, const double Density, const double Viscosity, const ProcessInfo &rCurrentProcessInfo)
Definition: three_step_updated_lagrangian_element.cpp:981
Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override
It creates a new element pointer and clones the previous element data.
Definition: three_step_updated_lagrangian_element.cpp:21
void ComputeBoundRHSVectorComplete(VectorType &BoundRHSVector, const double TimeStep, const double BoundRHSCoeffAcc, const double BoundRHSCoeffDev, const VectorType SpatialDefRate)
virtual void UpdateCauchyStress(unsigned int g, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:290
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: three_step_updated_lagrangian_element.h:77
std::size_t IndexType
Definition: three_step_updated_lagrangian_element.h:88
void ComputeLumpedMassMatrix(Matrix &rMassMatrix, const double Weight)
Definition: three_step_updated_lagrangian_element.cpp:269
void CalculatePSPGPressureSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: three_step_updated_lagrangian_element.cpp:341
std::string Info() const override
Turn back information as a string.
Definition: three_step_updated_lagrangian_element.h:322
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: three_step_updated_lagrangian_element.h:94
void CalculateTauFIC(double &TauOne, double ElemSize, const double Density, const double Viscosity, const ProcessInfo &rCurrentProcessInfo)
Definition: three_step_updated_lagrangian_element.cpp:766
BaseType::ElementalVariables ElementalVariables
Definition: three_step_updated_lagrangian_element.h:117
ThreeStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties)
Constuctor using geometry and properties.
Definition: three_step_updated_lagrangian_element.h:157
void ComputeBoundLHSMatrix(MatrixType &BoundLHSMatrix, const ShapeFunctionsType &rN, const double Weight)
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: three_step_updated_lagrangian_element.h:96
virtual void CalculateLocalMomentumEquations(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:358
ThreeStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using a geometry object.
Definition: three_step_updated_lagrangian_element.h:147
void CalculateMassMatrix(Matrix &rMassMatrix, const ProcessInfo &rCurrentProcessInfo) override
Add integration point contribution to the mass matrix.
Definition: three_step_updated_lagrangian_element.h:393
virtual double GetThetaMomentum() override
Definition: three_step_updated_lagrangian_element.h:381
virtual double GetThetaContinuity() override
Definition: three_step_updated_lagrangian_element.h:387
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: three_step_updated_lagrangian_element.h:99
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: three_step_updated_lagrangian_element.h:330
void AddMomentumRHSTerms(Vector &rRHSVector, const double Density, const array_1d< double, 3 > &rBodyForce, const double OldPressure, const ShapeFunctionsType &rN, const ShapeFunctionDerivativesType &rDN_DX, const double Weight)
Definition: three_step_updated_lagrangian_element.cpp:94
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Checks the input and that all required Kratos variables have been registered.
Definition: three_step_updated_lagrangian_element.h:311
void AddExternalForces(Vector &rRHSVector, const double Density, const ShapeFunctionsType &rN, const double Weight) override
Definition: three_step_updated_lagrangian_element.cpp:123
ConstitutiveLawType::Pointer ConstitutiveLawPointerType
Pointer type for constitutive laws.
Definition: three_step_updated_lagrangian_element.h:111
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: three_step_updated_lagrangian_element.h:105
BaseType::PropertiesType::Pointer pPropertiesType
Definition: three_step_updated_lagrangian_element.h:115
ThreeStepUpdatedLagrangianElement(ThreeStepUpdatedLagrangianElement const &rOther)
copy constructor
Definition: three_step_updated_lagrangian_element.h:163
void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:194
Node NodeType
Node type (default is: Node)
Definition: three_step_updated_lagrangian_element.h:74
void AddMomentumMassTerm(Matrix &rMassMatrix, const ShapeFunctionsType &rN, const double Weight)
Definition: three_step_updated_lagrangian_element.cpp:317
void ComputeStabLaplacianMatrix(MatrixType &StabLaplacianMatrix, const ShapeFunctionDerivativesType &rShapeDeriv, const double Weight)
Definition: three_step_updated_lagrangian_element.cpp:1021
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: three_step_updated_lagrangian_element.h:86
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:268
ThreeStepUpdatedLagrangianElement(IndexType NewId=0)
Default constuctor.
Definition: three_step_updated_lagrangian_element.h:129
void ComputeBoundRHSVector(VectorType &BoundRHSVector, const ShapeFunctionsType &rN, const double TimeStep, const double BoundRHSCoeffAcc, const double BoundRHSCoeffDev)
Element::PropertiesType PropertiesType
Definition: updated_lagrangian_element.h:142
void AddViscousTerm(MatrixType &rDampingMatrix, const ShapeFunctionDerivativesType &rShapeDeriv, const double Weight, const double theta_velocity)
virtual void CalcElasticPlasticCauchySplitted(ElementalVariables &rElementalVariables, const unsigned int g, const Vector &rN, const ProcessInfo &rCurrentProcessInfo, double &Density, double &DeviatoricCoeff, double &VolumetricCoeff) override
Definition: three_step_updated_lagrangian_element.h:404
std::size_t SizeType
Definition: three_step_updated_lagrangian_element.h:90
void InitializeElementalVariables(ElementalVariables &rElementalVariables) override
Definition: three_step_updated_lagrangian_element.cpp:1043
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(ThreeStepUpdatedLagrangianElement)
Pointer definition of ThreeStepUpdatedLagrangianElement.
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &rCurrentProcessInfo) const override
Provides the global indices for each one of this element's local rows.
Definition: three_step_updated_lagrangian_element.cpp:32
void InitializeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override
Initializes the element and all geometric information required for the problem.
Definition: three_step_updated_lagrangian_element.h:197
ConstitutiveLaw::Pointer mpConstitutiveLaw
Definition: three_step_updated_lagrangian_element.h:348
std::vector< std::size_t > EquationIdVectorType
Definition: three_step_updated_lagrangian_element.h:92
BaseType::PropertiesType PropertiesType
Definition: three_step_updated_lagrangian_element.h:113
void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_updated_lagrangian_element.h:199
virtual ~ThreeStepUpdatedLagrangianElement()
Destructor.
Definition: three_step_updated_lagrangian_element.h:166
ThreeStepUpdatedLagrangianElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: three_step_updated_lagrangian_element.h:138
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: updated_lagrangian_element.hpp:177
void ComputeBoundaryTermsForPressureSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ShapeFunctionsType &rN, const double lagMultiplier)
std::size_t IndexType
The definition of the index type.
Definition: key_hash.h:35
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
ModelPart::NodesContainerType NodesArrayType
Definition: gid_gauss_point_container.h:42
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
def load(f)
Definition: ode_solve.py:307
Definition: updated_lagrangian_element.h:74