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.
d_convection_diffusion_explicit.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: Riccardo Tosi
11 //
12 
13 #ifndef KRATOS_D_CONVECTION_DIFFUSION_EXPLICIT_H
14 #define KRATOS_D_CONVECTION_DIFFUSION_EXPLICIT_H
15 
16 // System includes
17 
18 
19 // External includes
20 
21 
22 // Project includes
24 
25 namespace Kratos
26 {
27 
30 
34 
38 
42 
46 
61 template< unsigned int TDim, unsigned int TNumNodes>
63 {
64 public:
67 
69  typedef typename BaseType::ElementData ElementData;
70  typedef Node NodeType;
73  typedef Vector VectorType;
74  typedef Matrix MatrixType;
75  typedef std::size_t IndexType;
76  typedef std::vector<std::size_t> EquationIdVectorType;
77  typedef std::vector< Dof<double>::Pointer > DofsVectorType;
79 
82 
86 
87  //Constructors.
89  IndexType NewId,
90  GeometryType::Pointer pGeometry);
91 
93  IndexType NewId,
94  GeometryType::Pointer pGeometry,
95  Properties::Pointer pProperties);
96 
98 
101 
105 
106 
110 
111  Element::Pointer Create(
112  IndexType NewId,
113  NodesArrayType const& ThisNodes,
114  Properties::Pointer pProperties) const override;
115 
116  Element::Pointer Create(
117  IndexType NewId,
118  GeometryType::Pointer pGeom,
119  Properties::Pointer pProperties) const override;
120 
122  MatrixType& rLeftHandSideMatrix,
123  VectorType& rRightHandSideVector,
124  const ProcessInfo& rCurrentProcessInfo) override;
125 
127  VectorType& rRightHandSideVector,
128  const ProcessInfo& rCurrentProcessInfo) override;
129 
130  void AddExplicitContribution(const ProcessInfo &rCurrentProcessInfo) override;
131 
132  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
133 
134  void FinalizeSolutionStep(const ProcessInfo& rCurrentProcessInfo) override;
135 
136  void Calculate(
137  const Variable<double>& rVariable,
138  double& Output,
139  const ProcessInfo& rCurrentProcessInfo) override;
140 
144 
148 
150  std::string Info() const override
151  {
152  return "DConvectionDiffusionExplicitElement #";
153  }
154 
156  void PrintInfo(std::ostream& rOStream) const override
157  {
158  rOStream << Info() << this->Id();
159  }
160 
162 
163 protected:
164 
167 
171 
175 
179 
180 
184 
188 
189  // Protected default constructor necessary for serialization
191  {
192  }
193 
195 
196 private:
197 
200 
201  BoundedVector<double, TNumNodes> mUnknownSubScale;
202 
206 
210 
211  friend class Serializer;
212 
213  void save(Serializer& rSerializer) const override
214  {
216  }
217 
218  void load(Serializer& rSerializer) override
219  {
221  }
222 
226 
231 
235 
236  void DCalculateRightHandSideInternal(
237  BoundedVector<double, TNumNodes>& rRightHandSideBoundedVector,
238  const ProcessInfo& rCurrentProcessInfo);
239 
240  void DCalculateOrthogonalSubgridScaleRHSInternal(
241  BoundedVector<double, TNumNodes>& rRightHandSideVector,
242  const ProcessInfo& rCurrentProcessInfo);
243 
244  void UpdateUnknownSubgridScaleGaussPoint(
245  ElementData& rData,
246  unsigned int g);
247 
248  void DCalculateTau(ElementData& rData);
249 
253 
254 
258 
259 
263 
265 
266 }; // Class DConvectionDiffusionExplicit
267 
269 
272 
273 
277 
278 
280 template< unsigned int TDim, unsigned int TNumNodes = TDim + 1>
281 inline std::istream& operator >>(
282  std::istream& rIStream,
284 {
285  return rIStream;
286 }
287 
289 template< unsigned int TDim, unsigned int TNumNodes = TDim + 1>
290 inline std::ostream& operator <<(
291  std::ostream& rOStream,
293 {
294  rThis.PrintInfo(rOStream);
295  rOStream << std::endl;
296  rThis.PrintData(rOStream);
297 
298  return rOStream;
299 }
301 
303 
304 } // namespace Kratos.
305 
306 #endif // KRATOS_D_CONVECTION_DIFFUSION_EXPLICIT_H
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
This element solves the convection-diffusion equation, stabilized with algebraic subgrid scale or ort...
Definition: d_convection_diffusion_explicit.h:63
Vector VectorType
Definition: d_convection_diffusion_explicit.h:73
virtual ~DConvectionDiffusionExplicit()
Default constuctor.
Definition: d_convection_diffusion_explicit.cpp:39
Node NodeType
Definition: d_convection_diffusion_explicit.h:70
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition: d_convection_diffusion_explicit.h:72
BaseType::ElementData ElementData
Definition: d_convection_diffusion_explicit.h:69
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: d_convection_diffusion_explicit.h:156
Geometry< NodeType > GeometryType
Definition: d_convection_diffusion_explicit.h:71
DConvectionDiffusionExplicit()
Definition: d_convection_diffusion_explicit.h:190
void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Definition: d_convection_diffusion_explicit.cpp:120
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: d_convection_diffusion_explicit.h:77
std::vector< std::size_t > EquationIdVectorType
Definition: d_convection_diffusion_explicit.h:76
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: d_convection_diffusion_explicit.cpp:68
std::size_t IndexType
Definition: d_convection_diffusion_explicit.h:75
GeometryData::IntegrationMethod IntegrationMethod
Definition: d_convection_diffusion_explicit.h:78
QSConvectionDiffusionExplicit< TDim, TNumNodes > BaseType
Definition: d_convection_diffusion_explicit.h:68
std::string Info() const override
Turn back information as a string.
Definition: d_convection_diffusion_explicit.h:150
void FinalizeSolutionStep(const ProcessInfo &rCurrentProcessInfo) override
Definition: d_convection_diffusion_explicit.cpp:135
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, Properties::Pointer pProperties) const override
Definition: d_convection_diffusion_explicit.cpp:45
Matrix MatrixType
Definition: d_convection_diffusion_explicit.h:74
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: d_convection_diffusion_explicit.cpp:82
void AddExplicitContribution(const ProcessInfo &rCurrentProcessInfo) override
Definition: d_convection_diffusion_explicit.cpp:95
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(DConvectionDiffusionExplicit)
Pointer definition of DConvectionDiffusionExplicit.
void Calculate(const Variable< double > &rVariable, double &Output, const ProcessInfo &rCurrentProcessInfo) override
Definition: d_convection_diffusion_explicit.cpp:167
Base class for all Elements.
Definition: element.h:60
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: element.h:1135
std::size_t IndexType
Definition: flags.h:74
IntegrationMethod
Definition: geometry_data.h:76
Geometry base class.
Definition: geometry.h:71
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
This element solves the convection-diffusion equation, stabilized with algebraic subgrid scale or ort...
Definition: qs_convection_diffusion_explicit.h:71
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
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
ProcessInfo
Definition: edgebased_PureConvection.py:116
def load(f)
Definition: ode_solve.py:307