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.
load_moment_director_5p_condition.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 
11 #if !defined(KRATOS_MOMENT_LOAD_DIRECTOR_5p_CONDITION_H_INCLUDED )
12 #define KRATOS_MOMENT_LOAD_DIRECTOR_5p_CONDITION_H_INCLUDED
13 
14 // System includes
15 
16 // External includes
17 
18 // Project includes
20 
21 #include "includes/condition.h"
22 
23 namespace Kratos
24 {
27  : public Condition
28  {
29  public:
32 
35 
37  typedef std::size_t SizeType;
38  typedef std::size_t IndexType;
39 
41  typedef typename GeometryType::Pointer GeometryPointerType;
42 
46 
49  IndexType NewId,
50  GeometryType::Pointer pGeometry)
51  : Condition(NewId, pGeometry)
52  {};
53 
56  IndexType NewId,
57  GeometryType::Pointer pGeometry,
58  PropertiesType::Pointer pProperties)
59  : Condition(NewId, pGeometry, pProperties)
60  {};
61 
64  {};
65 
68 
72 
74  Condition::Pointer Create(
75  IndexType NewId,
76  GeometryType::Pointer pGeom,
77  PropertiesType::Pointer pProperties
78  ) const final
79  {
80  return Kratos::make_intrusive<LoadMomentDirector5pCondition>(
81  NewId, pGeom, pProperties);
82  };
83 
85  Condition::Pointer Create(
86  IndexType NewId,
87  NodesArrayType const& ThisNodes,
88  PropertiesType::Pointer pProperties
89  ) const final
90  {
91  return Kratos::make_intrusive<LoadMomentDirector5pCondition>(
92  NewId, GetGeometry().Create(ThisNodes), pProperties);
93  };
94 
98 
106  VectorType& rRightHandSideVector,
107  const ProcessInfo& rCurrentProcessInfo) final
108  {
109  MatrixType left_hand_side_matrix = Matrix(0, 0);
110 
111  CalculateAll(left_hand_side_matrix, rRightHandSideVector,
112  rCurrentProcessInfo, false, true);
113  }
114 
122  MatrixType& rLeftHandSideMatrix,
123  const ProcessInfo& rCurrentProcessInfo) final
124  {
125  VectorType right_hand_side_vector = Vector(0);
126 
127  CalculateAll(rLeftHandSideMatrix, right_hand_side_vector,
128  rCurrentProcessInfo, true, false);
129  }
130 
140  MatrixType& rLeftHandSideMatrix,
141  VectorType& rRightHandSideVector,
142  const ProcessInfo& rCurrentProcessInfo) final
143  {
144  CalculateAll(rLeftHandSideMatrix, rRightHandSideVector,
145  rCurrentProcessInfo, true, true);
146  }
147 
153  void EquationIdVector(
154  EquationIdVectorType& rResult,
155  const ProcessInfo& rCurrentProcessInfo
156  ) const final;
157 
163  void GetDofList(
164  DofsVectorType& rElementalDofList,
165  const ProcessInfo& rCurrentProcessInfo
166  ) const final;
167 
176  void CalculateAll(
177  MatrixType& rLeftHandSideMatrix,
178  VectorType& rRightHandSideVector,
179  const ProcessInfo& rCurrentProcessInfo,
180  const bool CalculateStiffnessMatrixFlag,
181  const bool CalculateResidualVectorFlag
182  );
183 
185  const GeometryType& rGeometry,
186  Vector& rDeterminantOfJacobian);
187 
189  const GeometryType& rGeometry,
190  const Matrix& r_N,
191  const IndexType& point_number,
192  const array_1d<double, 3>& momentload);
193 
197 
199  std::string Info() const final
200  {
201  std::stringstream buffer;
202  buffer << "\"LoadCondition\" #" << Id();
203  return buffer.str();
204  }
205 
207  void PrintInfo(std::ostream& rOStream) const final
208  {
209  rOStream << "\"LoadCondition\" #" << Id();
210  }
211 
213  void PrintData(std::ostream& rOStream) const final
214  {
215  pGetGeometry()->PrintData(rOStream);
216  }
217 
219 
220  private:
223 
224  friend class Serializer;
225 
226  virtual void save(Serializer& rSerializer) const final
227  {
229  }
230 
231  virtual void load(Serializer& rSerializer) final
232  {
234  }
235 
237 
238  }; // Class LoadMomentDirector5pCondition
239 
240 } // namespace Kratos.
241 
242 #endif // KRATOS_MOMENT_LOAD_DIRECTOR_5p_CONDITION_H_INCLUDED defined
Base class for all Conditions.
Definition: condition.h:59
std::vector< std::size_t > EquationIdVectorType
Definition: condition.h:98
Matrix MatrixType
Definition: condition.h:90
std::vector< DofType::Pointer > DofsVectorType
Definition: condition.h:100
std::size_t IndexType
Definition: flags.h:74
GeometryType::Pointer pGetGeometry()
Returns the pointer to the geometry.
Definition: geometrical_object.h:140
GeometryType & GetGeometry()
Returns the reference of the geometry.
Definition: geometrical_object.h:158
Geometry base class.
Definition: geometry.h:71
IndexType Id() const
Definition: indexed_object.h:107
Condition for moment loads for the 5p shell based on directors.
Definition: load_moment_director_5p_condition.h:28
array_1d< double, 3 > calculateMomentLoadTimesDirectorTestFunction(const GeometryType &rGeometry, const Matrix &r_N, const IndexType &point_number, const array_1d< double, 3 > &momentload)
Definition: load_moment_director_5p_condition.cpp:94
LoadMomentDirector5pCondition()
Default constructor.
Definition: load_moment_director_5p_condition.h:63
void GetDofList(DofsVectorType &rElementalDofList, const ProcessInfo &rCurrentProcessInfo) const final
Sets on rConditionDofList the degrees of freedom of the considered element geometry.
Definition: load_moment_director_5p_condition.cpp:165
std::size_t IndexType
Definition: load_moment_director_5p_condition.h:38
std::size_t SizeType
Size types.
Definition: load_moment_director_5p_condition.h:37
std::string Info() const final
Turn back information as a string.
Definition: load_moment_director_5p_condition.h:199
void CalculateLeftHandSide(MatrixType &rLeftHandSideMatrix, const ProcessInfo &rCurrentProcessInfo) final
This is called during the assembling process in order to calculate the condition left hand side matri...
Definition: load_moment_director_5p_condition.h:121
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &rCurrentProcessInfo) const final
Sets on rResult the ID's of the element degrees of freedom.
Definition: load_moment_director_5p_condition.cpp:146
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) final
This is called during the assembling process in order to calculate the condition right hand side matr...
Definition: load_moment_director_5p_condition.h:105
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) final
This function provides a more general interface to the element.
Definition: load_moment_director_5p_condition.h:139
GeometryType::Pointer GeometryPointerType
Definition: load_moment_director_5p_condition.h:41
LoadMomentDirector5pCondition(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Constructor with Id, geometry and property.
Definition: load_moment_director_5p_condition.h:55
Condition::Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, PropertiesType::Pointer pProperties) const final
Create with Id, pointer to geometry and pointer to property.
Definition: load_moment_director_5p_condition.h:85
void PrintData(std::ostream &rOStream) const final
Print object's data.
Definition: load_moment_director_5p_condition.h:213
~LoadMomentDirector5pCondition()=default
Destructor.
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(LoadMomentDirector5pCondition)
Counted pointer definition of LoadCondition.
void PrintInfo(std::ostream &rOStream) const final
Print information about this object.
Definition: load_moment_director_5p_condition.h:207
void DeterminantOfJacobianInitial(const GeometryType &rGeometry, Vector &rDeterminantOfJacobian)
Definition: load_moment_director_5p_condition.cpp:114
Condition::Pointer Create(IndexType NewId, GeometryType::Pointer pGeom, PropertiesType::Pointer pProperties) const final
Create with Id, pointer to geometry and pointer to property.
Definition: load_moment_director_5p_condition.h:74
Geometry< Node > GeometryType
Definition: load_moment_director_5p_condition.h:40
void CalculateAll(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo, const bool CalculateStiffnessMatrixFlag, const bool CalculateResidualVectorFlag)
Definition: load_moment_director_5p_condition.cpp:24
LoadMomentDirector5pCondition(IndexType NewId, GeometryType::Pointer pGeometry)
Constructor with Id and geometry.
Definition: load_moment_director_5p_condition.h:48
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
#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
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
Internals::Matrix< double, AMatrix::dynamic, 1 > Vector
Definition: amatrix_interface.h:472
Internals::Matrix< double, AMatrix::dynamic, AMatrix::dynamic > Matrix
Definition: amatrix_interface.h:470
def load(f)
Definition: ode_solve.py:307
tuple const
Definition: ode_solve.py:403