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_fluid_navier_stokes_alpha_method.h
Go to the documentation of this file.
1 // | / |
2 // ' / __| _` | __| _ \ __|
3 // . \ | ( | | ( |\__ `
4 // _|\_\_| \__,_|\__|\___/ ____/
5 // Multi-Physics
6 //
7 // License: BSD License
8 // Kratos default license: kratos/license.txt
9 //
10 // Main author: Uxue Chasco
11 //
12 //
13 
14 #if !defined(KRATOS_TWO_FLUID_NAVIER_STOKES_ALPHA_METHOD)
15 #define KRATOS_TWO_FLUID_NAVIER_STOKES_ALPHA_METHOD
16 
17 // System includes
18 
19 // Project includes
20 #include "includes/define.h"
21 #include "includes/element.h"
22 #include "includes/variables.h"
23 #include "includes/serializer.h"
24 #include "includes/cfd_variables.h"
27 #include "utilities/geometry_utilities.h"
30 
31 namespace Kratos
32 {
33 
34 /*The "TwoFluidNavierStokesAlphaMethod" element is an element based on the Variational Multiscale Stabilization technique (VMS)
35 * which is designed for the solution of a two fluids (air and a non-air one) problems.
36 *
37 * A distinctive feature of the element is the use of 4 LOCAL enrichment functions, which allows to model
38 * a discontinuity in both the pressure field and in its gradient.
39 * The enrichment functions are obtained by duplicating all of the degrees of freedom of the element.
40 * Since the enrichment is performed elementwise, a purely local static condensation step is performed,
41 * meaning that no extra degrees of freedom are added..
42 *
43 * Since a jump in the pressure can be considered, the element shall be able to habdle moderate changes of the viscosity
44 * between the two fluids to be considered*/
45 
48 
52 
56 
60 
64 
65 template< class TElementData >
67 {
68 public:
69 
72 
75 
77 
78  typedef Node NodeType;
81  typedef Vector VectorType;
82  typedef Matrix MatrixType;
83  typedef std::size_t IndexType;
84  typedef std::size_t SizeType;
85  typedef std::vector<std::size_t> EquationIdVectorType;
86  typedef std::vector< Dof<double>::Pointer > DofsVectorType;
91  constexpr static unsigned int Dim = TwoFluidNavierStokes<TElementData>::Dim;
92  constexpr static unsigned int NumNodes = TwoFluidNavierStokes<TElementData>::NumNodes;
95  constexpr static unsigned int StrainSize = (Dim - 1) * 3;
96 
100 
102 
106 
108 
113 
115 
119  TwoFluidNavierStokesAlphaMethod(IndexType NewId, GeometryType::Pointer pGeometry);
120 
122 
127  TwoFluidNavierStokesAlphaMethod(IndexType NewId, GeometryType::Pointer pGeometry, Properties::Pointer pProperties);
128 
131 
135 
139 
141 
148  Element::Pointer Create(IndexType NewId,
149  NodesArrayType const& ThisNodes,
150  Properties::Pointer pProperties) const override;
151 
153 
160  Element::Pointer Create(IndexType NewId,
161  GeometryType::Pointer pGeom,
162  Properties::Pointer pProperties) const override;
163 
165  const Variable<double> &rVariable,
166  std::vector<double> &rOutput,
167  const ProcessInfo &rCurrentProcessInfo) override;
168 
169  void Calculate(
170  const Variable<double> &rVariable,
171  double &rOutput,
172  const ProcessInfo &rCurrentProcessInfo) override;
176 
180 
184 
188 
192 
196 
198 protected:
201 
205 
209 
213 
224  const TElementData &rData,
225  const Vector &rInterfaceWeights,
226  const Matrix &rEnrInterfaceShapeFunctionPos,
227  const Matrix &rEnrInterfaceShapeFunctionNeg,
228  const GeometryType::ShapeFunctionsGradientsType &rInterfaceShapeDerivatives,
229  MatrixType &rKeeTot,
230  VectorType &rRHSeeTot) override;
231 
239  TElementData& rData,
240  MatrixType& rLHS) override;
241 
249  TElementData& rData,
250  VectorType& rRHS) override;
251 
263  TElementData& rData,
264  MatrixType& rV,
265  MatrixType& rH,
266  MatrixType& rKee,
267  VectorType& rRHS_ee) override;
268 
274  void CalculateStrainRate(TElementData& rData) const override;
275 
282  double CalculateArtificialDynamicViscositySpecialization(TElementData &rData) const;
283 
287 
291 
295 
299 
303 
307 
308  friend class Serializer;
309 
310  void save(Serializer& rSerializer) const override;
311 
312  void load(Serializer& rSerializer) override;
313 
315 private:
318 
322 
326 
330 
333 
336 
338 
339 }; // Class TwoFluidNavierStokesAlphaMethod
342 
346 
348 template< class TElementData >
349 inline std::istream& operator >> (std::istream& rIStream,
351 {
352  return rIStream;
353 }
354 
356 template< class TElementData >
357 inline std::ostream& operator <<(std::ostream& rOStream,
359 {
360  rThis.PrintInfo(rOStream);
361  rOStream << std::endl;
362  rThis.PrintData(rOStream);
363 
364  return rOStream;
365 }
367 
368 } // namespace Kratos.
369 
370 #endif // KRATOS_TWO_FLUID_NAVIER_STOKES_ALPHA_METHOD
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
std::size_t IndexType
Definition: flags.h:74
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: fluid_element.hpp:584
This defines the geometrical object, base definition of the element and condition entities.
Definition: geometrical_object.h:58
Geometry base class.
Definition: geometry.h:71
This object defines an indexed object.
Definition: indexed_object.h:54
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
Definition: two_fluid_navier_stokes_alpha_method.h:67
Node NodeType
Definition: two_fluid_navier_stokes_alpha_method.h:78
virtual ~TwoFluidNavierStokesAlphaMethod()
Destructor.
Definition: two_fluid_navier_stokes_alpha_method.cpp:43
std::size_t SizeType
Definition: two_fluid_navier_stokes_alpha_method.h:84
TwoFluidNavierStokesAlphaMethod(IndexType NewId=0)
Default constuctor.
Definition: two_fluid_navier_stokes_alpha_method.cpp:24
void Calculate(const Variable< double > &rVariable, double &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: two_fluid_navier_stokes_alpha_method.cpp:100
Geometry< NodeType > GeometryType
Definition: two_fluid_navier_stokes_alpha_method.h:79
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(TwoFluidNavierStokesAlphaMethod)
Counted pointer of.
void ComputeGaussPointEnrichmentContributions(TElementData &rData, MatrixType &rV, MatrixType &rH, MatrixType &rKee, VectorType &rRHS_ee) override
Computes the pressure enrichment contributions This method computes the pressure enrichment contribut...
Matrix MatrixType
Definition: two_fluid_navier_stokes_alpha_method.h:82
constexpr static unsigned int StrainSize
Definition: two_fluid_navier_stokes_alpha_method.h:95
void CalculateStrainRate(TElementData &rData) const override
Calculate the strain rate In this function we calculate the strain rate at the mid step.
TwoFluidNavierStokes< TElementData >::ShapeFunctionDerivativesType ShapeFunctionDerivativesType
Definition: two_fluid_navier_stokes_alpha_method.h:89
constexpr static unsigned int BlockSize
Definition: two_fluid_navier_stokes_alpha_method.h:93
TwoFluidNavierStokes< TElementData >::ShapeFunctionDerivativesArrayType ShapeFunctionDerivativesArrayType
Definition: two_fluid_navier_stokes_alpha_method.h:90
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, Properties::Pointer pProperties) const override
Create a new element of this type.
Definition: two_fluid_navier_stokes_alpha_method.cpp:49
void ComputeGaussPointRHSContribution(TElementData &rData, VectorType &rRHS) override
Computes the RHS Gaus pt. contribution This method computes the contribution to the RHS of a Gauss pt...
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: two_fluid_navier_stokes_alpha_method.h:87
double CalculateArtificialDynamicViscositySpecialization(TElementData &rData) const
Calculate the artificial dynamic viscosity. In this function we calculate the artificial dynamic visc...
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: two_fluid_navier_stokes_alpha_method.h:86
void CalculateOnIntegrationPoints(const Variable< double > &rVariable, std::vector< double > &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Function to visualize the divergence field.
Definition: two_fluid_navier_stokes_alpha_method.cpp:67
void load(Serializer &rSerializer) override
Definition: two_fluid_navier_stokes_alpha_method.cpp:2024
void PressureGradientStabilization(const TElementData &rData, const Vector &rInterfaceWeights, const Matrix &rEnrInterfaceShapeFunctionPos, const Matrix &rEnrInterfaceShapeFunctionNeg, const GeometryType::ShapeFunctionsGradientsType &rInterfaceShapeDerivatives, MatrixType &rKeeTot, VectorType &rRHSeeTot) override
Computes the enriched LHS/RHS terms associated with the pressure stabilizations at the interface.
Definition: two_fluid_navier_stokes_alpha_method.cpp:1889
TwoFluidNavierStokes< TElementData >::ShapeFunctionsType ShapeFunctionsType
Definition: two_fluid_navier_stokes_alpha_method.h:88
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition: two_fluid_navier_stokes_alpha_method.h:80
constexpr static unsigned int Dim
Definition: two_fluid_navier_stokes_alpha_method.h:91
constexpr static unsigned int NumNodes
Definition: two_fluid_navier_stokes_alpha_method.h:92
constexpr static unsigned int LocalSize
Definition: two_fluid_navier_stokes_alpha_method.h:94
std::vector< std::size_t > EquationIdVectorType
Definition: two_fluid_navier_stokes_alpha_method.h:85
void save(Serializer &rSerializer) const override
Definition: two_fluid_navier_stokes_alpha_method.cpp:2017
void ComputeGaussPointLHSContribution(TElementData &rData, MatrixType &rLHS) override
Computes the LHS Gauss pt. contribution This method computes the contribution to the LHS of a Gauss p...
std::size_t IndexType
Definition: two_fluid_navier_stokes_alpha_method.h:83
Vector VectorType
Definition: two_fluid_navier_stokes_alpha_method.h:81
Definition: two_fluid_navier_stokes.h:67
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
FluidElement< TElementData >::ShapeFunctionsType ShapeFunctionsType
Definition: two_fluid_navier_stokes.h:86
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
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