www.mooseframework.org
Public Member Functions | Protected Attributes | List of all members
ContactLineSearchBase Class Reference

This class implements a custom line search for use with mechanical contact. More...

#include <ContactLineSearchBase.h>

Inheritance diagram for ContactLineSearchBase:
[legend]

Public Member Functions

 ContactLineSearchBase (const InputParameters &parameters)
 
void printContactInfo (const std::set< dof_id_type > &contact_set)
 Method for printing the contact information. More...
 
void insertSet (const std::set< dof_id_type > &mech_set)
 Unionize sets from different constraints. More...
 
virtual void reset ()
 Reset the line search data. More...
 

Protected Attributes

std::set< dof_id_type > _current_contact_state
 The current contact set. More...
 
std::set< dof_id_type > _old_contact_state
 The old contact set. More...
 
Real _user_ksp_rtol
 the linear tolerance set by the user in the input file More...
 
bool _user_ksp_rtol_set
 Whether the user linear tolerance has been set yet in this object. More...
 
Real _contact_lambda
 The multiplier of the newton step. More...
 
unsigned _allowed_lambda_cuts
 How many times the linsearch is allowed to cut lambda. More...
 
Real _contact_ltol
 What the linear tolerance should be while the contact state is changing. More...
 
bool _affect_ltol
 Whether to modify the linear tolerance. More...
 

Detailed Description

This class implements a custom line search for use with mechanical contact.

The line search is not fancy. It takes two parameters, set in the MOOSE Executioner block: contact_line_search_ltol and contact_line_search_allowed_lambda_cuts. The allowed_lambda_cuts parameter specifies the number of times the line search is allowed to cut lambda. If allowed to be cut, lambda will be reduced by half, and a new residual will be evaluated. If the residual is smaller with a smaller lambda, then cuts will continue until reaching allowed_lambda_cuts. If the residual is larger with a smaller lambda, then the line search is curtailed and the smaller residual is used. It's recommended that allowed_lambda_cuts be <= 3, with smaller values being used for smaller contact problems. This is to allow necessary residual increases when the transient problem requires significant changes in the contact state.

When the contact set is changing, the user may optionally use a looser linear tolerance set by the contact_line_search_ltol parameter. Then when the contact set is changing during the beginning of the Newton solve, unnecessary computational expense is avoided. Then when the contact set is resolved late in the Newton solve, the linear tolerance will return to the finer tolerance set through the traditional l_tol parameter.

Definition at line 39 of file ContactLineSearchBase.h.

Constructor & Destructor Documentation

◆ ContactLineSearchBase()

ContactLineSearchBase::ContactLineSearchBase ( const InputParameters &  parameters)

Definition at line 35 of file ContactLineSearchBase.C.

36  : LineSearch(parameters),
37  _user_ksp_rtol_set(false),
38  _allowed_lambda_cuts(getParam<unsigned>("allowed_lambda_cuts")),
39  _contact_ltol(getParam<Real>("contact_ltol")),
40  _affect_ltol(getParam<bool>("affect_ltol"))
41 {
42 }
bool _affect_ltol
Whether to modify the linear tolerance.
Real _contact_ltol
What the linear tolerance should be while the contact state is changing.
unsigned _allowed_lambda_cuts
How many times the linsearch is allowed to cut lambda.
bool _user_ksp_rtol_set
Whether the user linear tolerance has been set yet in this object.

Member Function Documentation

◆ insertSet()

void ContactLineSearchBase::insertSet ( const std::set< dof_id_type > &  mech_set)

Unionize sets from different constraints.

Definition at line 54 of file ContactLineSearchBase.C.

55 {
56  if (_current_contact_state.empty())
57  _current_contact_state = mech_set;
58  else
59  for (auto & node : mech_set)
60  _current_contact_state.insert(node);
61 }
std::set< dof_id_type > _current_contact_state
The current contact set.

◆ printContactInfo()

void ContactLineSearchBase::printContactInfo ( const std::set< dof_id_type > &  contact_set)

Method for printing the contact information.

Definition at line 45 of file ContactLineSearchBase.C.

Referenced by PetscContactLineSearch::lineSearch().

46 {
47  if (!contact_set.empty())
48  _console << contact_set.size() << " nodes in contact\n";
49  else
50  _console << "No nodes in contact\n";
51 }

◆ reset()

void ContactLineSearchBase::reset ( )
virtual

Reset the line search data.

Definition at line 64 of file ContactLineSearchBase.C.

65 {
66  _current_contact_state.clear();
67  zeroIts();
68 }
std::set< dof_id_type > _current_contact_state
The current contact set.

Member Data Documentation

◆ _affect_ltol

bool ContactLineSearchBase::_affect_ltol
protected

Whether to modify the linear tolerance.

Definition at line 80 of file ContactLineSearchBase.h.

Referenced by PetscContactLineSearch::lineSearch().

◆ _allowed_lambda_cuts

unsigned ContactLineSearchBase::_allowed_lambda_cuts
protected

How many times the linsearch is allowed to cut lambda.

Definition at line 74 of file ContactLineSearchBase.h.

Referenced by PetscContactLineSearch::lineSearch().

◆ _contact_lambda

Real ContactLineSearchBase::_contact_lambda
protected

The multiplier of the newton step.

Definition at line 71 of file ContactLineSearchBase.h.

Referenced by PetscContactLineSearch::lineSearch().

◆ _contact_ltol

Real ContactLineSearchBase::_contact_ltol
protected

What the linear tolerance should be while the contact state is changing.

Definition at line 77 of file ContactLineSearchBase.h.

Referenced by PetscContactLineSearch::lineSearch().

◆ _current_contact_state

std::set<dof_id_type> ContactLineSearchBase::_current_contact_state
protected

The current contact set.

Definition at line 61 of file ContactLineSearchBase.h.

Referenced by insertSet(), PetscContactLineSearch::lineSearch(), and reset().

◆ _old_contact_state

std::set<dof_id_type> ContactLineSearchBase::_old_contact_state
protected

The old contact set.

Definition at line 63 of file ContactLineSearchBase.h.

Referenced by PetscContactLineSearch::lineSearch().

◆ _user_ksp_rtol

Real ContactLineSearchBase::_user_ksp_rtol
protected

the linear tolerance set by the user in the input file

Definition at line 66 of file ContactLineSearchBase.h.

Referenced by PetscContactLineSearch::lineSearch().

◆ _user_ksp_rtol_set

bool ContactLineSearchBase::_user_ksp_rtol_set
protected

Whether the user linear tolerance has been set yet in this object.

Definition at line 68 of file ContactLineSearchBase.h.

Referenced by PetscContactLineSearch::lineSearch().


The documentation for this class was generated from the following files: