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.
two_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: April 2018 $
5 // Revision: $Revision: 0.0 $
6 //
7 //
8 
9 #if !defined(KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_ELEMENT_H_INCLUDED)
10 #define KRATOS_TWO_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 
65  public:
68 
71 
73 
75  typedef Node NodeType;
76 
79 
82 
84  typedef Vector VectorType;
85 
87  typedef Matrix MatrixType;
88 
89  typedef std::size_t IndexType;
90 
91  typedef std::size_t SizeType;
92 
93  typedef std::vector<std::size_t> EquationIdVectorType;
94 
95  typedef std::vector<Dof<double>::Pointer> DofsVectorType;
96 
98 
101 
104 
107 
110 
112  typedef ConstitutiveLawType::Pointer ConstitutiveLawPointerType;
113 
115 
116  typedef typename BaseType::PropertiesType::Pointer pPropertiesType;
117 
118  typedef typename BaseType::ElementalVariables ElementalVariables;
119 
123 
124  // Constructors.
125 
127 
131  {
132  }
133 
135 
139  TwoStepUpdatedLagrangianElement(IndexType NewId, const NodesArrayType &ThisNodes) : BaseType(NewId, ThisNodes)
140  {
141  }
142 
144 
148  TwoStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry) : BaseType(NewId, pGeometry)
149  {
150  }
151 
153 
158  TwoStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties) : BaseType(NewId, pGeometry, pProperties)
159  {
160  }
161 
163 
165 
168  {
169  }
170 
174 
178 
180 
187  Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes,
188  pPropertiesType pProperties) const override
189  {
190  return Element::Pointer(new TwoStepUpdatedLagrangianElement(NewId, this->GetGeometry().Create(ThisNodes), pProperties));
191  }
192 
193  Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override;
194 
195  void Initialize(const ProcessInfo &rCurrentProcessInfo) override{};
196 
198  void InitializeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override{};
199 
200  void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override{};
201 
203  void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix,
204  VectorType &rRightHandSideVector,
205  const ProcessInfo &rCurrentProcessInfo) override{};
206 
207  void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix,
208  const ProcessInfo &rCurrentProcessInfo) override
209  {
210  KRATOS_TRY;
211  KRATOS_THROW_ERROR(std::logic_error, "TwoStepUpdatedLagrangianElement::CalculateLeftHandSide not implemented", "");
212  KRATOS_CATCH("");
213  }
214 
215  void CalculateRightHandSide(VectorType &rRightHandSideVector,
216  const ProcessInfo &rCurrentProcessInfo) override{};
217 
218  // The following methods have different implementations depending on TDim
220 
227  const ProcessInfo &rCurrentProcessInfo) const override;
228 
230 
234  void GetDofList(DofsVectorType &rElementalDofList,
235  const ProcessInfo &rCurrentProcessInfo) const override;
236 
237  virtual void UpdateCauchyStress(unsigned int g, const ProcessInfo &rCurrentProcessInfo) override{};
238 
239  virtual void InitializeElementalVariables(ElementalVariables &rElementalVariables) override{};
240 
244 
248 
250 
258  int Check(const ProcessInfo &rCurrentProcessInfo) const override { return 0; };
259 
263 
267 
269  std::string Info() const override
270  {
271  std::stringstream buffer;
272  buffer << "TwoStepUpdatedLagrangianElement #" << this->Id();
273  return buffer.str();
274  }
275 
277  void PrintInfo(std::ostream &rOStream) const override {}
278 
279  // /// Print object's data.
280  // virtual void PrintData(std::ostream& rOStream) const;
281 
285 
287  protected:
290 
294 
295  ConstitutiveLaw::Pointer mpConstitutiveLaw = nullptr;
296 
300 
304 
305  virtual void CalculateLocalMomentumEquations(MatrixType &rLeftHandSideMatrix,
306  VectorType &rRightHandSideVector,
307  const ProcessInfo &rCurrentProcessInfo) override{};
308 
309  virtual void CalculateLocalContinuityEqForPressure(MatrixType &rLeftHandSideMatrix,
310  VectorType &rRightHandSideVector,
311  const ProcessInfo &rCurrentProcessInfo) override{};
312 
313  virtual void CalculateExplicitContinuityEquation(MatrixType &rLeftHandSideMatrix,
314  VectorType &rRightHandSideVector,
315  const ProcessInfo &rCurrentProcessInfo) override{};
316 
317  virtual double GetThetaMomentum() override
318  {
319  std::cout << "I SHOULD NOT ENTER HERE!" << std::endl;
320  return 0.0;
321  };
322 
323  virtual double GetThetaContinuity() override
324  {
325  std::cout << "I SHOULD NOT ENTER HERE!" << std::endl;
326  return 0.0;
327  };
328 
329  void CalculateMassMatrix(Matrix &rMassMatrix,
330  const ProcessInfo &rCurrentProcessInfo) override{};
331 
332  void ComputeMassMatrix(Matrix &rMassMatrix,
333  const ShapeFunctionsType &rN,
334  const double Weight,
335  double &MeanValue) override;
336 
337  void ComputeLumpedMassMatrix(Matrix &rMassMatrix,
338  const double Weight,
339  double &MeanValue) override;
340 
341  void AddExternalForces(Vector &rRHSVector,
342  const double Density,
343  const ShapeFunctionsType &rN,
344  const double Weight) override;
345 
346  void AddInternalForces(Vector &rRHSVector,
347  const ShapeFunctionDerivativesType &rDN_DX,
348  ElementalVariables &rElementalVariables,
349  const double Weight) override;
350 
351  void ComputeBulkMatrixLump(MatrixType &BulkMatrix,
352  const double Weight) override;
353 
354  virtual void ComputeBulkMatrixRHS(MatrixType &BulkMatrix,
355  const double Weight) override{};
356 
358  ElementalVariables &rElementalVariables,
359  const unsigned int g,
360  const Vector& rN,
361  const ProcessInfo &rCurrentProcessInfo,
362  double &Density,
363  double &DeviatoricCoeff,
364  double &VolumetricCoeff) override
365  {};
366 
370 
374 
378 
380 
381  private:
384 
388 
392 
393  friend class Serializer;
394 
395  void save(Serializer &rSerializer) const override
396  {
398  }
399 
400  void load(Serializer &rSerializer) override
401  {
403  }
404 
408 
412 
416 
420 
424 
427 
428  /* /// Copy constructor. */
429  /* TwoStepUpdatedLagrangianElement(TwoStepUpdatedLagrangianElement const& rOther); */
430 
432 
433  }; // Class TwoStepUpdatedLagrangianElement
434 
436 
439 
443 
445  template <unsigned int TDim>
446  inline std::istream &operator>>(std::istream &rIStream,
448  {
449  return rIStream;
450  }
451 
453  template <unsigned int TDim>
454  inline std::ostream &operator<<(std::ostream &rOStream,
456  {
457  rThis.PrintInfo(rOStream);
458  rOStream << std::endl;
459  rThis.PrintData(rOStream);
460 
461  return rOStream;
462  }
463 
464 } // namespace Kratos.
465 
466 #endif // KRATOS_TWO_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: two_step_updated_lagrangian_element.h:63
A stabilized element for the incompressible Navier-Stokes equations.
Definition: updated_lagrangian_element.h:64
#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
void AddExternalForces(Vector &rRHSVector, const double Density, const ShapeFunctionsType &rN, const double Weight) override
Definition: two_step_updated_lagrangian_element.cpp:171
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: two_step_updated_lagrangian_element.h:81
void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_element.h:200
TwoStepUpdatedLagrangianElement(TwoStepUpdatedLagrangianElement const &rOther)
copy constructor
Definition: two_step_updated_lagrangian_element.h:164
TwoStepUpdatedLagrangianElement(IndexType NewId=0)
Default constuctor.
Definition: two_step_updated_lagrangian_element.h:130
virtual void UpdateCauchyStress(unsigned int g, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_element.h:237
BaseType::PropertiesType::Pointer pPropertiesType
Definition: two_step_updated_lagrangian_element.h:116
void ComputeLumpedMassMatrix(Matrix &rMassMatrix, const double Weight, double &MeanValue) override
Definition: two_step_updated_lagrangian_element.cpp:84
virtual void ComputeBulkMatrixRHS(MatrixType &BulkMatrix, const double Weight) override
Definition: two_step_updated_lagrangian_element.h:354
Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override
It creates a new element pointer and clones the previous element data.
Definition: two_step_updated_lagrangian_element.cpp:21
void InitializeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override
Initializes the element and all geometric information required for the problem.
Definition: two_step_updated_lagrangian_element.h:198
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: two_step_updated_lagrangian_element.h:106
ConstitutiveLaw::Pointer mpConstitutiveLaw
Definition: two_step_updated_lagrangian_element.h:295
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: two_step_updated_lagrangian_element.h:97
virtual ~TwoStepUpdatedLagrangianElement()
Destructor.
Definition: two_step_updated_lagrangian_element.h:167
virtual double GetThetaContinuity() override
Definition: two_step_updated_lagrangian_element.h:323
UpdatedLagrangianElement< TDim > BaseType
Definition: two_step_updated_lagrangian_element.h:72
TwoStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using a geometry object.
Definition: two_step_updated_lagrangian_element.h:148
virtual void CalculateLocalMomentumEquations(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_element.h:305
ConstitutiveLawType::Pointer ConstitutiveLawPointerType
Pointer type for constitutive laws.
Definition: two_step_updated_lagrangian_element.h:112
BaseType::PropertiesType PropertiesType
Definition: two_step_updated_lagrangian_element.h:114
virtual void CalculateLocalContinuityEqForPressure(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_element.h:309
void ComputeMassMatrix(Matrix &rMassMatrix, const ShapeFunctionsType &rN, const double Weight, double &MeanValue) override
Definition: two_step_updated_lagrangian_element.cpp:139
virtual void CalcElasticPlasticCauchySplitted(ElementalVariables &rElementalVariables, const unsigned int g, const Vector &rN, const ProcessInfo &rCurrentProcessInfo, double &Density, double &DeviatoricCoeff, double &VolumetricCoeff) override
Definition: two_step_updated_lagrangian_element.h:357
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &rCurrentProcessInfo) const override
Provides the global indices for each one of this element's local rows.
Definition: two_step_updated_lagrangian_element.cpp:32
BaseType::ElementalVariables ElementalVariables
Definition: two_step_updated_lagrangian_element.h:118
Vector VectorType
Vector type for local contributions to the linear system.
Definition: two_step_updated_lagrangian_element.h:84
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: two_step_updated_lagrangian_element.h:95
ConstitutiveLaw ConstitutiveLawType
Reference type definition for constitutive laws.
Definition: two_step_updated_lagrangian_element.h:109
void GetDofList(DofsVectorType &rElementalDofList, const ProcessInfo &rCurrentProcessInfo) const override
Returns a list of the element's Dofs.
Definition: two_step_updated_lagrangian_element.cpp:58
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(TwoStepUpdatedLagrangianElement)
Pointer definition of TwoStepUpdatedLagrangianElement.
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_element.h:207
virtual double GetThetaMomentum() override
Definition: two_step_updated_lagrangian_element.h:317
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, pPropertiesType pProperties) const override
Create a new element of this type.
Definition: two_step_updated_lagrangian_element.h:187
virtual void InitializeElementalVariables(ElementalVariables &rElementalVariables) override
Definition: two_step_updated_lagrangian_element.h:239
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: two_step_updated_lagrangian_element.h:78
std::size_t SizeType
Definition: two_step_updated_lagrangian_element.h:91
Element::PropertiesType PropertiesType
Definition: updated_lagrangian_element.h:142
virtual void CalculateExplicitContinuityEquation(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_element.h:313
std::vector< std::size_t > EquationIdVectorType
Definition: two_step_updated_lagrangian_element.h:93
void ComputeBulkMatrixLump(MatrixType &BulkMatrix, const double Weight) override
Definition: two_step_updated_lagrangian_element.cpp:292
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: two_step_updated_lagrangian_element.h:100
void CalculateMassMatrix(Matrix &rMassMatrix, const ProcessInfo &rCurrentProcessInfo) override
Add integration point contribution to the mass matrix.
Definition: two_step_updated_lagrangian_element.h:329
void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_element.h:195
TwoStepUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties)
Constuctor using geometry and properties.
Definition: two_step_updated_lagrangian_element.h:158
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Calculate the element's local contribution to the system for the current step.
Definition: two_step_updated_lagrangian_element.h:203
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: two_step_updated_lagrangian_element.h:87
void AddInternalForces(Vector &rRHSVector, const ShapeFunctionDerivativesType &rDN_DX, ElementalVariables &rElementalVariables, const double Weight) override
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: two_step_updated_lagrangian_element.h:277
TwoStepUpdatedLagrangianElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: two_step_updated_lagrangian_element.h:139
Node NodeType
Node type (default is: Node)
Definition: two_step_updated_lagrangian_element.h:75
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: updated_lagrangian_element.hpp:177
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: two_step_updated_lagrangian_element.h:103
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Checks the input and that all required Kratos variables have been registered.
Definition: two_step_updated_lagrangian_element.h:258
std::size_t IndexType
Definition: two_step_updated_lagrangian_element.h:89
std::string Info() const override
Turn back information as a string.
Definition: two_step_updated_lagrangian_element.h:269
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_element.h:215
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