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.
calculate_gradient_Pouliot_2012.h
Go to the documentation of this file.
1 //
2 // Project Name: Kratos
3 // Last Modified by: $Author: gcasas $
4 // Date: $Date: 2016-03-12
5 //
6 
7 #if !defined(KRATOS_COMPUTE_GRADIENT_POULIOT_2012_H_INCLUDED )
8 #define KRATOS_COMPUTE_GRADIENT_POULIOT_2012_H_INCLUDED
9 
10 // System includes
11 #include <string>
12 #include <iostream>
13 
14 
15 // External includes
16 
17 
18 // Project includes
19 #include "containers/array_1d.h"
20 #include "includes/define.h"
21 #include "includes/element.h"
22 #include "includes/serializer.h"
23 #include "geometries/geometry.h"
24 #include "utilities/math_utils.h"
25 #include "utilities/geometry_utilities.h"
26 
27 // Application includes
28 #include "includes/variables.h"
31 
32 namespace Kratos
33 {
34 
37 
40 
44 
48 
52 
56 
58 
60 template< unsigned int TDim,
61  unsigned int TNumNodes = TDim + 1 >
62 class KRATOS_API(SWIMMING_DEM_APPLICATION) ComputeGradientPouliot2012 : public ComputeComponentGradientSimplex<TDim, TNumNodes>
63 {
64 public:
67 
70 
73  typedef Node NodeType;
74 
77 
80 
82  typedef Vector VectorType;
83 
85  typedef Matrix MatrixType;
86 
88 
89  typedef std::size_t IndexType;
90 
91  typedef std::size_t SizeType;
92 
93  typedef std::vector<std::size_t> EquationIdVectorType;
94 
95  typedef std::vector< Dof<double>::Pointer > DofsVectorType;
96 
98 
101 
104 
107 
111 
112  //Constructors.
113 
115 
119  BaseType(NewId)
120  {}
121 
123 
128  BaseType(NewId, ThisNodes)
129  {}
130 
132 
136  ComputeGradientPouliot2012(IndexType NewId, GeometryType::Pointer pGeometry) :
137  BaseType(NewId, pGeometry)
138  {}
139 
141 
146  ComputeGradientPouliot2012(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties) :
147  BaseType(NewId, pGeometry, pProperties)
148  {}
149 
152  {}
153 
154 
158 
159 
163 
165 
172  Element::Pointer Create(IndexType NewId, NodesArrayType const& ThisNodes,
173  PropertiesType::Pointer pProperties) const override
174  {
175  return Element::Pointer(new ComputeGradientPouliot2012(NewId, this->GetGeometry().Create(ThisNodes), pProperties));
176  }
177 
181  virtual void CalculateLocalSystem(MatrixType& rLeftHandSideMatrix, VectorType& rRightHandSideVector, const ProcessInfo& rCurrentProcessInfo) override;
182 
184 
190  virtual void EquationIdVector(EquationIdVectorType& rResult, const ProcessInfo& rCurrentProcessInfo) const override;
191 
193 
197  virtual void GetDofList(DofsVectorType& rElementalDofList, const ProcessInfo& rCurrentProcessInfo) const override;
198 
200 
208  virtual int Check(const ProcessInfo& rCurrentProcessInfo) const override;
209 
213 
214 
218 
220  virtual std::string Info() const override
221  {
222  std::stringstream buffer;
223  buffer << "ComputeGradientPouliot2012 #" << this->Id();
224  return buffer.str();
225  }
226 
228  virtual void PrintInfo(std::ostream& rOStream) const override
229  {
230  rOStream << "ComputeGradientPouliot2012" << TDim << "D";
231  }
232 
233 // /// Print object's data.
234 // virtual void PrintData(std::ostream& rOStream) const;
235 
239 
240 
242 
243 protected:
246 
247 
251 
252 
256 
257 
261 
262 
266 
267 
271 
272 
276 
277 
279 
280 private:
283 
287 
291 
292  friend class Serializer;
293 
297 
299 
307  virtual void AddPouliot2012LHS(MatrixType& rLeftHandSideMatrix, const ProcessInfo& rCurrentProcessInfo);
308 
309  virtual void AddPouliot2012StabilizationLHS(const double epsilon, MatrixType& rLeftHandSideMatrix, const ProcessInfo& rCurrentProcessInfo);
310 
311  virtual void AddFEMLaplacianStabilizationLHS(const double epsilon, MatrixType& rLeftHandSideMatrix, const ProcessInfo& rCurrentProcessInfo);
312 
313  void AssembleEdgeLHSContribution(const unsigned int edge[2], const array_1d<double, 3>& edge_normalized_vector, MatrixType& rLeftHandSideMatrix);
314 
315  virtual void AddPouliot2012RHS(VectorType& F, const ProcessInfo& rCurrentProcessInfo);
316 
317  virtual void AddStabilizationRHSContribution(VectorType& F,
318  const array_1d<double, TNumNodes>& rShapeFunc,
319  const double Weight);
320 
321  void CalculateStabilizationRHS(const double epsilon, VectorType& F, const ProcessInfo& rCurrentProcessInfo);
322 
323  void AssembleEdgeRHSContributionX(const unsigned int edge[2], const double h_edge_inv, const array_1d<double, 3>& edge_normalized_vector, VectorType& F);
324 
325  void AssembleEdgeRHSContributionY(const unsigned int edge[2], const double h_edge_inv, const array_1d<double, 3>& edge_normalized_vector, VectorType& F);
326 
327  void AssembleEdgeRHSContributionZ(const unsigned int edge[2], const double h_edge_inv, const array_1d<double, 3>& edge_normalized_vector, VectorType& F);
328 
329 
333 
334 
338 
339 
343 
344 
348 
351 
354 
356 
357 }; // Class ComputeGradientPouliot2012
358 
360 
363 
364 
368 
369 
371 template< unsigned int TDim >
372 inline std::istream& operator >>(std::istream& rIStream,
374 {
375  return rIStream;
376 }
377 
379 template< unsigned int TDim >
380 inline std::ostream& operator <<(std::ostream& rOStream,
382 {
383  rThis.PrintInfo(rOStream);
384  rOStream << std::endl;
385  rThis.PrintData(rOStream);
386 
387  return rOStream;
388 }
390 
392 
393 } // namespace Kratos.
394 
395 #endif // KRATOS_COMPUTE_GRADIENT_POULIOT_2012_H_INCLUDED defined
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
A post-processing element to recover the Laplacian from the velocity solution.
Definition: calculate_component_gradient_simplex_element.h:63
A post-processing element to recover the Laplacian from the velocity solution.
Definition: calculate_gradient_Pouliot_2012.h:63
std::vector< std::size_t > EquationIdVectorType
Definition: calculate_gradient_Pouliot_2012.h:93
ComputeGradientPouliot2012(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Constuctor using geometry and properties.
Definition: calculate_gradient_Pouliot_2012.h:146
Geometry< NodeType >::PointsArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: calculate_gradient_Pouliot_2012.h:79
Kratos::Vector ShapeFunctionsType
Type for shape function values container.
Definition: calculate_gradient_Pouliot_2012.h:100
ComputeGradientPouliot2012(IndexType NewId=0)
Default constuctor.
Definition: calculate_gradient_Pouliot_2012.h:118
ComputeGradientPouliot2012(IndexType NewId, const NodesArrayType &ThisNodes)
Constructor using an array of nodes.
Definition: calculate_gradient_Pouliot_2012.h:127
GeometryType::ShapeFunctionsGradientsType ShapeFunctionDerivativesArrayType
Type for an array of shape function gradient matrices.
Definition: calculate_gradient_Pouliot_2012.h:106
std::size_t SizeType
Definition: calculate_gradient_Pouliot_2012.h:91
Geometry< NodeType > GeometryType
Geometry type (using with given NodeType)
Definition: calculate_gradient_Pouliot_2012.h:76
virtual std::string Info() const override
Turn back information as a string.
Definition: calculate_gradient_Pouliot_2012.h:220
PointerVectorSet< Dof< double >, IndexedObject > DofsArrayType
Definition: calculate_gradient_Pouliot_2012.h:97
Node NodeType
Node type (default is: Node)
Definition: calculate_gradient_Pouliot_2012.h:73
Element::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, PropertiesType::Pointer pProperties) const override
Create a new element of this type.
Definition: calculate_gradient_Pouliot_2012.h:172
virtual ~ComputeGradientPouliot2012()
Destructor.
Definition: calculate_gradient_Pouliot_2012.h:151
Matrix MatrixType
Matrix type for local contributions to the linear system.
Definition: calculate_gradient_Pouliot_2012.h:85
Kratos::Matrix ShapeFunctionDerivativesType
Type for a matrix containing the shape function gradients.
Definition: calculate_gradient_Pouliot_2012.h:103
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(ComputeGradientPouliot2012)
Pointer definition of ComputeGradientPouliot2012.
ComputeComponentGradientSimplex< TDim, TNumNodes > BaseType
Definition: calculate_gradient_Pouliot_2012.h:71
Properties PropertiesType
Definition: calculate_gradient_Pouliot_2012.h:87
ComputeGradientPouliot2012(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor using a geometry object.
Definition: calculate_gradient_Pouliot_2012.h:136
Vector VectorType
Vector type for local contributions to the linear system.
Definition: calculate_gradient_Pouliot_2012.h:82
virtual void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: calculate_gradient_Pouliot_2012.h:228
std::size_t IndexType
Definition: calculate_gradient_Pouliot_2012.h:89
std::vector< Dof< double >::Pointer > DofsVectorType
Definition: calculate_gradient_Pouliot_2012.h:95
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: element.h:1135
std::size_t IndexType
Definition: flags.h:74
This defines the geometrical object, base definition of the element and condition entities.
Definition: geometrical_object.h:58
Geometry base class.
Definition: geometry.h:71
This object defines an indexed object.
Definition: indexed_object.h:54
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
A sorted associative container similar to an STL set, but uses a vector to store pointers to its data...
Definition: pointer_vector_set.h:72
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
Properties encapsulates data shared by different Elements or Conditions. It can store any type of dat...
Definition: properties.h:69
The serialization consists in storing the state of an object into a storage format like data file or ...
Definition: serializer.h:123
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
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
F
Definition: hinsberg_optimization.py:144