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.
compressible_perturbation_potential_flow_element.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 authors: Inigo Lopez, Marc Nuñez and Riccardo Rossi
11 //
12 
13 #if !defined(KRATOS_COMPRESSIBLE_PERTURBATION_POTENTIAL_FLOW_ELEMENT_H)
14 #define KRATOS_COMPRESSIBLE_PERTURBATION_POTENTIAL_FLOW_ELEMENT_H
15 
16 // System includes
17 
18 // External includes
19 
20 // Project includes
21 #include "includes/element.h"
22 #include "includes/kratos_flags.h"
23 #include "utilities/geometry_utilities.h"
25 namespace Kratos
26 {
29 
33 
37 
41 
45 
46 template <int Dim, int NumNodes>
48 {
49 public:
50  template <unsigned int TNumNodes, unsigned int TDim>
52  {
54  double vol;
55 
58  };
59 
62 
63  typedef Element BaseType;
64  static constexpr int TNumNodes = NumNodes;
65  static constexpr int TDim = Dim;
66 
71 
75 
76  // Constructors.
77 
79 
83  {
84  }
85 
90  : Element(NewId, ThisNodes)
91  {
92  }
93 
97  CompressiblePerturbationPotentialFlowElement(IndexType NewId, GeometryType::Pointer pGeometry)
98  : Element(NewId, pGeometry)
99  {
100  }
101 
106  GeometryType::Pointer pGeometry,
107  PropertiesType::Pointer pProperties)
108  : Element(NewId, pGeometry, pProperties)
109  {
110  }
111 
116 
121 
126  {
127  }
128 
132 
135 
138 
142 
143  Element::Pointer Create(IndexType NewId,
144  NodesArrayType const& ThisNodes,
145  PropertiesType::Pointer pProperties) const override;
146 
147  Element::Pointer Create(IndexType NewId,
148  GeometryType::Pointer pGeom,
149  PropertiesType::Pointer pProperties) const override;
150 
151  Element::Pointer Clone(IndexType NewId, NodesArrayType const& ThisNodes) const override;
152 
153  void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix,
154  VectorType& rRightHandSideVector,
155  const ProcessInfo& rCurrentProcessInfo) override;
156 
157  void CalculateRightHandSide(VectorType& rRightHandSideVector,
158  const ProcessInfo& rCurrentProcessInfo) override;
159 
160  void CalculateLeftHandSide(MatrixType& rLeftHandSideMatrix,
161  const ProcessInfo& rCurrentProcessInfo) override;
162 
163  void EquationIdVector(EquationIdVectorType& rResult, const ProcessInfo& CurrentProcessInfo) const override;
164 
165  void GetDofList(DofsVectorType& rElementalDofList, const ProcessInfo& rCurrentProcessInfo) const override;
166 
170 
171  void CalculateOnIntegrationPoints(const Variable<double>& rVariable,
172  std::vector<double>& rValues,
173  const ProcessInfo& rCurrentProcessInfo) override;
174 
175  void CalculateOnIntegrationPoints(const Variable<int>& rVariable,
176  std::vector<int>& rValues,
177  const ProcessInfo& rCurrentProcessInfo) override;
178 
180  std::vector<array_1d<double, 3>>& rValues,
181  const ProcessInfo& rCurrentProcessInfo) override;
182 
186 
187  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
188 
192 
194  std::string Info() const override;
195 
197  void PrintInfo(std::ostream& rOStream) const override;
198 
200  void PrintData(std::ostream& rOStream) const override;
201 
203 protected:
204 
205 private:
208 
209  void GetWakeDistances(array_1d<double, NumNodes>& distances) const;
210 
211  void GetEquationIdVectorNormalElement(EquationIdVectorType& rResult) const;
212 
213  void GetEquationIdVectorKuttaElement(EquationIdVectorType& rResult) const;
214 
215  void GetEquationIdVectorWakeElement(EquationIdVectorType& rResult) const;
216 
217  void GetDofListNormalElement(DofsVectorType& rElementalDofList) const;
218 
219  void GetDofListKuttaElement(DofsVectorType& rElementalDofList) const;
220 
221  void GetDofListWakeElement(DofsVectorType& rElementalDofList) const;
222 
223  void CalculateLeftHandSideNormalElement(MatrixType& rLeftHandSideMatrix,
224  const ProcessInfo& rCurrentProcessInfo);
225 
226  void CalculateRightHandSideNormalElement(VectorType& rRightHandSideVector,
227  const ProcessInfo& rCurrentProcessInfo);
228 
229  void CalculateLeftHandSideWakeElement(MatrixType& rLeftHandSideMatrix,
230  const ProcessInfo& rCurrentProcessInfo);
231 
232  BoundedMatrix<double, NumNodes, NumNodes> CalculateLeftHandSideWakeConditions(
233  const ElementalData<NumNodes, Dim>& rData,
234  const ProcessInfo& rCurrentProcessInfo);
235 
236  void CalculateRightHandSideWakeElement(VectorType& rRightHandSideVector,
237  const ProcessInfo& rCurrentProcessInfo);
238 
239  BoundedVector<double, NumNodes> CalculateRightHandSideWakeConditions(
240  const ElementalData<NumNodes, Dim>& rData,
241  const ProcessInfo& rCurrentProcessInfo,
242  const array_1d<double, Dim>& rDiff_velocity);
243 
244  void CalculateLeftHandSideContribution(BoundedMatrix<double, NumNodes, NumNodes>& rLhs_total,
245  const ProcessInfo& rCurrentProcessInfo,
246  const array_1d<double, Dim>& rVelocity,
247  const ElementalData<NumNodes, Dim>& rData);
248 
249  void CalculateRightHandSideContribution(BoundedVector<double, NumNodes>& rRhs_total,
250  const ProcessInfo& rCurrentProcessInfo,
251  const array_1d<double, Dim>& rVelocity,
252  const ElementalData<NumNodes, Dim>& rData);
253 
254  void CalculateLeftHandSideSubdividedElement(Matrix& lhs_positive,
255  Matrix& lhs_negative,
256  const ProcessInfo& rCurrentProcessInfo);
257  void CalculateVolumesSubdividedElement(double& rUpper_vol,
258  double& rLower_vol,
259  const ProcessInfo& rCurrentProcessInfo);
260 
261  void ComputeLHSGaussPointContribution(const double weight,
262  Matrix& lhs,
263  const ElementalData<NumNodes, Dim>& data) const;
264 
265  void AssignLeftHandSideSubdividedElement(Matrix& rLeftHandSideMatrix,
266  Matrix& lhs_positive,
267  Matrix& lhs_negative,
268  const BoundedMatrix<double, NumNodes, NumNodes>& rUpper_lhs_total,
269  const BoundedMatrix<double, NumNodes, NumNodes>& rLower_lhs_total,
270  const BoundedMatrix<double, NumNodes, NumNodes>& rLhs_wake_condition,
271  const ElementalData<NumNodes, Dim>& data) const;
272 
273  void AssignLeftHandSideWakeElement(MatrixType& rLeftHandSideMatrix,
274  const BoundedMatrix<double, NumNodes, NumNodes>& rUpper_lhs_total,
275  const BoundedMatrix<double, NumNodes, NumNodes>& rLower_lhs_total,
276  const BoundedMatrix<double, NumNodes, NumNodes>& rLhs_wake_condition,
277  const ElementalData<NumNodes, Dim>& rData) const;
278 
279  void AssignLeftHandSideWakeNode(MatrixType& rLeftHandSideMatrix,
280  const BoundedMatrix<double, NumNodes, NumNodes>& rUpper_lhs_total,
281  const BoundedMatrix<double, NumNodes, NumNodes>& rLower_lhs_total,
282  const BoundedMatrix<double, NumNodes, NumNodes>& rLhs_wake_condition,
283  const ElementalData<NumNodes, Dim>& rData,
284  unsigned int row) const;
285 
286  void AssignRightHandSideWakeNode(VectorType& rRightHandSideVector,
287  const BoundedVector<double, NumNodes>& rUpper_rhs,
288  const BoundedVector<double, NumNodes>& rLower_rhs,
289  const BoundedVector<double, NumNodes>& rWake_rhs,
290  const ElementalData<NumNodes, Dim>& rData,
291  unsigned int& rRow) const;
292 
296 
300 
301  friend class Serializer;
302 
303  void save(Serializer& rSerializer) const override;
304 
305  void load(Serializer& rSerializer) override;
306 
310 
314 
318 
320 
321 }; // Class CompressiblePerturbationPotentialFlowElement
322 
324 
327 
331 
333 
334 } // namespace Kratos.
335 
336 #endif // KRATOS_COMPRESSIBLE_PERTURBATION_POTENTIAL_FLOW_ELEMENT_H defined
Definition: compressible_perturbation_potential_flow_element.h:48
CompressiblePerturbationPotentialFlowElement(IndexType NewId, GeometryType::Pointer pGeometry)
Definition: compressible_perturbation_potential_flow_element.h:97
~CompressiblePerturbationPotentialFlowElement() override
Definition: compressible_perturbation_potential_flow_element.h:125
CompressiblePerturbationPotentialFlowElement(CompressiblePerturbationPotentialFlowElement &&rOther)=delete
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(CompressiblePerturbationPotentialFlowElement)
CompressiblePerturbationPotentialFlowElement(IndexType NewId=0)
Default constuctor.
Definition: compressible_perturbation_potential_flow_element.h:82
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Definition: compressible_perturbation_potential_flow_element.cpp:148
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: compressible_perturbation_potential_flow_element.cpp:266
CompressiblePerturbationPotentialFlowElement(CompressiblePerturbationPotentialFlowElement const &rOther)=delete
Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override
It creates a new element pointer and clones the previous element data.
Definition: compressible_perturbation_potential_flow_element.cpp:45
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &CurrentProcessInfo) const override
Definition: compressible_perturbation_potential_flow_element.cpp:89
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: compressible_perturbation_potential_flow_element.cpp:272
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, PropertiesType::Pointer pProperties) const override
It creates a new element pointer.
Definition: compressible_perturbation_potential_flow_element.cpp:26
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: compressible_perturbation_potential_flow_element.cpp:63
CompressiblePerturbationPotentialFlowElement(IndexType NewId, const NodesArrayType &ThisNodes)
Definition: compressible_perturbation_potential_flow_element.h:89
Element BaseType
Definition: compressible_perturbation_potential_flow_element.h:63
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: compressible_perturbation_potential_flow_element.cpp:55
CompressiblePerturbationPotentialFlowElement & operator=(CompressiblePerturbationPotentialFlowElement &&rOther)=delete
Move operator.
void GetDofList(DofsVectorType &rElementalDofList, const ProcessInfo &rCurrentProcessInfo) const override
Definition: compressible_perturbation_potential_flow_element.cpp:117
CompressiblePerturbationPotentialFlowElement & operator=(CompressiblePerturbationPotentialFlowElement const &rOther)=delete
Assignment operator.
std::string Info() const override
Turn back information as a string.
Definition: compressible_perturbation_potential_flow_element.cpp:258
CompressiblePerturbationPotentialFlowElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Definition: compressible_perturbation_potential_flow_element.h:105
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: compressible_perturbation_potential_flow_element.cpp:76
static constexpr int TDim
Definition: compressible_perturbation_potential_flow_element.h:65
static constexpr int TNumNodes
Definition: compressible_perturbation_potential_flow_element.h:64
void CalculateOnIntegrationPoints(const Variable< double > &rVariable, std::vector< double > &rValues, const ProcessInfo &rCurrentProcessInfo) override
Definition: compressible_perturbation_potential_flow_element.cpp:175
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
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
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
AMatrix::MatrixRow< const TExpressionType > row(AMatrix::MatrixExpression< TExpressionType, TCategory > const &TheExpression, std::size_t RowIndex)
Definition: amatrix_interface.h:649
lhs
Definition: generate_frictional_mortar_condition.py:297
data
Definition: mesh_to_mdpa_converter.py:59
def load(f)
Definition: ode_solve.py:307
Definition: compressible_perturbation_potential_flow_element.h:52
double vol
Definition: compressible_perturbation_potential_flow_element.h:54
BoundedMatrix< double, TNumNodes, TDim > DN_DX
Definition: compressible_perturbation_potential_flow_element.h:56
array_1d< double, TNumNodes > distances
Definition: compressible_perturbation_potential_flow_element.h:53
array_1d< double, TNumNodes > N
Definition: compressible_perturbation_potential_flow_element.h:57
array_1d< double, TNumNodes > potentials
Definition: compressible_perturbation_potential_flow_element.h:53