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.
transonic_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: Eloisa Baez Jones, Inigo Lopez, Marc Nuñez and Riccardo Rossi
11 //
12 
13 #if !defined(KRATOS_TRANSONIC_PERTURBATION_POTENTIAL_FLOW_ELEMENT_H)
14 #define KRATOS_TRANSONIC_PERTURBATION_POTENTIAL_FLOW_ELEMENT_H
15 
16 // System includes
17 
18 // External includes
19 
20 // Project includes
21 #include "includes/element.h"
23 
24 namespace Kratos
25 {
28 
32 
36 
40 
44 
45 template <int TDim, int TNumNodes>
47 {
48 public:
51 
52  typedef Element BaseType;
53 
56 
61 
65 
66  // Constructors.
67 
69 
73  {
74  }
75 
80  const NodesArrayType& ThisNodes)
81  : Element(NewId, ThisNodes)
82  {
83  }
84 
89  GeometryType::Pointer pGeometry)
90  : Element(NewId, pGeometry)
91  {
92  }
93 
98  GeometryType::Pointer pGeometry,
99  PropertiesType::Pointer pProperties)
100  : Element(NewId, pGeometry, pProperties)
101  {
102  }
103 
108 
113 
118  {
119  }
120 
124 
127 
130 
134 
135  Element::Pointer Create(IndexType NewId,
136  NodesArrayType const& ThisNodes,
137  PropertiesType::Pointer pProperties) const override;
138 
139  Element::Pointer Create(IndexType NewId,
140  GeometryType::Pointer pGeom,
141  PropertiesType::Pointer pProperties) const override;
142 
143  Element::Pointer Clone(IndexType NewId,
144  NodesArrayType const& ThisNodes) const override;
145 
146  void Initialize(const ProcessInfo& rCurrentProcessInfo) override;
147 
148  void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix,
149  VectorType& rRightHandSideVector,
150  const ProcessInfo& rCurrentProcessInfo) override;
151 
152  void CalculateRightHandSide(VectorType& rRightHandSideVector,
153  const ProcessInfo& rCurrentProcessInfo) override;
154 
155  void CalculateLeftHandSide(MatrixType& rLeftHandSideMatrix,
156  const ProcessInfo& rCurrentProcessInfo) override;
157 
159  const ProcessInfo& CurrentProcessInfo) const override;
160 
161  void GetDofList(DofsVectorType& rElementalDofList,
162  const ProcessInfo& rCurrentProcessInfo) const override;
163 
167 
168  void CalculateOnIntegrationPoints(const Variable<double>& rVariable,
169  std::vector<double>& rValues,
170  const ProcessInfo& rCurrentProcessInfo) override;
171 
172  void CalculateOnIntegrationPoints(const Variable<int>& rVariable,
173  std::vector<int>& rValues,
174  const ProcessInfo& rCurrentProcessInfo) override;
175 
177  std::vector<array_1d<double, 3>>& rValues,
178  const ProcessInfo& rCurrentProcessInfo) override;
179 
183 
184  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
185 
189 
191  std::string Info() const override;
192 
194  void PrintInfo(std::ostream& rOStream) const override;
195 
197  void PrintData(std::ostream& rOStream) const override;
198 
200 protected:
201 
202  BoundedVector<double, TNumNodes + 1> AssembleDensityDerivativeAndShapeFunctions(const double densityDerivativeWRTVelocitySquared, const double densityDerivativeWRTUpwindVelocitySquared, const array_1d<double, TDim>& velocity, const array_1d<double, TDim>& upwindVelocity,const ProcessInfo& rCurrentProcessInfo);
203 
204  void CalculateRightHandSideNormalElement(VectorType& rRightHandSideVector,
205  const ProcessInfo& rCurrentProcessInfo);
206 
207  void CalculateLeftHandSideNormalElement(MatrixType& rLeftHandSideMatrix,
208  const ProcessInfo& rCurrentProcessInfo);
209 
210  void CalculateLeftHandSideWakeElement(MatrixType& rLeftHandSideMatrix,
211  const ProcessInfo& rCurrentProcessInfo);
212 
213  virtual void CalculateRightHandSideWakeElement(VectorType& rRightHandSideVector,
214  const ProcessInfo& rCurrentProcessInfo);
215 
216 
217  virtual void AssembleSupersonicLeftHandSide(MatrixType& rLeftHandSideMatrix,
218  const double densityDerivativeWRTVelocity,
219  const double densityDerivativeWRTUpwindVelocity,
221  const array_1d<double, TDim> upwindVelocity,
222  const ProcessInfo& rCurrentProcessInfo);
223 
225  const ProcessInfo& rCurrentProcessInfo,
226  const array_1d<double, TDim>& rVelocity,
227  const ElementalData& rData);
228 
230  const double rDensity,
231  const array_1d<double, TDim>& rVelocity);
232 
234 
235  bool CheckUpwindElement();
236 
237  void pSetUpwindElement(GlobalPointer<Element> pUpwindElement);
238 
240 
241  virtual void FindUpwindElement(const ProcessInfo& rCurrentProcessInfo);
242 
243 private:
247 
248  GlobalPointer<Element> mpUpwindElement;
249 
252 
253  void GetWakeDistances(array_1d<double,
254  TNumNodes>& distances) const;
255 
256  void GetEquationIdVectorExtendedElement(EquationIdVectorType& rResult) const;
257 
258  void AddUpwindEquationId(EquationIdVectorType& rResult) const;
259 
260  void GetEquationIdVectorNormalElement(EquationIdVectorType& rResult) const;
261 
262  void GetEquationIdVectorKuttaElement(EquationIdVectorType& rResult) const;
263 
264  void GetEquationIdVectorWakeElement(EquationIdVectorType& rResult) const;
265 
266  void GetDofListNormalElement(DofsVectorType& rElementalDofList) const;
267 
268  void GetDofListKuttaElement(DofsVectorType& rElementalDofList) const;
269 
270  void GetDofListWakeElement(DofsVectorType& rElementalDofList) const;
271 
272  void CalculateLeftHandSideSubsonicElement(MatrixType& rLeftHandSideMatrix,
273  const ProcessInfo& rCurrentProcessInfo);
274 
275  // void CalculateRightHandSideSupersonicElement(VectorType& rRightHandSideVector,
276  // const ProcessInfo& rCurrentProcessInfo);
277 
278  BoundedMatrix<double, TNumNodes, TNumNodes> CalculateLeftHandSideWakeConditions(
279  const ElementalData& rData,
280  const ProcessInfo& rCurrentProcessInfo);
281 
282  BoundedVector<double, TNumNodes> CalculateRightHandSideWakeConditions(
283  const ElementalData& rData,
284  const ProcessInfo& rCurrentProcessInfo,
285  const array_1d<double, TDim>& rDiff_velocity);
286 
287  void CalculateLeftHandSideSubdividedElement(Matrix& lhs_positive,
288  Matrix& lhs_negative,
289  const ProcessInfo& rCurrentProcessInfo);
290  void CalculateVolumesSubdividedElement(double& rUpper_vol,
291  double& rLower_vol,
292  const ProcessInfo& rCurrentProcessInfo);
293 
294  void ComputeLHSGaussPointContribution(const double weight,
295  Matrix& lhs,
296  const ElementalData& data) const;
297 
298  void AssignLeftHandSideSubdividedElement(
299  Matrix& rLeftHandSideMatrix,
300  Matrix& lhs_positive,
301  Matrix& lhs_negative,
302  const BoundedMatrix<double, TNumNodes, TNumNodes>& rUpper_lhs_total,
303  const BoundedMatrix<double, TNumNodes, TNumNodes>& rLower_lhs_total,
304  const BoundedMatrix<double, TNumNodes, TNumNodes>& rLhs_wake_condition,
305  const ElementalData& data) const;
306 
307  void AssignLeftHandSideWakeElement(MatrixType& rLeftHandSideMatrix,
308  const BoundedMatrix<double, TNumNodes, TNumNodes>& rUpper_lhs_total,
309  const BoundedMatrix<double, TNumNodes, TNumNodes>& rLower_lhs_total,
310  const BoundedMatrix<double, TNumNodes, TNumNodes>& rLhs_wake_condition,
311  const ElementalData& rData) const;
312 
313  void AssignLeftHandSideWakeNode(MatrixType& rLeftHandSideMatrix,
314  const BoundedMatrix<double, TNumNodes, TNumNodes>& rUpper_lhs_total,
315  const BoundedMatrix<double, TNumNodes, TNumNodes>& rLower_lhs_total,
316  const BoundedMatrix<double, TNumNodes, TNumNodes>& rLhs_wake_condition,
317  const ElementalData& rData,
318  unsigned int row) const;
319 
320  void AssignRightHandSideWakeNode(VectorType& rRightHandSideVector,
321  const BoundedVector<double, TNumNodes>& rUpper_rhs,
322  const BoundedVector<double, TNumNodes>& rLower_rhs,
323  const BoundedVector<double, TNumNodes>& rWake_rhs,
324  const ElementalData& rData,
325  unsigned int& rRow) const;
326 
327  array_1d<size_t, TNumNodes> GetAssemblyKey(const GeometryType& rGeom, const GeometryType& rUpwindGeom, const ProcessInfo& rCurrentProcessInfo);
328 
329  void FindUpwindEdge(GeometryType& rUpwindEdge,
330  const ProcessInfo& rCurrentProcessInfo);
331 
332  void GetElementGeometryBoundary(GeometriesArrayType& rElementGeometryBoundary);
333 
334 
335  void SelectUpwindElement(std::vector<IndexType>& rUpwindElementNodesIds,
336  GlobalPointersVector<Element>& rUpwindElementCandidates);
337 
338  int GetAdditionalUpwindNodeIndex() const;
339 
343 
347 
348  friend class Serializer;
349 
350  void save(Serializer& rSerializer) const override;
351 
352  void load(Serializer& rSerializer) override;
353 
357 
361 
365 
367 
368 }; // Class TransonicPerturbationPotentialFlowElement
369 
371 
374 
378 
380 
381 } // namespace Kratos.
382 
383 #endif // KRATOS_TRANSONIC_PERTURBATION_POTENTIAL_FLOW_ELEMENT_H defined
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
This class is a wrapper for a pointer to a data that is located in a different rank.
Definition: global_pointer.h:44
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
Definition: transonic_perturbation_potential_flow_element.h:47
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Definition: transonic_perturbation_potential_flow_element.cpp:198
Element::Pointer Clone(IndexType NewId, NodesArrayType const &ThisNodes) const override
It creates a new element pointer and clones the previous element data.
Definition: transonic_perturbation_potential_flow_element.cpp:49
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: transonic_perturbation_potential_flow_element.cpp:351
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: transonic_perturbation_potential_flow_element.cpp:345
BoundedVector< double, TNumNodes+1 > AssembleDensityDerivativeAndShapeFunctions(const double densityDerivativeWRTVelocitySquared, const double densityDerivativeWRTUpwindVelocitySquared, const array_1d< double, TDim > &velocity, const array_1d< double, TDim > &upwindVelocity, const ProcessInfo &rCurrentProcessInfo)
Definition: transonic_perturbation_potential_flow_element.cpp:1163
TransonicPerturbationPotentialFlowElement(TransonicPerturbationPotentialFlowElement const &rOther)=delete
Element BaseType
Definition: transonic_perturbation_potential_flow_element.h:52
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: transonic_perturbation_potential_flow_element.cpp:66
~TransonicPerturbationPotentialFlowElement() override
Definition: transonic_perturbation_potential_flow_element.h:117
virtual void CalculateLeftHandSideContribution(BoundedMatrix< double, TNumNodes, TNumNodes > &rLhs_total, const ProcessInfo &rCurrentProcessInfo, const array_1d< double, TDim > &rVelocity, const ElementalData &rData)
Definition: transonic_perturbation_potential_flow_element.cpp:870
virtual void CalculateRightHandSideWakeElement(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: transonic_perturbation_potential_flow_element.cpp:765
TransonicPerturbationPotentialFlowElement(IndexType NewId, const NodesArrayType &ThisNodes)
Definition: transonic_perturbation_potential_flow_element.h:79
TransonicPerturbationPotentialFlowElement(TransonicPerturbationPotentialFlowElement &&rOther)=delete
TransonicPerturbationPotentialFlowElement & operator=(TransonicPerturbationPotentialFlowElement const &rOther)=delete
Assignment operator.
virtual void AssembleSupersonicLeftHandSide(MatrixType &rLeftHandSideMatrix, const double densityDerivativeWRTVelocity, const double densityDerivativeWRTUpwindVelocity, const array_1d< double, TDim > velocity, const array_1d< double, TDim > upwindVelocity, const ProcessInfo &rCurrentProcessInfo)
Definition: transonic_perturbation_potential_flow_element.cpp:1126
TransonicPerturbationPotentialFlowElement(IndexType NewId, GeometryType::Pointer pGeometry)
Definition: transonic_perturbation_potential_flow_element.h:88
virtual void FindUpwindElement(const ProcessInfo &rCurrentProcessInfo)
Definition: transonic_perturbation_potential_flow_element.cpp:1218
PointerVector< GeometryType > GeometriesArrayType
Definition: transonic_perturbation_potential_flow_element.h:54
void Initialize(const ProcessInfo &rCurrentProcessInfo) override
Definition: transonic_perturbation_potential_flow_element.cpp:60
array_1d< double, 3 > GetEdgeNormal(const GeometryType &rEdge)
Definition: transonic_perturbation_potential_flow_element.cpp:1278
void GetDofList(DofsVectorType &rElementalDofList, const ProcessInfo &rCurrentProcessInfo) const override
Definition: transonic_perturbation_potential_flow_element.cpp:158
TransonicPerturbationPotentialFlowElement(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Definition: transonic_perturbation_potential_flow_element.h:97
GlobalPointer< Element > pGetUpwindElement() const
Definition: transonic_perturbation_potential_flow_element.cpp:361
std::string Info() const override
Turn back information as a string.
Definition: transonic_perturbation_potential_flow_element.cpp:337
void CalculateLeftHandSideNormalElement(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo)
Definition: transonic_perturbation_potential_flow_element.cpp:631
virtual void CalculateRightHandSideContribution(BoundedVector< double, TNumNodes > &rRhs_total, const double rDensity, const array_1d< double, TDim > &rVelocity)
Definition: transonic_perturbation_potential_flow_element.cpp:619
void CalculateLeftHandSideWakeElement(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo)
Definition: transonic_perturbation_potential_flow_element.cpp:685
bool CheckUpwindElement()
Definition: transonic_perturbation_potential_flow_element.cpp:375
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: transonic_perturbation_potential_flow_element.cpp:92
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(TransonicPerturbationPotentialFlowElement)
void CalculateRightHandSideNormalElement(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo)
Definition: transonic_perturbation_potential_flow_element.cpp:580
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: transonic_perturbation_potential_flow_element.cpp:76
void pSetUpwindElement(GlobalPointer< Element > pUpwindElement)
Definition: transonic_perturbation_potential_flow_element.cpp:369
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, PropertiesType::Pointer pProperties) const override
It creates a new element pointer.
Definition: transonic_perturbation_potential_flow_element.cpp:26
TransonicPerturbationPotentialFlowElement & operator=(TransonicPerturbationPotentialFlowElement &&rOther)=delete
Move operator.
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &CurrentProcessInfo) const override
Definition: transonic_perturbation_potential_flow_element.cpp:125
void CalculateOnIntegrationPoints(const Variable< double > &rVariable, std::vector< double > &rValues, const ProcessInfo &rCurrentProcessInfo) override
Definition: transonic_perturbation_potential_flow_element.cpp:225
TransonicPerturbationPotentialFlowElement(IndexType NewId=0)
Default constuctor.
Definition: transonic_perturbation_potential_flow_element.h:72
PotentialFlowUtilities::ElementalData< TNumNodes, TDim > ElementalData
Definition: transonic_perturbation_potential_flow_element.h:55
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
float velocity
Definition: PecletTest.py:54
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: potential_flow_utilities.h:39