ViscoelasticPM: unconsistent computation of contact stiffness and uncoherent effective restitution coefficient with a fixed particle

Dear all,

Noé Rebourcier is starting his PhD in our lab and he found out that the Ip2_ViscElMat_ViscElMat_ViscElPhys associated to the Law2_ScGeom_ViscElPhys_Basic has some unconsistencies. This ip2 has different behaviors when prescribing either (young, poisson, en, es) of the ViscElMat material, or directly (kn, ks, cn, cs). In the first case, the effective contact parameters calculated for two identitical particles of a given material leads (I would say logically) to a contact property identical to the particle's one. On the contrary, if you prescribe an identical parameter for the two particles when prescribing (kn, ks, cn, cs), for example kn1 = kn2 = kn in the ViscElMat, the parameters of the contact will be half of the prescribed one: Kn = kn1 kn2/(kn1+kn2) = kn/2. I wrote a small note that I attach to the question if this is not clear enough. descriptionViscoElastic.pdf

I believe that this is not intuitive and that it should be modified in the law, by writting Kn = 2 kn1 kn2/(kn1+kn2). The consequence for anyone using the contact law with (kn, ks, cn, cs) and not anticipating the change will be to have them doubled in their setup. Are you ok if I change this in the source code ?

Related to my first question, we found out with Noé a problem in the Law2_ScGeom_ViscElPhys_Basic when considering an impact between a moving and a fixed particle. In this case, when considering the same material for the free particle and the fixed one, the effective restitution coefficient observed is different from the one prescribed. There is actually a problem in Ip2_ViscElMat_ViscElMat_ViscElPhys in the reduced mass calculation when one particle is fixed. When one particle is fixed, it is as if its mass is infinite. Therefore, in order to recover the prescribed restitution coefficient, it would be necessary to set the mass of the fixed particle to infinity, so that the reduced mass: mred = m1m2/(m1 + m2) ~ m1m2/m1 = m2. This is not done in the contact law as the dynamic or fixed character of the particles are never checked. Two options are possible:

  1. That each user defines a different material for the fixed particles in the simulations and set effectively an almost infinite density for the fixed particles.
  2. To put the reduced mass to m2 in the ip2 if particle 1 is fixed, and conversely.

I believe that the second option is more intuitive, and I would like to know if you would agree and if it is something that I could technically do in the contact law.

Raphael

Edited by Raphael Maurin