Fix issue #184 (aka question #694548) (!628) · Merge requests · yade-dev / trunk

This merge request corresponds to the fix I presented in my e-mail of Tue, 02 Mar 2021 in Yade-dev for the issue #184 (closed) (aka question #694548). This issue concerns the unphysical collision between two PFacet objects.

In short, I believe issue #184 (closed) was caused by a permutation in the p1, p2, p3 weights that are used to calculate the forces on each individual gridNode from the force of the collision.

In more details:

The motion of PFacet is dictated by that of the 3 gridNode objects it connects. Each gridNode recieves a different force value calculated from the force of the interaction. The function Law2_ScGridCoGeom_FrictPhys_CundallStrack::go (pkg/common/Grid.cpp:537) does that using weights:

scene->forces.addForce(geom->id3, geom->weight[0] * -force);
scene->forces.addForce(geom->id4, geom->weight[1] * -force);
scene->forces.addForce(geom->id5, geom->weight[2] * -force);

Those weights are calculated in pkg/common/PFacet.cpp (Ig2_Sphere_PFacet_ScGridCoGeom::projection()) as

     ‖v₁‖²(v₀·v₂) - (v₁·v₂)(v₀·v₁)
p₁ = —————————————————————————————        (1a)
         ‖v₀‖²‖v₁‖² - (v₀·v₁)²
     ‖v₀‖²(v₁·v₂) - (v₀·v₂)(v₀·v₁)
p₂ = —————————————————————————————        (1b)
         ‖v₀‖²‖v₁‖² - (v₀·v₁)²
p₃ = 1 - p₂ - p₃                          (1c)

But actually, it should be:

     ‖v₁‖²(v₀·v₂) - (v₁·v₂)(v₀·v₁)
p₂ = —————————————————————————————        (2a)
         ‖v₀‖²‖v₁‖² - (v₀·v₁)²
     ‖v₀‖²(v₁·v₂) - (v₀·v₂)(v₀·v₁)
p₃ = —————————————————————————————        (2b)
         ‖v₀‖²‖v₁‖² - (v₀·v₁)²
p₁ = 1 - p₂ - p₃                          (2c)

Details of the calculations here. Some more details of the fix are in the e-mail, in particular link to calculations (which will soon be dead: should I make a LaTeX version of it?).

@klausthoeni I hope we have not been colliding…

Fix issue #184 (aka question #694548)

Merge request reports