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_dem_coupled.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: Joaquin Gonzalez-Usua
11 //
12 
13 #ifndef KRATOS_QS_VMS_DEM_COUPLED_H
14 #define KRATOS_QS_VMS_DEM_COUPLED_H
15 
16 // System includes
17 #include <string>
18 #include <iostream>
19 
20 // External includes
21 
22 // Project includes
23 #include "includes/checks.h"
24 #include "includes/element.h"
25 #include "includes/serializer.h"
26 #include "geometries/geometry.h"
27 #include "includes/cfd_variables.h"
28 
29 // Application includes
31 #include "custom_elements/qs_vms.h"
33 
34 namespace Kratos
35 {
36 
39 
42 
46 
50 
54 
58 
59 template< class TElementData >
60 class QSVMSDEMCoupled : public QSVMS<TElementData>
61 {
62 public:
65 
68 
70  typedef Node NodeType;
71 
74 
77 
79  typedef Vector VectorType;
80 
82  typedef Matrix MatrixType;
83 
84  typedef std::size_t IndexType;
85 
86  typedef std::size_t SizeType;
87 
88  typedef std::vector<std::size_t> EquationIdVectorType;
89 
90  typedef std::vector< Dof<double>::Pointer > DofsVectorType;
91 
93 
96 
99 
102 
104 
105  constexpr static unsigned int Dim = QSVMS<TElementData>::Dim;
106  constexpr static unsigned int NumNodes = QSVMS<TElementData>::NumNodes;
107  constexpr static unsigned int BlockSize = QSVMS<TElementData>::BlockSize;
108  constexpr static unsigned int LocalSize = QSVMS<TElementData>::LocalSize;
109  constexpr static unsigned int StrainSize = QSVMS<TElementData>::StrainSize;
110 
114 
115  //Constructors.
116 
118 
121  QSVMSDEMCoupled(IndexType NewId = 0);
122 
124 
128  QSVMSDEMCoupled(IndexType NewId, const NodesArrayType& ThisNodes);
129 
131 
135  QSVMSDEMCoupled(IndexType NewId, GeometryType::Pointer pGeometry);
136 
138 
143  QSVMSDEMCoupled(IndexType NewId, GeometryType::Pointer pGeometry, Properties::Pointer pProperties);
144 
146  ~QSVMSDEMCoupled() override;
147 
151 
152 
156 
157 
159 
166  Element::Pointer Create(
167  IndexType NewId,
168  NodesArrayType const& ThisNodes,
169  Properties::Pointer pProperties) const override;
170 
172 
179  Element::Pointer Create(
180  IndexType NewId,
181  GeometryType::Pointer pGeom,
182  Properties::Pointer pProperties) const override;
183 
184  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
185 
187 
188  void Calculate(
189  const Variable<double>& rVariable,
190  double& rOutput,
191  const ProcessInfo& rCurrentProcessInfo) override;
192 
193  void Calculate(
194  const Variable<array_1d<double, 3>>& rVariable,
195  array_1d<double, 3>& rOutput,
196  const ProcessInfo& rCurrentProcessInfo) override;
197 
199  const Variable<array_1d<double, 3>>& rVariable,
200  std::vector<array_1d<double, 3>>& rOutput,
201  const ProcessInfo& rCurrentProcessInfo) override;
202 
204  const Variable<double>& rVariable,
205  std::vector<double>& rOutput,
206  const ProcessInfo& rCurrentProcessInfo) override;
207 
209  Variable<Matrix> const& rVariable,
210  std::vector<Matrix>& rValues,
211  ProcessInfo const& rCurrentProcessInfo) override;
212 
213  void InitializeNonLinearIteration(const ProcessInfo& rCurrentProcessInfo) override;
214 
215  void FinalizeNonLinearIteration(const ProcessInfo& rCurrentProcessInfo) override;
216 
217  void UpdateSubscaleVelocity(const TElementData& rData);
218 
219  void EquationIdVector(
220  EquationIdVectorType& rResult,
221  const ProcessInfo& rCurrentProcessInfo) const override;
222 
224  VectorType& rRightHandSideVector,
225  const ProcessInfo& rCurrentProcessInfo) override;
226 
230 
234 
235  int Check(const ProcessInfo &rCurrentProcessInfo) const override;
236 
240 
242  std::string Info() const override;
243 
244 
246  void PrintInfo(std::ostream& rOStream) const override;
247 
248 
252 
253 
255 
256 protected:
257 
260 
261 
267  // Velocity subscale history, stored at integration points
270 
274 
275 
279 
280  // Protected interface of FluidElement ////////////////////////////////////
281 
283  const TElementData& rData,
284  const array_1d<double,3> &rConvectionVelocity,
285  array_1d<double,3>& rResidual) const override;
286 
287  void MomentumProjTerm(
288  const TElementData& rData,
289  const array_1d<double,3>& rConvectionVelocity,
290  array_1d<double,3> &rMomentumRHS) const override;
291 
293  TElementData& rData,
294  BoundedMatrix<double,NumNodes*(Dim+1),NumNodes*(Dim+1)>& rLHS,
295  VectorType& rLocalRHS);
296 
298  TElementData& rData,
299  MatrixType &rMassMatrix);
300 
302  void CalculateTau(
303  const TElementData& rData,
304  const array_1d<double,3> &Velocity,
306  double &TauTwo) const;
307 
309  const ProcessInfo &rCurrentProcessInfo) override;
310 
312  TElementData& rData,
313  unsigned int IntegrationPointIndex,
314  double Weight,
315  const typename TElementData::MatrixRowType& rN,
316  const typename TElementData::ShapeDerivativesType& rDN_DX,
317  const typename TElementData::ShapeFunctionsSecondDerivativesType& rDDN_DDX) const;
318 
319  void AddVelocitySystem(
320  TElementData& rData,
321  MatrixType &rLocalLHS,
322  VectorType &rLocalRHS) override;
323 
324  void CalculateMassMatrix(MatrixType& rMassMatrix,
325  const ProcessInfo& rCurrentProcessInfo) override;
326 
328  VectorType& rRightHandSideVector,
329  const ProcessInfo& rCurrentProcessInfo) override;
330 
332  const TElementData& rData);
333 
334  void AddMassLHS(
335  TElementData& rData,
336  MatrixType& rMassMatrix) override;
337 
338  void MassProjTerm(
339  const TElementData& rData,
340  double &rMassRHS) const override;
341 
342  void SubscaleVelocity(
343  const TElementData& rData,
344  array_1d<double,3> &rVelocitySubscale) const override;
345 
346  void SubscalePressure(
347  const TElementData& rData,
348  double &rPressureSubscale) const override;
349 
353 
354 
358 
359 
363 
365 
366 private:
367 
370 
374 
378 
379  friend class Serializer;
380 
381  void save(Serializer& rSerializer) const override;
382 
383  void load(Serializer& rSerializer) override;
384 
388 
389 
393 
394 
398 
399 
403 
404 
408 
410  QSVMSDEMCoupled& operator=(QSVMSDEMCoupled const& rOther);
411 
413  QSVMSDEMCoupled(QSVMSDEMCoupled const& rOther);
415 
416 
417 }; // Class QSVMSDEMCoupled
418 
420 
423 
424 
428 
429 
431 template< class TElementData >
432 inline std::istream& operator >>(std::istream& rIStream,
434 {
435  return rIStream;
436 }
437 
439 template< class TElementData >
440 inline std::ostream& operator <<(std::ostream& rOStream,
441  const QSVMSDEMCoupled<TElementData>& rThis)
442 {
443  rThis.PrintInfo(rOStream);
444  rOStream << std::endl;
445  rThis.PrintData(rOStream);
446 
447  return rOStream;
448 }
450 
452 
453 } // namespace Kratos.
454 
455 #endif // KRATOS_QS_VMS_DEM_COUPLED_H
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
std::vector< std::size_t > EquationIdVectorType
Definition: element.h:98
std::size_t IndexType
Definition: flags.h:74
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: fluid_element.hpp:584
IntegrationMethod
Definition: geometry_data.h:76
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
Definition: qs_vms_dem_coupled.h:61
constexpr static unsigned int NumNodes
Definition: qs_vms_dem_coupled.h:106
constexpr static unsigned int Dim
Definition: qs_vms_dem_coupled.h:105
constexpr static unsigned int BlockSize
Definition: qs_vms_dem_coupled.h:107
DenseVector< array_1d< double, Dim > > mPreviousVelocity
Definition: qs_vms_dem_coupled.h:269
void CalculateOnIntegrationPoints(const Variable< array_1d< double, 3 >> &rVariable, std::vector< array_1d< double, 3 >> &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: qs_vms_dem_coupled.cpp:136
void CalculateLocalVelocityContribution(MatrixType &rDampMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
CalculateLocalVelocityContribution Calculate the local contribution in terms of velocity and pressure...
Definition: qs_vms_dem_coupled.cpp:807
void AddReactionStabilization(TElementData &rData, BoundedMatrix< double, NumNodes *(Dim+1), NumNodes *(Dim+1)> &rLHS, VectorType &rLocalRHS)
Definition: qs_vms_dem_coupled.cpp:452
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &rCurrentProcessInfo) const override
Definition: qs_vms_dem_coupled.cpp:280
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(QSVMSDEMCoupled)
Pointer definition of QSVMSDEMCoupled.
constexpr static unsigned int StrainSize
Definition: qs_vms_dem_coupled.h:109
void CalculateTau(const TElementData &rData, const array_1d< double, 3 > &Velocity, BoundedMatrix< double, Dim, Dim > &TauOne, double &TauTwo) const
Definition: qs_vms_dem_coupled.cpp:912
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: qs_vms_dem_coupled.h:101
void FinalizeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override
Definition: qs_vms_dem_coupled.cpp:326
GeometryData::IntegrationMethod GetIntegrationMethod() const override
GetIntegrationMethod Return the integration order to be used.
Definition: qs_vms_dem_coupled.cpp:227
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: qs_vms_dem_coupled.h:82
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: qs_vms_dem_coupled.h:73
void AddVelocitySystem(TElementData &rData, MatrixType &rLocalLHS, VectorType &rLocalRHS) override
Definition: qs_vms_dem_coupled.cpp:595
~QSVMSDEMCoupled() override
Destructor.
Definition: qs_vms_dem_coupled.cpp:54
void Calculate(const Variable< double > &rVariable, double &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: qs_vms_dem_coupled.cpp:258
void CalculateResistanceTensor(const TElementData &rData)
Definition: qs_vms_dem_coupled.cpp:849
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: qs_vms_dem_coupled.h:98
void AddMassStabilization(TElementData &rData, MatrixType &rMassMatrix)
Definition: qs_vms_dem_coupled.cpp:536
std::size_t SizeType
Definition: qs_vms_dem_coupled.h:86
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, Properties::Pointer pProperties) const override
Create a new element of this type.
Definition: qs_vms_dem_coupled.cpp:58
void SubscalePressure(const TElementData &rData, double &rPressureSubscale) const override
Definition: qs_vms_dem_coupled.cpp:1036
std::string Info() const override
Turn back information as a string.
Definition: qs_vms_dem_coupled.cpp:288
void UpdateIntegrationPointDataSecondDerivatives(TElementData &rData, unsigned int IntegrationPointIndex, double Weight, const typename TElementData::MatrixRowType &rN, const typename TElementData::ShapeDerivativesType &rDN_DX, const typename TElementData::ShapeFunctionsSecondDerivativesType &rDDN_DDX) const
Definition: qs_vms_dem_coupled.cpp:348
DenseVector< BoundedMatrix< double, Dim, Dim > > mViscousResistanceTensor
Definition: qs_vms_dem_coupled.h:266
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Definition: qs_vms_dem_coupled.cpp:239
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: qs_vms_dem_coupled.h:90
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: qs_vms_dem_coupled.cpp:297
DenseVector< array_1d< double, Dim > > mPredictedSubscaleVelocity
Definition: qs_vms_dem_coupled.h:268
int mInterpolationOrder
Definition: qs_vms_dem_coupled.h:265
void AlgebraicMomentumResidual(const TElementData &rData, const array_1d< double, 3 > &rConvectionVelocity, array_1d< double, 3 > &rResidual) const override
Definition: qs_vms_dem_coupled.cpp:373
void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Set up the element for solution.
Definition: qs_vms_dem_coupled.cpp:180
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: qs_vms_dem_coupled.h:95
constexpr static unsigned int LocalSize
Definition: qs_vms_dem_coupled.h:108
void MomentumProjTerm(const TElementData &rData, const array_1d< double, 3 > &rConvectionVelocity, array_1d< double, 3 > &rMomentumRHS) const override
Definition: qs_vms_dem_coupled.cpp:412
QSVMSDEMCoupled(IndexType NewId=0)
Default constuctor.
Definition: qs_vms_dem_coupled.cpp:30
std::size_t IndexType
Definition: qs_vms_dem_coupled.h:84
void AddMassLHS(TElementData &rData, MatrixType &rMassMatrix) override
Definition: qs_vms_dem_coupled.cpp:860
void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override
Definition: qs_vms_dem_coupled.cpp:303
Vector VectorType
Vector type for local contributions to the linear system.
Definition: qs_vms_dem_coupled.h:79
void UpdateSubscaleVelocity(const TElementData &rData)
Definition: qs_vms_dem_coupled.cpp:1058
void MassProjTerm(const TElementData &rData, double &rMassRHS) const override
Definition: qs_vms_dem_coupled.cpp:890
void CalculateMassMatrix(MatrixType &rMassMatrix, const ProcessInfo &rCurrentProcessInfo) override
MassMatrix Calculate the local mass matrix.
Definition: qs_vms_dem_coupled.cpp:771
void CalculateProjections(const ProcessInfo &rCurrentProcessInfo) override
Definition: qs_vms_dem_coupled.cpp:956
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
CalculateRightHandSide Return an empty matrix of appropriate size. This element does not have a local...
Definition: qs_vms_dem_coupled.cpp:361
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: qs_vms_dem_coupled.h:92
Node NodeType
Node type (default is: Node)
Definition: qs_vms_dem_coupled.h:70
GeometryType::ShapeFunctionsSecondDerivativesType ShapeFunctionsSecondDerivativesType
Definition: qs_vms_dem_coupled.h:103
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: qs_vms_dem_coupled.h:76
std::vector< std::size_t > EquationIdVectorType
Definition: qs_vms_dem_coupled.h:88
void SubscaleVelocity(const TElementData &rData, array_1d< double, 3 > &rVelocitySubscale) const override
Definition: qs_vms_dem_coupled.cpp:1015
Definition: qs_vms.h:52
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
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