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.
symbolic_stokes.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: Ruben Zorrilla
11 //
12 
13 #ifndef KRATOS_SYMBOLIC_STOKES_H
14 #define KRATOS_SYMBOLIC_STOKES_H
15 
16 // External includes
17 
18 // System includes
19 
20 // Project includes
21 #include "includes/define.h"
22 #include "includes/element.h"
23 #include "includes/serializer.h"
24 #include "geometries/geometry.h"
25 
26 // Application includes
27 #include "includes/cfd_variables.h"
29 
30 namespace Kratos
31 {
32 
35 
38 
42 
46 
50 
54 
55 template< class TElementData >
56 class SymbolicStokes : public FluidElement<TElementData>
57 {
58 public:
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 = FluidElement<TElementData>::Dim;
100  constexpr static unsigned int NumNodes = FluidElement<TElementData>::NumNodes;
101  constexpr static unsigned int BlockSize = FluidElement<TElementData>::BlockSize;
102  constexpr static unsigned int LocalSize = FluidElement<TElementData>::LocalSize;
103 
104  constexpr static unsigned int StrainSize = (Dim-1)*3;
105 
109 
110  //Constructors.
111 
113 
116  SymbolicStokes(IndexType NewId = 0);
117 
119 
124  IndexType NewId,
125  const NodesArrayType& ThisNodes);
126 
128 
133  IndexType NewId,
134  GeometryType::Pointer pGeometry);
135 
137 
143  IndexType NewId,
144  GeometryType::Pointer pGeometry,
145  Properties::Pointer pProperties);
146 
148  virtual ~SymbolicStokes();
149 
153 
154 
158 
159 
161 
168  Element::Pointer Create(
169  IndexType NewId,
170  NodesArrayType const& ThisNodes,
171  Properties::Pointer pProperties) const override;
172 
174 
181  Element::Pointer Create(
182  IndexType NewId,
183  GeometryType::Pointer pGeom,
184  Properties::Pointer pProperties) const override;
185 
189 
190  int Check(const ProcessInfo &rCurrentProcessInfo) const override;
191 
195 
196  const Parameters GetSpecifications() const override;
197 
199  std::string Info() const override;
200 
201 
203  void PrintInfo(std::ostream& rOStream) const override;
204 
206 protected:
209 
210 
214 
215 
219 
221  TElementData& rData,
222  MatrixType& rLHS,
223  VectorType& rRHS) override;
224 
226  TElementData& rData,
227  MatrixType& rLHS) override;
228 
230  TElementData& rData,
231  VectorType& rRHS) override;
232 
233  void AddBoundaryTraction(
234  TElementData& rData,
235  const Vector& rUnitNormal,
236  MatrixType& rLHS,
237  VectorType& rRHS) override;
238 
239  void Calculate(
240  const Variable<double>& rVariable,
241  double& rOutput,
242  const ProcessInfo& rCurrentProcessInfo) override;
243 
245  TElementData& rData,
246  MatrixType& rLHS);
247 
249  TElementData& rData,
250  VectorType& rRHS);
251 
255 
256 
260 
261 
265 
266 
268 
269 private:
270 
273 
277 
281 
282  friend class Serializer;
283 
284  void save(Serializer& rSerializer) const override;
285 
286  void load(Serializer& rSerializer) override;
287 
291 
292 
296 
297 
301 
302 
306 
307 
311 
313  SymbolicStokes& operator=(SymbolicStokes const& rOther);
314 
316  SymbolicStokes(SymbolicStokes const& rOther);
317 
319 }; // Class SymbolicStokes
320 
322 
325 
326 
330 
332 template< class TElementData >
333 inline std::istream& operator >>(
334  std::istream& rIStream,
336 {
337  return rIStream;
338 }
339 
341 template< class TElementData >
342 inline std::ostream& operator <<(
343  std::ostream& rOStream,
344  const SymbolicStokes<TElementData>& rThis)
345 {
346  rThis.PrintInfo(rOStream);
347  rOStream << std::endl;
348  rThis.PrintData(rOStream);
349 
350  return rOStream;
351 }
353 
355 } // namespace Kratos.
356 
357 #endif // KRATOS_SYMBOLIC_STOKES_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
MatrixRow< Matrix > ShapeFunctionsType
Type for shape function values container.
Definition: fluid_element.h:95
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: fluid_element.hpp:584
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
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
Definition: symbolic_stokes.h:57
constexpr static unsigned int LocalSize
Definition: symbolic_stokes.h:102
std::vector< std::size_t > EquationIdVectorType
Definition: symbolic_stokes.h:84
void AddTimeIntegratedSystem(TElementData &rData, MatrixType &rLHS, VectorType &rRHS) override
Definition: symbolic_stokes.cpp:167
constexpr static unsigned int Dim
Definition: symbolic_stokes.h:99
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Definition: symbolic_stokes.cpp:87
const Parameters GetSpecifications() const override
This method provides the specifications/requirements of the element.
Definition: symbolic_stokes.cpp:104
void AddTimeIntegratedRHS(TElementData &rData, VectorType &rRHS) override
Definition: symbolic_stokes.cpp:185
FluidElement< TElementData >::ShapeFunctionsType ShapeFunctionsType
Type for shape function values container.
Definition: symbolic_stokes.h:91
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: symbolic_stokes.cpp:193
std::string Info() const override
Turn back information as a string.
Definition: symbolic_stokes.cpp:145
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: symbolic_stokes.h:78
constexpr static unsigned int BlockSize
Definition: symbolic_stokes.h:101
FluidElement< TElementData >::ShapeFunctionDerivativesType ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: symbolic_stokes.h:94
constexpr static unsigned int NumNodes
Definition: symbolic_stokes.h:100
Node NodeType
Node type (default is: Node)
Definition: symbolic_stokes.h:66
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: symbolic_stokes.h:86
void Calculate(const Variable< double > &rVariable, double &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: symbolic_stokes.cpp:243
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(SymbolicStokes)
Pointer definition of SymbolicStokes.
virtual ~SymbolicStokes()
Destructor.
Definition: symbolic_stokes.cpp:58
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: symbolic_stokes.h:88
FluidElement< TElementData >::ShapeFunctionDerivativesArrayType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: symbolic_stokes.h:97
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: symbolic_stokes.h:72
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: symbolic_stokes.h:69
SymbolicStokes(IndexType NewId=0)
Default constuctor.
Definition: symbolic_stokes.cpp:31
Vector VectorType
Vector type for local contributions to the linear system.
Definition: symbolic_stokes.h:75
std::size_t SizeType
Definition: symbolic_stokes.h:82
constexpr static unsigned int StrainSize
Definition: symbolic_stokes.h:104
void ComputeGaussPointLHSContribution(TElementData &rData, MatrixType &rLHS)
std::size_t IndexType
Definition: symbolic_stokes.h:80
void AddTimeIntegratedLHS(TElementData &rData, MatrixType &rLHS) override
Definition: symbolic_stokes.cpp:177
void ComputeGaussPointRHSContribution(TElementData &rData, VectorType &rRHS)
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, Properties::Pointer pProperties) const override
Create a new element of this type.
Definition: symbolic_stokes.cpp:65
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: symbolic_stokes.cpp:153
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