Apply a HCSITS (Hard-contact semi-implicit timestepping solver) relaxation to all contacts.
Call this kernel on all bodies to perform on relaxation step of the solver. There exist different friction and relaxation models. See comments on their respective methods for an extended documentation. Use setRelaxationModel() to set it. getDeltaMax() can be used to query the maximum change in the impulses applied wrt to the last iteration.
Use setMaxSubIterations() to set the maximum number of iterations performed to solve optimization problems arising during the relaxation (in InelasticGeneralizedMaximumDissipationContact and InelasticCoulombContactByDecoupling).
Use setCor() to set the coefficient of restitution between 0 and 1 for the InelasticProjectedGaussSeidel model. All other models describe inelastic contacts.
Example of Using the HCSITS kernels in a simulation:
mesa_pd::mpi::ReduceProperty reductionKernel;
mesa_pd::mpi::BroadcastProperty broadcastKernel;
kernel::IntegrateBodiesHCSITS integration;
kernel::InitContactsForHCSITS initContacts(1);
initContacts.setFriction(0,0,
real_t(0.2));
kernel::InitBodiesForHCSITS initBodies;
kernel::HCSITSRelaxationStep relaxationStep;
relaxationStep.setCor(
real_t(0.1));
cs.forEachContact(false, kernel::SelectAll(), ca, initContacts, ca, pa);
ps.forEachParticle(false, kernel::SelectAll(), pa, initBodies, pa, dt);
reductionKernel.operator()<VelocityCorrectionNotification>(*ps);
broadcastKernel.operator()<VelocityUpdateNotification>(*ps);
for(int j = 0; j < 10; j++){
cs.forEachContact(false, kernel::SelectAll(), ca, relaxationStep, ca, pa, dt);
reductionKernel.operator()<VelocityCorrectionNotification>(*ps);
broadcastKernel.operator()<VelocityUpdateNotification>(*ps);
}
ps.forEachParticle(false, kernel::SelectAll(), pa, integration, pa, dt);
|
| | HCSITSRelaxationStep () |
| |
| const size_t & | getMaxSubIterations () const |
| |
| void | setMaxSubIterations (size_t v) |
| |
| const RelaxationModel & | getRelaxationModel () const |
| |
| void | setRelaxationModel (RelaxationModel v) |
| |
| const real_t & | getDeltaMax () const |
| |
| void | setDeltaMax (real_t v) |
| |
| const real_t & | getCor () const |
| |
| void | setCor (real_t v) |
| |
| template<typename CAccessor , typename PAccessor > |
| void | operator() (size_t cid, CAccessor &ca, PAccessor &pa, real_t dt) |
| | Perform relaxation for a single contact. More...
|
| |
| template<typename CAccessor , typename PAccessor > |
| real_t | relaxInelasticFrictionlessContacts (size_t cid, real_t dtinv, CAccessor &ca, PAccessor &pa) |
| | Relaxes The contact with ID cid once. More...
|
| |
| template<typename CAccessor , typename PAccessor > |
| real_t | relaxApproximateInelasticCoulombContactsByDecoupling (size_t cid, real_t dtinv, CAccessor &ca, PAccessor &pa) |
| | Relaxes all contacts once. More...
|
| |
| template<typename CAccessor , typename PAccessor > |
| real_t | relaxInelasticCoulombContactsByDecoupling (size_t cid, real_t dtinv, CAccessor &ca, PAccessor &pa) |
| | Relaxes all contacts once. More...
|
| |
| template<typename CAccessor , typename PAccessor > |
| real_t | relaxInelasticCoulombContactsByOrthogonalProjections (size_t cid, real_t dtinv, bool approximate, CAccessor &ca, PAccessor &pa) |
| | Relaxes all contacts once. More...
|
| |
| template<typename CAccessor , typename PAccessor > |
| real_t | relaxInelasticGeneralizedMaximumDissipationContacts (size_t cid, real_t dtinv, CAccessor &ca, PAccessor &pa) |
| | Relaxes all contacts once. More...
|
| |
| template<typename CAccessor , typename PAccessor > |
| real_t | relaxInelasticContactsByProjectedGaussSeidel (size_t cid, real_t dtinv, CAccessor &ca, PAccessor &pa) |
| | Relaxes all contacts once. More...
|
| |
template<typename CAccessor , typename PAccessor >
| real_t walberla::mesa_pd::kernel::HCSITSRelaxationStep::relaxInelasticCoulombContactsByOrthogonalProjections |
( |
size_t |
cid, |
|
|
real_t |
dtinv, |
|
|
bool |
approximate, |
|
|
CAccessor & |
ca, |
|
|
PAccessor & |
pa |
|
) |
| |
|
inline |
Relaxes all contacts once.
The contact model is for inelastic unilateral contacts with Coulomb friction.
- Parameters
-
| dtinv | The inverse of the current time step. |
| approximate | Use the approximate model showing bouncing. |
- Returns
- The largest variation of contact impulses in the L-infinity norm.
This function is to be called from resolveContacts(). The iterative method to solve the contact problem is e.g. described in the article "A matrix-free cone complementarity approach for
solving large-scale, nonsmooth, rigid body dynamics" by A. Tasora and M. Anitescu in Computer Methods in Applied Mechanics and Engineering (Volume 200, Issues 5–8, 15 January 2011, Pages 439-453). The contact model is purely quadratic and convergence should be good but depends on a parameter. The one-contact problem has a unique solution. The frictional reactions for a dynamic contact converge to those that directly oppose slip. However, the contact is not perfectly inelastic for dynamic contacts but bounces. These vertical motions tend to go to zero for smaller time steps and can be interpreted as exaggerated vertical motions coming from micro asperities (see "Optimization-based simulation of nonsmooth rigid multibody
dynamics" by M. Anitescu in Mathematical Programming (Volume 105, Issue 1, January 2006, Pages 113-143). These motions can be prevented by a small change in the iteration proposed in "The bipotential method: a constructive approach to design the complete contact law with friction and improved numerical algorithms" by G. De Saxce and Z-Q. Feng in Mathematical and Computer Modelling (Volume 28, Issue 4, 1998, Pages 225-245). Which iteration is used is controlled with the approximate parameter.