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.
helmholtz_element.h
Go to the documentation of this file.
1 // | / |
2 // ' / __| _` | __| _ \ __|
3 // . \ | ( | | ( |\__ `
4 // _|\_\_| \__,_|\__|\___/ ____/
5 // Multi-Physics
6 //
7 // License: BSD License
8 // license: OptimizationApplication/license.txt
9 //
10 // Main authors: Reza Najian Asl
11 // Suneth Warnakulasuriya
12 //
13 
14 #pragma once
15 
16 // System includes
17 
18 // External includes
19 
20 // Project includes
21 #include "includes/define.h"
22 #include "includes/element.h"
24 
25 // Application includes
26 
27 namespace Kratos {
28 
31 
33 
40 template<class TDataContainer>
41 class KRATOS_API(OPTIMIZATION_APPLICATION) HelmholtzElement : public Element
42 {
43 public:
46 
47  using IndexType = std::size_t;
48 
49  using PointType = Node;
50 
51  using PointPtrType = Node::Pointer;
52 
54 
55  static constexpr IndexType NumberOfNodes = TDataContainer::NumberOfNodes;
56 
57  static constexpr IndexType DataDimension = TDataContainer::TargetVariablesList.size();
58 
59  static constexpr IndexType LocalSize = NumberOfNodes * DataDimension;
60 
63 
67 
70  IndexType NewId,
71  GeometryType::Pointer pGeometry);
72 
74  IndexType NewId,
75  GeometryType::Pointer pGeometry,
76  PropertiesType::Pointer pProperties);
77 
79  virtual ~HelmholtzElement() = default;
80 
84 
92  Element::Pointer Create(
93  IndexType NewId,
94  GeometryType::Pointer pGeom,
95  PropertiesType::Pointer pProperties) const override;
96 
104  Element::Pointer Create(
105  IndexType NewId,
106  NodesArrayType const& ThisNodes,
107  PropertiesType::Pointer pProperties) const override;
108 
116  Element::Pointer Clone(
117  IndexType NewId,
118  const NodesArrayType& rThisNodes) const override;
119 
125  void EquationIdVector(
126  EquationIdVectorType& rResult,
127  const ProcessInfo& rCurrentProcessInfo) const override;
128 
134  void GetDofList(
135  DofsVectorType& rElementalDofList,
136  const ProcessInfo& rCurrentProcessInfo) const override;
137 
143  void GetValuesVector(
144  Vector& rValues,
145  int Step = 0) const override;
146 
154  void CalculateLocalSystem(
155  MatrixType& rLeftHandSideMatrix,
156  VectorType& rRightHandSideVector,
157  const ProcessInfo& rCurrentProcessInfo) override;
158 
164  void CalculateLeftHandSide(
165  MatrixType& rLeftHandSideMatrix,
166  const ProcessInfo& rCurrentProcessInfo) override;
167 
173  void CalculateRightHandSide(
174  VectorType& rRightHandSideVector,
175  const ProcessInfo& rCurrentProcessInfo) override;
176 
182  void CalculateMassMatrix(
183  MatrixType& rMassMatrix,
184  const ProcessInfo& rCurrentProcessInfo) override;
185 
191  void Calculate(
192  const Variable<double>& rVariable,
193  double& rOutput,
194  const ProcessInfo& rCurrentProcessInfo) override;
195 
203  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
204 
206 
207 protected:
208  // Protected default constructor necessary for serialization
209  HelmholtzElement() : Element(), mDataContainer(this->GetGeometry())
210  {
211  }
212 
213 private:
216 
217  TDataContainer mDataContainer;
218 
222  friend class Serializer;
223 
224  void save(Serializer& rSerializer) const override
225  {
227  }
228 
229  void load(Serializer& rSerializer) override
230  {
232  }
233 
237 
242  void CalculateStiffnessMatrix(
243  MatrixType& rStiffnessMatrix,
244  const ProcessInfo& rCurrentProcessInfo) const;
245 
247 
248 }; // Class HelmholtzElement
249 
251 
252 } // namespace Kratos.
Base class for all Elements.
Definition: element.h:60
std::vector< DofType::Pointer > DofsVectorType
Definition: element.h:100
std::vector< std::size_t > EquationIdVectorType
Definition: element.h:98
std::size_t IndexType
Definition: flags.h:74
Geometry base class.
Definition: geometry.h:71
SizeType size() const
Definition: geometry.h:518
Short class definition.
Definition: helmholtz_element.h:42
virtual ~HelmholtzElement()=default
Destructor.
HelmholtzElement()
Definition: helmholtz_element.h:209
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(HelmholtzElement)
Counted pointer of HelmholtzElement.
Node::Pointer PointPtrType
Definition: helmholtz_element.h:51
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
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
TDataType Calculate(GeometryType &dummy, const Variable< TDataType > &rVariable)
Definition: add_geometries_to_python.cpp:103
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
ProcessInfo
Definition: edgebased_PureConvection.py:116
def load(f)
Definition: ode_solve.py:307