Energy conservation in YADE

There is probably different possible ways answering your question.. Hoping that others won't be too much different, here is a first one.

As for the time dependancy, I guess you re speaking on timestep, isn't it ? You may give a look to Figure 1.19 of [Duriez2009a] (my PhD thesis, sorry..) on https://tel.archives-ouvertes.fr/tel-00462072/document

No need to understand french, it shows the mechanical (kinetic + elastic) energy of a mass-spring system. This system is simulated using a Finite Difference scheme quite similar to "Yade's one", and energy data are plotted according to timesteps.
As you can see on the Figure, the discrepancy with respect to the constant initial value corresponding to theory (thick "Theorique" line) increases with the timestep.

As for myself, I consider this as related with an approximate nature for computations done within DEM (as with any numerical method), with velocity and acting forces being constant throughout time intervals corresponding to one timestep. Leading to approximated energy balances. (dW = f . v dt is not equivalent to W = ExtremeValueOfF . ExtremeValueOfV . DeltaT)

You may have another point of view, though.. Considering that the discrepancy arises from using values that do not correspond to the same time, with e.g. velocity values being associated to the middle of a time interval (at t = t_i + DeltaT / 2), while positions correspond to t = t_i or t_{i+1}. Which may allow to be consistent with theory.
Other people could probably elaborate this point of view better than me ;-)

As for me, I do not feel so much at ease with this point of view, since I think the half of the timestep (DeltaT / 2) never enters in Yade computations (as in any other DEM code ?)

Thanks Duriez. So according to you the integration scheme is the culprit. If I use a more precise integration scheme(higher order), can I control the energy growth.

Yes, provided you - like me - are reluctant to consider velocities as corresponding to the half of the time steps, it comes directly from the integration scheme, from my point of view.

However, if I understood (and remember) well discussions with Bruno, other (higher order) finite difference schemes would not necessarily improve the situation. But I can not elaborate on this point..

Sir,
I have send the Energy curve to you through mail. Please just have a look. My doubt is, why in my case enegy grows and reaches a constant and not behave in the way shown in your thesis.

Does your system stop to evolve (e.g. the impacting particle is at rest) when the energy levels off ?
If yes, it is a major difference with "my" example, where the mass keeps oscillating.

Hi, the main source of energy sink is usually due to the fact that there is no accurate creation-deletion of contacts. When a new contact is created during a timestep it is only reflected at the next timestep.
It makes more accurate integration schemes useless usually. Even if they are nominally of order 3, 4 or more, the dominant error is due to this failure to precisely acknowledge the contact duration, then everything is of order 2 ultimately (unless adaptative timestepping is used, a nightmare with thousands of particles). Worst that that, some higher order and/or implicit methods have the property to not conserve energy even for simple systems.
If your initial impact is the fastest of all collisions in the problem it explain why the energy change is most noticable at this time.
As always with numerical methods the result will be more accurate if the timestep is decreased.
Bruno