yade-users team mailing list archive

[Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>]

> I tried the sample example with the two balls and the time step 
> estimated by PWave. I attach the figure.
Mmmm... It is not "that" better. I expected better matching with 
corrected kinetic (as in the simple mass-spring test).

> Now I suppose that if the kinetic energy is computed correctly (is 
> it?) the problem relies on the elastic energy. You said yourself that 
> the total energy is "approximately" constant.
ret+=.5*(b->state->mass*(b->state->vel+b->state->accel*dt/2).SquaredLength()+(b->state->angVel+b->state->angAccel*dt/2).Dot(diagMult(b->state->inertia,(b->state->angVel+b->state->angAccel*dt/2))));

Ah! I see. At time "t" we have contact forces F(t-dt) and velocities 
V(t-dt/2).
So, you need in fact to __substract__ accel*dt/2, and you will get at 
time "t" the energy of time "t-dt". Sorry if I wrote something different 
in previous mails.
It explains why it match better with non-shifted curves, you shifted the 
wrong way. ;)

Now replying to another question :

> plasticDissipation += 
(shearDisp+(1/currentContactPhysics->ks)*(shearForce-prevForce)).dot(shearForce)
> Is shearDisp the total shear displacement, right? I would compute the 
dissipation incrementally, something like
> plasticDissipation += [delta_us + (Fs_max_current - Fs_max_prev)/ks] 
* Fs_max_current

The equations are almost the same. shearDisp is the increment.
The problem in your equation is it doesn't handle the transition 
plastic->elastic or elastic->plastic (when the increment of displacement 
is partly elastic and partly plastic), it is the reason to use 
shearForce/prevForce instead of max_current/max_previous. If the 
increment is fully plastic, it makes no difference : in such case 
shearForce=max_current and prevForce=max_previous.

Bruno