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.
qs_vms.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 authors: Jordi Cotela
11 //
12 
13 #ifndef KRATOS_QS_VMS_H
14 #define KRATOS_QS_VMS_H
15 
16 #include "includes/define.h"
17 #include "includes/element.h"
18 #include "includes/serializer.h"
19 #include "geometries/geometry.h"
20 
21 #include "includes/cfd_variables.h"
24 
25 namespace Kratos
26 {
27 
30 
33 
37 
41 
45 
49 
50 template< class TElementData >
51 class QSVMS : public FluidElement<TElementData>
52 {
53 public:
56 
59 
62 
64  typedef Node NodeType;
65 
68 
71 
73  typedef Vector VectorType;
74 
76  typedef Matrix MatrixType;
77 
78  typedef std::size_t IndexType;
79 
80  typedef std::size_t SizeType;
81 
82  typedef std::vector<std::size_t> EquationIdVectorType;
83 
84  typedef std::vector< Dof<double>::Pointer > DofsVectorType;
85 
87 
90 
93 
96 
97  constexpr static unsigned int Dim = BaseType::Dim;
98  constexpr static unsigned int NumNodes = BaseType::NumNodes;
99  constexpr static unsigned int BlockSize = BaseType::BlockSize;
100  constexpr static unsigned int LocalSize = BaseType::LocalSize;
101  constexpr static unsigned int StrainSize = BaseType::StrainSize;
102 
106 
107  //Constructors.
108 
110 
113  QSVMS(IndexType NewId = 0);
114 
116 
120  QSVMS(IndexType NewId, const NodesArrayType& ThisNodes);
121 
123 
127  QSVMS(IndexType NewId, GeometryType::Pointer pGeometry);
128 
130 
135  QSVMS(IndexType NewId, GeometryType::Pointer pGeometry, Properties::Pointer pProperties);
136 
138  ~QSVMS() override;
139 
143 
144 
148 
149 
151 
158  Element::Pointer Create(IndexType NewId,
159  NodesArrayType const& ThisNodes,
160  Properties::Pointer pProperties) const override;
161 
163 
170  Element::Pointer Create(IndexType NewId,
171  GeometryType::Pointer pGeom,
172  Properties::Pointer pProperties) const override;
173 
174 
175  void Calculate(
176  const Variable<double>& rVariable,
177  double& rOutput,
178  const ProcessInfo& rCurrentProcessInfo) override;
179 
180 
181  void Calculate(
182  const Variable<array_1d<double, 3 > >& rVariable,
183  array_1d<double, 3 > & rOutput,
184  const ProcessInfo& rCurrentProcessInfo) override;
185 
186 
187  void Calculate(
188  const Variable<Vector >& rVariable,
189  Vector& Output,
190  const ProcessInfo& rCurrentProcessInfo) override;
191 
192  void Calculate(
193  const Variable<Matrix >& rVariable,
194  Matrix& Output,
195  const ProcessInfo& rCurrentProcessInfo) override;
196 
200 
202  Variable<array_1d<double, 3>> const& rVariable,
203  std::vector<array_1d<double, 3>>& rValues,
204  ProcessInfo const& rCurrentProcessInfo) override;
205 
207  Variable<double> const& rVariable,
208  std::vector<double>& rValues,
209  ProcessInfo const& rCurrentProcessInfo) override;
210 
212  Variable<array_1d<double, 6>> const& rVariable,
213  std::vector<array_1d<double, 6>>& rValues,
214  ProcessInfo const& rCurrentProcessInfo) override;
215 
217  Variable<Vector> const& rVariable,
218  std::vector<Vector>& rValues,
219  ProcessInfo const& rCurrentProcessInfo) override;
220 
222  Variable<Matrix> const& rVariable,
223  std::vector<Matrix>& rValues,
224  ProcessInfo const& rCurrentProcessInfo) override;
225 
229 
230  int Check(const ProcessInfo &rCurrentProcessInfo) const override;
231 
235 
236  const Parameters GetSpecifications() const override;
237 
239  std::string Info() const override;
240 
241 
243  void PrintInfo(std::ostream& rOStream) const override;
244 
245 
249 
250 
252 
253 protected:
254 
257 
258 
262 
263 
267 
268 
272 
273  // Protected interface of FluidElement ////////////////////////////////////
274 
276  TElementData& rData,
277  MatrixType& rLHS,
278  VectorType& rRHS) override;
279 
281  TElementData& rData,
282  MatrixType& rLHS) override;
283 
285  TElementData& rData,
286  VectorType& rRHS) override;
287 
288  void AddVelocitySystem(
289  TElementData& rData,
290  MatrixType& rLocalLHS,
291  VectorType& rLocalRHS) override;
292 
293  void AddMassLHS(
294  TElementData& rData,
295  MatrixType& rMassMatrix) override;
296 
297  // This function integrates the traction over a cut. It is only required to implement embedded formulations
298  void AddBoundaryTraction(
299  TElementData& rData,
300  const Vector& rUnitNormal,
301  MatrixType& rLHS,
302  VectorType& rRHS) override;
303 
304  // Implementation details of QSVMS ////////////////////////////////////////
305 
307  TElementData& rData,
308  MatrixType& rMassMatrix);
309 
310  virtual void AddViscousTerm(
311  const TElementData& rData,
313  VectorType& rRHS);
314 
322  KRATOS_DEPRECATED virtual double EffectiveViscosity(
323  TElementData& rData,
324  double ElementSize);
325 
326  virtual void CalculateTau(
327  const TElementData& rData,
328  const array_1d<double,3> &Velocity,
329  double &TauOne,
330  double &TauTwo) const;
331 
332  virtual void CalculateProjections(const ProcessInfo &rCurrentProcessInfo);
333 
334  virtual void MomentumProjTerm(
335  const TElementData& rData,
336  const array_1d<double,3>& rConvectionVelocity,
337  array_1d<double,3>& rMomentumRHS) const;
338 
339  virtual void MassProjTerm(
340  const TElementData& rData,
341  double& rMassRHS) const;
342 
343  virtual void SubscaleVelocity(
344  const TElementData& rData,
345  array_1d<double,3>& rVelocitySubscale) const;
346 
347  virtual void SubscalePressure(
348  const TElementData& rData,
349  double &rPressureSubscale) const;
350 
351  virtual void AlgebraicMomentumResidual(
352  const TElementData& rData,
353  const array_1d<double,3> &rConvectionVelocity,
354  array_1d<double,3>& rResidual) const;
355 
356  virtual void AlgebraicMassResidual(
357  const TElementData& rData,
358  double& rMomentumRes) const;
359 
360  virtual void OrthogonalMomentumResidual(
361  const TElementData& rData,
362  const array_1d<double,3> &rConvectionVelocity,
363  array_1d<double,3>& rResidual) const;
364 
365  virtual void OrthogonalMassResidual(
366  const TElementData& rData,
367  double& rMassRes) const;
368 
372 
373 
377 
378 
382 
384 
385 private:
386 
389 
393 
397 
398  friend class Serializer;
399 
400  void save(Serializer& rSerializer) const override;
401 
402  void load(Serializer& rSerializer) override;
403 
407 
408 
412 
413 
417 
418 
422 
423 
427 
429  QSVMS& operator=(QSVMS const& rOther);
430 
432  QSVMS(QSVMS const& rOther);
433 
435 
436 
437 }; // Class QSVMS
438 
440 
443 
444 
448 
449 
451 template< class TElementData >
452 inline std::istream& operator >>(std::istream& rIStream,
453  QSVMS<TElementData>& rThis)
454 {
455  return rIStream;
456 }
457 
459 template< class TElementData >
460 inline std::ostream& operator <<(std::ostream& rOStream,
461  const QSVMS<TElementData>& rThis)
462 {
463  rThis.PrintInfo(rOStream);
464  rOStream << std::endl;
465  rThis.PrintData(rOStream);
466 
467  return rOStream;
468 }
470 
472 
473 } // namespace Kratos.
474 
475 #endif // KRATOS_QS_VMS_H
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
std::size_t IndexType
Definition: flags.h:74
Large Displacement Lagrangian Element for 3D and 2D geometries. (base class)
Definition: fluid_element.h:61
static constexpr unsigned int Dim
Definition: fluid_element.h:105
static constexpr unsigned int NumNodes
Definition: fluid_element.h:107
static constexpr unsigned int LocalSize
Definition: fluid_element.h:111
static constexpr unsigned int BlockSize
Definition: fluid_element.h:109
static constexpr unsigned int StrainSize
Definition: fluid_element.h:113
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
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
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
Definition: qs_vms.h:52
void AddTimeIntegratedLHS(TElementData &rData, MatrixType &rLHS) override
Definition: qs_vms.cpp:379
void Calculate(const Variable< double > &rVariable, double &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: qs_vms.cpp:73
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: qs_vms.h:84
constexpr static unsigned int BlockSize
Definition: qs_vms.h:99
constexpr static unsigned int LocalSize
Definition: qs_vms.h:100
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: qs_vms.h:86
std::size_t IndexType
Definition: qs_vms.h:78
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: qs_vms.h:76
QSVMS(IndexType NewId=0)
Default constuctor.
Definition: qs_vms.cpp:30
virtual void MassProjTerm(const TElementData &rData, double &rMassRHS) const
Definition: qs_vms.cpp:355
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Definition: qs_vms.cpp:115
virtual void SubscalePressure(const TElementData &rData, double &rPressureSubscale) const
Definition: qs_vms.cpp:783
constexpr static unsigned int Dim
Definition: qs_vms.h:97
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: qs_vms.h:89
virtual void SubscaleVelocity(const TElementData &rData, array_1d< double, 3 > &rVelocitySubscale) const
Definition: qs_vms.cpp:763
virtual void CalculateProjections(const ProcessInfo &rCurrentProcessInfo)
Definition: qs_vms.cpp:708
std::string Info() const override
Turn back information as a string.
Definition: qs_vms.cpp:264
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: qs_vms.h:95
virtual void MomentumProjTerm(const TElementData &rData, const array_1d< double, 3 > &rConvectionVelocity, array_1d< double, 3 > &rMomentumRHS) const
Definition: qs_vms.cpp:337
virtual void CalculateTau(const TElementData &rData, const array_1d< double, 3 > &Velocity, double &TauOne, double &TauTwo) const
Definition: qs_vms.cpp:681
std::size_t SizeType
Definition: qs_vms.h:80
void AddTimeIntegratedRHS(TElementData &rData, VectorType &rRHS) override
Definition: qs_vms.cpp:386
void AddMassStabilization(TElementData &rData, MatrixType &rMassMatrix)
Definition: qs_vms.cpp:535
virtual void AddViscousTerm(const TElementData &rData, BoundedMatrix< double, LocalSize, LocalSize > &rLHS, VectorType &rRHS)
Definition: qs_vms.cpp:620
virtual void AlgebraicMassResidual(const TElementData &rData, double &rMomentumRes) const
Definition: qs_vms.cpp:306
void CalculateOnIntegrationPoints(Variable< array_1d< double, 3 >> const &rVariable, std::vector< array_1d< double, 3 >> &rValues, ProcessInfo const &rCurrentProcessInfo) override
const Parameters GetSpecifications() const override
This method provides the specifications/requirements of the element.
Definition: qs_vms.cpp:230
void AddTimeIntegratedSystem(TElementData &rData, MatrixType &rLHS, VectorType &rRHS) override
Definition: qs_vms.cpp:369
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: qs_vms.h:67
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(QSVMS)
Pointer definition of QSVMS.
virtual KRATOS_DEPRECATED double EffectiveViscosity(TElementData &rData, double ElementSize)
EffectiveViscosity Evaluate the total kinematic viscosity at a given integration point....
Definition: qs_vms.cpp:642
Node NodeType
Node type (default is: Node)
Definition: qs_vms.h:64
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: qs_vms.h:70
constexpr static unsigned int NumNodes
Definition: qs_vms.h:98
std::vector< std::size_t > EquationIdVectorType
Definition: qs_vms.h:82
void Calculate(const Variable< array_1d< double, 3 > > &rVariable, array_1d< double, 3 > &rOutput, const ProcessInfo &rCurrentProcessInfo) override
void AddBoundaryTraction(TElementData &rData, const Vector &rUnitNormal, MatrixType &rLHS, VectorType &rRHS) override
Adds the boundary traction component along a cut plane for embedded formulations. This method adds th...
Definition: qs_vms.cpp:577
constexpr static unsigned int StrainSize
Definition: qs_vms.h:101
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: qs_vms.cpp:273
virtual void OrthogonalMassResidual(const TElementData &rData, double &rMassRes) const
Definition: qs_vms.cpp:326
Vector VectorType
Vector type for local contributions to the linear system.
Definition: qs_vms.h:73
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, Properties::Pointer pProperties) const override
Create a new element of this type.
Definition: qs_vms.cpp:60
void AddMassLHS(TElementData &rData, MatrixType &rMassMatrix) override
Definition: qs_vms.cpp:502
void AddVelocitySystem(TElementData &rData, MatrixType &rLocalLHS, VectorType &rLocalRHS) override
Definition: qs_vms.cpp:393
virtual void OrthogonalMomentumResidual(const TElementData &rData, const array_1d< double, 3 > &rConvectionVelocity, array_1d< double, 3 > &rResidual) const
Definition: qs_vms.cpp:314
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: qs_vms.h:92
virtual void AlgebraicMomentumResidual(const TElementData &rData, const array_1d< double, 3 > &rConvectionVelocity, array_1d< double, 3 > &rResidual) const
Definition: qs_vms.cpp:282
~QSVMS() override
Destructor.
Definition: qs_vms.cpp:53
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
#define KRATOS_DEPRECATED
Definition: define.h:738
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
def load(f)
Definition: ode_solve.py:307