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.
fic.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 // Ignasi de Pouplana
12 //
13 
14 #ifndef KRATOS_FIC_FLUID_ELEMENT_H
15 #define KRATOS_FIC_FLUID_ELEMENT_H
16 
17 #include "includes/define.h"
18 #include "includes/element.h"
19 #include "includes/serializer.h"
20 #include "geometries/geometry.h"
21 
22 #include "includes/cfd_variables.h"
25 
26 namespace Kratos
27 {
28 
31 
34 
38 
42 
46 
50 
54 
55 template< class TElementData >
56 class FIC : public FluidElement<TElementData>
57 {
58 public:
61 
64 
67 
69  typedef Node NodeType;
70 
73 
76 
78  typedef Vector VectorType;
79 
81  typedef Matrix MatrixType;
82 
83  typedef std::size_t IndexType;
84 
85  typedef std::size_t SizeType;
86 
87  typedef std::vector<std::size_t> EquationIdVectorType;
88 
89  typedef std::vector< Dof<double>::Pointer > DofsVectorType;
90 
92 
95 
98 
101 
102  constexpr static unsigned int Dim = BaseType::Dim;
103  constexpr static unsigned int NumNodes = BaseType::NumNodes;
104  constexpr static unsigned int BlockSize = BaseType::BlockSize;
105  constexpr static unsigned int LocalSize = BaseType::LocalSize;
106  constexpr static unsigned int StrainSize = BaseType::StrainSize;
107 
111 
112  //Constructors.
113 
115 
118  FIC(IndexType NewId = 0);
119 
121 
125  FIC(IndexType NewId, const NodesArrayType& ThisNodes);
126 
128 
132  FIC(IndexType NewId, GeometryType::Pointer pGeometry);
133 
135 
140  FIC(IndexType NewId, GeometryType::Pointer pGeometry, Properties::Pointer pProperties);
141 
143  ~FIC() override;
144 
148 
149 
153 
154 
156 
163  Element::Pointer Create(IndexType NewId,
164  NodesArrayType const& ThisNodes,
165  Properties::Pointer pProperties) const override;
166 
168 
175  Element::Pointer Create(IndexType NewId,
176  GeometryType::Pointer pGeom,
177  Properties::Pointer pProperties) const override;
178 
179  void Calculate(
180  const Variable<double>& rVariable,
181  double& rOutput,
182  const ProcessInfo& rCurrentProcessInfo) override;
183 
184 
185  void Calculate(
186  const Variable<array_1d<double, 3 > >& rVariable,
187  array_1d<double, 3 > & rOutput,
188  const ProcessInfo& rCurrentProcessInfo) override;
189 
190 
191  void Calculate(
192  const Variable<Vector >& rVariable,
193  Vector& Output,
194  const ProcessInfo& rCurrentProcessInfo) override;
195 
196  void Calculate(
197  const Variable<Matrix >& rVariable,
198  Matrix& Output,
199  const ProcessInfo& rCurrentProcessInfo) override;
200 
204 
206  Variable<array_1d<double, 3>> const& rVariable,
207  std::vector<array_1d<double, 3>>& rValues,
208  ProcessInfo const& rCurrentProcessInfo) override;
209 
211  Variable<double> const& rVariable,
212  std::vector<double>& rValues,
213  ProcessInfo const& rCurrentProcessInfo) override;
214 
216  Variable<array_1d<double, 6>> const& rVariable,
217  std::vector<array_1d<double, 6>>& rValues,
218  ProcessInfo const& rCurrentProcessInfo) override;
219 
221  Variable<Vector> const& rVariable,
222  std::vector<Vector>& rValues,
223  ProcessInfo const& rCurrentProcessInfo) override;
224 
226  Variable<Matrix> const& rVariable,
227  std::vector<Matrix>& rValues,
228  ProcessInfo const& rCurrentProcessInfo) override;
229 
233 
234  int Check(const ProcessInfo &rCurrentProcessInfo) const override;
235 
239 
240  const Parameters GetSpecifications() const override;
241 
243  std::string Info() const override;
244 
245 
247  void PrintInfo(std::ostream& rOStream) const override;
248 
249 
253 
254 
256 
257 protected:
258 
261 
262 
266 
267 
271 
272 
276 
278  TElementData& rData,
279  MatrixType& rLHS,
280  VectorType& rRHS) override;
281 
283  TElementData& rData,
284  MatrixType& rLHS) override;
285 
287  TElementData& rData,
288  VectorType& rRHS) override;
289 
290  void AddVelocitySystem(
291  TElementData& rData,
292  MatrixType& rLocalLHS,
293  VectorType& rLocalRHS) override;
294 
295  void AddMassLHS(
296  TElementData& rData,
297  MatrixType& rMassMatrix) override;
298 
299  // This function integrates the traction over a cut. It is only required to implement embedded formulations
300  void AddBoundaryTraction(
301  TElementData& rData,
302  const Vector& rUnitNormal,
303  MatrixType& rLHS,
304  VectorType& rRHS) override;
305 
306  // Implementation details of FIC ////////////////////////////////////////
307 
309  TElementData& rData,
310  MatrixType& rMassMatrix);
311 
312  void AddViscousTerm(
313  const TElementData& rData,
315  VectorType& rRHS);
316 
317  virtual void CalculateTau(
318  const TElementData& rData,
319  const array_1d<double,3> &Velocity,
320  double &TauIncompr,
321  double &TauMomentum,
322  array_1d<double,3> &TauGrad) const;
323 
324  virtual void CalculateTauGrad(
325  const TElementData& rData,
326  array_1d<double,3> &TauGrad) const;
327 
328  virtual void AlgebraicMomentumResidual(
329  const TElementData& rData,
330  const Vector& rConvection,
331  array_1d<double,3>& rResidual) const;
332 
336 
337 
341 
342 
346 
348 
349 private:
350 
353 
357 
361 
362  friend class Serializer;
363 
364  void save(Serializer& rSerializer) const override;
365 
366  void load(Serializer& rSerializer) override;
367 
371 
372 
376 
377 
381 
382 
386 
387 
391 
393  FIC& operator=(FIC const& rOther);
394 
396  FIC(FIC const& rOther);
397 
399 
400 
401 }; // Class FIC
402 
404 
407 
408 
412 
413 
415 template< class TElementData >
416 inline std::istream& operator >>(std::istream& rIStream,
417  FIC<TElementData>& rThis)
418 {
419  return rIStream;
420 }
421 
423 template< class TElementData >
424 inline std::ostream& operator <<(std::ostream& rOStream,
425  const FIC<TElementData>& rThis)
426 {
427  rThis.PrintInfo(rOStream);
428  rOStream << std::endl;
429  rThis.PrintData(rOStream);
430 
431  return rOStream;
432 }
434 
435 
436 namespace Internals {
437 
438 template <class TElementData, bool TDataKnowsAboutTimeIntegration>
440  public:
441  static void AddSystem(FIC<TElementData>* pElement,
442  TElementData& rData, Kratos::Matrix& rLHS, Vector& rRHS);
443 };
444 
445 template <class TElementData>
446 class FICSpecializedAddTimeIntegratedSystem<TElementData, true> {
447  public:
448  static void AddSystem(FIC<TElementData>* pElement,
449  TElementData& rData, Kratos::Matrix& rLHS, Vector& rRHS);
450 };
451 
452 template <class TElementData>
453 class FICSpecializedAddTimeIntegratedSystem<TElementData, false> {
454  public:
455  static void AddSystem(FIC<TElementData>* pElement,
456  TElementData& rData, Kratos::Matrix& rLHS, Vector& rRHS);
457 };
458 
460 
461 } // namespace Internals
462 
463 } // namespace Kratos.
464 
465 #endif // KRATOS_FIC_FLUID_ELEMENT_H
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Definition: fic.h:57
Large Displacement Lagrangian Element for 3D and 2D geometries. (base class)
Definition: fluid_element.h:61
static constexpr unsigned int Dim
Definition: fluid_element.h:105
static constexpr unsigned int NumNodes
Definition: fluid_element.h:107
static constexpr unsigned int LocalSize
Definition: fluid_element.h:111
static constexpr unsigned int BlockSize
Definition: fluid_element.h:109
static constexpr unsigned int StrainSize
Definition: fluid_element.h:113
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
static void AddSystem(FIC< TElementData > *pElement, TElementData &rData, Kratos::Matrix &rLHS, Vector &rRHS)
static void AddSystem(FIC< TElementData > *pElement, TElementData &rData, Kratos::Matrix &rLHS, Vector &rRHS)
static void AddSystem(FIC< TElementData > *pElement, TElementData &rData, Kratos::Matrix &rLHS, Vector &rRHS)
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
std::string Info() const override
Turn back information as a string.
Definition: fic.cpp:215
void AddMassLHS(TElementData &rData, MatrixType &rMassMatrix) override
Definition: fic.cpp:389
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: fic.h:81
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: fic.h:97
constexpr static unsigned int NumNodes
Definition: fic.h:103
constexpr static unsigned int BlockSize
Definition: fic.h:104
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: fic.h:100
void Calculate(const Variable< array_1d< double, 3 > > &rVariable, array_1d< double, 3 > &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Vector VectorType
Vector type for local contributions to the linear system.
Definition: fic.h:78
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: fic.h:89
std::size_t SizeType
Definition: fic.h:85
constexpr static unsigned int LocalSize
Definition: fic.h:105
void Calculate(const Variable< double > &rVariable, double &rOutput, const ProcessInfo &rCurrentProcessInfo) override
Definition: fic.cpp:74
void CalculateOnIntegrationPoints(Variable< array_1d< double, 3 >> const &rVariable, std::vector< array_1d< double, 3 >> &rValues, ProcessInfo const &rCurrentProcessInfo) override
std::size_t IndexType
Definition: fic.h:83
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: fic.h:94
virtual void CalculateTau(const TElementData &rData, const array_1d< double, 3 > &Velocity, double &TauIncompr, double &TauMomentum, array_1d< double, 3 > &TauGrad) const
Definition: fic.cpp:525
std::vector< std::size_t > EquationIdVectorType
Definition: fic.h:87
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: fic.cpp:224
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: fic.h:75
~FIC() override
Destructor.
Definition: fic.cpp:54
void AddMassStabilization(TElementData &rData, MatrixType &rMassMatrix)
Definition: fic.cpp:417
void AddTimeIntegratedSystem(TElementData &rData, MatrixType &rLHS, VectorType &rRHS) override
Definition: fic.cpp:257
virtual void CalculateTauGrad(const TElementData &rData, array_1d< double, 3 > &TauGrad) const
Definition: fic.cpp:581
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: fic.h:72
void AddVelocitySystem(TElementData &rData, MatrixType &rLocalLHS, VectorType &rLocalRHS) override
Definition: fic.cpp:279
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, Properties::Pointer pProperties) const override
Create a new element of this type.
Definition: fic.cpp:61
void AddViscousTerm(const TElementData &rData, BoundedMatrix< double, LocalSize, LocalSize > &rLHS, VectorType &rRHS)
Definition: fic.cpp:504
const Parameters GetSpecifications() const override
This method provides the specifications/requirements of the element.
Definition: fic.cpp:181
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Definition: fic.cpp:113
virtual void AlgebraicMomentumResidual(const TElementData &rData, const Vector &rConvection, array_1d< double, 3 > &rResidual) const
Definition: fic.cpp:233
Node NodeType
Node type (default is: Node)
Definition: fic.h:69
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(FIC)
Pointer definition of FIC.
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: fic.h:91
constexpr static unsigned int Dim
Definition: fic.h:102
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: fic.cpp:461
void AddTimeIntegratedLHS(TElementData &rData, MatrixType &rLHS) override
Definition: fic.cpp:267
FIC(IndexType NewId=0)
Default constuctor.
Definition: fic.cpp:31
void AddTimeIntegratedRHS(TElementData &rData, VectorType &rRHS) override
Definition: fic.cpp:273
constexpr static unsigned int StrainSize
Definition: fic.h:106
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