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_condition.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 #if !defined(KRATOS_HELMHOLTZ_SURF_SHAPE_CONDITION_3D_H_INCLUDED )
13 #define KRATOS_HELMHOLTZ_SURF_SHAPE_CONDITION_3D_H_INCLUDED
14 
15 // System includes
16 
17 // External includes
18 
19 // Project includes
20 #include "includes/condition.h"
21 #include "includes/model_part.h"
23 #include "utilities/geometry_utilities.h"
25 
26 namespace Kratos
27 {
28 
31 
35 
39 
43 
47 
54 class KRATOS_API(OPTIMIZATION_APPLICATION) HelmholtzSurfShapeCondition
55  : public Condition
56 {
57 public:
58 
61 
64 
67 
70 
73 
76 
79 
82 
83  // Counted pointer of HelmholtzSurfShapeCondition
85 
89 
90  // Constructor void
92  {};
93 
94  // Constructor using an array of nodes
95  HelmholtzSurfShapeCondition( IndexType NewId, GeometryType::Pointer pGeometry ):Condition(NewId,pGeometry)
96  {};
97 
98  // Constructor using an array of nodes with properties
99  HelmholtzSurfShapeCondition( IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties ):Condition(NewId,pGeometry,pProperties)
100  {};
101 
104 
105  // Destructor
107  {};
108 
112 
115 
119 
127  Condition::Pointer Create(
128  IndexType NewId,
129  NodesArrayType const& ThisNodes,
130  PropertiesType::Pointer pProperties
131  ) const override;
132 
140  Condition::Pointer Create(
141  IndexType NewId,
142  GeometryType::Pointer pGeom,
143  PropertiesType::Pointer pProperties
144  ) const override;
145 
152  Condition::Pointer Clone (
153  IndexType NewId,
154  NodesArrayType const& ThisNodes
155  ) const override;
156 
162  void EquationIdVector(
163  EquationIdVectorType& rResult,
164  const ProcessInfo& rCurrentProcessInfo
165  ) const override;
166 
172  void GetDofList(
173  DofsVectorType& ElementalDofList,
174  const ProcessInfo& rCurrentProcessInfo
175  ) const override;
176 
182  void GetValuesVector(
183  Vector& rValues,
184  int Step = 0
185  ) const override;
186 
187  void Calculate(const Variable<Matrix>& rVariable,
188  Matrix& rOutput, const ProcessInfo& rCurrentProcessInfo) override;
189 
190  void Calculate(const Variable<double>& rVariable,
191  double& rOutput, const ProcessInfo& rCurrentProcessInfo) override;
192 
201  void CalculateLocalSystem(
202  MatrixType& rLeftHandSideMatrix,
203  VectorType& rRightHandSideVector,
204  const ProcessInfo& rCurrentProcessInfo
205  ) override;
206 
212  void CalculateRightHandSide(
213  VectorType& rRightHandSideVector,
214  const ProcessInfo& rCurrentProcessInfo
215  ) override;
216 
222  int Check( const ProcessInfo& rCurrentProcessInfo ) const override;
223 
224 
228 
232 
236 
238  std::string Info() const override
239  {
240  std::stringstream buffer;
241  buffer << "Base surface filter Condition #" << Id();
242  return buffer.str();
243  }
244 
246 
247  void PrintInfo(std::ostream& rOStream) const override
248  {
249  rOStream << "Base surface filter Condition #" << Id();
250  }
251 
253  void PrintData(std::ostream& rOStream) const override
254  {
255  pGetGeometry()->PrintData(rOStream);
256  }
257 
261 
262 protected:
263 
266 
270 
274 
278 
282 
286 
290 
291 private:
294 
295 
299 
300 
301 
305 
309 
310  void CalculateNormal(VectorType & r_n) const;
311  void GetParentElementShapeFunctionsValues(MatrixType& rNMatrix,const IntegrationMethod& rIntegrationMethod, const ProcessInfo& rCurrentProcessInfo) const;
312  void GetParentElementShapeFunctionsGlobalGradients(MatrixType& rDN_DX,const IndexType PointNumber,const IntegrationMethod& rIntegrationMethod, const ProcessInfo& rCurrentProcessInfo) const;
313  void CalculateSurfaceMassMatrix(MatrixType& rMassMatrix,const ProcessInfo& rCurrentProcessInfo) const;
314  void CalculateSurfaceStiffnessMatrix(MatrixType& rStiffnessMatrix,const ProcessInfo& rCurrentProcessInfo) const;
315 
319 
320 
324 
328 
329  friend class Serializer;
330 
331  void save( Serializer& rSerializer ) const override;
332 
333  void load( Serializer& rSerializer ) override;
334 
335 }; // class HelmholtzSurfShapeCondition.
336 
340 
341 
345 
346 } // namespace Kratos.
347 
348 #endif // KRATOS_HELMHOLTZ_SURF_SHAPE_CONDITION_3D_H_INCLUDED defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Base class for all Conditions.
Definition: condition.h:59
std::size_t IndexType
Definition: flags.h:74
std::size_t IndexType
Defines the index type.
Definition: geometrical_object.h:73
IntegrationMethod
Definition: geometry_data.h:76
Geometry base class.
Definition: geometry.h:71
This is the class of surface PDE/Helmholtz-based surface filtering.
Definition: helmholtz_surf_shape_condition.h:56
BaseType::GeometryType GeometryType
Definition of the geometry type with given NodeType.
Definition: helmholtz_surf_shape_condition.h:78
HelmholtzSurfShapeCondition(IndexType NewId, GeometryType::Pointer pGeometry)
Definition: helmholtz_surf_shape_condition.h:95
BaseType::NodesArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: helmholtz_surf_shape_condition.h:81
~HelmholtzSurfShapeCondition() override
Definition: helmholtz_surf_shape_condition.h:106
Condition BaseType
We define the base class Condition.
Definition: helmholtz_surf_shape_condition.h:63
std::string Info() const override
Turn back information as a string.
Definition: helmholtz_surf_shape_condition.h:238
BaseType::SizeType SizeType
Definition of the size type.
Definition: helmholtz_surf_shape_condition.h:69
BaseType::IndexType IndexType
Dfinition of the index type.
Definition: helmholtz_surf_shape_condition.h:66
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: helmholtz_surf_shape_condition.h:253
BaseType::NodeType NodeType
Definition of the node type.
Definition: helmholtz_surf_shape_condition.h:72
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(HelmholtzSurfShapeCondition)
HelmholtzSurfShapeCondition(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Definition: helmholtz_surf_shape_condition.h:99
HelmholtzSurfShapeCondition()
Definition: helmholtz_surf_shape_condition.h:91
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: helmholtz_surf_shape_condition.h:247
BaseType::PropertiesType PropertiesType
Definition of the properties type.
Definition: helmholtz_surf_shape_condition.h:75
This class defines the node.
Definition: node.h:65
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
Variable class contains all information needed to store and retrive data from a data container.
Definition: variable.h:63
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
std::size_t SizeType
Definition: nurbs_utilities.h:41
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
Properties PropertiesType
Definition: regenerate_pfem_pressure_conditions_process.h:26
ModelPart::NodesContainerType NodesArrayType
Definition: gid_gauss_point_container.h:42
def load(f)
Definition: ode_solve.py:307