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.
d_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_D_VMS_H
14 #define KRATOS_D_VMS_H
15 
16 //#define KRATOS_D_VMS_SUBSCALE_ERROR_INSTRUMENTATION
17 
18 #include "includes/define.h"
19 #include "includes/element.h"
20 #include "includes/serializer.h"
21 #include "geometries/geometry.h"
22 
23 #include "includes/cfd_variables.h"
24 #include "custom_elements/qs_vms.h"
26 
27 namespace Kratos
28 {
29 
32 
35 
39 
43 
47 
51 
52 template< class TElementData >
53 class DVMS : public QSVMS<TElementData>
54 {
55 public:
58 
61 
64 
66  typedef Node NodeType;
67 
70 
73 
75  typedef Vector VectorType;
76 
78  typedef Matrix MatrixType;
79 
80  typedef std::size_t IndexType;
81 
82  typedef std::size_t SizeType;
83 
84  typedef std::vector<std::size_t> EquationIdVectorType;
85 
86  typedef std::vector< Dof<double>::Pointer > DofsVectorType;
87 
89 
92 
95 
98 
99  constexpr static unsigned int Dim = BaseType::Dim;
100  constexpr static unsigned int NumNodes = BaseType::NumNodes;
101  constexpr static unsigned int BlockSize = BaseType::BlockSize;
102  constexpr static unsigned int LocalSize = BaseType::LocalSize;
103  constexpr static unsigned int StrainSize = BaseType::StrainSize;
104 
108 
109  //Constructors.
110 
112 
115  DVMS(IndexType NewId = 0);
116 
118 
122  DVMS(IndexType NewId, const NodesArrayType& ThisNodes);
123 
125 
129  DVMS(IndexType NewId, GeometryType::Pointer pGeometry);
130 
132 
137  DVMS(IndexType NewId, GeometryType::Pointer pGeometry, Properties::Pointer pProperties);
138 
140  ~DVMS() override;
141 
145 
146 
150 
151 
153 
160  Element::Pointer Create(IndexType NewId,
161  NodesArrayType const& ThisNodes,
162  Properties::Pointer pProperties) const override;
163 
165 
172  Element::Pointer Create(IndexType NewId,
173  GeometryType::Pointer pGeom,
174  Properties::Pointer pProperties) const override;
175 
176 
177  void Calculate(
178  const Variable<double>& rVariable,
179  double& rOutput,
180  const ProcessInfo& rCurrentProcessInfo) override;
181 
182 
183  void Calculate(
184  const Variable<array_1d<double, 3 > >& rVariable,
185  array_1d<double, 3 > & rOutput,
186  const ProcessInfo& rCurrentProcessInfo) override;
187 
188 
189  void Calculate(
190  const Variable<Vector >& rVariable,
191  Vector& Output,
192  const ProcessInfo& rCurrentProcessInfo) override;
193 
194  void Calculate(
195  const Variable<Matrix >& rVariable,
196  Matrix& Output,
197  const ProcessInfo& rCurrentProcessInfo) override;
198 
200 
201  virtual void Initialize(const ProcessInfo &rCurrentProcessInfo) override;
202 
204  virtual void FinalizeSolutionStep(const ProcessInfo& rCurrentProcessInfo) override;
205 
207  virtual void InitializeNonLinearIteration(const ProcessInfo& rCurrentProcessInfo) override;
208 
212 
214  Variable<array_1d<double, 3>> const& rVariable,
215  std::vector<array_1d<double, 3>>& rValues,
216  ProcessInfo const& rCurrentProcessInfo) override;
217 
219  Variable<double> const& rVariable,
220  std::vector<double>& rValues,
221  ProcessInfo const& rCurrentProcessInfo) override;
222 
224  Variable<array_1d<double, 6>> const& rVariable,
225  std::vector<array_1d<double, 6>>& rValues,
226  ProcessInfo const& rCurrentProcessInfo) override;
227 
229  Variable<Vector> const& rVariable,
230  std::vector<Vector>& rValues,
231  ProcessInfo const& rCurrentProcessInfo) override;
232 
234  Variable<Matrix> const& rVariable,
235  std::vector<Matrix>& rValues,
236  ProcessInfo const& rCurrentProcessInfo) override;
237 
241 
242  int Check(const ProcessInfo &rCurrentProcessInfo) const override;
243 
247 
248  const Parameters GetSpecifications() const override;
249 
251  std::string Info() const override;
252 
253 
255  void PrintInfo(std::ostream& rOStream) const override;
256 
257 
261 
262 
264 
265 protected:
266 
269 
270  constexpr static double mTauC1 = 8.0;
271  constexpr static double mTauC2 = 2.0;
272  constexpr static double mSubscalePredictionVelocityTolerance = 1e-14;
273  constexpr static double mSubscalePredictionResidualTolerance = 1e-14;
274  constexpr static unsigned int mSubscalePredictionMaxIterations = 10;
275 
279 
280  // Velocity subscale history, stored at integration points
283 
284  #ifdef KRATOS_D_VMS_SUBSCALE_ERROR_INSTRUMENTATION
285  std::vector< double > mSubscaleIterationError;
286  std::vector< unsigned int > mSubscaleIterationCount;
287  #endif
288 
289 
293 
294 
298 
299  // Protected interface of FluidElement ////////////////////////////////////
300 
301  void AddVelocitySystem(
302  TElementData& rData,
303  MatrixType& rLocalLHS,
304  VectorType& rLocalRHS) override;
305 
306  virtual void AddMassLHS(
307  TElementData& rData,
308  MatrixType& rMassMatrix) override;
309 
310  // Implementation details of DVMS /////////////////////////////////////////
311 
312  virtual void AddMassStabilization(
313  TElementData& rData,
314  MatrixType& rMassMatrix);
315 
316  void CalculateProjections(const ProcessInfo &rCurrentProcessInfo) override;
317 
319  const TElementData& rData,
320  const array_1d<double,3> &Velocity,
321  double &TauOne,
322  double &TauTwo,
323  double &TauP) const;
324 
325  virtual void SubscaleVelocity(
326  const TElementData& rData,
327  array_1d<double,3>& rVelocitySubscale) const override;
328 
329  virtual void SubscalePressure(
330  const TElementData& rData,
331  double &rPressureSubscale) const override;
332 
334  const TElementData& rData) const;
335 
337  const TElementData& rData);
341 
342 
346 
347 
351 
353 
354 private:
355 
358 
362 
366 
367  friend class Serializer;
368 
369  void save(Serializer& rSerializer) const override;
370 
371  void load(Serializer& rSerializer) override;
372 
376 
377 
381 
385 
386 
390 
391 
395 
397  DVMS& operator=(DVMS const& rOther);
398 
400  DVMS(DVMS const& rOther);
401 
403 
404 
405 }; // Class DVMS
406 
408 
411 
412 
416 
417 
419 template< class TElementData >
420 inline std::istream& operator >>(std::istream& rIStream,
421  DVMS<TElementData>& rThis)
422 {
423  return rIStream;
424 }
425 
427 template< class TElementData >
428 inline std::ostream& operator <<(std::ostream& rOStream,
429  const DVMS<TElementData>& rThis)
430 {
431  rThis.PrintInfo(rOStream);
432  rOStream << std::endl;
433  rThis.PrintData(rOStream);
434 
435  return rOStream;
436 }
438 
440 
441 } // namespace Kratos.
442 
443 #endif // KRATOS_D_VMS_H
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Definition: d_vms.h:54
constexpr static unsigned int StrainSize
Definition: d_vms.h:103
Node NodeType
Node type (default is: Node)
Definition: d_vms.h:66
virtual void UpdateSubscaleVelocityPrediction(const TElementData &rData)
Definition: d_vms.cpp:792
void CalculateOnIntegrationPoints(Variable< array_1d< double, 3 >> const &rVariable, std::vector< array_1d< double, 3 >> &rValues, ProcessInfo const &rCurrentProcessInfo) override
virtual void AddMassStabilization(TElementData &rData, MatrixType &rMassMatrix)
Definition: d_vms.cpp:578
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Definition: d_vms.cpp:200
std::size_t SizeType
Definition: d_vms.h:82
~DVMS() override
Destructor.
Definition: d_vms.cpp:60
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: d_vms.h:91
std::vector< std::size_t > EquationIdVectorType
Definition: d_vms.h:84
constexpr static unsigned int Dim
Definition: d_vms.h:99
void CalculateStabilizationParameters(const TElementData &rData, const array_1d< double, 3 > &Velocity, double &TauOne, double &TauTwo, double &TauP) const
Definition: d_vms.cpp:688
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: d_vms.h:69
const Parameters GetSpecifications() const override
This method provides the specifications/requirements of the element.
Definition: d_vms.cpp:352
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: d_vms.h:94
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: d_vms.h:86
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: d_vms.h:97
constexpr static double mSubscalePredictionResidualTolerance
Definition: d_vms.h:273
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: d_vms.h:78
void Calculate(const Variable< array_1d< double, 3 > > &rVariable, array_1d< double, 3 > &rOutput, const ProcessInfo &rCurrentProcessInfo) override
void CalculateProjections(const ProcessInfo &rCurrentProcessInfo) override
Definition: d_vms.cpp:624
DenseVector< array_1d< double, Dim > > mPredictedSubscaleVelocity
Definition: d_vms.h:281
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: d_vms.h:72
Vector VectorType
Vector type for local contributions to the linear system.
Definition: d_vms.h:75
virtual void SubscalePressure(const TElementData &rData, double &rPressureSubscale) const override
Definition: d_vms.cpp:744
std::size_t IndexType
Definition: d_vms.h:80
DVMS(IndexType NewId=0)
Default constuctor.
Definition: d_vms.cpp:29
constexpr static unsigned int mSubscalePredictionMaxIterations
Definition: d_vms.h:274
constexpr static double mTauC2
Definition: d_vms.h:271
virtual void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Set up the element.
Definition: d_vms.cpp:121
virtual void InitializeNonLinearIteration(const ProcessInfo &rCurrentProcessInfo) override
Predict the value of the small scale velocity for the current iteration.
Definition: d_vms.cpp:177
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, Properties::Pointer pProperties) const override
Create a new element of this type.
Definition: d_vms.cpp:67
std::string Info() const override
Turn back information as a string.
Definition: d_vms.cpp:386
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: d_vms.cpp:395
void Calculate(const Variable< double > &rVariable, double &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: d_vms.cpp:80
constexpr static unsigned int LocalSize
Definition: d_vms.h:102
DenseVector< array_1d< double, Dim > > mOldSubscaleVelocity
Definition: d_vms.h:282
constexpr static unsigned int BlockSize
Definition: d_vms.h:101
virtual array_1d< double, 3 > FullConvectiveVelocity(const TElementData &rData) const
Definition: d_vms.cpp:778
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: d_vms.h:88
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(DVMS)
Pointer definition of DVMS.
constexpr static double mTauC1
Definition: d_vms.h:270
constexpr static unsigned int NumNodes
Definition: d_vms.h:100
virtual void SubscaleVelocity(const TElementData &rData, array_1d< double, 3 > &rVelocitySubscale) const override
Definition: d_vms.cpp:713
virtual void AddMassLHS(TElementData &rData, MatrixType &rMassMatrix) override
Definition: d_vms.cpp:545
constexpr static double mSubscalePredictionVelocityTolerance
Definition: d_vms.h:272
void AddVelocitySystem(TElementData &rData, MatrixType &rLocalLHS, VectorType &rLocalRHS) override
Definition: d_vms.cpp:407
virtual void FinalizeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override
Update the values of tracked small scale quantities.
Definition: d_vms.cpp:150
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
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
constexpr static unsigned int BlockSize
Definition: qs_vms.h:99
constexpr static unsigned int LocalSize
Definition: qs_vms.h:100
constexpr static unsigned int Dim
Definition: qs_vms.h:97
constexpr static unsigned int NumNodes
Definition: qs_vms.h:98
constexpr static unsigned int StrainSize
Definition: qs_vms.h:101
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
e
Definition: run_cpp_mpi_tests.py:31