Difference in unbalancedForce() between yade and yadedaily

Asked by Sacha Duverger on 2020-02-27

Hello,

While I was performing triaxial tests I noticed that the same script gave different values for unbalancedForce() depending on which version of yade I was using: yade 2018.02b or yadedaily 20191105-2658~35899c7~bionic1. It seems that the problem only occurs with clumps.

Here is a MWE:

##################

from yade import plot
import numpy as np

O.bodies.appendClumped([sphere((0,0,0), 1), sphere((1.5,0,0), 1)])

wallIds = O.bodies.append(aabbWalls(thickness=1e-10))

O.engines=[ForceResetter()
     ,InsertionSortCollider([Bo1_Sphere_Aabb(), Bo1_Box_Aabb()])
     ,InteractionLoop(
         [Ig2_Sphere_Sphere_ScGeom(), Ig2_Box_Sphere_ScGeom6D()],
         [Ip2_FrictMat_FrictMat_FrictPhys()],
         [Law2_ScGeom_FrictPhys_CundallStrack()])
     ,GlobalStiffnessTimeStepper()
     ,NewtonIntegrator(damping=0, gravity=(0,-10,0))
     ,PyRunner(command='yade.plot.addData(uF=unbalancedForce())', iterPeriod=1)
    ]

O.run(1000, wait=True)

print(np.mean(plot.data['uF']))

##################

yade 2018.02b outputs 0.11062286226025234 .
yadedaily 20191105-2658~35899c7~bionic1 outputs 75.5814592275 .

Which version of Yade should I trust?

Thanks in advance

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Jérôme Duriez
Solved:
2020-03-03
Last query:
2020-03-03
Last reply:
2020-03-02
Best Jérôme Duriez (jduriez) said : #1

What do you observe for O.interactions.has(0,1) ?

True for yade 2018.02b and False for the more recent one ?

If yes, that is another illustration/confirmation of the fact that more recent versions should hopefully be more correct.

Clump members are not supposed to interact with each other (at least in modern design [*] but I do not think this has ever been otherwise ?) and yade 2018.02b looks wrong in this sense.

[*] https://gitlab.com/yade-dev/trunk/-/blob/master/pkg/common/Collider.cpp#L27

Chareyre (bruno-chareyre-9) said : #2

My dear Sacha, thanks a lot for revealing two issues in a row. Do you know
that you can report them on gitlab directly (gitlab issue = bug report)?
I hope that one will be as fast as previous one but I have doubts. If you
have candidate scripts that could help keeping clumps constantly working
they are welcome, we apparently miss a few of them. You can check
«checkScripts» in sources for examples.
Bruno

Le jeu. 27 févr. 2020 15:37, Sacha Duverger <
<email address hidden>> a écrit :

> New question #689025 on Yade:
> https://answers.launchpad.net/yade/+question/689025
>
> Hello,
>
> While I was performing triaxial tests I noticed that the same script gave
> different values for unbalancedForce() depending on which version of yade I
> was using: yade 2018.02b or yadedaily 20191105-2658~35899c7~bionic1. It
> seems that the problem only occurs with clumps.
>
> Here is a MWE:
>
> ##################
>
> from yade import plot
> import numpy as np
>
> O.bodies.appendClumped([sphere((0,0,0), 1), sphere((1.5,0,0), 1)])
>
> wallIds = O.bodies.append(aabbWalls(thickness=1e-10))
>
> O.engines=[ForceResetter()
> ,InsertionSortCollider([Bo1_Sphere_Aabb(),
> Bo1_Box_Aabb()])
> ,InteractionLoop(
> [Ig2_Sphere_Sphere_ScGeom(),
> Ig2_Box_Sphere_ScGeom6D()],
> [Ip2_FrictMat_FrictMat_FrictPhys()],
> [Law2_ScGeom_FrictPhys_CundallStrack()])
> ,GlobalStiffnessTimeStepper()
> ,NewtonIntegrator(damping=0, gravity=(0,-10,0))
>
> ,PyRunner(command='yade.plot.addData(uF=unbalancedForce())', iterPeriod=1)
> ]
>
> O.run(1000, wait=True)
>
> print(np.mean(plot.data['uF']))
>
> ##################
>
> yade 2018.02b outputs 0.11062286226025234 .
> yadedaily 20191105-2658~35899c7~bionic1 outputs 75.5814592275 .
>
> Which version of Yade should I trust?
>
> Thanks in advance
>
> --
> You received this question notification because your team yade-users is
> an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>
>
>

Sacha Duverger (schmxprr) said : #3

Yes, O.interactions.has(0,1) returns True with yade 2018.02b and False with yadedaily.
No, I didn't know that anyone could report an issue. I'll think of it next time I face a problem that looks like a bug.

Thank you for your answers, I'll use the most recent version of yade I can from now on.

Sacha Duverger (schmxprr) said : #4

Thanks Jérôme Duriez, that solved my question.

So you both believe there is no issue?...
Beside unbalanced force, are there differences in the trajectories of the clump between 2018.02b and now?
Bruno

Ah... I see.
If you plot the evolution of the unbalanced force it looks weird, with huge values (~1e3) appearing at a fixed frequency.
The issue is that when there is no contact force (which is necessarily the case in the initial position, and also periodically when the bouncing pass through the initial position again), the unbalanced force is virtually infinite (weight/zero).
The huge values you get are thus simply the result of numerical noise, in principle you should get infinity.

No surprise if different (even slightly different) versions give very different results in such case. They can still be both correct (based on your example at least).

I'm unsure this is related to O.interactions.has(0,1), if you think so please elaborate.

Bruno

Sacha Duverger (schmxprr) said : #7

I took off the box and printed the clump's final position in the MWE of the first post (cf new MWE). Yade 2018.02b and yadedaily 20200301-3527~67223fc~bionic1 give the exact same final position, so it looks like the trajectories are the same.

As for O.interactions.has(0,1), I understood with Jérôme's message that it should be "False" since no contact is created among clump members of a same clump ([*]). If it is "True" there is an interaction in O.interactions that is not supposed to be there, the computation of unbalancedForce() could then be wrong.

[*] : https://yade-dem.org/doc/formulation.html#clumps-rigid-aggregates

New MWE:

########################

c_id = O.bodies.appendClumped([sphere((0,0,0), 1), sphere((1.5,0,0), 1)])[0]

O.engines=[ForceResetter()
     ,InsertionSortCollider([Bo1_Sphere_Aabb(), Bo1_Box_Aabb()])
     ,InteractionLoop(
         [Ig2_Sphere_Sphere_ScGeom(), Ig2_Box_Sphere_ScGeom6D()],
         [Ip2_FrictMat_FrictMat_FrictPhys()],
         [Law2_ScGeom_FrictPhys_CundallStrack()])
     ,NewtonIntegrator(damping=0, gravity=(0,-10,0))
    ]
O.dt = 1e-6
O.run(1000, wait=True)

print(O.bodies[c_id].state.pos)

Chareyre (bruno-chareyre-9) said : #8

So, this issue seems to be solved indeed.
Thanks for feedback.
B

On Wed, 4 Mar 2020 at 15:47, Sacha Duverger <
<email address hidden>> wrote:

> Question #689025 on Yade changed:
> https://answers.launchpad.net/yade/+question/689025
>
> Sacha Duverger posted a new comment:
> I took off the box and printed the clump's final position in the MWE of
> the first post (cf new MWE). Yade 2018.02b and yadedaily
> 20200301-3527~67223fc~bionic1 give the exact same final position, so it
> looks like the trajectories are the same.
>
> As for O.interactions.has(0,1), I understood with Jérôme's message that
> it should be "False" since no contact is created among clump members of
> a same clump ([*]). If it is "True" there is an interaction in
> O.interactions that is not supposed to be there, the computation of
> unbalancedForce() could then be wrong.
>
> [*] : https://yade-dem.org/doc/formulation.html#clumps-rigid-aggregates
>
> New MWE:
>
> ########################
>
> c_id = O.bodies.appendClumped([sphere((0,0,0), 1), sphere((1.5,0,0),
> 1)])[0]
>
> O.engines=[ForceResetter()
> ,InsertionSortCollider([Bo1_Sphere_Aabb(),
> Bo1_Box_Aabb()])
> ,InteractionLoop(
> [Ig2_Sphere_Sphere_ScGeom(),
> Ig2_Box_Sphere_ScGeom6D()],
> [Ip2_FrictMat_FrictMat_FrictPhys()],
> [Law2_ScGeom_FrictPhys_CundallStrack()])
> ,NewtonIntegrator(damping=0, gravity=(0,-10,0))
> ]
> O.dt = 1e-6
> O.run(1000, wait=True)
>
> print(O.bodies[c_id].state.pos)
>
> --
> You received this question notification because your team yade-users is
> an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

--
--
_______________
Bruno Chareyre
Associate Professor
ENSE³ - Grenoble INP
Lab. 3SR
BP 53
38041 Grenoble cedex 9
Tél : +33 4 56 52 86 21
________________

Email too brief?
Here's why: email charter
<https://marcuselliott.co.uk/wp-content/uploads/2017/04/emailCharter.jpg>