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.
mpm_grid_base_load_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 // Main authors: Bodhinanda Chandra
11 //
12 
13 
14 #if !defined(KRATOS_MPM_GRID_BASE_LOAD_CONDITION_3D_H_INCLUDED )
15 #define KRATOS_MPM_GRID_BASE_LOAD_CONDITION_3D_H_INCLUDED
16 
17 // System includes
18 
19 // External includes
20 
21 // Project includes
22 #include "includes/define.h"
23 #include "includes/condition.h"
25 
26 namespace Kratos
27 {
28 
31 
35 
39 
43 
47 
49  : public Condition
50 {
51 public:
52 
54  typedef std::size_t SizeType;
56 
57  // Counted pointer of MPMGridBaseLoadCondition
59 
63 
64  // Constructor void
66  {};
67 
68  // Constructor using an array of nodes
69  MPMGridBaseLoadCondition( IndexType NewId, GeometryType::Pointer pGeometry ):Condition(NewId,pGeometry)
70  {};
71 
72  // Constructor using an array of nodes with properties
73  MPMGridBaseLoadCondition( IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties ):Condition(NewId,pGeometry,pProperties)
74  {};
75 
76  // Destructor
78  {};
79 
83 
84 
88 
94  void EquationIdVector(
95  EquationIdVectorType& rResult,
96  const ProcessInfo& rCurrentProcessInfo
97  ) const override;
98 
104  void GetDofList(
105  DofsVectorType& ElementalDofList,
106  const ProcessInfo& rCurrentProcessInfo
107  ) const override;
108 
114  void GetValuesVector(
115  Vector& rValues,
116  int Step = 0
117  ) const override;
118 
125  Vector& rValues,
126  int Step = 0
127  ) const override;
128 
135  Vector& rValues,
136  int Step = 0
137  ) const override;
138 
148  MatrixType& rLeftHandSideMatrix,
149  VectorType& rRightHandSideVector,
150  const ProcessInfo& rCurrentProcessInfo
151  ) override;
152 
159  VectorType& rRightHandSideVector,
160  const ProcessInfo& rCurrentProcessInfo
161  ) override;
162 
168  void CalculateMassMatrix(
169  MatrixType& rMassMatrix,
170  const ProcessInfo& rCurrentProcessInfo
171  ) override;
172 
180  MatrixType& rDampingMatrix,
181  const ProcessInfo& rCurrentProcessInfo
182  ) override;
183 
193  void AddExplicitContribution(const VectorType& rRHS,
194  const Variable<VectorType>& rRHSVariable,
195  const Variable<array_1d<double,3> >& rDestinationVariable,
196  const ProcessInfo& rCurrentProcessInfo) override;
197 
205  int Check( const ProcessInfo& rCurrentProcessInfo ) const override;
206 
210  bool HasRotDof(){return (GetGeometry()[0].HasDofFor(ROTATION_X) && GetGeometry().size() == 2);};
211 
212  unsigned int GetBlockSize()
213  {
214  unsigned int dimension = GetGeometry().WorkingSpaceDimension();
215  if( HasRotDof() ) // if it has rotations
216  {
217  if(dimension == 2)
218  return 3;
219  else if(dimension == 3)
220  return 6;
221  else
222  KRATOS_ERROR << "the conditions only works for 2D and 3D elements";
223  }
224  else
225  {
226  return dimension;
227  }
228  }
229 
233 
234 
238 
239 
243 
247 
248 protected:
249 
252 
256 
260 
264 
273  virtual void CalculateAll(
274  MatrixType& rLeftHandSideMatrix,
275  VectorType& rRightHandSideVector,
276  const ProcessInfo& rCurrentProcessInfo,
277  const bool CalculateStiffnessMatrixFlag,
278  const bool CalculateResidualVectorFlag
279  );
280 
287  virtual double GetIntegrationWeight(
288  const GeometryType::IntegrationPointsArrayType& IntegrationPoints,
289  const unsigned int PointNumber,
290  const double detJ
291  );
292 
296 
300 
304 
305 private:
308 
309 
313 
314 
315 
319 
323 
324 
328 
329 
333 
337 
338  friend class Serializer;
339 
340  void save( Serializer& rSerializer ) const override
341  {
343  }
344 
345  void load( Serializer& rSerializer ) override
346  {
348  }
349 
350 }; // class MPMGridBaseLoadCondition.
351 
355 
356 
360 
361 } // namespace Kratos.
362 
363 #endif // KRATOS_MPM_GRID_BASE_LOAD_CONDITION_3D_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 & GetGeometry()
Returns the reference of the geometry.
Definition: geometrical_object.h:158
std::vector< IntegrationPointType > IntegrationPointsArrayType
Definition: geometry.h:161
SizeType WorkingSpaceDimension() const
Definition: geometry.h:1287
Definition: mpm_grid_base_load_condition.h:50
virtual double GetIntegrationWeight(const GeometryType::IntegrationPointsArrayType &IntegrationPoints, const unsigned int PointNumber, const double detJ)
Definition: mpm_grid_base_load_condition.cpp:275
MPMGridBaseLoadCondition(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Definition: mpm_grid_base_load_condition.h:73
~MPMGridBaseLoadCondition() override
Definition: mpm_grid_base_load_condition.h:77
void GetFirstDerivativesVector(Vector &rValues, int Step=0) const override
Definition: mpm_grid_base_load_condition.cpp:132
void CalculateDampingMatrix(MatrixType &rDampingMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: mpm_grid_base_load_condition.cpp:228
virtual void CalculateAll(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo, const bool CalculateStiffnessMatrixFlag, const bool CalculateResidualVectorFlag)
Definition: mpm_grid_base_load_condition.cpp:242
void GetSecondDerivativesVector(Vector &rValues, int Step=0) const override
Definition: mpm_grid_base_load_condition.cpp:161
MPMGridBaseLoadCondition(IndexType NewId, GeometryType::Pointer pGeometry)
Definition: mpm_grid_base_load_condition.h:69
void GetValuesVector(Vector &rValues, int Step=0) const override
Definition: mpm_grid_base_load_condition.cpp:103
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(MPMGridBaseLoadCondition)
int Check(const ProcessInfo &rCurrentProcessInfo) const override
Definition: mpm_grid_base_load_condition.cpp:255
void CalculateLocalSystem(MatrixType &rLeftHandSideMatrix, VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: mpm_grid_base_load_condition.cpp:202
unsigned int GetBlockSize()
Definition: mpm_grid_base_load_condition.h:212
void AddExplicitContribution(const VectorType &rRHS, const Variable< VectorType > &rRHSVariable, const Variable< array_1d< double, 3 > > &rDestinationVariable, const ProcessInfo &rCurrentProcessInfo) override
Definition: mpm_grid_base_load_condition.cpp:287
std::size_t SizeType
Definition: mpm_grid_base_load_condition.h:54
MPMGridBaseLoadCondition()
Definition: mpm_grid_base_load_condition.h:65
void CalculateMassMatrix(MatrixType &rMassMatrix, const ProcessInfo &rCurrentProcessInfo) override
Definition: mpm_grid_base_load_condition.cpp:214
void CalculateRightHandSide(VectorType &rRightHandSideVector, const ProcessInfo &rCurrentProcessInfo) override
Definition: mpm_grid_base_load_condition.cpp:190
void GetDofList(DofsVectorType &ElementalDofList, const ProcessInfo &rCurrentProcessInfo) const override
Definition: mpm_grid_base_load_condition.cpp:67
void EquationIdVector(EquationIdVectorType &rResult, const ProcessInfo &rCurrentProcessInfo) const override
Definition: mpm_grid_base_load_condition.cpp:27
bool HasRotDof()
Definition: mpm_grid_base_load_condition.h:210
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
#define KRATOS_SERIALIZE_SAVE_BASE_CLASS(Serializer, BaseType)
Definition: define.h:812
#define KRATOS_SERIALIZE_LOAD_BASE_CLASS(Serializer, BaseType)
Definition: define.h:815
#define KRATOS_ERROR
Definition: exception.h:161
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
int dimension
Definition: isotropic_damage_automatic_differentiation.py:123
def load(f)
Definition: ode_solve.py:307