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

Petsc implementation of the contact line search (based on the Petsc LineSearchShell) More...

#include <PetscContactLineSearch.h>

Inheritance diagram for PetscContactLineSearch:
[legend]

Public Member Functions

 PetscContactLineSearch (const InputParameters &parameters)
 
virtual void lineSearch () override
 
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

PetscNonlinearSolver< Real > * _solver
 
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

Petsc implementation of the contact line search (based on the Petsc LineSearchShell)

Definition at line 37 of file PetscContactLineSearch.h.

Constructor & Destructor Documentation

◆ PetscContactLineSearch()

PetscContactLineSearch::PetscContactLineSearch ( const InputParameters &  parameters)

Definition at line 29 of file PetscContactLineSearch.C.

30  : ContactLineSearchBase(parameters)
31 {
32  _solver = dynamic_cast<PetscNonlinearSolver<Real> *>(
33  _fe_problem.getNonlinearSystem().nonlinearSolver());
34  if (!_solver)
35  mooseError(
36  "This line search operates only with Petsc, so Petsc must be your nonlinear solver.");
37 }
ContactLineSearchBase(const InputParameters &parameters)
PetscNonlinearSolver< Real > * _solver

Member Function Documentation

◆ insertSet()

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

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.

◆ lineSearch()

void PetscContactLineSearch::lineSearch ( )
overridevirtual

Definition at line 40 of file PetscContactLineSearch.C.

41 {
42  PetscBool changed_y = PETSC_FALSE, changed_w = PETSC_FALSE;
43  PetscErrorCode ierr;
44  Vec X, F, Y, W, G, W1;
45  SNESLineSearch line_search;
46  PetscReal fnorm, xnorm, ynorm, gnorm;
47  PetscBool domainerror;
48  PetscReal ksp_rtol, ksp_abstol, ksp_dtol;
49  PetscInt ksp_maxits;
50  KSP ksp;
51  SNES snes = _solver->snes();
52 
53  ierr = SNESGetLineSearch(snes, &line_search);
54  LIBMESH_CHKERR(ierr);
55  ierr = SNESLineSearchGetVecs(line_search, &X, &F, &Y, &W, &G);
56  LIBMESH_CHKERR(ierr);
57  ierr = SNESLineSearchGetNorms(line_search, &xnorm, &fnorm, &ynorm);
58  LIBMESH_CHKERR(ierr);
59  ierr = SNESLineSearchGetSNES(line_search, &snes);
60  LIBMESH_CHKERR(ierr);
61  ierr = SNESLineSearchSetReason(line_search, SNES_LINESEARCH_SUCCEEDED);
62  LIBMESH_CHKERR(ierr);
63  ierr = SNESGetKSP(snes, &ksp);
64  LIBMESH_CHKERR(ierr);
65  ierr = KSPGetTolerances(ksp, &ksp_rtol, &ksp_abstol, &ksp_dtol, &ksp_maxits);
66  LIBMESH_CHKERR(ierr);
67  ierr = VecDuplicate(W, &W1);
68  LIBMESH_CHKERR(ierr);
69 
70  if (!_user_ksp_rtol_set)
71  {
72  _user_ksp_rtol = ksp_rtol;
73  _user_ksp_rtol_set = true;
74  }
75 
76  ++_nl_its;
77 
78  /* precheck */
79  ierr = SNESLineSearchPreCheck(line_search, X, Y, &changed_y);
80  LIBMESH_CHKERR(ierr);
81 
82  /* temporary update */
83  _contact_lambda = 1.;
84  ierr = VecWAXPY(W, -_contact_lambda, Y, X);
85  LIBMESH_CHKERR(ierr);
86 
87  /* compute residual to determine whether contact state has changed since the last non-linear
88  * residual evaluation */
89  _current_contact_state.clear();
90  ierr = SNESComputeFunction(snes, W, F);
91  LIBMESH_CHKERR(ierr);
92  ierr = SNESGetFunctionDomainError(snes, &domainerror);
93  LIBMESH_CHKERR(ierr);
94  if (domainerror)
95  {
96  ierr = SNESLineSearchSetReason(line_search, SNES_LINESEARCH_FAILED_DOMAIN);
97  LIBMESH_CHKERR(ierr);
98  }
99  ierr = VecNorm(F, NORM_2, &fnorm);
100  LIBMESH_CHKERR(ierr);
101  std::set<dof_id_type> contact_state_stored = _current_contact_state;
102  _current_contact_state.clear();
103  printContactInfo(contact_state_stored);
104 
105  if (_affect_ltol)
106  {
107  if (contact_state_stored != _old_contact_state)
108  {
109  KSPSetTolerances(ksp, _contact_ltol, ksp_abstol, ksp_dtol, ksp_maxits);
110  _console << "Contact set changed since previous non-linear iteration!\n";
111  }
112  else
113  KSPSetTolerances(ksp, _user_ksp_rtol, ksp_abstol, ksp_dtol, ksp_maxits);
114  }
115 
116  size_t ls_its = 0;
117  while (ls_its < _allowed_lambda_cuts)
118  {
119  _contact_lambda *= 0.5;
120  /* update */
121  ierr = VecWAXPY(W1, -_contact_lambda, Y, X);
122  LIBMESH_CHKERR(ierr);
123 
124  ierr = SNESComputeFunction(snes, W1, G);
125  LIBMESH_CHKERR(ierr);
126  ierr = SNESGetFunctionDomainError(snes, &domainerror);
127  LIBMESH_CHKERR(ierr);
128  if (domainerror)
129  {
130  ierr = SNESLineSearchSetReason(line_search, SNES_LINESEARCH_FAILED_DOMAIN);
131  LIBMESH_CHKERR(ierr);
132  }
133  ierr = VecNorm(G, NORM_2, &gnorm);
134  LIBMESH_CHKERR(ierr);
135  if (gnorm < fnorm)
136  {
137  VecCopy(G, F);
138  LIBMESH_CHKERR(ierr);
139  VecCopy(W1, W);
140  LIBMESH_CHKERR(ierr);
141  fnorm = gnorm;
142  contact_state_stored.swap(_current_contact_state);
143  _current_contact_state.clear();
144  printContactInfo(contact_state_stored);
145  ++ls_its;
146  }
147  else
148  break;
149  }
150 
151  ierr = VecScale(Y, _contact_lambda);
152  LIBMESH_CHKERR(ierr);
153  /* postcheck */
154  ierr = SNESLineSearchPostCheck(line_search, X, Y, W, &changed_y, &changed_w);
155  LIBMESH_CHKERR(ierr);
156 
157  if (changed_y)
158  {
159  ierr = VecWAXPY(W, -1., Y, X);
160  LIBMESH_CHKERR(ierr);
161  }
162 
163  if (changed_w || changed_y)
164  {
165  ierr = SNESComputeFunction(snes, W, F);
166  LIBMESH_CHKERR(ierr);
167  ierr = SNESGetFunctionDomainError(snes, &domainerror);
168  LIBMESH_CHKERR(ierr);
169  if (domainerror)
170  {
171  ierr = SNESLineSearchSetReason(line_search, SNES_LINESEARCH_FAILED_DOMAIN);
172  LIBMESH_CHKERR(ierr);
173  }
174  contact_state_stored.swap(_current_contact_state);
175  _current_contact_state.clear();
176  printContactInfo(contact_state_stored);
177  }
178 
179  /* copy the solution over */
180  ierr = VecCopy(W, X);
181  LIBMESH_CHKERR(ierr);
182 
183  ierr = SNESLineSearchComputeNorms(line_search);
184  LIBMESH_CHKERR(ierr);
185 
186  ierr = VecDestroy(&W1);
187  LIBMESH_CHKERR(ierr);
188 
189  _old_contact_state = std::move(contact_state_stored);
190 }
bool _affect_ltol
Whether to modify the linear tolerance.
std::set< dof_id_type > _old_contact_state
The old contact set.
Real _contact_ltol
What the linear tolerance should be while the contact state is changing.
std::set< dof_id_type > _current_contact_state
The current contact set.
Real _contact_lambda
The multiplier of the newton step.
unsigned _allowed_lambda_cuts
How many times the linsearch is allowed to cut lambda.
void printContactInfo(const std::set< dof_id_type > &contact_set)
Method for printing the contact information.
Real _user_ksp_rtol
the linear tolerance set by the user in the input file
PetscNonlinearSolver< Real > * _solver
bool _user_ksp_rtol_set
Whether the user linear tolerance has been set yet in this object.

◆ printContactInfo()

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

Method for printing the contact information.

Definition at line 45 of file ContactLineSearchBase.C.

Referenced by 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 ( )
virtualinherited

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
protectedinherited

Whether to modify the linear tolerance.

Definition at line 80 of file ContactLineSearchBase.h.

Referenced by lineSearch().

◆ _allowed_lambda_cuts

unsigned ContactLineSearchBase::_allowed_lambda_cuts
protectedinherited

How many times the linsearch is allowed to cut lambda.

Definition at line 74 of file ContactLineSearchBase.h.

Referenced by lineSearch().

◆ _contact_lambda

Real ContactLineSearchBase::_contact_lambda
protectedinherited

The multiplier of the newton step.

Definition at line 71 of file ContactLineSearchBase.h.

Referenced by lineSearch().

◆ _contact_ltol

Real ContactLineSearchBase::_contact_ltol
protectedinherited

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

Definition at line 77 of file ContactLineSearchBase.h.

Referenced by lineSearch().

◆ _current_contact_state

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

The current contact set.

Definition at line 61 of file ContactLineSearchBase.h.

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

◆ _old_contact_state

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

The old contact set.

Definition at line 63 of file ContactLineSearchBase.h.

Referenced by lineSearch().

◆ _solver

PetscNonlinearSolver<Real>* PetscContactLineSearch::_solver
protected

Definition at line 45 of file PetscContactLineSearch.h.

Referenced by lineSearch(), and PetscContactLineSearch().

◆ _user_ksp_rtol

Real ContactLineSearchBase::_user_ksp_rtol
protectedinherited

the linear tolerance set by the user in the input file

Definition at line 66 of file ContactLineSearchBase.h.

Referenced by lineSearch().

◆ _user_ksp_rtol_set

bool ContactLineSearchBase::_user_ksp_rtol_set
protectedinherited

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

Definition at line 68 of file ContactLineSearchBase.h.

Referenced by lineSearch().


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