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_fluid_element.h
Go to the documentation of this file.
1 //
2 // Project Name: KratosFluidDynamicsApplication $
3 // Last modified by: $Author: AFranci $
4 // Date: $Date: February 2016 $
5 // Revision: $Revision: 0.0 $
6 //
7 //
8 
9 #if !defined(KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_V_P_IMPLICIT_FLUID_ELEMENT_H_INCLUDED)
10 #define KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_V_P_IMPLICIT_FLUID_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"
26 
28 
29 namespace Kratos
30 {
31 
34 
37 
41 
45 
49 
53 
55 
57  template <unsigned int TDim>
59  {
60 
61  public:
64 
67 
70 
72  typedef Node NodeType;
73 
76 
79 
81  typedef Vector VectorType;
82 
84  typedef Matrix MatrixType;
85 
86  typedef std::size_t IndexType;
87 
88  typedef std::size_t SizeType;
89 
90  typedef std::vector<std::size_t> EquationIdVectorType;
91 
92  typedef std::vector<Dof<double>::Pointer> DofsVectorType;
93 
95 
98 
101 
104 
106 
107  typedef typename BaseType::PropertiesType::Pointer pPropertiesType;
108 
110 
112 
115 
117  typedef ConstitutiveLawType::Pointer ConstitutiveLawPointerType;
118 
122 
123  // Constructors.
124 
126 
130  {
131  }
132 
134 
139  {
140  }
141 
143 
147  TwoStepUpdatedLagrangianVPImplicitFluidElement(IndexType NewId, GeometryType::Pointer pGeometry) : BaseType(NewId, pGeometry)
148  {
149  }
150 
152 
157  TwoStepUpdatedLagrangianVPImplicitFluidElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties) : BaseType(NewId, pGeometry, pProperties)
158  {
159  }
160 
162 
164  {
165  }
166 
169  {
170  }
171 
175 
179 
181 
188  Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes,
189  pPropertiesType pProperties) const override
190  {
191  return Element::Pointer(new TwoStepUpdatedLagrangianVPImplicitFluidElement(NewId, BaseType::GetGeometry().Create(ThisNodes), pProperties));
192  }
193 
194  Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override;
195 
196  void Initialize(const ProcessInfo &rCurrentProcessInfo) override;
197 
199  void InitializeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override{};
200 
201  void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override{};
202 
203  void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix,
204  const ProcessInfo &rCurrentProcessInfo) override
205  {
206  KRATOS_TRY;
207  KRATOS_THROW_ERROR(std::logic_error, "TwoStepUpdatedLagrangianVPImplicitFluidElement::CalculateLeftHandSide not implemented", "");
208  KRATOS_CATCH("");
209  }
210 
211  void CalculateRightHandSide(VectorType &rRightHandSideVector,
212  const ProcessInfo &rCurrentProcessInfo) override
213  {
214  KRATOS_TRY;
215  KRATOS_THROW_ERROR(std::logic_error, "TwoStepUpdatedLagrangianVPImplicitFluidElement::CalculateRightHandSide not implemented", "");
216  KRATOS_CATCH("");
217  }
218 
219  // The following methods have different implementations depending on TDim
221 
229 
234  void InitializeElementalVariables(ElementalVariables &rElementalVariables) override;
235 
239 
243 
245 
253  int Check(const ProcessInfo &rCurrentProcessInfo) const override;
254 
258 
262 
264  std::string Info() const override
265  {
266  std::stringstream buffer;
267  buffer << "TwoStepUpdatedLagrangianVPImplicitFluidElement #" << BaseType::Id();
268  return buffer.str();
269  }
270 
272  void PrintInfo(std::ostream &rOStream) const override
273  {
274  rOStream << "TwoStepUpdatedLagrangianVPImplicitFluidElement" << TDim << "D";
275  }
276 
277  // /// Print object's data.
278  // virtual void PrintData(std::ostream& rOStream) const;
279 
283 
285 
286  protected:
289 
293 
294  ConstitutiveLaw::Pointer mpConstitutiveLaw = nullptr;
295 
299 
303 
312  double &MeanValue,
313  const ShapeFunctionDerivativesType &rShapeDeriv,
314  const double secondLame,
315  double &bulkModulus,
316  const double Weight,
317  double &MeanValueMass,
318  const double TimeStep);
319 
321  MatrixType StiffnessMatrix,
322  double &meanValueStiff,
323  double &bulkCoefficient,
324  double timeStep) override;
325 
326  void ComputeBulkMatrixLump(MatrixType &BulkMatrix,
327  const double Weight) override;
328 
330  ElementalVariables &rElementalVariables,
331  const unsigned int g,
332  const Vector& rN,
333  const ProcessInfo &rCurrentProcessInfo,
334  double &Density,
335  double &DeviatoricCoeff,
336  double &VolumetricCoeff) override;
337 
338  double GetThetaMomentum() override { return 0.5; };
339 
340  double GetThetaContinuity() override { return 1.0; };
341 
342  void UpdateStressTensor(ElementalVariables &rElementalVariables);
343 
344  void SetYieldedElements(ElementalVariables &rElementalVariables);
345 
349 
353 
357 
359 
360  private:
363 
367 
371 
372  friend class Serializer;
373 
374  void save(Serializer &rSerializer) const override
375  {
377  }
378 
379  void load(Serializer &rSerializer) override
380  {
382  }
383 
387 
391 
395 
399 
403 
406 
407  /* /// Copy constructor. */
408  /* TwoStepUpdatedLagrangianVPImplicitFluidElement(TwoStepUpdatedLagrangianVPImplicitFluidElement const& rOther); */
409 
411 
412  }; // Class TwoStepUpdatedLagrangianVPImplicitFluidElement
413 
415 
418 
422 
424  template <unsigned int TDim>
425  inline std::istream &operator>>(std::istream &rIStream,
427  {
428  return rIStream;
429  }
430 
432  template <unsigned int TDim>
433  inline std::ostream &operator<<(std::ostream &rOStream,
435  {
436  rThis.PrintInfo(rOStream);
437  rOStream << std::endl;
438  rThis.PrintData(rOStream);
439 
440  return rOStream;
441  }
442 
443 } // namespace Kratos.
444 
445 #endif // KRATOS_TWO_STEP_UPDATED_LAGRANGIAN_V_P_FLUID_ELEMENT defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Definition: constitutive_law.h:47
Base class for all Elements.
Definition: element.h:60
virtual void MassMatrix(MatrixType &rMassMatrix, const ProcessInfo &rCurrentProcessInfo)
Definition: element.h:926
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 class is a vector which stores global pointers.
Definition: global_pointers_vector.h:61
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_fluid_element.h:59
#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
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:103
GlobalPointersVector< Node > NodeWeakPtrVectorType
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:111
void UpdateStressTensor(ElementalVariables &rElementalVariables)
BaseType::PropertiesType::Pointer pPropertiesType
Definition: two_step_updated_lagrangian_element.h:116
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_fluid_element.cpp:56
BaseType::PropertiesType::Pointer pPropertiesType
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:107
void CalcElasticPlasticCauchySplitted(ElementalVariables &rElementalVariables, const unsigned int g, const Vector &rN, const ProcessInfo &rCurrentProcessInfo, double &Density, double &DeviatoricCoeff, double &VolumetricCoeff) override
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_fluid_element.h:199
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:272
Vector VectorType
Vector type for local contributions to the linear system.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:81
void SetYieldedElements(ElementalVariables &rElementalVariables)
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:211
virtual ~TwoStepUpdatedLagrangianVPImplicitFluidElement()
Destructor.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:168
TwoStepUpdatedLagrangianVPImplicitFluidElement(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using a geometry object.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:147
ConstitutiveLaw::Pointer mpConstitutiveLaw
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:294
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:75
TwoStepUpdatedLagrangianVPImplicitElement< TDim > BaseType
base type:
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:69
double GetThetaMomentum() override
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:338
double GetThetaContinuity() override
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:340
Node NodeType
Node type (default is: Node)
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:72
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:84
void ComputeBulkMatrixLump(MatrixType &BulkMatrix, const double Weight) override
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.cpp:385
void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:201
TwoStepUpdatedLagrangianVPImplicitFluidElement(IndexType NewId, GeometryType::Pointer pGeometry, pPropertiesType pProperties)
Constuctor using geometry and properties.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:157
TwoStepUpdatedLagrangianVPImplicitFluidElement(TwoStepUpdatedLagrangianVPImplicitFluidElement const &rOther)
copy constructor
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:163
void ComputeMeanValueMaterialTangentMatrix(ElementalVariables &rElementalVariables, double &MeanValue, const ShapeFunctionDerivativesType &rShapeDeriv, const double secondLame, double &bulkModulus, const double Weight, double &MeanValueMass, const double TimeStep)
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:100
void InitializeElementalVariables(ElementalVariables &rElementalVariables) override
Provides the global indices for each one of this element's local rows.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.cpp:407
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:97
BaseType::PropertiesType PropertiesType
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:105
TwoStepUpdatedLagrangianVPImplicitFluidElement(IndexType NewId=0)
Default constuctor.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:129
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_fluid_element.h:188
TwoStepUpdatedLagrangianVPImplicitFluidElement(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:138
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:78
BaseType::PropertiesType PropertiesType
Definition: two_step_updated_lagrangian_V_P_implicit_element.h:111
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:94
void ComputeBulkReductionCoefficient(MatrixType MassMatrix, MatrixType StiffnessMatrix, double &meanValueStiff, double &bulkCoefficient, double timeStep) override
std::string Info() const override
Turn back information as a string.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:264
std::size_t SizeType
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:88
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:203
void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.cpp:38
ConstitutiveLawType::Pointer ConstitutiveLawPointerType
Pointer type for constitutive laws.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:117
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:92
ConstitutiveLaw ConstitutiveLawType
Reference type definition for constitutive laws.
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:114
std::size_t IndexType
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:86
BaseType::ElementalVariables ElementalVariables
Definition: two_step_updated_lagrangian_V_P_implicit_element.h:115
BaseType::ElementalVariables ElementalVariables
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:109
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(TwoStepUpdatedLagrangianVPImplicitFluidElement)
Pointer definition of TwoStepUpdatedLagrangianVPImplicitFluidElement.
std::vector< std::size_t > EquationIdVectorType
Definition: two_step_updated_lagrangian_V_P_implicit_fluid_element.h:90
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: updated_lagrangian_element.hpp:177
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_fluid_element.cpp:26
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