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_V_P_implicit_nodally_integrated_element.h
Go to the documentation of this file.
1 //
2 // Project Name: KratosFluidDynamicsApplication $
3 // Last modified by: $Author: AFranci $
4 // Date: $Date: January 2016 $
5 // Revision: $Revision: 0.0 $
6 //
7 //
8 
9 #if !defined(KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_V_P_IMPLICIT_NODALLY_INTEGRATED_ELEMENT_H_INCLUDED)
10 #define KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_V_P_IMPLICIT_NODALLY_INTEGRATED_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/element.h" */
22 #include "includes/serializer.h"
23 #include "geometries/geometry.h"
24 #include "utilities/math_utils.h"
25 #include "utilities/geometry_utilities.h"
26 
28 
30 
31 #include "includes/model_part.h"
32 /* #include "includes/node.h" */
33 
34 namespace Kratos
35 {
36 
39 
42 
46 
47  typedef GlobalPointersVector<Element> ElementWeakPtrVectorType;
48 
52 
56 
60 
62 
65  template <unsigned int TDim>
67  {
68  public:
71 
74 
77 
79  typedef Node NodeType;
80 
83 
86 
89 
92 
93  typedef std::size_t IndexType;
94 
95  typedef std::size_t SizeType;
96 
97  typedef std::vector<std::size_t> EquationIdVectorType;
98 
99  typedef std::vector<Dof<double>::Pointer> DofsVectorType;
100 
102 
105 
108 
111 
112  /* typedef Element::PropertiesType::Pointer PropertiesType::Pointer; */
113 
115 
116  typedef typename BaseType::PropertiesType::Pointer pPropertiesType;
117 
119 
123 
124  //Constructors.
125 
127 
131  {
132  }
133 
135 
140  {
141  }
142 
144 
148  TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement(IndexType NewId, GeometryType::Pointer pGeometry) : BaseType(NewId, pGeometry)
149  {
150  }
151 
153 
158  TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties) : BaseType(NewId, pGeometry, pProperties)
159  {
160  }
161 
163 
165  {
166  }
167 
170  {
171  }
172 
176 
180 
182 
189  Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes,
190  pPropertiesType pProperties) const override
191  {
192  return Element::Pointer(new TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement(NewId, this->GetGeometry().Create(ThisNodes), pProperties));
193  }
194 
195  Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override;
196 
197  void Initialize(const ProcessInfo &rCurrentProcessInfo) override{};
198 
200  void InitializeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override{};
201 
202  void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override;
203 
204  void CalcElementalStrains(ElementalVariables &rElementalVariables,
205  const ProcessInfo &rCurrentProcessInfo,
206  const ShapeFunctionDerivativesType &rDN_DX);
207 
208  void GetNodesPosition(Vector &rValues,
209  const ProcessInfo &rCurrentProcessInfo,
210  double theta);
211 
213  Matrix &rNContainer,
214  Vector &rGaussWeights) override;
215 
216  void InitializeElementalVariables(ElementalVariables &rElementalVariables) override
217  {
218  KRATOS_TRY;
219  KRATOS_THROW_ERROR(std::logic_error, "InitializeElementalVariables", "");
220  KRATOS_CATCH("");
221  };
222 
223  void CalculateLocalMomentumEquations(MatrixType &rLeftHandSideMatrix,
224  VectorType &rRightHandSideVector,
225  const ProcessInfo &rCurrentProcessInfo) override{};
226 
228 
236  int Check(const ProcessInfo &rCurrentProcessInfo) const override;
237 
238  void CalculateElementalLaplacian(MatrixType &rLeftHandSideMatrix,
239  VectorType &rRightHandSideVector,
240  const ProcessInfo &rCurrentProcessInfo);
241 
242  void CalculateElementalLaplacianAndTau(MatrixType &rLeftHandSideMatrix,
243  VectorType &rRightHandSideVector,
244  const ProcessInfo &rCurrentProcessInfo);
245 
246  void CalculateVolumetricStabilizedTerms(MatrixType &rLeftHandSideMatrix,
247  VectorType &rRightHandSideVector,
248  const ProcessInfo &rCurrentProcessInfo);
249 
250  void CalculateElementalContinuityEqForPressure(MatrixType &rLeftHandSideMatrix,
251  VectorType &rRightHandSideVector,
252  const ProcessInfo &rCurrentProcessInfo);
253 
254  void CalculateLocalContinuityEqForPressure(MatrixType &rLeftHandSideMatrix,
255  VectorType &rRightHandSideVector,
256  const ProcessInfo &rCurrentProcessInfo) override;
257 
259  VectorType &rRightHandSideVector,
260  const ProcessInfo &rCurrentProcessInfo);
261 
262  void ComputeBoundLHSMatrix(MatrixType &BoundLHSMatrix,
263  const ShapeFunctionsType &rN,
264  const double Weight) override;
265 
267  const double TimeStep,
268  const double BoundRHSCoeffAcc,
269  const double BoundRHSCoeffDev,
270  const VectorType SpatialDefRate);
271 
273  const double TimeStep,
274  const double BoundRHSCoeffAcc,
275  const double BoundRHSCoeffDev);
276 
277  void ComputeStabLaplacianMatrix(MatrixType &StabLaplacianMatrix,
278  const ShapeFunctionDerivativesType &rShapeDeriv,
279  const double Weight) override;
280 
281  void CalculateTauFIC(double &TauOne,
282  double ElemSize,
283  const double Density,
284  const double Viscosity,
285  const ProcessInfo &rCurrentProcessInfo) override;
286 
287  void AddStabilizationNodalTermsRHS(VectorType &rRightHandSideVector,
288  const double Tau,
289  const double Density,
290  const double Weight,
291  const ShapeFunctionDerivativesType &rDN_DX,
292  const SizeType i) override;
293 
297 
301 
303  std::string Info() const override
304  {
305  std::stringstream buffer;
306  buffer << "TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement #" << BaseType::Id();
307  return buffer.str();
308  }
309 
311  void PrintInfo(std::ostream &rOStream) const override
312  {
313  rOStream << "TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement" << TDim << "D";
314  }
315 
316  protected:
317  void NodalFreeSurfaceLength(unsigned int nodeIndex);
318 
319  private:
320  friend class Serializer;
321 
322  void save(Serializer &rSerializer) const override
323  {
325  }
326 
327  void load(Serializer &rSerializer) override
328  {
330  }
331 
335 
339 
343 
347 
351 
354 
355  /* /// Copy constructor. */
356  /* TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement(TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement const& rOther); */
357 
359 
360  }; // Class TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement
361 
363 
366 
370 
372  template <unsigned int TDim>
373  inline std::istream &operator>>(std::istream &rIStream,
375  {
376  return rIStream;
377  }
378 
380  template <unsigned int TDim>
381  inline std::ostream &operator<<(std::ostream &rOStream,
383  {
384  rThis.PrintInfo(rOStream);
385  rOStream << std::endl;
386  rThis.PrintData(rOStream);
387 
388  return rOStream;
389  }
390 
391 } // namespace Kratos.
392 
393 #endif // KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_V_P_IMPLICIT_NODALLY_INTEGRATED_ELEMENT defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Base class for all Elements.
Definition: element.h:60
std::size_t SizeType
Definition: element.h:94
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_V_P_implicit_element.h:63
A stabilized element for the incompressible Navier-Stokes equations.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:67
#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
BaseType::PropertiesType::Pointer pPropertiesType
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:116
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:311
void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:32
void CalculateVolumetricStabilizedTerms(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:882
virtual ~TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement()
Destructor.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:169
TwoStepUpdatedLagrangianVPImplicitElement< TDim > BaseType
base type:
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:76
std::size_t IndexType
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:93
BaseType::PropertiesType::Pointer pPropertiesType
Definition: two_step_updated_lagrangian_element.h:116
std::string Info() const override
Turn back information as a string.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:303
void CalculateGeometryData(ShapeFunctionDerivativesArrayType &rDN_DX, Matrix &rNContainer, Vector &rGaussWeights) override
Determine integration point weights and shape funcition derivatives from the element's geometry.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:686
Kratos::Vector VectorType
Vector type for local contributions to the linear system.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:88
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement)
Pointer definition of TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement.
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:99
void CalculateTauFIC(double &TauOne, double ElemSize, const double Density, const double Viscosity, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:1269
void InitializeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override
Initializes the element and all geometric information required for the problem.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:200
void CalculateElementalContinuityEqForPressure(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:951
void CalculateElementalLaplacian(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:770
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:101
BaseType::PropertiesType PropertiesType
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:114
void ComputeElementalBoundRHSVector(VectorType &BoundRHSVector, const double TimeStep, const double BoundRHSCoeffAcc, const double BoundRHSCoeffDev)
TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement(IndexType NewId=0)
Default constuctor.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:130
void ComputeBoundRHSVectorComplete(VectorType &BoundRHSVector, const double TimeStep, const double BoundRHSCoeffAcc, const double BoundRHSCoeffDev, const VectorType SpatialDefRate)
void ComputeStabLaplacianMatrix(MatrixType &StabLaplacianMatrix, const ShapeFunctionDerivativesType &rShapeDeriv, const double Weight) override
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:2059
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:110
std::size_t SizeType
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:95
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:104
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:85
void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:197
BaseType::PropertiesType PropertiesType
Definition: two_step_updated_lagrangian_V_P_implicit_element.h:111
void CalcElementalStrains(ElementalVariables &rElementalVariables, const ProcessInfo &rCurrentProcessInfo, const ShapeFunctionDerivativesType &rDN_DX)
void CalculateElementalLaplacianAndTau(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:817
Node NodeType
Node type (default is: Node)
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:79
void ComputeBoundLHSMatrix(MatrixType &BoundLHSMatrix, const ShapeFunctionsType &rN, const double Weight) override
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Checks the input and that all required Kratos variables have been registered.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:708
void CalculateLocalMomentumEquations(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:223
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:107
TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using a geometry object.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:148
void InitializeElementalVariables(ElementalVariables &rElementalVariables) override
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:216
BaseType::ElementalVariables ElementalVariables
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:118
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:82
Kratos::Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:91
TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties)
Constuctor using geometry and properties.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:158
void AddStabilizationNodalTermsRHS(VectorType &rRightHandSideVector, const double Tau, const double Density, const double Weight, const ShapeFunctionDerivativesType &rDN_DX, const SizeType i) override
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:2081
BaseType::ElementalVariables ElementalVariables
Definition: two_step_updated_lagrangian_V_P_implicit_element.h:115
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_V_P_implicit_nodally_integrated_element.cpp:21
void CalculateStabilizingTermsContinuityEqForPressure(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:1186
void GetNodesPosition(Vector &rValues, const ProcessInfo &rCurrentProcessInfo, double theta)
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, pPropertiesType pProperties) const override
Create a new element of this type.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:189
std::vector< std::size_t > EquationIdVectorType
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:97
void CalculateLocalContinuityEqForPressure(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.cpp:1103
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: updated_lagrangian_element.hpp:177
GlobalPointersVector< Element > ElementWeakPtrVectorType
Definition: build_model_part_boundary_process.hpp:60
TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:139
TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement(TwoStepUpdatedLagrangianVPImplicitNodallyIntegratedElement const &rOther)
copy constructor
Definition: two_step_updated_lagrangian_V_P_implicit_nodally_integrated_element.h:164
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
integer i
Definition: TensorModule.f:17
Definition: updated_lagrangian_element.h:74