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.
fractional_step_discontinuous.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: Riccardo Rossi
11 //
12 
13 #if !defined(KRATOS_FRACTIONAL_STEP_DISCONTINUOUS_H_INCLUDED )
14 #define KRATOS_FRACTIONAL_STEP_DISCONTINUOUS_H_INCLUDED
15 
16 // System includes
17 #include <string>
18 #include <iostream>
19 
20 
21 // External includes
22 
23 
24 // Project includes
25 #include "containers/array_1d.h"
26 #include "includes/define.h"
27 #include "fractional_step.h"
28 #include "includes/serializer.h"
29 #include "geometries/geometry.h"
30 #include "utilities/math_utils.h"
32 
33 // Application includes
35 
36 namespace Kratos
37 {
38 
41 
44 
48 
52 
56 
60 
62 
64 template< unsigned int TDim >
66 {
67 public:
70 
73 
75  typedef Node NodeType;
76 
79 
82 
84  typedef Vector VectorType;
85 
87  typedef Matrix MatrixType;
88 
89  typedef std::size_t IndexType;
90 
91  typedef std::size_t SizeType;
92 
93  typedef std::vector<std::size_t> EquationIdVectorType;
94 
95  typedef std::vector< Dof<double>::Pointer > DofsVectorType;
96 
98 
101 
104 
107 
111 
112  //Constructors.
113 
115 
119  FractionalStep<TDim>(NewId),
120  medge_areas(array_1d<double,(TDim-1)*3 >( (TDim-1)*3, 0.0))
121  {}
122 
124 
129  FractionalStep<TDim>(NewId, ThisNodes),
130  medge_areas(array_1d<double,(TDim-1)*3 >( (TDim-1)*3, 0.0))
131  {}
132 
134 
138  FractionalStepDiscontinuous(IndexType NewId, GeometryType::Pointer pGeometry) :
139  FractionalStep<TDim>(NewId, pGeometry),
140  medge_areas(array_1d<double,(TDim-1)*3 >( (TDim-1)*3, 0.0))
141  {}
142 
144 
149  FractionalStepDiscontinuous(IndexType NewId, GeometryType::Pointer pGeometry, Element::PropertiesType::Pointer pProperties) :
150  FractionalStep<TDim>(NewId, pGeometry, pProperties),
151  medge_areas(array_1d<double,(TDim-1)*3 >( (TDim-1)*3, 0.0))
152  {}
153 
156  {}
157 
158 
166  void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix,
167  VectorType& rRightHandSideVector,
168  const ProcessInfo& rCurrentProcessInfo) override;
169 
173 
175 
182  Element::Pointer Create(IndexType NewId, NodesArrayType const& ThisNodes,
183  Element::PropertiesType::Pointer pProperties) const override
184  {
185  return Kratos::make_intrusive< FractionalStepDiscontinuous<TDim> >(NewId, this->GetGeometry().Create(ThisNodes), pProperties);
186  }
187 
188  Element::Pointer Create(IndexType NewId, Element::GeometryType::Pointer pGeom, Element::PropertiesType::Pointer pProperties) const override
189  {
190  return Kratos::make_intrusive< FractionalStepDiscontinuous<TDim> >(NewId, pGeom, pProperties);
191  }
192 
196 
200 
202 
210  /* virtual int Check(const ProcessInfo& rCurrentProcessInfo);*/
211 
212 
218  void Calculate(const Variable<double>& rVariable,
219  double& rOutput,
220  const ProcessInfo& rCurrentProcessInfo) override;
226  void Calculate(const Variable<array_1d<double, 3 > >& rVariable,
227  array_1d<double, 3 > & rOutput,
228  const ProcessInfo& rCurrentProcessInfo) override;
229 
233 
234 
238 
240  std::string Info() const override
241  {
242  std::stringstream buffer;
243  buffer << "FractionalStepDiscontinuous #" << this->Id();
244  return buffer.str();
245  }
246 
248  void PrintInfo(std::ostream& rOStream) const override
249  {
250  rOStream << "FractionalStepDiscontinuous" << TDim << "D";
251  }
252 
253 // /// Print object's data.
254 // virtual void PrintData(std::ostream& rOStream) const;
255 
259 
260 
262 
263 protected:
266 
267 
272 
273 
277  void CalculateLocalPressureSystem(MatrixType& rLeftHandSideMatrix,
278  VectorType& rRightHandSideVector,
279  const ProcessInfo& rCurrentProcessInfo) override;
280 
281  void AddMomentumSystemTerms(Matrix& rLHSMatrix,
282  Vector& rRHSVector,
283  const double Density,
284  const Vector& rConvOperator,
285  const array_1d<double,3>& rBodyForce,
286  const double OldPressure,
287  const double TauOne,
288  const double TauTwo,
289  const array_1d<double,3>& rMomentumProjection,
290  const double MassProjection,
291  const ShapeFunctionsType& rN,
292  const ShapeFunctionDerivativesType& rDN_DX,
293  const double Weight) override;
297 
300  Matrix& rNContainer,
301  Vector& rGaussWeights) override;
302 
303 
304 
305 
306 
310 
311 
315 
316 
320 
321 
323 
324 private:
327 
331 
335 
336  friend class Serializer;
337 
338  void save(Serializer& rSerializer) const override
339  {
340  typedef FractionalStep<TDim> basetype;
341  KRATOS_SERIALIZE_SAVE_BASE_CLASS(rSerializer, basetype );
342  }
343 
344  void load(Serializer& rSerializer) override
345  {
346  typedef FractionalStep<TDim> basetype;
347  KRATOS_SERIALIZE_LOAD_BASE_CLASS(rSerializer, basetype );
348  }
349 
353 
354 
358 
359 
363 
364 
368 
369 
373 
376 
379 
381 
382 }; // Class FractionalStepDiscontinuous
383 
385 
388 
389 
393 
394 
396 template< unsigned int TDim >
397 inline std::istream& operator >>(std::istream& rIStream,
399 {
400  return rIStream;
401 }
402 
404 template< unsigned int TDim >
405 inline std::ostream& operator <<(std::ostream& rOStream,
407 {
408  rThis.PrintInfo(rOStream);
409  rOStream << std::endl;
410  rThis.PrintData(rOStream);
411 
412  return rOStream;
413 }
415 
417 
418 } // namespace Kratos.
419 
420 #endif // KRATOS_FRACTIONAL_STEP_DISCONTINUOUS_H_INCLUDED defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: element.h:1135
std::size_t IndexType
Definition: flags.h:74
A stabilized element for the incompressible Navier-Stokes equations.
Definition: fractional_step_discontinuous.h:66
void CalculateGeometryData(ShapeFunctionDerivativesArrayType &rDN_DX, Matrix &rNContainer, Vector &rGaussWeights) override
Determine integration point weights and shape funcition derivatives from the element's geometry.
Definition: fractional_step_discontinuous.cpp:511
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: fractional_step_discontinuous.h:87
array_1d< double,(TDim-1) *3 > medge_areas
Definition: fractional_step_discontinuous.h:271
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: fractional_step_discontinuous.h:97
FractionalStepDiscontinuous(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using a geometry object.
Definition: fractional_step_discontinuous.h:138
FractionalStepDiscontinuous(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: fractional_step_discontinuous.h:128
Element::Pointer Create(IndexType NewId, Element::GeometryType::Pointer pGeom, Element::PropertiesType::Pointer pProperties) const override
Definition: fractional_step_discontinuous.h:188
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Calculate the element's local contribution to the system for the current step.
Definition: fractional_step_discontinuous.cpp:628
void Calculate(const Variable< double > &rVariable, double &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Checks the input and that all required Kratos variables have been registered.
Definition: fractional_step_discontinuous.cpp:458
Node NodeType
Node type (default is: Node)
Definition: fractional_step_discontinuous.h:75
std::string Info() const override
Turn back information as a string.
Definition: fractional_step_discontinuous.h:240
void AddMomentumSystemTerms(Matrix &rLHSMatrix, Vector &rRHSVector, const double Density, const Vector &rConvOperator, const array_1d< double, 3 > &rBodyForce, const double OldPressure, const double TauOne, const double TauTwo, const array_1d< double, 3 > &rMomentumProjection, const double MassProjection, const ShapeFunctionsType &rN, const ShapeFunctionDerivativesType &rDN_DX, const double Weight) override
Definition: fractional_step_discontinuous.cpp:58
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: fractional_step_discontinuous.h:106
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: fractional_step_discontinuous.h:248
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, Element::PropertiesType::Pointer pProperties) const override
Create a new element of this type.
Definition: fractional_step_discontinuous.h:182
FractionalStepDiscontinuous(IndexType NewId, GeometryType::Pointer pGeometry, Element::PropertiesType::Pointer pProperties)
Constuctor using geometry and properties.
Definition: fractional_step_discontinuous.h:149
std::vector< std::size_t > EquationIdVectorType
Definition: fractional_step_discontinuous.h:93
~FractionalStepDiscontinuous() override
Destructor.
Definition: fractional_step_discontinuous.h:155
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: fractional_step_discontinuous.h:81
std::size_t SizeType
Definition: fractional_step_discontinuous.h:91
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(FractionalStepDiscontinuous)
Pointer definition of FractionalStepDiscontinuous.
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: fractional_step_discontinuous.h:103
void CalculateLocalPressureSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: fractional_step_discontinuous.cpp:130
Vector VectorType
Vector type for local contributions to the linear system.
Definition: fractional_step_discontinuous.h:84
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: fractional_step_discontinuous.h:95
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: fractional_step_discontinuous.h:78
FractionalStepDiscontinuous(IndexType NewId=0)
Default constuctor.
Definition: fractional_step_discontinuous.h:118
std::size_t IndexType
Definition: fractional_step_discontinuous.h:89
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: fractional_step_discontinuous.h:100
A stabilized element for the incompressible Navier-Stokes equations.
Definition: fractional_step.h:65
GeometryType & GetGeometry()
Returns the reference of the geometry.
Definition: geometrical_object.h:158
Geometry base class.
Definition: geometry.h:71
virtual Pointer Create(PointsArrayType const &rThisPoints) const
Creates a new geometry pointer.
Definition: geometry.h:813
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
Short class definition.
Definition: array_1d.h:61
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
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
TABLE_NUMBER_ANGULAR_VELOCITY TABLE_NUMBER_MOMENT I33 BEAM_INERTIA_ROT_UNIT_LENGHT_Y KRATOS_DEFINE_APPLICATION_VARIABLE(DEM_APPLICATION, double, BEAM_INERTIA_ROT_UNIT_LENGHT_Z) typedef std double
Definition: DEM_application_variables.h:182
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