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_surf_shape_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, https://github.com/RezaNajian
11 //
12 
13 #if !defined(KRATOS_HELMHOLTZ_SURF_SHAPE_ELEMENT_H_INCLUDED )
14 #define KRATOS_HELMHOLTZ_SURF_SHAPE_ELEMENT_H_INCLUDED
15 
16 // System includes
17 
18 // External includes
19 
20 // Project includes
21 #include "includes/define.h"
22 #include "includes/element.h"
24 #include "utilities/geometry_utilities.h"
28 #include "includes/variables.h"
30 
31 
32 namespace Kratos
33 {
34 
37 
41 
45 
49 
53 
55 
57 class KRATOS_API(OPTIMIZATION_APPLICATION) HelmholtzSurfShapeElement
58  : public Element
59 {
60 public:
63 
64  typedef Node PointType;
65  typedef Node::Pointer PointPtrType;
69 
72 
76 
78  HelmholtzSurfShapeElement(IndexType NewId, GeometryType::Pointer pGeometry);
79  HelmholtzSurfShapeElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties);
80 
82  virtual ~HelmholtzSurfShapeElement();
83 
84 
88 
89 
93 
94  Element::Pointer Create(IndexType NewId, NodesArrayType const& ThisNodes, PropertiesType::Pointer pProperties) const override;
95 
96  Element::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, PropertiesType::Pointer pProperties) const override;
97 
98  void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix, VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
99 
100  void CalculateLeftHandSide(MatrixType& rLeftHandSideMatrix, const ProcessInfo& rCurrentProcessInfo) override;
101 
102  void CalculateRightHandSide(VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
103 
104  void EquationIdVector(EquationIdVectorType& rResult, const ProcessInfo& rCurrentProcessInfo) const override;
105 
106  void GetDofList(DofsVectorType& ElementalDofList, const ProcessInfo& CurrentProcessInfo) const override;
107 
108  void GetValuesVector(VectorType &rValues, int Step = 0) const override;
109 
110  void Calculate(const Variable<Matrix>& rVariable, Matrix& rOutput, const ProcessInfo& rCurrentProcessInfo) override;
111 
112  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
113 
117 
118 
122 
126 
130 
131 
133 
134 protected:
137 
138 
142 
143 
147 
148 
152 
153 
157 
158 
162 
163 
167 
168  // Protected default constructor necessary for serialization
170  {
171  }
172 
174 
175 private:
178 
179 
183 
184 
188  friend class Serializer;
189 
190  void save(Serializer& rSerializer) const override
191  {
193  }
194 
195  void load(Serializer& rSerializer) override
196  {
198  }
199 
203 
204 
208 
209  void CalculateSurfaceMassMatrix(MatrixType& rMassMatrix,const ProcessInfo& rCurrentProcessInfo) const;
210  void CalculateSurfaceStiffnessMatrix(MatrixType& rStiffnessMatrix,const ProcessInfo& rCurrentProcessInfo) const;
211  void GetPseudoBulkSurfaceShapeFunctionsValues(MatrixType& rNMatrix,const IntegrationMethod& rIntegrationMethod, const ProcessInfo& rCurrentProcessInfo) const;
212  void CalculatePseudoBulkSurfaceDN_DXMatrix(MatrixType& rDN_DX, const IntegrationMethod& rIntegrationMethod, const IndexType PointNumber, const ProcessInfo& rCurrentProcessInfo) const;
213  void CalculateAvgSurfUnitNormal(VectorType & rNormal) const;
214  void CalculateCMatrix(MatrixType& rCMatrix, const IntegrationMethod& rIntegrationMethod, const IndexType PointNumber) const;
215  void CalculateBMatrix(MatrixType& rBMatrix, const MatrixType& rDN_DX_tMatrix, const IntegrationMethod& rIntegrationMethod, const IndexType PointNumber) const;
216 
220 
221 
225 
226 
230 
232  //HelmholtzSurfShapeElement& operator=(const HelmholtzSurfShapeElement& rOther);
233 
235  //HelmholtzSurfShapeElement(const HelmholtzSurfShapeElement& rOther);
236 
237 
239 
240 }; // Class HelmholtzSurfShapeElement
241 
243 
246 
247 
251 
252 
254 /* inline std::istream& operator >> (std::istream& rIStream,
255  HelmholtzSurfShapeElement& rThis);
256 */
258 /* inline std::ostream& operator << (std::ostream& rOStream,
259  const HelmholtzSurfShapeElement& rThis)
260  {
261  rThis.PrintInfo(rOStream);
262  rOStream << std::endl;
263  rThis.PrintData(rOStream);
264 
265  return rOStream;
266  }*/
268 
269 } // namespace Kratos.
270 
271 #endif // KRATOS_HELMHOLTZ_SURF_SHAPE_ELEMENT_H_INCLUDED defined
272 
273 
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
Short class definition.
Definition: helmholtz_surf_shape_element.h:59
Geometry< PointType > GeometryType
Definition: helmholtz_surf_shape_element.h:66
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(HelmholtzSurfShapeElement)
Counted pointer of HelmholtzSurfShapeElement.
HelmholtzSurfShapeElement()
Definition: helmholtz_surf_shape_element.h:169
Tetrahedra3D4< PointType > TetrahedraGeometryType
Definition: helmholtz_surf_shape_element.h:68
Node::Pointer PointPtrType
Definition: helmholtz_surf_shape_element.h:65
Node PointType
Definition: helmholtz_surf_shape_element.h:64
Pyramid3D5< PointType > PyramidGeometryType
Definition: helmholtz_surf_shape_element.h:67
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
A five node pyramid geometry with linear shape functions.
Definition: pyramid_3d_5.h:59
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
A four node tetrahedra geometry with linear shape functions.
Definition: tetrahedra_3d_4.h:59
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
std::size_t IndexType
The definition of the index type.
Definition: key_hash.h:35
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
GeometryData::IntegrationMethod IntegrationMethod
Type definition for integration methods.
Definition: exact_mortar_segmentation_utility.h:57
array_1d< double, 3 > VectorType
Definition: solid_mechanics_application_variables.cpp:19
ProcessInfo
Definition: edgebased_PureConvection.py:116
def load(f)
Definition: ode_solve.py:307