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.
convection_diffusion_reaction_cross_wind_stabilized_element.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: Suneth Warnakulasuriya
11 //
12 
13 #if !defined(KRATOS_CONVECTION_DIFFUSION_REACTION_CROSS_WIND_STABILIZED_ELEMENT_H_INCLUDED)
14 #define KRATOS_CONVECTION_DIFFUSION_REACTION_CROSS_WIND_STABILIZED_ELEMENT_H_INCLUDED
15 
16 // System includes
17 #include <cmath>
18 
19 // External includes
20 
21 // Project includes
22 #include "includes/define.h"
23 #include "includes/element.h"
24 
25 // Application includes
27 
28 namespace Kratos
29 {
32 
33 template <unsigned int TDim,
34  unsigned int TNumNodes,
35  class TConvectionDiffusionReactionData>
37 : public ConvectionDiffusionReactionElement<TDim, TNumNodes, TConvectionDiffusionReactionData>
38 {
39 public:
42 
44 
46  using NodeType = Node;
47 
50 
53 
55  using VectorType = Vector;
56 
58  using MatrixType = Matrix;
59 
60  using IndexType = std::size_t;
61 
62  using EquationIdVectorType = std::vector<IndexType>;
63 
64  using DofsVectorType = std::vector<Dof<double>::Pointer>;
65 
67 
70 
71  using ConvectionDiffusionReactionDataType = TConvectionDiffusionReactionData;
72 
75 
80 
84 
89  IndexType NewId = 0)
90  : BaseType(NewId)
91  {
92  }
93 
98  IndexType NewId,
99  const NodesArrayType& ThisNodes)
100  : BaseType(NewId, ThisNodes)
101  {
102  }
103 
108  IndexType NewId,
109  GeometryType::Pointer pGeometry)
110  : BaseType(NewId, pGeometry)
111  {
112  }
113 
118  IndexType NewId,
119  GeometryType::Pointer pGeometry,
120  typename PropertiesType::Pointer pProperties)
121  : BaseType(NewId, pGeometry, pProperties)
122  {
123  }
124 
130  : BaseType(rOther)
131  {
132  }
133 
138 
142 
155  Element::Pointer Create(
156  IndexType NewId,
157  NodesArrayType const& ThisNodes,
158  typename PropertiesType::Pointer pProperties) const override
159  {
160  KRATOS_TRY
161  return Kratos::make_intrusive<CurrentElementType>(
162  NewId, this->GetGeometry().Create(ThisNodes), pProperties);
163  KRATOS_CATCH("");
164  }
165 
173  Element::Pointer Create(
174  IndexType NewId,
175  GeometryType::Pointer pGeom,
176  typename PropertiesType::Pointer pProperties) const override
177  {
178  KRATOS_TRY
179  return Kratos::make_intrusive<CurrentElementType>(NewId, pGeom, pProperties);
180  KRATOS_CATCH("");
181  }
182 
190  Element::Pointer Clone(
191  IndexType NewId,
192  NodesArrayType const& ThisNodes) const override
193  {
194  KRATOS_TRY
195  return Kratos::make_intrusive<CurrentElementType>(
196  NewId, this->GetGeometry().Create(ThisNodes), this->pGetProperties());
197  KRATOS_CATCH("");
198  }
199 
214  VectorType& rRightHandSideVector,
215  const ProcessInfo& rCurrentProcessInfo) override;
216 
224  MatrixType& rDampingMatrix,
225  VectorType& rRightHandSideVector,
226  const ProcessInfo& rCurrentProcessInfo) override;
227 
243  void CalculateMassMatrix(
244  MatrixType& rMassMatrix,
245  const ProcessInfo& rCurrentProcessInfo) override;
246 
254  MatrixType& rDampingMatrix,
255  const ProcessInfo& rCurrentProcessInfo) override;
256 
260 
262  std::string Info() const override
263  {
264  std::stringstream buffer;
265  buffer << "ConvectionDiffusionReactionCrossWindStabilizedElement #" << this->Id();
266  return buffer.str();
267  }
268 
269  void PrintInfo(std::ostream& rOStream) const override
270  {
271  rOStream << "CDRCrossWind" << TConvectionDiffusionReactionData::GetName();
272  }
273 
277 
279 
280 protected:
283 
286  const Matrix& rParameterDerivatives) const;
287 
288  virtual double GetDeltaTime(
289  const ProcessInfo& rProcessInfo) const;
290 
304 
306 
307 private:
310 
311  friend class Serializer;
312 
313  void save(Serializer& rSerializer) const override
314  {
315  KRATOS_TRY
316 
318 
319  KRATOS_CATCH("");
320  }
321  void load(Serializer& rSerializer) override
322  {
323  KRATOS_TRY
324 
326 
327  KRATOS_CATCH("");
328  }
329 
331 }; // Class ConvectionDiffusionReactionCrossWindStabilizedElement
332 
336 
337 template <unsigned int TDim, unsigned int TNumNodes, class TConvectionDiffusionReactionData>
338 inline std::istream& operator>>(
339  std::istream& rIStream,
341 
343 template <unsigned int TDim, unsigned int TNumNodes, class TConvectionDiffusionReactionData>
344 inline std::ostream& operator<<(
345  std::ostream& rOStream,
347 {
348  rThis.PrintInfo(rOStream);
349  rOStream << " : " << std::endl;
350  rThis.PrintData(rOStream);
351  return rOStream;
352 }
353 
355 
356 } // namespace Kratos.
357 
358 #endif // KRATOS_CONVECTION_DIFFUSION_REACTION_CROSS_WIND_STABILIZED_ELEMENT_H_INCLUDED defined
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:38
ShapeFunctionDerivativesArrayType GetGeometryParameterDerivatives() const
Get the Geometry Parameter Derivatives object.
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.cpp:358
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:269
void CalculateMassMatrix(MatrixType &rMassMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.cpp:136
void CalculateDampingMatrix(MatrixType &rDampingMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.cpp:205
virtual double GetDeltaTime(const ProcessInfo &rProcessInfo) const
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.cpp:348
std::string Info() const override
Turn back information as a string.
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:262
ConvectionDiffusionReactionCrossWindStabilizedElement(IndexType NewId, GeometryType::Pointer pGeometry, typename PropertiesType::Pointer pProperties)
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:117
ConvectionDiffusionReactionCrossWindStabilizedElement(IndexType NewId=0)
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:88
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(ConvectionDiffusionReactionCrossWindStabilizedElement)
Element::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, typename PropertiesType::Pointer pProperties) const override
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:173
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, typename PropertiesType::Pointer pProperties) const override
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:155
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:58
ConvectionDiffusionReactionCrossWindStabilizedElement(IndexType NewId, const NodesArrayType &ThisNodes)
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:97
void CalculateContravariantMetricTensor(BoundedMatrix< double, TDim, TDim > &rOutput, const Matrix &rParameterDerivatives) const
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.cpp:340
ConvectionDiffusionReactionCrossWindStabilizedElement(IndexType NewId, GeometryType::Pointer pGeometry)
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:107
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:69
Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:190
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.cpp:46
void CalculateLocalVelocityContribution(MatrixType &rDampingMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
CalculateLocalVelocityContribution Calculate the local contribution in terms of velocity and pressure...
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.cpp:122
TConvectionDiffusionReactionData ConvectionDiffusionReactionDataType
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:71
ConvectionDiffusionReactionCrossWindStabilizedElement(ConvectionDiffusionReactionCrossWindStabilizedElement const &rOther)
Definition: convection_diffusion_reaction_cross_wind_stabilized_element.h:128
Definition: convection_diffusion_reaction_element.h:36
typename BaseType::PropertiesType PropertiesType
Definition: convection_diffusion_reaction_element.h:67
Base class for all Elements.
Definition: element.h:60
std::vector< DofType::Pointer > DofsVectorType
Definition: element.h:100
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: element.h:1135
std::vector< std::size_t > EquationIdVectorType
Definition: element.h:98
PropertiesType::Pointer pGetProperties()
returns the pointer to the property of the element. Does not throw an error, to allow copying of elem...
Definition: element.h:1014
std::size_t IndexType
Definition: flags.h:74
This defines the geometrical object, base definition of the element and condition entities.
Definition: geometrical_object.h:58
GeometryType & GetGeometry()
Returns the reference of the geometry.
Definition: geometrical_object.h:158
Geometry base class.
Definition: geometry.h:71
GeometryData::ShapeFunctionsGradientsType ShapeFunctionsGradientsType
Definition: geometry.h:189
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
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
#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
Internals::Matrix< double, AMatrix::dynamic, 1 > Vector
Definition: amatrix_interface.h:472
Internals::Matrix< double, AMatrix::dynamic, AMatrix::dynamic > Matrix
Definition: amatrix_interface.h:470
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