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.
displacement_contact_criteria.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
23 
24 namespace Kratos
25 {
28 
31 
35 
39 
43 
46 
54 template< class TSparseSpace,
55  class TDenseSpace >
57  : public ConvergenceCriteria< TSparseSpace, TDenseSpace >
58 {
59 public:
60 
63 
66 
68  KRATOS_DEFINE_LOCAL_FLAG( PRINTING_OUTPUT );
69  KRATOS_DEFINE_LOCAL_FLAG( TABLE_IS_INITIALIZED );
70  KRATOS_DEFINE_LOCAL_FLAG( ROTATION_DOF_IS_CONSIDERED );
71 
74 
77 
80 
83 
86 
88  using SparseSpaceType = TSparseSpace;
89 
91  using TablePrinterPointerType = TableStreamUtility::Pointer;
92 
94  using IndexType = std::size_t;
95 
99 
104  : BaseType()
105  {
106  }
107 
113  : BaseType()
114  {
115  // Validate and assign defaults
116  ThisParameters = this->ValidateAndAssignParameters(ThisParameters, this->GetDefaultParameters());
117  this->AssignSettings(ThisParameters);
118  }
119 
130  const double DispRatioTolerance,
131  const double DispAbsTolerance,
132  const double RotRatioTolerance,
133  const double RotAbsTolerance,
134  const bool PrintingOutput = false
135  )
136  : BaseType()
137  {
138  // Set local flags
139  mOptions.Set(DisplacementContactCriteria::PRINTING_OUTPUT, PrintingOutput);
140  mOptions.Set(DisplacementContactCriteria::TABLE_IS_INITIALIZED, false);
141  mOptions.Set(DisplacementContactCriteria::ROTATION_DOF_IS_CONSIDERED, false);
142 
143  // The displacement solution
144  mDispRatioTolerance = DispRatioTolerance;
145  mDispAbsTolerance = DispAbsTolerance;
146 
147  // The rotation solution
148  mRotRatioTolerance = RotRatioTolerance;
149  mRotAbsTolerance = RotAbsTolerance;
150  }
151 
152  // Copy constructor.
154  :BaseType(rOther)
155  ,mOptions(rOther.mOptions)
156  ,mDispRatioTolerance(rOther.mDispRatioTolerance)
157  ,mDispAbsTolerance(rOther.mDispAbsTolerance)
158  ,mRotRatioTolerance(rOther.mRotRatioTolerance)
159  ,mRotAbsTolerance(rOther.mRotAbsTolerance)
160  {
161  }
162 
164  ~DisplacementContactCriteria() override = default;
165 
169 
173 
178  typename BaseType::Pointer Create(Parameters ThisParameters) const override
179  {
180  return Kratos::make_shared<ClassType>(ThisParameters);
181  }
182 
193  ModelPart& rModelPart,
194  DofsArrayType& rDofSet,
195  const TSystemMatrixType& rA,
196  const TSystemVectorType& rDx,
197  const TSystemVectorType& rb
198  ) override
199  {
200  if (SparseSpaceType::Size(rDx) != 0) { //if we are solving for something
201  // Initialize
202  double disp_solution_norm = 0.0, disp_increase_norm = 0.0;
203  IndexType disp_dof_num(0);
204  double rot_solution_norm = 0.0, rot_increase_norm = 0.0;
205  IndexType rot_dof_num(0);
206 
207  // Auxiliary values
208  struct AuxValues {
209  std::size_t dof_id = 0;
210  double dof_value = 0.0, dof_incr = 0.0;
211  };
212 
213  // Auxiliary displacement DoF check
214  const std::function<bool(const VariableData&)> check_without_rot =
215  [](const VariableData& rCurrVar) -> bool {return true;};
216  const std::function<bool(const VariableData&)> check_with_rot =
217  [](const VariableData& rCurrVar) -> bool {return ((rCurrVar == DISPLACEMENT_X) || (rCurrVar == DISPLACEMENT_Y) || (rCurrVar == DISPLACEMENT_Z));};
218  const auto* p_check_disp = (mOptions.Is(DisplacementContactCriteria::ROTATION_DOF_IS_CONSIDERED)) ? &check_with_rot : &check_without_rot;
219 
220  // Loop over Dofs
222  std::tie(disp_solution_norm,disp_increase_norm,disp_dof_num,rot_solution_norm,rot_increase_norm,rot_dof_num) = block_for_each<SixReduction>(rDofSet, AuxValues(), [p_check_disp,&rDx](Dof<double>& rDof, AuxValues& aux_values) {
223  if (rDof.IsFree()) {
224  aux_values.dof_id = rDof.EquationId();
225  aux_values.dof_value = rDof.GetSolutionStepValue(0);
226  aux_values.dof_incr = rDx[aux_values.dof_id];
227 
228  const auto& r_curr_var = rDof.GetVariable();
229  if ((*p_check_disp)(r_curr_var)) {
230  return std::make_tuple(std::pow(aux_values.dof_value, 2),std::pow(aux_values.dof_incr, 2),1,0.0,0.0,0);
231  } else {
232  KRATOS_DEBUG_ERROR_IF_NOT((r_curr_var == ROTATION_X) || (r_curr_var == ROTATION_Y) || (r_curr_var == ROTATION_Z)) << "Variable must be a ROTATION and it is: " << r_curr_var.Name() << std::endl;
233  return std::make_tuple(0.0,0.0,0,std::pow(aux_values.dof_value, 2),std::pow(aux_values.dof_incr, 2),1);
234  }
235  }
236  return std::make_tuple(0.0,0.0,0,0.0,0.0,0);
237  });
238 
239  if(disp_increase_norm == 0.0) disp_increase_norm = 1.0;
240  if(disp_solution_norm == 0.0) disp_solution_norm = 1.0;
241  if(rot_increase_norm == 0.0) rot_increase_norm = 1.0;
242  if(rot_solution_norm == 0.0) rot_solution_norm = 1.0;
243 
244  const double disp_ratio = std::sqrt(disp_increase_norm/disp_solution_norm);
245  const double disp_abs = std::sqrt(disp_increase_norm)/ static_cast<double>(disp_dof_num);
246  const double rot_ratio = std::sqrt(rot_increase_norm/rot_solution_norm);
247  const double rot_abs = std::sqrt(rot_increase_norm)/ static_cast<double>(rot_dof_num);
248 
249  // The process info of the model part
250  ProcessInfo& r_process_info = rModelPart.GetProcessInfo();
251 
252  // We print the results // TODO: Replace for the new log
253  if (rModelPart.GetCommunicator().MyPID() == 0 && this->GetEchoLevel() > 0) {
254  if (r_process_info.Has(TABLE_UTILITY)) {
255  std::cout.precision(4);
256  TablePrinterPointerType p_table = r_process_info[TABLE_UTILITY];
257  auto& Table = p_table->GetTable();
258  if (mOptions.Is(DisplacementContactCriteria::ROTATION_DOF_IS_CONSIDERED)) {
259  Table << disp_ratio << mDispRatioTolerance << disp_abs << mDispAbsTolerance << rot_ratio << mRotRatioTolerance << rot_abs << mRotAbsTolerance;
260  } else {
261  Table << disp_ratio << mDispRatioTolerance << disp_abs << mDispAbsTolerance;
262  }
263  } else {
264  std::cout.precision(4);
265  if (mOptions.IsNot(DisplacementContactCriteria::PRINTING_OUTPUT)) {
266  KRATOS_INFO("DisplacementContactCriteria") << BOLDFONT("DoF ONVERGENCE CHECK") << "\tSTEP: " << r_process_info[STEP] << "\tNL ITERATION: " << r_process_info[NL_ITERATION_NUMBER] << std::endl;
267  KRATOS_INFO("DisplacementContactCriteria") << BOLDFONT("\tDISPLACEMENT: RATIO = ") << disp_ratio << BOLDFONT(" EXP.RATIO = ") << mDispRatioTolerance << BOLDFONT(" ABS = ") << disp_abs << BOLDFONT(" EXP.ABS = ") << mDispAbsTolerance << std::endl;
268  if (mOptions.Is(DisplacementContactCriteria::ROTATION_DOF_IS_CONSIDERED)) {
269  KRATOS_INFO("DisplacementContactCriteria") << BOLDFONT("\tROTATION: RATIO = ") << rot_ratio << BOLDFONT(" EXP.RATIO = ") << mRotRatioTolerance << BOLDFONT(" ABS = ") << rot_abs << BOLDFONT(" EXP.ABS = ") << mRotAbsTolerance << std::endl;
270  }
271  } else {
272  KRATOS_INFO("DisplacementContactCriteria") << "DoF ONVERGENCE CHECK" << "\tSTEP: " << r_process_info[STEP] << "\tNL ITERATION: " << r_process_info[NL_ITERATION_NUMBER] << std::endl;
273  KRATOS_INFO("DisplacementContactCriteria") << "\tDISPLACEMENT: RATIO = " << disp_ratio << " EXP.RATIO = " << mDispRatioTolerance << " ABS = " << disp_abs << " EXP.ABS = " << mDispAbsTolerance << std::endl;
274  if (mOptions.Is(DisplacementContactCriteria::ROTATION_DOF_IS_CONSIDERED)) {
275  KRATOS_INFO("DisplacementContactCriteria") << "\tROTATION: RATIO = " << rot_ratio << " EXP.RATIO = " << mRotRatioTolerance << " ABS = " << rot_abs << " EXP.ABS = " << mRotAbsTolerance << std::endl;
276  }
277  }
278  }
279  }
280 
281  // We check if converged
282  const bool disp_converged = (disp_ratio <= mDispRatioTolerance || disp_abs <= mDispAbsTolerance);
283  const bool rot_converged = mOptions.Is(DisplacementContactCriteria::ROTATION_DOF_IS_CONSIDERED) ? (rot_ratio <= mRotRatioTolerance || rot_abs <= mRotAbsTolerance) : true;
284 
285  if (disp_converged && rot_converged) {
286  if (rModelPart.GetCommunicator().MyPID() == 0 && this->GetEchoLevel() > 0) {
287  if (r_process_info.Has(TABLE_UTILITY)) {
288  TablePrinterPointerType p_table = r_process_info[TABLE_UTILITY];
289  auto& table = p_table->GetTable();
290  if (mOptions.IsNot(DisplacementContactCriteria::PRINTING_OUTPUT))
291  table << BOLDFONT(FGRN(" Achieved"));
292  else
293  table << "Achieved";
294  } else {
295  if (mOptions.IsNot(DisplacementContactCriteria::PRINTING_OUTPUT))
296  KRATOS_INFO("DisplacementContactCriteria") << BOLDFONT("\tDoF") << " convergence is " << BOLDFONT(FGRN("achieved")) << std::endl;
297  else
298  KRATOS_INFO("DisplacementContactCriteria") << "\tDoF convergence is achieved" << std::endl;
299  }
300  }
301  return true;
302  } else {
303  if (rModelPart.GetCommunicator().MyPID() == 0 && this->GetEchoLevel() > 0) {
304  if (r_process_info.Has(TABLE_UTILITY)) {
305  TablePrinterPointerType p_table = r_process_info[TABLE_UTILITY];
306  auto& table = p_table->GetTable();
307  if (mOptions.IsNot(DisplacementContactCriteria::PRINTING_OUTPUT))
308  table << BOLDFONT(FRED(" Not achieved"));
309  else
310  table << "Not achieved";
311  } else {
312  if (mOptions.IsNot(DisplacementContactCriteria::PRINTING_OUTPUT))
313  KRATOS_INFO("DisplacementContactCriteria") << BOLDFONT("\tDoF") << " convergence is " << BOLDFONT(FRED(" not achieved")) << std::endl;
314  else
315  KRATOS_INFO("DisplacementContactCriteria") << "\tDoF convergence is not achieved" << std::endl;
316  }
317  }
318  return false;
319  }
320  }
321  else // In this case all the displacements are imposed!
322  return true;
323  }
324 
329  void Initialize( ModelPart& rModelPart ) override
330  {
331  // Initialize
332  BaseType::mConvergenceCriteriaIsInitialized = true;
333 
334  // Check rotation dof
335  mOptions.Set(DisplacementContactCriteria::ROTATION_DOF_IS_CONSIDERED, ContactUtilities::CheckModelPartHasRotationDoF(rModelPart));
336 
337  // Initialize header
338  ProcessInfo& r_process_info = rModelPart.GetProcessInfo();
339  if (r_process_info.Has(TABLE_UTILITY) && mOptions.IsNot(DisplacementContactCriteria::TABLE_IS_INITIALIZED)) {
340  TablePrinterPointerType p_table = r_process_info[TABLE_UTILITY];
341  auto& r_table = p_table->GetTable();
342  r_table.AddColumn("DP RATIO", 10);
343  r_table.AddColumn("EXP. RAT", 10);
344  r_table.AddColumn("ABS", 10);
345  r_table.AddColumn("EXP. ABS", 10);
346  if (mOptions.Is(DisplacementContactCriteria::ROTATION_DOF_IS_CONSIDERED)) {
347  r_table.AddColumn("RT RATIO", 10);
348  r_table.AddColumn("EXP. RAT", 10);
349  r_table.AddColumn("ABS", 10);
350  r_table.AddColumn("EXP. ABS", 10);
351  }
352  r_table.AddColumn("CONVERGENCE", 15);
353  mOptions.Set(DisplacementContactCriteria::TABLE_IS_INITIALIZED, true);
354  }
355  }
356 
362  {
363  Parameters default_parameters = Parameters(R"(
364  {
365  "name" : "displacement_contact_criteria",
366  "ensure_contact" : false,
367  "print_convergence_criterion" : false,
368  "displacement_relative_tolerance" : 1.0e-4,
369  "displacement_absolute_tolerance" : 1.0e-9,
370  "rotation_relative_tolerance" : 1.0e-4,
371  "rotation_absolute_tolerance" : 1.0e-9
372  })");
373 
374  // Getting base class default parameters
375  const Parameters base_default_parameters = BaseType::GetDefaultParameters();
376  default_parameters.RecursivelyAddMissingParameters(base_default_parameters);
377  return default_parameters;
378  }
379 
384  static std::string Name()
385  {
386  return "displacement_contact_criteria";
387  }
388 
392 
396 
400 
402  std::string Info() const override
403  {
404  return "DisplacementContactCriteria";
405  }
406 
408  void PrintInfo(std::ostream& rOStream) const override
409  {
410  rOStream << Info();
411  }
412 
414  void PrintData(std::ostream& rOStream) const override
415  {
416  rOStream << Info();
417  }
418 
420 protected:
423 
427 
431 
435 
440  void AssignSettings(const Parameters ThisParameters) override
441  {
442  BaseType::AssignSettings(ThisParameters);
443 
444  // The displacement solution
445  mDispRatioTolerance = ThisParameters["displacement_relative_tolerance"].GetDouble();
446  mDispAbsTolerance = ThisParameters["displacement_absolute_tolerance"].GetDouble();
447 
448  // The rotation solution
449  mRotRatioTolerance = ThisParameters["rotation_relative_tolerance"].GetDouble();
450  mRotAbsTolerance = ThisParameters["rotation_absolute_tolerance"].GetDouble();
451 
452  // Set local flags
453  mOptions.Set(DisplacementContactCriteria::PRINTING_OUTPUT, ThisParameters["print_convergence_criterion"].GetBool());
454  mOptions.Set(DisplacementContactCriteria::TABLE_IS_INITIALIZED, false);
455  mOptions.Set(DisplacementContactCriteria::ROTATION_DOF_IS_CONSIDERED, false);
456  }
457 
459 private:
462 
466 
467  Flags mOptions;
468 
469  double mDispRatioTolerance;
470  double mDispAbsTolerance;
471 
472  double mRotRatioTolerance;
473  double mRotAbsTolerance;
474 
476 }; // Kratos DisplacementContactCriteria
477 
480 
482 template<class TSparseSpace, class TDenseSpace>
483 const Kratos::Flags DisplacementContactCriteria<TSparseSpace, TDenseSpace>::PRINTING_OUTPUT(Kratos::Flags::Create(1));
484 template<class TSparseSpace, class TDenseSpace>
485 const Kratos::Flags DisplacementContactCriteria<TSparseSpace, TDenseSpace>::TABLE_IS_INITIALIZED(Kratos::Flags::Create(2));
486 template<class TSparseSpace, class TDenseSpace>
487 const Kratos::Flags DisplacementContactCriteria<TSparseSpace, TDenseSpace>::ROTATION_DOF_IS_CONSIDERED(Kratos::Flags::Create(3));
488 }
std::string Info() const override
Turn back information as a string.
Definition: periodic_interface_process.hpp:93
virtual int MyPID() const
Definition: communicator.cpp:91
static bool CheckModelPartHasRotationDoF(ModelPart &rModelPart)
This method checks that the modelpart has a rotation DoF.
Definition: contact_utilities.cpp:147
This is the base class to define the different convergence criterion considered.
Definition: convergence_criteria.h:58
int GetEchoLevel()
This returns the level of echo for the solving strategy.
Definition: convergence_criteria.h:209
virtual Parameters ValidateAndAssignParameters(Parameters ThisParameters, const Parameters DefaultParameters) const
This method validate and assign default parameters.
Definition: convergence_criteria.h:466
TSparseSpace::MatrixType TSystemMatrixType
Matrix type definition.
Definition: convergence_criteria.h:72
ModelPart::DofsArrayType DofsArrayType
DoF array type definition.
Definition: convergence_criteria.h:81
TSparseSpace::VectorType TSystemVectorType
Vector type definition.
Definition: convergence_criteria.h:74
bool Has(const Variable< TDataType > &rThisVariable) const
Checks if the data container has a value associated with a given variable.
Definition: data_value_container.h:382
Convergence criteria for contact problems.
Definition: displacement_contact_criteria.h:58
Dof represents a degree of freedom (DoF).
Definition: dof.h:86
bool IsFree() const
Definition: dof.h:382
Definition: flags.h:58
void Set(const Flags ThisFlag)
Definition: flags.cpp:33
bool Is(Flags const &rOther) const
Definition: flags.h:274
static Flags Create(IndexType ThisPosition, bool Value=true)
Definition: flags.h:138
bool IsNot(Flags const &rOther) const
Definition: flags.h:291
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
Communicator & GetCommunicator()
Definition: model_part.h:1821
ProcessInfo & GetProcessInfo()
Definition: model_part.h:1746
This class provides to Kratos a data structure for I/O based on the standard of JSON.
Definition: kratos_parameters.h:59
double GetDouble() const
This method returns the double contained in the current Parameter.
Definition: kratos_parameters.cpp:657
void RecursivelyAddMissingParameters(const Parameters &rDefaultParameters)
This function is designed to verify that the parameters under testing contain at least all parameters...
Definition: kratos_parameters.cpp:1457
ProcessInfo holds the current value of different solution parameters.
Definition: process_info.h:59
utility function to do a sum reduction
Definition: reduction_utilities.h:68
This class represents the value of its variable depending to other variable.
Definition: piecewize_linear_table.h:63
static IndexType Size(VectorType const &rV)
return size of vector rV
Definition: ublas_space.h:190
This class is the base of variables and variable's components which contains their common data.
Definition: variable_data.h:49
#define FGRN(x)
Definition: color_utilities.h:27
#define FRED(x)
Definition: color_utilities.h:26
#define BOLDFONT(x)
Definition: color_utilities.h:34
TableStreamUtility::Pointer TablePrinterPointerType
The table stream definition TODO: Replace by logger.
Definition: displacement_contact_criteria.h:91
bool PostCriteria(ModelPart &rModelPart, DofsArrayType &rDofSet, const TSystemMatrixType &rA, const TSystemVectorType &rDx, const TSystemVectorType &rb) override
Compute relative and absolute error.
Definition: displacement_contact_criteria.h:192
void PrintData(std::ostream &rOStream) const override
Print object's data.
Definition: displacement_contact_criteria.h:414
typename BaseType::DofsArrayType DofsArrayType
The dofs array type.
Definition: displacement_contact_criteria.h:79
KRATOS_DEFINE_LOCAL_FLAG(ROTATION_DOF_IS_CONSIDERED)
TSparseSpace SparseSpaceType
The sparse space used.
Definition: displacement_contact_criteria.h:88
KRATOS_CLASS_POINTER_DEFINITION(DisplacementContactCriteria)
Pointer definition of DisplacementContactCriteria.
DisplacementContactCriteria(Kratos::Parameters ThisParameters)
Default constructor. (with parameters)
Definition: displacement_contact_criteria.h:112
KRATOS_DEFINE_LOCAL_FLAG(PRINTING_OUTPUT)
Local Flags.
BaseType::Pointer Create(Parameters ThisParameters) const override
Create method.
Definition: displacement_contact_criteria.h:178
typename BaseType::TSystemMatrixType TSystemMatrixType
The sparse matrix type.
Definition: displacement_contact_criteria.h:82
void PrintInfo(std::ostream &rOStream) const override
Print information about this object.
Definition: displacement_contact_criteria.h:408
std::size_t IndexType
The index type definition.
Definition: displacement_contact_criteria.h:94
KRATOS_DEFINE_LOCAL_FLAG(TABLE_IS_INITIALIZED)
Parameters GetDefaultParameters() const override
This method provides the defaults parameters to avoid conflicts between the different constructors.
Definition: displacement_contact_criteria.h:361
void Initialize(ModelPart &rModelPart) override
This function initialize the convergence criteria.
Definition: displacement_contact_criteria.h:329
static std::string Name()
Returns the name of the class as used in the settings (snake_case format)
Definition: displacement_contact_criteria.h:384
DisplacementContactCriteria()
Default constructor.
Definition: displacement_contact_criteria.h:103
DisplacementContactCriteria(const double DispRatioTolerance, const double DispAbsTolerance, const double RotRatioTolerance, const double RotAbsTolerance, const bool PrintingOutput=false)
Default constructor.
Definition: displacement_contact_criteria.h:129
void AssignSettings(const Parameters ThisParameters) override
This method assigns settings to member variables.
Definition: displacement_contact_criteria.h:440
~DisplacementContactCriteria() override=default
Destructor.
std::string Info() const override
Turn back information as a string.
Definition: displacement_contact_criteria.h:402
typename BaseType::TSystemVectorType TSystemVectorType
The dense vector type.
Definition: displacement_contact_criteria.h:85
DisplacementContactCriteria(DisplacementContactCriteria const &rOther)
Definition: displacement_contact_criteria.h:153
#define KRATOS_INFO(label)
Definition: logger.h:250
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
REACTION_CHECK_STIFFNESS_FACTOR INNER_LOOP_ITERATION DISTANCE_THRESHOLD ACTIVE_CHECK_FACTOR AUXILIAR_COORDINATES NORMAL_GAP WEIGHTED_GAP WEIGHTED_SCALAR_RESIDUAL bool
Definition: contact_structural_mechanics_application_variables.h:93
Definition: reduction_utilities.h:310