www.mooseframework.org
ContactSlipDamper.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
8 #include "ContactSlipDamper.h"
9 #include "FEProblem.h"
10 #include "DisplacedProblem.h"
11 #include "AuxiliarySystem.h"
12 #include "PenetrationLocator.h"
13 #include "NearestNodeLocator.h"
14 
15 template <>
16 InputParameters
18 {
19  InputParameters params = validParams<GeneralDamper>();
20  params.addRequiredParam<std::vector<int>>(
21  "master", "IDs of the master surfaces for which slip reversals should be damped");
22  params.addRequiredParam<std::vector<int>>(
23  "slave", "IDs of the slave surfaces for which slip reversals should be damped");
24  params.addParam<Real>(
25  "max_iterative_slip", std::numeric_limits<Real>::max(), "Maximum iterative slip");
26  params.addRangeCheckedParam<Real>("min_damping_factor",
27  0.0,
28  "min_damping_factor < 1.0",
29  "Minimum permissible value for damping factor");
30  params.addParam<Real>("damping_threshold_factor",
31  1.0e3,
32  "If previous iterations's slip is below "
33  "the slip tolerance, only damp a slip "
34  "reversal if the slip magnitude is "
35  "greater than than this factor times "
36  "the old slip.");
37  params.addParam<bool>("debug_output", false, "Output detailed debugging information");
38  return params;
39 }
40 
41 ContactSlipDamper::ContactSlipDamper(const InputParameters & parameters)
42  : GeneralDamper(parameters),
43  _aux_sys(parameters.get<FEProblemBase *>("_fe_problem_base")->getAuxiliarySystem()),
44  _displaced_problem(parameters.get<FEProblemBase *>("_fe_problem_base")->getDisplacedProblem()),
45  _num_contact_nodes(0),
46  _num_sticking(0),
47  _num_slipping(0),
48  _num_slipping_friction(0),
49  _num_stick_locked(0),
50  _num_slip_reversed(0),
51  _max_iterative_slip(parameters.get<Real>("max_iterative_slip")),
52  _min_damping_factor(parameters.get<Real>("min_damping_factor")),
53  _damping_threshold_factor(parameters.get<Real>("damping_threshold_factor")),
54  _debug_output(parameters.get<bool>("debug_output"))
55 {
56  if (!_displaced_problem)
57  mooseError("Must have displaced problem to use ContactSlipDamper");
58 
59  std::vector<int> master = parameters.get<std::vector<int>>("master");
60  std::vector<int> slave = parameters.get<std::vector<int>>("slave");
61 
62  unsigned int num_interactions = master.size();
63  if (num_interactions != slave.size())
64  mooseError("Sizes of master surface and slave surface lists must match in ContactSlipDamper");
65  if (num_interactions == 0)
66  mooseError("Must define at least one master/slave pair in ContactSlipDamper");
67 
68  for (unsigned int i = 0; i < master.size(); ++i)
69  {
70  std::pair<int, int> ms_pair(master[i], slave[i]);
71  _interactions.insert(ms_pair);
72  }
73 }
74 
75 void
77 {
78  GeometricSearchData & displaced_geom_search_data = _displaced_problem->geomSearchData();
79  std::map<std::pair<unsigned int, unsigned int>, PenetrationLocator *> * penetration_locators =
80  &displaced_geom_search_data._penetration_locators;
81 
82  for (pl_iterator plit = penetration_locators->begin(); plit != penetration_locators->end();
83  ++plit)
84  {
85  PenetrationLocator & pen_loc = *plit->second;
86 
87  if (operateOnThisInteraction(pen_loc))
88  {
89  std::vector<dof_id_type> & slave_nodes = pen_loc._nearest_node._slave_nodes;
90 
91  for (unsigned int i = 0; i < slave_nodes.size(); i++)
92  {
93  dof_id_type slave_node_num = slave_nodes[i];
94 
95  if (pen_loc._penetration_info[slave_node_num])
96  {
97  PenetrationInfo & info = *pen_loc._penetration_info[slave_node_num];
98  const Node * node = info._node;
99 
100  if (node->processor_id() == processor_id())
101  // && info.isCaptured()) //TODO maybe just set this
102  // everywhere?
103  info._slip_reversed = false;
104  }
105  }
106  }
107  }
108 }
109 
110 Real
111 ContactSlipDamper::computeDamping(const NumericVector<Number> & solution,
112  const NumericVector<Number> & /*update*/)
113 {
114  // Do new contact search to update positions of slipped nodes
115  _displaced_problem->updateMesh(solution, *_aux_sys.currentSolution());
116 
117  Real damping = 1.0;
118 
119  _num_contact_nodes = 0;
120  _num_sticking = 0;
121  _num_slipping = 0;
123  _num_stick_locked = 0;
124  _num_slip_reversed = 0;
125 
126  GeometricSearchData & displaced_geom_search_data = _displaced_problem->geomSearchData();
127  std::map<std::pair<unsigned int, unsigned int>, PenetrationLocator *> * penetration_locators =
128  &displaced_geom_search_data._penetration_locators;
129 
130  for (pl_iterator plit = penetration_locators->begin(); plit != penetration_locators->end();
131  ++plit)
132  {
133  PenetrationLocator & pen_loc = *plit->second;
134 
135  if (operateOnThisInteraction(pen_loc))
136  {
137  std::vector<dof_id_type> & slave_nodes = pen_loc._nearest_node._slave_nodes;
138 
139  for (unsigned int i = 0; i < slave_nodes.size(); i++)
140  {
141  dof_id_type slave_node_num = slave_nodes[i];
142 
143  if (pen_loc._penetration_info[slave_node_num])
144  {
145  PenetrationInfo & info = *pen_loc._penetration_info[slave_node_num];
146  const Node * node = info._node;
147 
148  if (node->processor_id() == processor_id())
149  {
150  if (info.isCaptured())
151  {
153  if (info._mech_status == PenetrationInfo::MS_STICKING)
154  _num_sticking++;
155  else if (info._mech_status == PenetrationInfo::MS_SLIPPING)
156  _num_slipping++;
157  else if (info._mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
159  if (info._stick_locked_this_step >= 2) // TODO get from contact interaction
161 
162  RealVectorValue tangential_inc_slip_prev_iter =
163  info._incremental_slip_prev_iter -
164  (info._incremental_slip_prev_iter * info._normal) * info._normal;
165  RealVectorValue tangential_inc_slip =
166  info._incremental_slip - (info._incremental_slip * info._normal) * info._normal;
167 
168  RealVectorValue tangential_it_slip =
169  tangential_inc_slip - tangential_inc_slip_prev_iter;
170  Real node_damping_factor = 1.0;
171  if ((tangential_inc_slip_prev_iter * tangential_inc_slip < 0.0) &&
172  info._mech_status == PenetrationInfo::MS_SLIPPING_FRICTION)
173  {
174  info._slip_reversed = true;
176  Real prev_iter_slip_mag = tangential_inc_slip_prev_iter.norm();
177  RealVectorValue prev_iter_slip_dir =
178  tangential_inc_slip_prev_iter / prev_iter_slip_mag;
179  Real cur_it_slip_in_old_dir = tangential_it_slip * prev_iter_slip_dir;
180 
181  if (prev_iter_slip_mag > info._slip_tol ||
182  cur_it_slip_in_old_dir > -_damping_threshold_factor * prev_iter_slip_mag)
183  node_damping_factor =
184  1.0 - (cur_it_slip_in_old_dir + prev_iter_slip_mag) / cur_it_slip_in_old_dir;
185 
186  if (node_damping_factor < 0.0)
187  mooseError("Damping factor can't be negative");
188 
189  if (node_damping_factor < _min_damping_factor)
190  node_damping_factor = _min_damping_factor;
191  }
192 
193  if (tangential_it_slip.norm() > _max_iterative_slip)
194  node_damping_factor =
195  (tangential_it_slip.norm() - _max_iterative_slip) / tangential_it_slip.norm();
196 
197  if (_debug_output && node_damping_factor < 1.0)
198  _console << "Damping node: " << node->id()
199  << " prev iter slip: " << info._incremental_slip_prev_iter
200  << " curr iter slip: " << info._incremental_slip
201  << " slip_tol: " << info._slip_tol
202  << " damping factor: " << node_damping_factor << "\n";
203 
204  if (node_damping_factor < damping)
205  damping = node_damping_factor;
206  }
207  }
208  }
209  }
210  }
211  }
212  _console << std::flush;
213  _communicator.sum(_num_contact_nodes);
214  _communicator.sum(_num_sticking);
215  _communicator.sum(_num_slipping);
216  _communicator.sum(_num_slipping_friction);
217  _communicator.sum(_num_stick_locked);
218  _communicator.sum(_num_slip_reversed);
219  _communicator.min(damping);
220 
221  _console << " ContactSlipDamper: Damping #Cont #Stick #Slip #SlipFric #StickLock "
222  "#SlipRev\n";
223 
224  _console << std::right << std::setw(29) << damping << std::setw(10) << _num_contact_nodes
225  << std::setw(10) << _num_sticking << std::setw(10) << _num_slipping << std::setw(10)
226  << _num_slipping_friction << std::setw(11) << _num_stick_locked << std::setw(10)
227  << _num_slip_reversed << "\n\n";
228  _console << std::flush;
229 
230  return damping;
231 }
232 
233 bool
234 ContactSlipDamper::operateOnThisInteraction(const PenetrationLocator & pen_loc)
235 {
236  bool operate_on_this_interaction = false;
237  std::set<std::pair<int, int>>::iterator ipit;
238  std::pair<int, int> ms_pair(pen_loc._master_boundary, pen_loc._slave_boundary);
239  ipit = _interactions.find(ms_pair);
240  if (ipit != _interactions.end())
241  operate_on_this_interaction = true;
242  return operate_on_this_interaction;
243 }
virtual Real computeDamping(const NumericVector< Number > &solution, const NumericVector< Number > &update)
Compute the amount of damping.
std::set< std::pair< int, int > > _interactions
virtual void timestepSetup()
MooseSharedPointer< DisplacedProblem > _displaced_problem
InputParameters validParams< ContactSlipDamper >()
std::map< std::pair< unsigned int, unsigned int >, PenetrationLocator * >::iterator pl_iterator
Convenient typedef for frequently used iterator.
ContactSlipDamper(const InputParameters &parameters)
AuxiliarySystem & _aux_sys
bool operateOnThisInteraction(const PenetrationLocator &pen_loc)
Determine whether the damper should operate on the interaction corresponding to the supplied Penetrat...