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.
base_load_condition.h
Go to the documentation of this file.
1 // KRATOS ___| | | |
2 // \___ \ __| __| | | __| __| | | __| _` | |
3 // | | | | | ( | | | | ( | |
4 // _____/ \__|_| \__,_|\___|\__|\__,_|_| \__,_|_| MECHANICS
5 //
6 // License: BSD License
7 // license: StructuralMechanicsApplication/license.txt
8 //
9 // Main authors: Riccardo Rossi
10 // Vicente Mataix Ferrandiz
11 //
12 
13 #pragma once
14 
15 // System includes
16 
17 // External includes
18 
19 // Project includes
20 #include "includes/condition.h"
22 
23 namespace Kratos
24 {
25 
28 
32 
36 
40 
44 
52 class KRATOS_API(STRUCTURAL_MECHANICS_APPLICATION) BaseLoadCondition
53  : public Condition
54 {
55 public:
56 
59 
62 
65 
68 
71 
74 
77 
80 
81  // Counted pointer of BaseLoadCondition
83 
87 
88  // Constructor void
90  {};
91 
92  // Constructor using an array of nodes
93  BaseLoadCondition( IndexType NewId, GeometryType::Pointer pGeometry ):Condition(NewId,pGeometry)
94  {};
95 
96  // Constructor using an array of nodes with properties
97  BaseLoadCondition( IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties ):Condition(NewId,pGeometry,pProperties)
98  {};
99 
101  BaseLoadCondition(BaseLoadCondition const& rOther);
102 
103  // Destructor
105  {};
106 
110 
113 
117 
125  Condition::Pointer Create(
126  IndexType NewId,
127  NodesArrayType const& ThisNodes,
128  PropertiesType::Pointer pProperties
129  ) const override;
130 
138  Condition::Pointer Create(
139  IndexType NewId,
140  GeometryType::Pointer pGeom,
141  PropertiesType::Pointer pProperties
142  ) const override;
143 
150  Condition::Pointer Clone (
151  IndexType NewId,
152  NodesArrayType const& ThisNodes
153  ) const override;
154 
160  void EquationIdVector(
161  EquationIdVectorType& rResult,
162  const ProcessInfo& rCurrentProcessInfo
163  ) const override;
164 
170  void GetDofList(
171  DofsVectorType& ElementalDofList,
172  const ProcessInfo& rCurrentProcessInfo
173  ) const override;
174 
180  void GetValuesVector(
181  Vector& rValues,
182  int Step = 0
183  ) const override;
184 
190  void GetFirstDerivativesVector(
191  Vector& rValues,
192  int Step = 0
193  ) const override;
194 
200  void GetSecondDerivativesVector(
201  Vector& rValues,
202  int Step = 0
203  ) const override;
204 
213  void CalculateLocalSystem(
214  MatrixType& rLeftHandSideMatrix,
215  VectorType& rRightHandSideVector,
216  const ProcessInfo& rCurrentProcessInfo
217  ) override;
218 
224  void CalculateRightHandSide(
225  VectorType& rRightHandSideVector,
226  const ProcessInfo& rCurrentProcessInfo
227  ) override;
228 
234  void CalculateMassMatrix(
235  MatrixType& rMassMatrix,
236  const ProcessInfo& rCurrentProcessInfo
237  ) override;
238 
244  void CalculateDampingMatrix(
245  MatrixType& rDampingMatrix,
246  const ProcessInfo& rCurrentProcessInfo
247  ) override;
248 
256  void AddExplicitContribution(const VectorType& rRHS,
257  const Variable<VectorType>& rRHSVariable,
258  const Variable<array_1d<double,3> >& rDestinationVariable,
259  const ProcessInfo& rCurrentProcessInfo
260  ) override;
261 
267  int Check( const ProcessInfo& rCurrentProcessInfo ) const override;
268 
273  virtual bool HasRotDof() const;
274 
279  unsigned int GetBlockSize() const
280  {
281  unsigned int dim = GetGeometry().WorkingSpaceDimension();
282  if( HasRotDof() ) { // if it has rotations
283  if(dim == 2)
284  return 3;
285  else if(dim == 3)
286  return 6;
287  else
288  KRATOS_ERROR << "The conditions only works for 2D and 3D elements";
289  } else {
290  return dim;
291  }
292  }
293 
297 
301 
305 
311  const Parameters GetSpecifications() const override;
312 
314  std::string Info() const override
315  {
316  std::stringstream buffer;
317  buffer << "Base load Condition #" << Id();
318  return buffer.str();
319  }
320 
322 
323  void PrintInfo(std::ostream& rOStream) const override
324  {
325  rOStream << "Base load Condition #" << Id();
326  }
327 
329  void PrintData(std::ostream& rOStream) const override
330  {
331  pGetGeometry()->PrintData(rOStream);
332  }
333 
337 
338 protected:
339 
342 
346 
350 
354 
363  virtual void CalculateAll(
364  MatrixType& rLeftHandSideMatrix,
365  VectorType& rRightHandSideVector,
366  const ProcessInfo& rCurrentProcessInfo,
367  const bool CalculateStiffnessMatrixFlag,
368  const bool CalculateResidualVectorFlag
369  );
370 
377  virtual double GetIntegrationWeight(
378  const GeometryType::IntegrationPointsArrayType& IntegrationPoints,
379  const SizeType PointNumber,
380  const double detJ
381  ) const;
382 
386 
390 
394 
395 private:
398 
399 
403 
404 
405 
409 
413 
414 
418 
419 
423 
427 
428  friend class Serializer;
429 
430  void save( Serializer& rSerializer ) const override;
431 
432  void load( Serializer& rSerializer ) override;
433 
434 }; // class BaseLoadCondition.
435 
439 
440 
444 
445 } // namespace Kratos.
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
This is the base class of all the load conditions on StructuralMechanicsApplication.
Definition: base_load_condition.h:54
BaseType::GeometryType GeometryType
Definition of the geometry type with given NodeType.
Definition: base_load_condition.h:76
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: base_load_condition.h:329
std::string Info() const override
Turn back information as a string.
Definition: base_load_condition.h:314
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: base_load_condition.h:323
BaseType::NodesArrayType NodesArrayType
Definition of nodes container type, redefined from GeometryType.
Definition: base_load_condition.h:79
BaseType::SizeType SizeType
Definition of the size type.
Definition: base_load_condition.h:67
BaseType::PropertiesType PropertiesType
Definition of the properties type.
Definition: base_load_condition.h:73
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(BaseLoadCondition)
unsigned int GetBlockSize() const
This method computes the DoF block size.
Definition: base_load_condition.h:279
BaseType::NodeType NodeType
Definition of the node type.
Definition: base_load_condition.h:70
BaseLoadCondition(IndexType NewId, GeometryType::Pointer pGeometry, PropertiesType::Pointer pProperties)
Definition: base_load_condition.h:97
BaseLoadCondition(IndexType NewId, GeometryType::Pointer pGeometry)
Definition: base_load_condition.h:93
~BaseLoadCondition() override
Definition: base_load_condition.h:104
Condition BaseType
We define the base class Condition.
Definition: base_load_condition.h:61
BaseLoadCondition()
Definition: base_load_condition.h:89
BaseType::IndexType IndexType
Dfinition of the index type.
Definition: base_load_condition.h:64
Base class for all Conditions.
Definition: condition.h:59
std::size_t IndexType
Definition: flags.h:74
std::size_t IndexType
Defines the index type.
Definition: geometrical_object.h:73
Geometry base class.
Definition: geometry.h:71
std::vector< IntegrationPointType > IntegrationPointsArrayType
Definition: geometry.h:161
This class defines the node.
Definition: node.h:65
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
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
std::size_t IndexType
The definition of the index type.
Definition: key_hash.h:35
#define KRATOS_ERROR
Definition: exception.h:161
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
std::size_t SizeType
Definition: nurbs_utilities.h:41
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
Properties PropertiesType
Definition: regenerate_pfem_pressure_conditions_process.h:26
std::size_t SizeType
The definition of the size type.
Definition: mortar_classes.h:43
ModelPart::NodesContainerType NodesArrayType
Definition: gid_gauss_point_container.h:42
def load(f)
Definition: ode_solve.py:307
int dim
Definition: sensitivityMatrix.py:25