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.
ALM_frictionless_mortar_contact_condition.h
Go to the documentation of this file.
1 // KRATOS ______ __ __ _____ __ __ __
2 // / ____/___ ____ / /_____ ______/ /_/ ___// /________ _______/ /___ ___________ _/ /
3 // / / / __ \/ __ \/ __/ __ `/ ___/ __/\__ \/ __/ ___/ / / / ___/ __/ / / / ___/ __ `/ /
4 // / /___/ /_/ / / / / /_/ /_/ / /__/ /_ ___/ / /_/ / / /_/ / /__/ /_/ /_/ / / / /_/ / /
5 // \____/\____/_/ /_/\__/\__,_/\___/\__//____/\__/_/ \__,_/\___/\__/\__,_/_/ \__,_/_/ MECHANICS
6 //
7 // License: BSD License
8 // license: ContactStructuralMechanicsApplication/license.txt
9 //
10 // Main authors: Vicente Mataix Ferrandiz
11 //
12 
13 #pragma once
14 
15 // System includes
16 
17 // External includes
18 
19 // Project includes
21 
22 namespace Kratos
23 {
24 
27 
31 
33  using SizeType = std::size_t;
34 
38 
42 
46 
60 template< SizeType TDim, SizeType TNumNodes, bool TNormalVariation, const SizeType TNumNodesMaster = TNumNodes >
61 class KRATOS_API(CONTACT_STRUCTURAL_MECHANICS_APPLICATION) AugmentedLagrangianMethodFrictionlessMortarContactCondition
62  : public MortarContactCondition<TDim, TNumNodes, FrictionalCase::FRICTIONLESS, TNormalVariation, TNumNodesMaster>
63 {
64 public:
67 
70 
73 
75  using MortarConditionMatrices = typename BaseType::MortarConditionMatrices;
76 
79 
82 
85 
88 
91 
93  using GeometryPointerType = typename ConditionBaseType::GeometryType::Pointer;
94 
97 
100 
102  using PropertiesPointerType = typename PropertiesType::Pointer;
103 
106 
109 
111  using PointType = Point;
112 
115 
118 
120  using ConditionArrayListType = std::vector<array_1d<PointType, TDim>>;
121 
124 
127 
130 
133 
135  static constexpr IndexType MatrixSize = TDim * (TNumNodes + TNumNodesMaster) + TNumNodes;
136 
140 
143  : BaseType()
144  {
145  }
146 
147  // Constructor 1
149  IndexType NewId,
150  GeometryPointerType pGeometry
151  ):BaseType(NewId, pGeometry)
152  {
153  }
154 
155  // Constructor 2
157  IndexType NewId,
158  GeometryPointerType pGeometry,
159  PropertiesPointerType pProperties
160  ):BaseType( NewId, pGeometry, pProperties )
161  {
162  }
163 
164  // Constructor 3
166  IndexType NewId,
167  GeometryPointerType pGeometry,
168  PropertiesPointerType pProperties,
169  GeometryType::Pointer pMasterGeometry
170  ):BaseType( NewId, pGeometry, pProperties, pMasterGeometry )
171  {
172  }
173 
176  {
177  }
178 
181 
185 
186 
190 
198  Condition::Pointer Create(
199  IndexType NewId,
200  NodesArrayType const& rThisNodes,
201  PropertiesPointerType pProperties
202  ) const override;
203 
211  Condition::Pointer Create(
212  IndexType NewId,
213  GeometryPointerType pGeom,
214  PropertiesPointerType pProperties
215  ) const override;
216 
225  Condition::Pointer Create(
226  IndexType NewId,
227  GeometryPointerType pGeom,
228  PropertiesPointerType pProperties,
229  GeometryPointerType pMasterGeom
230  ) const override;
231 
232  /********************************************************************************/
233  /**************** METHODS TO CALCULATE MORTAR CONDITION MATRICES ****************/
234  /********************************************************************************/
235 
245  PairedCondition* pCondition,
246  Vector& rLocalRHS,
247  const MortarConditionMatrices& rMortarConditionMatrices,
248  const DerivativeDataType& rDerivativeData,
249  const IndexType rActiveInactive,
250  const ProcessInfo& rCurrentProcessInfo
251  );
252 
253  /******************************************************************/
254  /********** AUXILLIARY METHODS FOR GENERAL CALCULATIONS ***********/
255  /******************************************************************/
256 
262  void EquationIdVector(
263  EquationIdVectorType& rResult,
264  const ProcessInfo& rCurrentProcessInfo
265  ) const override;
266 
272  void GetDofList(
273  DofsVectorType& rConditionalDofList,
274  const ProcessInfo& rCurrentProcessInfo
275  ) const override;
276 
284  int Check(const ProcessInfo& rCurrentProcessInfo) const override;
285 
289 
293 
297 
299  std::string Info() const override
300  {
301  std::stringstream buffer;
302  buffer << "AugmentedLagrangianMethodFrictionlessMortarContactCondition #" << this->Id();
303  return buffer.str();
304  }
305 
307  void PrintInfo(std::ostream& rOStream) const override
308  {
309  rOStream << "AugmentedLagrangianMethodFrictionlessMortarContactCondition #" << this->Id();
310  }
311 
313  void PrintData(std::ostream& rOStream) const override
314  {
315  PrintInfo(rOStream);
316  this->GetParentGeometry().PrintData(rOStream);
317  this->GetPairedGeometry().PrintData(rOStream);
318  }
319 
323 
325 
326 protected:
329 
333 
337 
341 
342  /********************************************************************************/
343  /**************** METHODS TO CALCULATE MORTAR CONDITION MATRICES ****************/
344  /********************************************************************************/
345 
354  Matrix& rLocalLHS,
355  const MortarConditionMatrices& rMortarConditionMatrices,
356  const DerivativeDataType& rDerivativeData,
357  const IndexType rActiveInactive,
358  const ProcessInfo& rCurrentProcessInfo
359  ) override;
360 
368  void CalculateLocalRHS(
369  Vector& rLocalRHS,
370  const MortarConditionMatrices& rMortarConditionMatrices,
371  const DerivativeDataType& rDerivativeData,
372  const IndexType rActiveInactive,
373  const ProcessInfo& rCurrentProcessInfo
374  ) override;
375 
376  /******************************************************************/
377  /********** AUXILLIARY METHODS FOR GENERAL CALCULATIONS ***********/
378  /******************************************************************/
379 
385  IndexType GetActiveInactiveValue(const GeometryType& CurrentGeometry) const override
386  {
387  IndexType value = 0;
388  for (IndexType i_node = 0; i_node < TNumNodes; ++i_node)
389  if (CurrentGeometry[i_node].Is(ACTIVE) == true)
390  value += 1 << i_node;
391 
392  return value;
393  }
394 
398 
402 
406 
408 private:
411 
415 
419 
423 
427 
431 
435 
436  // Serialization
437 
438  friend class Serializer;
439 
440  void save(Serializer& rSerializer) const override
441  {
443  }
444 
445  void load(Serializer& rSerializer) override
446  {
447  KRATOS_SERIALIZE_LOAD_BASE_CLASS( rSerializer, BaseType );
448  }
449 
451 
452 }; // Class AugmentedLagrangianMethodFrictionlessMortarContactCondition
453 
455 
458 
462 
464 
465 }// namespace Kratos.
AugmentedLagrangianMethodFrictionlessMortarContactCondition.
Definition: ALM_frictionless_mortar_contact_condition.h:63
AugmentedLagrangianMethodFrictionlessMortarContactCondition(IndexType NewId, GeometryPointerType pGeometry)
Definition: ALM_frictionless_mortar_contact_condition.h:148
AugmentedLagrangianMethodFrictionlessMortarContactCondition()
Default constructor.
Definition: ALM_frictionless_mortar_contact_condition.h:142
typename ConditionBaseType::GeometryType::Pointer GeometryPointerType
Pointer type for the geometry used in the condition.
Definition: ALM_frictionless_mortar_contact_condition.h:93
AugmentedLagrangianMethodFrictionlessMortarContactCondition(AugmentedLagrangianMethodFrictionlessMortarContactCondition const &rOther)
Copy constructor.
Definition: ALM_frictionless_mortar_contact_condition.h:175
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: ALM_frictionless_mortar_contact_condition.h:307
std::vector< array_1d< PointType, TDim > > ConditionArrayListType
Type for the array list of conditions with points.
Definition: ALM_frictionless_mortar_contact_condition.h:120
typename std::conditional< TDim==2, LineType, TriangleType >::type DecompositionType
Type definition for the decomposition based on dimension.
Definition: ALM_frictionless_mortar_contact_condition.h:129
typename BaseType::MortarConditionMatrices MortarConditionMatrices
Type for the matrices used in the mortar contact condition.
Definition: ALM_frictionless_mortar_contact_condition.h:75
AugmentedLagrangianMethodFrictionlessMortarContactCondition(IndexType NewId, GeometryPointerType pGeometry, PropertiesPointerType pProperties)
Definition: ALM_frictionless_mortar_contact_condition.h:156
Condition::Pointer Create(IndexType NewId, GeometryPointerType pGeom, PropertiesPointerType pProperties, GeometryPointerType pMasterGeom) const override
Creates a new element pointer from an existing geometry.
typename ConditionBaseType::IndexType IndexType
Type for the index used in the condition.
Definition: ALM_frictionless_mortar_contact_condition.h:90
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: ALM_frictionless_mortar_contact_condition.h:313
std::string Info() const override
Turn back information as a string.
Definition: ALM_frictionless_mortar_contact_condition.h:299
KRATOS_CLASS_INTRUSIVE_POINTER_DEFINITION(AugmentedLagrangianMethodFrictionlessMortarContactCondition)
Counted pointer of AugmentedLagrangianMethodFrictionlessMortarContactCondition.
AugmentedLagrangianMethodFrictionlessMortarContactCondition(IndexType NewId, GeometryPointerType pGeometry, PropertiesPointerType pProperties, GeometryType::Pointer pMasterGeometry)
Definition: ALM_frictionless_mortar_contact_condition.h:165
typename GeometryType::IntegrationPointsArrayType IntegrationPointsType
Type for the integration points in the geometry.
Definition: ALM_frictionless_mortar_contact_condition.h:117
typename PropertiesType::Pointer PropertiesPointerType
Pointer type for the properties used in the condition.
Definition: ALM_frictionless_mortar_contact_condition.h:102
static void StaticCalculateLocalRHS(PairedCondition *pCondition, Vector &rLocalRHS, const MortarConditionMatrices &rMortarConditionMatrices, const DerivativeDataType &rDerivativeData, const IndexType rActiveInactive, const ProcessInfo &rCurrentProcessInfo)
Calculates the local contibution of the RHS.
IndexType GetActiveInactiveValue(const GeometryType &CurrentGeometry) const override
Returns a value depending of the active/inactive set.
Definition: ALM_frictionless_mortar_contact_condition.h:385
void CalculateLocalLHS(Matrix &rLocalLHS, const MortarConditionMatrices &rMortarConditionMatrices, const DerivativeDataType &rDerivativeData, const IndexType rActiveInactive, const ProcessInfo &rCurrentProcessInfo) override
Calculates the local contibution of the LHS.
Base class for all Conditions.
Definition: condition.h:59
std::vector< std::size_t > EquationIdVectorType
Definition: condition.h:98
std::vector< DofType::Pointer > DofsVectorType
Definition: condition.h:100
This data will be used to compute the derivatives.
Definition: mortar_classes.h:638
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
std::vector< IntegrationPointType > IntegrationPointsArrayType
Definition: geometry.h:161
An two node 2D line geometry with linear shape functions.
Definition: line_2d_2.h:65
MortarContactCondition.
Definition: mortar_contact_condition.h:78
This is a base class for the conditions paired.
Definition: paired_condition.h:53
Point class.
Definition: point.h:59
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
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
A three node 3D triangle geometry with linear shape functions.
Definition: triangle_3d_3.h:77
#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
Definition: binary_expression.cpp:25
Vector VectorType
Definition: geometrical_transformation_utilities.h:56
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
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
type
Definition: generate_gid_list_file.py:35
def load(f)
Definition: ode_solve.py:307