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_first_order_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_FIRST_ORDER_UPDATED_LAGRANGIAN_ELEMENT_H_INCLUDED)
10 #define KRATOS_THREE_STEP_FIRST_ORDER_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 
118 
122 
123  //Constructors.
124 
126 
130  {
131  }
132 
134 
139  {
140  }
141 
143 
147  ThreeStepFirstOrderUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry) : BaseType(NewId, pGeometry)
148  {
149  }
150 
152 
157  ThreeStepFirstOrderUpdatedLagrangianElement(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 ThreeStepFirstOrderUpdatedLagrangianElement(NewId, this->GetGeometry().Create(ThisNodes), pProperties));
190  }
191 
192  Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override;
193 
195  void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix,
196  VectorType &rRightHandSideVector,
197  const ProcessInfo &rCurrentProcessInfo) override;
198 
199  void CalculateFirstVelocitySystem(MatrixType &rLeftHandSideMatrix,
200  VectorType &rRightHandSideVector,
201  const ProcessInfo &rCurrentProcessInfo);
202 
208  void Calculate(const Variable<array_1d<double, 3>> &rVariable,
209  array_1d<double, 3> &rOutput,
210  const ProcessInfo &rCurrentProcessInfo) override;
211 
215 
219 
223 
227 
229  std::string Info() const override
230  {
231  std::stringstream buffer;
232  buffer << "ThreeStepFirstOrderUpdatedLagrangianElement #" << this->Id();
233  return buffer.str();
234  }
235 
236  // /// Print object's data.
237  // virtual void PrintData(std::ostream& rOStream) const;
238 
242 
244  protected:
247 
251 
252  ConstitutiveLaw::Pointer mpConstitutiveLaw = nullptr;
253 
257 
261 
262  void CalculateStandardFSPressureSystem(MatrixType &rLeftHandSideMatrix,
263  VectorType &rRightHandSideVector,
264  const ProcessInfo &rCurrentProcessInfo);
265 
269 
273 
277 
279 
280  private:
283 
287 
291 
292  friend class Serializer;
293 
294  void save(Serializer &rSerializer) const override
295  {
297  }
298 
299  void load(Serializer &rSerializer) override
300  {
302  }
303 
307 
311 
315 
319 
323 
326 
327  /* /// Copy constructor. */
328  /* ThreeStepFirstOrderUpdatedLagrangianElement(ThreeStepFirstOrderUpdatedLagrangianElement const& rOther); */
329 
331 
332  }; // Class ThreeStepFirstOrderUpdatedLagrangianElement
333 
335 
338 
342 
344  template <unsigned int TDim>
345  inline std::istream &operator>>(std::istream &rIStream,
347  {
348  return rIStream;
349  }
350 
352  template <unsigned int TDim>
353  inline std::ostream &operator<<(std::ostream &rOStream,
355  {
356  rThis.PrintInfo(rOStream);
357  rOStream << std::endl;
358  rThis.PrintData(rOStream);
359 
360  return rOStream;
361  }
362 
363 } // namespace Kratos.
364 
365 #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_first_order_updated_lagrangian_element.h:63
A stabilized element for the incompressible Navier-Stokes equations.
Definition: three_step_updated_lagrangian_element.h:63
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
Node NodeType
Node type (default is: Node)
Definition: three_step_first_order_updated_lagrangian_element.h:74
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: three_step_first_order_updated_lagrangian_element.h:80
ThreeStepFirstOrderUpdatedLagrangianElement(ThreeStepFirstOrderUpdatedLagrangianElement const &rOther)
copy constructor
Definition: three_step_first_order_updated_lagrangian_element.h:163
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: three_step_first_order_updated_lagrangian_element.h:102
ThreeStepUpdatedLagrangianElement< TDim > BaseType
Definition: three_step_first_order_updated_lagrangian_element.h:71
ThreeStepFirstOrderUpdatedLagrangianElement(IndexType NewId=0)
Default constuctor.
Definition: three_step_first_order_updated_lagrangian_element.h:129
ThreeStepFirstOrderUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using a geometry object.
Definition: three_step_first_order_updated_lagrangian_element.h:147
ConstitutiveLawType::Pointer ConstitutiveLawPointerType
Pointer type for constitutive laws.
Definition: three_step_first_order_updated_lagrangian_element.h:111
void CalculateStandardFSPressureSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: three_step_first_order_updated_lagrangian_element.cpp:223
std::size_t IndexType
Definition: three_step_first_order_updated_lagrangian_element.h:88
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: three_step_first_order_updated_lagrangian_element.h:86
std::size_t SizeType
Definition: three_step_first_order_updated_lagrangian_element.h:90
ConstitutiveLaw ConstitutiveLawType
Reference type definition for constitutive laws.
Definition: three_step_first_order_updated_lagrangian_element.h:108
virtual ~ThreeStepFirstOrderUpdatedLagrangianElement()
Destructor.
Definition: three_step_first_order_updated_lagrangian_element.h:166
BaseType::ElementalVariables ElementalVariables
Definition: three_step_updated_lagrangian_element.h:117
void Calculate(const Variable< array_1d< double, 3 >> &rVariable, array_1d< double, 3 > &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: three_step_first_order_updated_lagrangian_element.cpp:160
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_first_order_updated_lagrangian_element.cpp:21
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_first_order_updated_lagrangian_element.cpp:36
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: three_step_first_order_updated_lagrangian_element.h:96
ConstitutiveLaw::Pointer mpConstitutiveLaw
Definition: three_step_first_order_updated_lagrangian_element.h:252
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: three_step_updated_lagrangian_element.h:330
std::string Info() const override
Turn back information as a string.
Definition: three_step_first_order_updated_lagrangian_element.h:229
ThreeStepFirstOrderUpdatedLagrangianElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties)
Constuctor using geometry and properties.
Definition: three_step_first_order_updated_lagrangian_element.h:157
void CalculateFirstVelocitySystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: three_step_first_order_updated_lagrangian_element.cpp:65
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: three_step_first_order_updated_lagrangian_element.h:94
BaseType::PropertiesType PropertiesType
Definition: three_step_first_order_updated_lagrangian_element.h:113
BaseType::ElementalVariables ElementalVariables
Definition: three_step_first_order_updated_lagrangian_element.h:117
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(ThreeStepFirstOrderUpdatedLagrangianElement)
Pointer definition of ThreeStepFirstOrderUpdatedLagrangianElement.
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, pPropertiesType pProperties) const override
Create a new element of this type.
Definition: three_step_first_order_updated_lagrangian_element.h:186
BaseType::PropertiesType::Pointer pPropertiesType
Definition: three_step_first_order_updated_lagrangian_element.h:115
BaseType::PropertiesType PropertiesType
Definition: three_step_updated_lagrangian_element.h:113
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: three_step_first_order_updated_lagrangian_element.h:77
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: three_step_first_order_updated_lagrangian_element.h:105
ThreeStepFirstOrderUpdatedLagrangianElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: three_step_first_order_updated_lagrangian_element.h:138
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: updated_lagrangian_element.hpp:177
std::vector< std::size_t > EquationIdVectorType
Definition: three_step_first_order_updated_lagrangian_element.h:92
Vector VectorType
Vector type for local contributions to the linear system.
Definition: three_step_first_order_updated_lagrangian_element.h:83
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: three_step_first_order_updated_lagrangian_element.h:99
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