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.
eulerian_conv_diff.h
Go to the documentation of this file.
1 // KRATOS ___ ___ _ ___ __ ___ ___ ___ ___
2 // / __/ _ \| \| \ \ / /__| \_ _| __| __|
3 // | (_| (_) | .` |\ V /___| |) | || _|| _|
4 // \___\___/|_|\_| \_/ |___/___|_| |_| APPLICATION
5 //
6 // License: BSD License
7 // Kratos default license: kratos/license.txt
8 //
9 // Main authors: Riccardo Rossi
10 //
11 
12 #if !defined(KRATOS_EULERIAN_CONVECTION_DIFFUSION_ELEMENT_INCLUDED )
13 #define KRATOS_EULERIAN_CONVECTION_DIFFUSION_ELEMENT_INCLUDED
14 
15 
16 // System includes
17 
18 
19 // External includes
20 
21 
22 // Project includes
23 #include "includes/define.h"
24 #include "includes/element.h"
26 #include "includes/variables.h"
27 #include "includes/cfd_variables.h"
28 #include "includes/serializer.h"
29 #include "utilities/math_utils.h"
30 #include "utilities/geometry_utilities.h"
31 
32 
33 
34 namespace Kratos
35 {
36 
38 template< unsigned int TDim, unsigned int TNumNodes>
39 class KRATOS_API(CONVECTION_DIFFUSION_APPLICATION) EulerianConvectionDiffusionElement
40  : public Element
41 {
42 public:
45 
46 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
47 
49 
51  {
52  }
53 
54  EulerianConvectionDiffusionElement(IndexType NewId, GeometryType::Pointer pGeometry)
55  : Element(NewId, pGeometry)
56  {}
57 
58  EulerianConvectionDiffusionElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
59  : Element(NewId, pGeometry, pProperties)
60  {}
61 
64 
65 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
66 
67  Element::Pointer Create(
68  IndexType NewId,
69  NodesArrayType const& ThisNodes,
70  PropertiesType::Pointer pProperties
71  ) const override
72  {
73  return Kratos::make_intrusive<EulerianConvectionDiffusionElement>(NewId, GetGeometry().Create(ThisNodes), pProperties);
74  }
75 
76  Element::Pointer Create(
77  IndexType NewId,
78  GeometryType::Pointer pGeom,
79  PropertiesType::Pointer pProperties
80  ) const override
81  {
82  return Kratos::make_intrusive<EulerianConvectionDiffusionElement>(NewId, pGeom, pProperties);
83  }
84 
85  void EquationIdVector(EquationIdVectorType& rResult, const ProcessInfo& rCurrentProcessInfo) const override;
86 
87  void GetDofList(DofsVectorType& ElementalDofList, const ProcessInfo& rCurrentProcessInfo) const override;
88 
89  void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix, VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
90 
91  void CalculateRightHandSide(VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
92 
93 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
94 
95  std::string Info() const override
96  {
97  return "EulerianConvectionDiffusionElement #";
98  }
99 
101 
102  void PrintInfo(std::ostream& rOStream) const override
103  {
104  rOStream << Info() << Id();
105  }
106 
107 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
108 
109 protected:
110 
112  {
113  double theta;
114  double dyn_st_beta;
115  double dt_inv;
117  double conductivity;
119  double density;
120  double beta;
121  double div_v;
122 
128  };
129 
130  void InitializeEulerianElement(ElementVariables& rVariables, const ProcessInfo& rCurrentProcessInfo);
131 
132  void CalculateGeometry(BoundedMatrix<double,TNumNodes,TDim>& rDN_DX, double& rVolume);
133 
134  double ComputeH(BoundedMatrix<double,TNumNodes,TDim>& rDN_DX);
135 
136  void GetNodalValues(ElementVariables& rVariables, const ProcessInfo& rCurrentProcessInfo) const;
137 
138  double CalculateTau(const ElementVariables& rVariables, double norm_vel, double h);
139 
140  // Member Variables
141 
142 
143 //----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
144 
145 private:
146 
147  // Serialization
148 
149  friend class Serializer;
150 
151  void save(Serializer& rSerializer) const override
152  {
154  }
155 
156  void load(Serializer& rSerializer) override
157  {
159  }
160 
161 
162 }; // Class EulerianConvectionDiffusionElement
163 
164 } // namespace Kratos.
165 
166 #endif // KRATOS_EULERIAN_CONVECTION_DIFFUSION_ELEMENT_INCLUDED defined
std::string Info() const override
Turn back information as a string.
Definition: periodic_interface_process.hpp:93
Base class for all Elements.
Definition: element.h:60
formulation described in https://docs.google.com/document/d/13a_zGLj6xORDuLgoOG5LwHI6BwShvfO166opZ815...
Definition: eulerian_conv_diff.h:41
EulerianConvectionDiffusionElement()
Default constructor.
Definition: eulerian_conv_diff.h:50
virtual ~EulerianConvectionDiffusionElement()
Destructor.
Definition: eulerian_conv_diff.h:63
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(EulerianConvectionDiffusionElement)
Counted pointer of.
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, PropertiesType::Pointer pProperties) const override
It creates a new element pointer.
Definition: eulerian_conv_diff.h:67
Element::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, PropertiesType::Pointer pProperties) const override
It creates a new element pointer.
Definition: eulerian_conv_diff.h:76
EulerianConvectionDiffusionElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Definition: eulerian_conv_diff.h:58
EulerianConvectionDiffusionElement(IndexType NewId, GeometryType::Pointer pGeometry)
Definition: eulerian_conv_diff.h:54
std::string Info() const override
Turn back information as a string.
Definition: eulerian_conv_diff.h:95
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: eulerian_conv_diff.h:102
std::size_t IndexType
Definition: flags.h:74
Definition: amatrix_interface.h:41
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
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_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
Matrix MatrixType
Definition: geometrical_transformation_utilities.h:55
Modeler::Pointer Create(const std::string &ModelerName, Model &rModel, const Parameters ModelParameters)
Checks if the modeler is registered.
Definition: modeler_factory.cpp:30
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
h
Definition: generate_droplet_dynamics.py:91
def load(f)
Definition: ode_solve.py:307
double div_v
Definition: eulerian_conv_diff.h:121
double lumping_factor
Definition: eulerian_conv_diff.h:116
double conductivity
Definition: eulerian_conv_diff.h:117
double specific_heat
Definition: eulerian_conv_diff.h:118
double theta
Definition: eulerian_conv_diff.h:113
double dyn_st_beta
Definition: eulerian_conv_diff.h:114
array_1d< double, TNumNodes > volumetric_source
Definition: eulerian_conv_diff.h:125
array_1d< array_1d< double, 3 >, TNumNodes > v
Definition: eulerian_conv_diff.h:126
array_1d< double, TNumNodes > phi
Definition: eulerian_conv_diff.h:123
double dt_inv
Definition: eulerian_conv_diff.h:115
double beta
Definition: eulerian_conv_diff.h:120
array_1d< array_1d< double, 3 >, TNumNodes > vold
Definition: eulerian_conv_diff.h:127
array_1d< double, TNumNodes > phi_old
Definition: eulerian_conv_diff.h:124
double density
Definition: eulerian_conv_diff.h:119