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.
stokes_wall_condition.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 
14 #ifndef KRATOS_STOKES_WALL_CONDITION_H
15 #define KRATOS_STOKES_WALL_CONDITION_H
16 
17 // System includes
18 #include <string>
19 #include <iostream>
20 
21 
22 // External includes
23 
24 
25 // Project includes
26 #include "includes/define.h"
27 #include "includes/serializer.h"
28 #include "includes/condition.h"
29 #include "includes/process_info.h"
30 
31 // Application includes
34 #include "includes/cfd_variables.h"
35 #include "includes/kratos_flags.h"
36 
37 namespace Kratos
38 {
41 
44 
48 
52 
56 
60 
62 
66 template< unsigned int TDim, unsigned int TNumNodes = TDim >
68 {
69 public:
72 
75 
76  typedef Node NodeType;
77 
79 
81 
83 
84  typedef Vector VectorType;
85 
86  typedef Matrix MatrixType;
87 
88  typedef std::size_t IndexType;
89 
90  typedef std::size_t SizeType;
91 
92  typedef std::vector<std::size_t> EquationIdVectorType;
93 
94  typedef std::vector< Dof<double>::Pointer > DofsVectorType;
95 
97 
101 
103 
107  Condition(NewId)
108  {
109  }
110 
112 
117  const NodesArrayType& ThisNodes):
118  Condition(NewId,ThisNodes)
119  {
120  }
121 
123 
128  GeometryType::Pointer pGeometry):
129  Condition(NewId,pGeometry)
130  {
131  }
132 
134 
140  GeometryType::Pointer pGeometry,
141  PropertiesType::Pointer pProperties):
142  Condition(NewId,pGeometry,pProperties)
143  {
144  }
145 
148  Condition(rOther)
149  {
150  }
151 
153  ~StokesWallCondition() override {}
154 
155 
159 
162  {
163  Condition::operator=(rOther);
164 
165  return *this;
166  }
167 
171 
173 
178  Condition::Pointer Create(IndexType NewId, NodesArrayType const& ThisNodes, PropertiesType::Pointer pProperties) const override
179  {
180  return Kratos::make_intrusive<StokesWallCondition>(NewId, GetGeometry().Create(ThisNodes), pProperties);
181  }
182 
183 
185 
190  Condition::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, PropertiesType::Pointer pProperties) const override {
191  return Kratos::make_intrusive< StokesWallCondition >(NewId, pGeom, pProperties);
192  }
193 
194 
196 
200  MatrixType& rLeftHandSideMatrix,
201  VectorType& rRightHandSideVector,
202  const ProcessInfo& rCurrentProcessInfo) override
203  {
204  const SizeType BlockSize = TDim + 1;
205  const SizeType LocalSize = BlockSize * TNumNodes;
206 
207  if (rLeftHandSideMatrix.size1() != LocalSize)
208  rLeftHandSideMatrix.resize(LocalSize,LocalSize);
209 
210  if (rRightHandSideVector.size() != LocalSize)
211  rRightHandSideVector.resize(LocalSize);
212 
213  noalias(rLeftHandSideMatrix) = ZeroMatrix(LocalSize,LocalSize);
214  noalias(rRightHandSideVector) = ZeroVector(LocalSize);
215 
216  if( this->Is(OUTLET) )
217  {
218  ApplyNeumannCondition(rLeftHandSideMatrix,rRightHandSideVector);
219  }
220  }
221 
223 
227  MatrixType& rLeftHandSideMatrix,
228  const ProcessInfo& rCurrentProcessInfo) override
229  {
230  const SizeType BlockSize = TDim + 1;
231  const SizeType LocalSize = BlockSize * TNumNodes;
232 
233  if (rLeftHandSideMatrix.size1() != LocalSize)
234  rLeftHandSideMatrix.resize(LocalSize,LocalSize);
235 
236  noalias(rLeftHandSideMatrix) = ZeroMatrix(LocalSize,LocalSize);
237  }
238 
240 
244  VectorType& rRightHandSideVector,
245  const ProcessInfo& rCurrentProcessInfo) override
246  {
247  const SizeType BlockSize = TDim + 1;
248  const SizeType LocalSize = BlockSize * TNumNodes;
249 
250  if (rRightHandSideVector.size() != LocalSize)
251  rRightHandSideVector.resize(LocalSize);
252 
253  noalias(rRightHandSideVector) = ZeroVector(LocalSize);
254 
255  if( this->Is(OUTLET) )
256  {
257  Matrix tmp;
258  ApplyNeumannCondition(tmp,rRightHandSideVector);
259  }
260  }
261 
262 
263 
265  int Check(const ProcessInfo& rCurrentProcessInfo) const override
266  {
267  KRATOS_TRY;
268 
269  int Check = Condition::Check(rCurrentProcessInfo); // Checks id > 0 and area > 0
270 
271  if (Check != 0)
272  {
273  return Check;
274  }
275  else
276  {
277  // Checks on nodes
278  // Check that the element's nodes contain all required SolutionStepData and Degrees of freedom
279  for(unsigned int i=0; i<this->GetGeometry().size(); ++i)
280  {
281  if(this->GetGeometry()[i].SolutionStepsDataHas(VELOCITY) == false)
282  KRATOS_THROW_ERROR(std::invalid_argument,"missing VELOCITY variable on solution step data for node ",this->GetGeometry()[i].Id());
283  if(this->GetGeometry()[i].SolutionStepsDataHas(PRESSURE) == false)
284  KRATOS_THROW_ERROR(std::invalid_argument,"missing PRESSURE variable on solution step data for node ",this->GetGeometry()[i].Id());
285  if(this->GetGeometry()[i].SolutionStepsDataHas(EXTERNAL_PRESSURE) == false)
286  KRATOS_THROW_ERROR(std::invalid_argument,"missing EXTERNAL_PRESSURE variable on solution step data for node ",this->GetGeometry()[i].Id());
287  if(this->GetGeometry()[i].HasDofFor(VELOCITY_X) == false ||
288  this->GetGeometry()[i].HasDofFor(VELOCITY_Y) == false ||
289  this->GetGeometry()[i].HasDofFor(VELOCITY_Z) == false)
290  KRATOS_THROW_ERROR(std::invalid_argument,"missing VELOCITY component degree of freedom on node ",this->GetGeometry()[i].Id());
291  if(this->GetGeometry()[i].HasDofFor(PRESSURE) == false)
292  KRATOS_THROW_ERROR(std::invalid_argument,"missing PRESSURE component degree of freedom on node ",this->GetGeometry()[i].Id());
293  }
294 
295  return Check;
296  }
297 
298  KRATOS_CATCH("");
299  }
300 
301 
303 
308  EquationIdVectorType& rResult,
309  const ProcessInfo& rCurrentProcessInfo) const override;
310 
311 
313 
318  DofsVectorType& ConditionDofList,
319  const ProcessInfo& CurrentProcessInfo) const override;
320 
321 
325 
326 
330 
331 
335 
337  std::string Info() const override
338  {
339  std::stringstream buffer;
340  buffer << "StokesWallCondition" << TDim << "D";
341  return buffer.str();
342  }
343 
345  void PrintInfo(std::ostream& rOStream) const override
346  {
347  rOStream << "StokesWallCondition";
348  }
349 
351  void PrintData(std::ostream& rOStream) const override {}
352 
353 
357 
358 
360 
361 protected:
364 
365 
369 
370 
374 
375 
379 
380 
381  void ApplyNeumannCondition(MatrixType &rLocalMatrix, VectorType &rLocalVector);
382 
383 
387 
388 
392 
393 
397 
398 
400 
401 private:
404 
405 
409 
410 
414 
415  friend class Serializer;
416 
417  void save(Serializer& rSerializer) const override
418  {
420  }
421 
422  void load(Serializer& rSerializer) override
423  {
425  }
426 
430 
431 
435 
436 
440 
441 
445 
446 
450 
451 
453 
454 }; // Class StokesWallCondition
455 
456 
458 
461 
462 
466 
467 
469 template< unsigned int TDim, unsigned int TNumNodes >
470 inline std::istream& operator >> (std::istream& rIStream,
472 {
473  return rIStream;
474 }
475 
477 template< unsigned int TDim, unsigned int TNumNodes >
478 inline std::ostream& operator << (std::ostream& rOStream,
480 {
481  rThis.PrintInfo(rOStream);
482  rOStream << std::endl;
483  rThis.PrintData(rOStream);
484 
485  return rOStream;
486 }
487 
489 
491 
492 
493 } // namespace Kratos.
494 
495 #endif // KRATOS_STOKES_WALL_CONDITION_H
Base class for all Conditions.
Definition: condition.h:59
std::size_t SizeType
Definition: condition.h:94
std::vector< DofType::Pointer > DofsVectorType
Definition: condition.h:100
virtual int Check(const ProcessInfo &rCurrentProcessInfo) const
Definition: condition.h:854
Condition & operator=(Condition const &rOther)
Assignment operator.
Definition: condition.h:181
std::size_t IndexType
Definition: flags.h:74
bool Is(Flags const &rOther) const
Definition: flags.h:274
GeometryType & GetGeometry()
Returns the reference of the geometry.
Definition: geometrical_object.h:158
Geometry base class.
Definition: geometry.h:71
SizeType size() const
Definition: geometry.h:518
IndexType const & Id() const
Id of this Geometry.
Definition: geometry.h:964
This object defines an indexed object.
Definition: indexed_object.h:54
void resize(std::size_t NewSize1, std::size_t NewSize2, bool preserve=0)
Definition: amatrix_interface.h:224
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
Properties encapsulates data shared by different Elements or Conditions. It can store any type of dat...
Definition: properties.h:69
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
Implements a wall condition for the monolithic stokes formulation - based on BDF2.
Definition: stokes_wall_condition.h:68
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: stokes_wall_condition.h:96
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: stokes_wall_condition.h:345
~StokesWallCondition() override
Destructor.
Definition: stokes_wall_condition.h:153
Vector VectorType
Definition: stokes_wall_condition.h:84
Geometry< NodeType > GeometryType
Definition: stokes_wall_condition.h:80
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: stokes_wall_condition.h:94
StokesWallCondition(IndexType NewId=0)
Default constructor.
Definition: stokes_wall_condition.h:106
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) override
Return a matrix of the correct size, filled with zeros (for compatibility with time schemes).
Definition: stokes_wall_condition.h:226
void GetDofList(DofsVectorType &ConditionDofList, const ProcessInfo &CurrentProcessInfo) const override
Returns a list of the element's Dofs.
StokesWallCondition(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Constructor using Properties.
Definition: stokes_wall_condition.h:139
Properties PropertiesType
Definition: stokes_wall_condition.h:78
std::vector< std::size_t > EquationIdVectorType
Definition: stokes_wall_condition.h:92
std::string Info() const override
Turn back information as a string.
Definition: stokes_wall_condition.h:337
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Check that all data required by this condition is available and reasonable.
Definition: stokes_wall_condition.h:265
Condition::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, PropertiesType::Pointer pProperties) const override
Create a new StokesWallCondition object.
Definition: stokes_wall_condition.h:178
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition: stokes_wall_condition.h:82
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(StokesWallCondition)
Pointer definition of StokesWallCondition.
Matrix MatrixType
Definition: stokes_wall_condition.h:86
StokesWallCondition(StokesWallCondition const &rOther)
Copy constructor.
Definition: stokes_wall_condition.h:147
StokesWallCondition(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: stokes_wall_condition.h:116
Node NodeType
Definition: stokes_wall_condition.h:76
void ApplyNeumannCondition(MatrixType &rLocalMatrix, VectorType &rLocalVector)
Definition: stokes_wall_condition.cpp:156
StokesWallCondition & operator=(StokesWallCondition const &rOther)
Assignment operator.
Definition: stokes_wall_condition.h:161
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Return local right hand side of the correct size, filled with zeros (for compatibility with time sche...
Definition: stokes_wall_condition.h:243
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: stokes_wall_condition.h:351
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &rCurrentProcessInfo) const override
Provides the global indices for each one of this element's local rows.
std::size_t SizeType
Definition: stokes_wall_condition.h:90
Condition::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, PropertiesType::Pointer pProperties) const override
Create a new StokesWallCondition object.
Definition: stokes_wall_condition.h:190
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Return local contributions of the correct size, filled with zeros (for compatibility with time scheme...
Definition: stokes_wall_condition.h:199
StokesWallCondition(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using Geometry.
Definition: stokes_wall_condition.h:127
std::size_t IndexType
Definition: stokes_wall_condition.h:88
#define KRATOS_THROW_ERROR(ExceptionType, ErrorMessage, MoreInfo)
Definition: define.h:77
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_TRY
Definition: define.h:109
#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
KratosZeroVector< double > ZeroVector
Definition: amatrix_interface.h:561
KratosZeroMatrix< double > ZeroMatrix
Definition: amatrix_interface.h:559
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
T & noalias(T &TheMatrix)
Definition: amatrix_interface.h:484
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
tuple tmp
Definition: generate_total_lagrangian_mixed_volumetric_strain_element.py:98
def load(f)
Definition: ode_solve.py:307
integer i
Definition: TensorModule.f:17