Viscoelastic Force Model

Hi Andy,

2014-05-20 20:11 GMT+02:00 Andy Torres <email address hidden>:
> I have been working in a toy model, trying to figure out how all the viscoelastic forces works, and I think there is a bug here, because there are forces acting even when the particles are not touching each other. This happens if the particles separate after a collision and then meet again.
If penetration depth < 0, then the interaction should be removed, so
no forces at all: line 70-71, 45 [1]

[1] https://github.com/yade/trunk/blob/master/pkg/dem/ViscoelasticPM.cpp#L45

> The simulation I made is just two aligned spheres ( one right above the other, so only normal force is involved) one fixed and the other free to fall under gravity force and bounce on the first fixed sphere.
> In the first bounce, the behavior of the system is the right one, the Penetration Depth (PD) starts at 0, the normal Force (Fn) has a large initial value because of the velocity of the falling particle at time of contact. The interaction doesn't disappear after the bounce, but the code has an "if" to set the force to zero and prevent an attractive force due to a viscous component:
> (line 104-106 of /pkg/dem/ViscoelasticPM.cpp)
> const Real normForceReal = phys.kn * geom.penetrationDepth + phys.cn * normalVelocity;
> if (normForceReal < 0) {
> phys.normalForce = Vector3r::Zero();

Please provide minimal working example, so we can test it. Better with
the proposal, how can it be fixed.

Regards

Anton

Hi Anton,
thanks for the quick anser. In the code, there is an "if"
l.35 if (computeForceTorqueViscEl(_geom, _phys, I, force, torque1, torque2) and (I->isActive)) {...}
l.44 else {scene->interactions->requestErase(I); return;}
should this lines remove the interaction? (if so, seems like its not working)
why YADE keep the interaction when the particles are not in contact?

i just try this code:(with yadedaily Yade 1.07.0-116-c00c469~jessie)

#!/usr/bin/env python
# -*- coding: utf-8 -*-
from ya#!/usr/bin/env python
# -*- coding: utf-8 -*-
from yade import pack, plot
import numpy as np

Corr='-test'

Arch1='energ'+Corr+'.out'
deltaT=5e-5
grav=-9.81
zi=10.0
sumado=Vector3(0.0,0.0,.0)
# "PARTÍCULAS"
idmatVisEl=O.materials.append(ViscElMat(kn=1e9,cn=300e4,ks=1.4e4,cs=200e4,frictionAngle=0.36,label='matVisEl'))

EsfMo = O.bodies.append([utils.sphere(center=(0,0,zi),material='matVisEl', radius=1.0,fixed=False)])
EsFix = O.bodies.append([utils.sphere(center=(0,0,0), material='matVisEl',radius=1.0,fixed=True)])

O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
      [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
      [Law2_ScGeom_ViscElPhys_Basic()]
   ),
 NewtonIntegrator(gravity=(0,0,grav),damping=0.0),
 PyRunner(command='dissip(sumado)',iterPeriod=1,label='potDisip'),#Cuando comienzo a calcular las interacciones
]
def dissip(sumado):
 g = open (Arch1,'a')
 v=O.bodies[0].state.vel[2]
 z=O.bodies[0].state.pos[2]
 m=O.bodies[0].state.mass
 Fn=0.0
 for i in O.interactions:
  PD=i.geom.penetrationDepth
  YFn=i.phys.normalForce
  Fn=YFn[2]
  g.write(str(O.iter)+'\t'+str(z)+'\t'+str(v)+'\t'+'\t'+str(Fn)+'\t'+str(PD)+'\n')
 g.closed
O.dt=deltaT
O.saveTmp()de import pack, plot
import numpy as np

Corr='-test'

Arch1='energ'+Corr+'.out'
deltaT=5e-5
grav=-9.81
zi=10.0
sumado=Vector3(0.0,0.0,.0)
# "PARTÍCULAS"
idmatVisEl=O.materials.append(ViscElMat(kn=1e9,cn=300e4,ks=1.4e4,cs=200e4,frictionAngle=0.36,label='matVisEl'))

EsfMo = O.bodies.append([utils.sphere(center=(0,0,zi),material='matVisEl', radius=1.0,fixed=False)])
EsFix = O.bodies.append([utils.sphere(center=(0,0,0), material='matVisEl',radius=1.0,fixed=True)])

O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
      [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
      [Law2_ScGeom_ViscElPhys_Basic()]
   ),
 NewtonIntegrator(gravity=(0,0,grav),damping=0.0),
 PyRunner(command='dissip(sumado)',iterPeriod=1,label='potDisip'),#Cuando comienzo a calcular las interacciones
]
def dissip(sumado):
 g = open (Arch1,'a')
 v=O.bodies[0].state.vel[2]
 z=O.bodies[0].state.pos[2]
 m=O.bodies[0].state.mass
 Fn=0.0
 for i in O.interactions:
  PD=i.geom.penetrationDepth
  YFn=i.phys.normalForce
  Fn=YFn[2]
  g.write(str(O.iter)+'\t'+str(z)+'\t'+str(v)+'\t'+'\t'+str(Fn)+'\t'+str(PD)+'\n')
 g.closed
O.dt=deltaT
O.saveTmp()

and using gnuplot, i plot:

Kn=5e8
cn=1.5e6
dt=5e-5
m=4188.79
g=9.81
set term wxt 0
#set terminal eps
#set output "NormalForce.eps"
set title "Penetration Depth, normal Velcity and Force Vs time"
set xrange[25000:45000]
set grid
plot \
 "energ-test.out" every 15 u 1:($5) ps 0.5 pt 1 title "PD",\
 "energ-test.out" every 15 u 1:($3/100) ps 0.5 title "normVel/100",\
 "energ-test.out" every 25 u 1:($4/1e8) ps 0.5 pt 1 title "Fn from yade/1e8",\
 "energ-test.out" every 25 u 1:((-$5*Kn+$3*cn)/1e8) ps 0.5 pt 4 title "(calculated Fn=-PD*Kn+normVel*Cn)/1e8"

Regards,
Andy

Hi Andy,

thanks for raising the question. Sometimes it is important to
check the correctness of basic functionality.

See an attached, slightly modified, script and diagram. I put the first
particle (upper) not at the height 10m, but directly on the
border of the second one, adding an initial velocity, which
is almost the same, like you have before the impact. Also I have
added the plot function, so you can see all parameters during simulation.

As you can see, there are no attractive forces, if there is no penetration.
The normal force is 0.

Regards

Anton

2014-05-20 21:26 GMT+02:00 Andy Torres <email address hidden>:
> Question #248996 on Yade changed:
> https://answers.launchpad.net/yade/+question/248996
>
> Andy Torres posted a new comment:
> Hi Anton,
> thanks for the quick anser. In the code, there is an "if"
> l.35 if (computeForceTorqueViscEl(_geom, _phys, I, force, torque1, torque2) and (I->isActive)) {...}
> l.44 else {scene->interactions->requestErase(I); return;}
> should this lines remove the interaction? (if so, seems like its not working)
> why YADE keep the interaction when the particles are not in contact?
>
>
> i just try this code:(with yadedaily Yade 1.07.0-116-c00c469~jessie)
>
> #!/usr/bin/env python
> # -*- coding: utf-8 -*-
> from ya#!/usr/bin/env python
> # -*- coding: utf-8 -*-
> from yade import pack, plot
> import numpy as np
>
> Corr='-test'
>
> Arch1='energ'+Corr+'.out'
> deltaT=5e-5
> grav=-9.81
> zi=10.0
> sumado=Vector3(0.0,0.0,.0)
> # "PARTÍCULAS"
> idmatVisEl=O.materials.append(ViscElMat(kn=1e9,cn=300e4,ks=1.4e4,cs=200e4,frictionAngle=0.36,label='matVisEl'))
>
> EsfMo = O.bodies.append([utils.sphere(center=(0,0,zi),material='matVisEl', radius=1.0,fixed=False)])
> EsFix = O.bodies.append([utils.sphere(center=(0,0,0), material='matVisEl',radius=1.0,fixed=True)])
>
> O.engines=[
> ForceResetter(),
> InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
> InteractionLoop(
> [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
> [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
> [Law2_ScGeom_ViscElPhys_Basic()]
> ),
> NewtonIntegrator(gravity=(0,0,grav),damping=0.0),
> PyRunner(command='dissip(sumado)',iterPeriod=1,label='potDisip'),#Cuando comienzo a calcular las interacciones
> ]
> def dissip(sumado):
> g = open (Arch1,'a')
> v=O.bodies[0].state.vel[2]
> z=O.bodies[0].state.pos[2]
> m=O.bodies[0].state.mass
> Fn=0.0
> for i in O.interactions:
> PD=i.geom.penetrationDepth
> YFn=i.phys.normalForce
> Fn=YFn[2]
> g.write(str(O.iter)+'\t'+str(z)+'\t'+str(v)+'\t'+'\t'+str(Fn)+'\t'+str(PD)+'\n')
> g.closed
> O.dt=deltaT
> O.saveTmp()de import pack, plot
> import numpy as np
>
> Corr='-test'
>
> Arch1='energ'+Corr+'.out'
> deltaT=5e-5
> grav=-9.81
> zi=10.0
> sumado=Vector3(0.0,0.0,.0)
> # "PARTÍCULAS"
> idmatVisEl=O.materials.append(ViscElMat(kn=1e9,cn=300e4,ks=1.4e4,cs=200e4,frictionAngle=0.36,label='matVisEl'))
>
> EsfMo = O.bodies.append([utils.sphere(center=(0,0,zi),material='matVisEl', radius=1.0,fixed=False)])
> EsFix = O.bodies.append([utils.sphere(center=(0,0,0), material='matVisEl',radius=1.0,fixed=True)])
>
> O.engines=[
> ForceResetter(),
> InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
> InteractionLoop(
> [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
> [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
> [Law2_ScGeom_ViscElPhys_Basic()]
> ),
> NewtonIntegrator(gravity=(0,0,grav),damping=0.0),
> PyRunner(command='dissip(sumado)',iterPeriod=1,label='potDisip'),#Cuando comienzo a calcular las interacciones
> ]
> def dissip(sumado):
> g = open (Arch1,'a')
> v=O.bodies[0].state.vel[2]
> z=O.bodies[0].state.pos[2]
> m=O.bodies[0].state.mass
> Fn=0.0
> for i in O.interactions:
> PD=i.geom.penetrationDepth
> YFn=i.phys.normalForce
> Fn=YFn[2]
> g.write(str(O.iter)+'\t'+str(z)+'\t'+str(v)+'\t'+'\t'+str(Fn)+'\t'+str(PD)+'\n')
> g.closed
> O.dt=deltaT
> O.saveTmp()
>
> and using gnuplot, i plot:
>
> Kn=5e8
> cn=1.5e6
> dt=5e-5
> m=4188.79
> g=9.81
> set term wxt 0
> #set terminal eps
> #set output "NormalForce.eps"
> set title "Penetration Depth, normal Velcity and Force Vs time"
> set xrange[25000:45000]
> set grid
> plot \
> "energ-test.out" every 15 u 1:($5) ps 0.5 pt 1 title "PD",\
> "energ-test.out" every 15 u 1:($3/100) ps 0.5 title "normVel/100",\
> "energ-test.out" every 25 u 1:($4/1e8) ps 0.5 pt 1 title "Fn from yade/1e8",\
> "energ-test.out" every 25 u 1:((-$5*Kn+$3*cn)/1e8) ps 0.5 pt 4 title "(calculated Fn=-PD*Kn+normVel*Cn)/1e8"
>
> Regards,
> Andy
>
> --
> You received this question notification because you are a member of
> yade-users, which 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

I forgot, that this Question-interface does not support attachments.

#!/usr/bin/env python
# -*- coding: utf-8 -*-
from yade import pack, plot
import numpy as np

Corr='-test'

Arch1='energ'+Corr+'.out'
deltaT=5e-5
grav=-9.81
zi=10.0
sumado=Vector3(0.0,0.0,.0)
# "PARTÍCULAS"
idmatVisEl=O.materials.append(ViscElMat(kn=1e9,cn=300e4,ks=1.4e4,cs=200e4,frictionAngle=0.36,label='matVisEl'))

EsfMo = O.bodies.append(utils.sphere(center=(0,0,2.0),material='matVisEl', radius=1.0,fixed=False))
EsFix = O.bodies.append(utils.sphere(center=(0,0,0), material='matVisEl',radius=1.0,fixed=True))

O.bodies[EsfMo].state.vel=Vector3(0,0,-12.528367)

O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
      [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
      [Law2_ScGeom_ViscElPhys_Basic()]
   ),
        NewtonIntegrator(gravity=(0,0,grav),damping=0.0),
        PyRunner(command='dissip(sumado)',iterPeriod=1,label='potDisip'),#Cuando comienzo a calcular las interacciones
        PyRunner(command='addPlotData()',iterPeriod=1),
]

from yade import qt

qt.View()

# Function to add data to plot
def addPlotData():
  posZ = O.bodies[EsfMo].state.pos[2]
  velZ = O.bodies[EsfMo].state.vel[2]
  PDcalc = 2.0 - posZ
  PDgeom = 0.0

  try:
    PDgeom = O.interactions[0,1].geom.penetrationDepth
  except:
    pass

  Fn = 0.0
  try:
    Fn = O.interactions[0,1].phys.normalForce[2]
  except:
    pass

  FnCalculated = 5e8*PDgeom - 1.5e6*velZ

  plot.addData(it1=O.iter, it2=O.iter, it3=O.iter, PDcalc = PDcalc, PDgeom = PDgeom, velZ = velZ, Fn = Fn, FnYade = -Fn, FnCalculated = FnCalculated)

plot.plots={'it1':('PDcalc','PDgeom'), 'it2':('velZ'), 'it3':('Fn'), 'PDgeom':('FnCalculated', 'FnYade')}; plot.plot()

def dissip(sumado):
        g = open (Arch1,'a')
        v=O.bodies[0].state.vel[2]
        z=O.bodies[0].state.pos[2]
        m=O.bodies[0].state.mass
        Fn=0.0
        for i in O.interactions:
                PD=i.geom.penetrationDepth
                YFn=i.phys.normalForce
                Fn=YFn[2]
                g.write(str(O.iter)+'\t'+str(z)+'\t'+str(v)+'\t'+'\t'+str(Fn)+'\t'+str(PD)+'\n')
        g.closed
O.dt=deltaT

O.run(200, True)
O.wait
plot.saveGnuplot('sim-data_Sphere')

Hi Anton,
Yes, actually I had prepeared some attachments to make the thing easier, but... haha
ok, I checked your script but with 5e4 time steps, and when I plot PD, Vel/100 and Fn/1e8 ("1/100" and "1/1e8" factors are just to have all in the same scale) I have exactly the same problem.
I mean, the first bounse its all right, like you said, but, in the second one, if you plot all together (PD,Vel/100 and Fn/1e8 Vs iterN) and make some zoom to that (between iteration numbers 14400-14800), you will see that the Fn != 0 even when PD < 0 (meaning that there is forces acting even if ther is no contact, right?). And with this plot, its posible to see that the force is actually acting over the particle, because the velocity also change before the particles make contact.
Its seems like the "if" of line 44 is not removing the interaction like it shoud, right?
please, tell me if I'm not beeing clear enough, i have the files done and ready to send.
again, tanks .
Regards,
Andy

Sorry Andy, I do not get your idea. There is an updated script, which
lets you to analyze all parameters after the second impact e.g. between
interactions 14400-14800.

If you look at the blue line (normal force), it is always 0, when
there is no penetration.
And we never get an attractive force during penetration, what is
restricted by lines
105-106 in ViscoelasticPM.cpp. From my point of view, all is right.

> its posible to see that the force is actually acting over the particle, because the velocity
> also change before the particles make contact.

Sure, you have a gravity there.

Anton

================================================

#!/usr/bin/env python
# -*- coding: utf-8 -*-
from yade import pack, plot
import numpy as np

Corr='-test'

Arch1='energ'+Corr+'.out'
deltaT=5e-5
grav=-9.81
zi=10.0
sumado=Vector3(0.0,0.0,.0)
# "PARTÍCULAS"
idmatVisEl=O.materials.append(ViscElMat(kn=1e9,cn=300e4,ks=1.4e4,cs=200e4,frictionAngle=0.36,label='matVisEl'))

EsfMo = O.bodies.append(utils.sphere(center=(0,0,2.0),material='matVisEl',
radius=1.0,fixed=False))
EsFix = O.bodies.append(utils.sphere(center=(0,0,0),
material='matVisEl',radius=1.0,fixed=True))

O.bodies[EsfMo].state.vel=Vector3(0,0,-12.528367)

O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
      [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
      [Law2_ScGeom_ViscElPhys_Basic()]
   ),
        NewtonIntegrator(gravity=(0,0,grav),damping=0.0),
        PyRunner(command='dissip(sumado)',iterPeriod=1,label='potDisip'),#Cuando
comienzo a calcular las interacciones
        PyRunner(command='addPlotData()',iterPeriod=1,dead=True,
label='PlotAdd'),
]

from yade import qt

qt.View()

# Function to add data to plot
def addPlotData():
  posZ = O.bodies[EsfMo].state.pos[2]
  velZ = O.bodies[EsfMo].state.vel[2]
  PDcalc = 2.0 - posZ
  PDgeom = 0.0

  try:
    PDgeom = O.interactions[0,1].geom.penetrationDepth
  except:
    pass

  Fn = 0.0
  try:
    Fn = O.interactions[0,1].phys.normalForce[2]
  except:
    pass

  FnCalculated = 5e8*PDgeom - 1.5e6*velZ

  plot.addData(it1=O.iter, PDgeom = PDgeom, velZ = velZ/100.0, FnYade
= -Fn/1e8, FnCalculated = FnCalculated/1e8)

plot.plots={'it1':('PDgeom','velZ','FnYade','FnCalculated')}; plot.plot()

def dissip(sumado):
        g = open (Arch1,'a')
        v=O.bodies[0].state.vel[2]
        z=O.bodies[0].state.pos[2]
        m=O.bodies[0].state.mass
        Fn=0.0
        for i in O.interactions:
                PD=i.geom.penetrationDepth
                YFn=i.phys.normalForce
                Fn=YFn[2]

g.write(str(O.iter)+'\t'+str(z)+'\t'+str(v)+'\t'+'\t'+str(Fn)+'\t'+str(PD)+'\n')
        g.closed
O.dt=deltaT

O.run(14400, True)
O.wait
PlotAdd.dead=False
O.run(400, True)
O.wait
plot.saveGnuplot('sim-data_Sphere')

OK! maybe im missing something here, but I think we almost have it!

>If you look at the blue line (normal force), it is always 0, when
>there is no penetration.
is that right? if you check the plot you gave me, Fn start growing from iterN 14520, but the PD is negative until iterN 14600. Shouldn't be start changing all together?

>> its posible to see that the force is actually acting over the particle, because the velocity
>> also change before the particles make contact.
>Sure, you have a gravity there.
yes, true, i was not clear here, the Vel under gravity, should be linear, but in this case is not.

I do not know, what is wrong, here is my data-file.
Before iteration 14580 there is no penetration/force,
they are both 0.0.

Anton

==============================================
==============================================

# FnCalculated FnYade PDgeom it1 velZ
0.0528854137881 -0.0 0.0 14577 -0.0352569425254
0.0528927712881 -0.0 0.0 14578 -0.0352618475254
0.0529001287881 -0.0 0.0 14579 -0.0352667525254
0.0529074862881 -0.0 0.0 14580 -0.0352716575254
0.0527418052399 0.0536958685668 0.000157676455738 14581
-0.0346356153075
0.0526552110261 0.0536076956226 0.000330854532275 14582
-0.0340006255765
0.0525545767687 0.0535052266655 0.000500857660158 14583
-0.0333668589786
0.0524401838845 0.0533887482432 0.000667691955051 14584
-0.0327344827395
...
...
...
0.000843280448325 0.000851162962522 0.00320248714932 14717
0.0101127701989
0.000587246526371 0.000590461193354 0.00315192329833 14718
0.0101149133102
0.000335744255672 0.000334373693617 0.00310134873178 14719
0.0101139996021
8.87675495041e-05 8.28942656184e-05 0.00305077873376 14720
0.0101100840795
-0.000156627052485 -0.0 0.00300022831337 14721 0.0101051790795
-0.000401899029474 -0.0 0.00294970241797 14722 0.0101002740795
-0.000647048381463 -0.0 0.00289920104757 14723 0.0100953690795
-0.000892075108451 -0.0 0.00284872420217 14724 0.0100904640795

Anton

2014-05-21 17:26 GMT+02:00 Andy Torres <email address hidden>:
> Question #248996 on Yade changed:
> https://answers.launchpad.net/yade/+question/248996
>
> Andy Torres posted a new comment:
> OK! maybe im missing something here, but I think we almost have it!
>
>>If you look at the blue line (normal force), it is always 0, when
>>there is no penetration.
> is that right? if you check the plot you gave me, Fn start growing from iterN 14520, but the PD is negative until iterN 14600. Shouldn't be start changing all together?
>
>>> its posible to see that the force is actually acting over the particle, because the velocity
>>> also change before the particles make contact.
>>Sure, you have a gravity there.
> yes, true, i was not clear here, the Vel under gravity, should be linear, but in this case is not.
>
> --
> You received this question notification because you are a member of
> yade-users, which 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

oh! i see, i dont know which is the problem then... you have the right jump i was expecting in Fn, but when i run the same script here, the FnYade and the FnCalculated are equal... (and not 0 like the FnCalculaded-column1 in your data) any idea about this problem? im using yadedaily, and the Yade version is 1.07.0-116-c00c469~jessie

Sorry to keep with this, but i cant figure out why i have different behavior.
should I try compiling version instead daily version? (i had some pyton problems with that)

Here my data.. Fn starts growing at iteration 14521, when PD < 0 when it should start at iter 14524, when PD > 0.
any idea

# FnCalculated FnYade PDgeom it1 velZ
-0.0091216254158 -0.0 -0.0123028173408 14510 -0.0349283075254
-0.00824106022767 -0.0 -0.0121281758032 14511 -0.0349332125254
-0.00736037241453 -0.0 -0.0119535097405 14512 -0.0349381175254
-0.0064795619764 -0.0 -0.0117788191529 14513 -0.0349430225254
-0.00559862891326 -0.0 -0.0116041040403 14514 -0.0349479275254
-0.00471757322513 -0.0 -0.0114293644026 14515 -0.0349528325254
-0.00383639491199 -0.0 -0.01125460024 14516 -0.0349577375254
-0.00295509397385 -0.0 -0.0110798115524 14517 -0.0349626425254
-0.00207367041072 -0.0 -0.0109049983398 14518 -0.0349675475254
-0.00119212422258 -0.0 -0.0107301606021 14519 -0.0349724525254
-0.000310455409448 -0.0 -0.0105552983395 14520 -0.0349773575254
0.000561238031991 0.000563978528688 -0.0103804115519 14521 -0.0349755305276
0.00141727900528 0.00143562629518 -0.0102055338992 14522 -0.034963299001
0.00225769231088 0.00229136148031 -0.0100307174042 14523 -0.0349408528881
0.00308250696873 0.00313121363309 -0.0098560131398 14524 -0.0349083817785
0.00389175613405 0.0039552165132 -0.00968147123091 14525 -0.0348660748591
0.0046854770134 0.00476340800553 -0.00950714085662 14526 -0.0348141208643
0.00546371078104 0.005555830035 -0.00933307025229 14527 -0.0347527080283
0.0062265024956 0.00633252848174 -0.00915930671215 14528 -0.0346820240376
0.00697390101711 0.00709355309654 -0.00898589659196 14529 -0.0346022559846
0.0077059589243 0.00783895741672 -0.00881288531204 14530 -0.034513590323
.
.
.
0.0271560590089 0.0276436593244 -0.000558670074082 14588 -0.0199662729195
0.0271674085979 0.0276552158319 -0.000458838709485 14589 -0.0196410680969
0.027170570422 0.0276584353003 -0.000360633369 14590 -0.0193158248447
0.027165690139 0.0276534660431 -0.000264054244777 14591 -0.0189906409086
0.0271529131986 0.0276404561617 -0.000169101040234 14592 -0.0186656122665
0.0271323848028 0.0276195535053 -7.57729789012e-05 14593 -0.0183408331316

0.0271042498669 0.0275909056311 1.59311867567e-05 14594 -0.0180163959554

0.0270686529812 0.0275546597657 0.000106013166534 14595 -0.0176923914324
0.0270257383742 0.027510962767 0.000194475123695 14596 -0.0173689085038
0.0269756498756 0.0274599610868 0.000281319666215 14597 -0.017046034363
0.0269185308806 0.0274018007347 0.00036654983803 14598 -0.0167238544603

Thanks Anton!
I updated to yadedaily 1.09 version and the FnYADE makes what we all expect.
I think that is important to emphasize that for the 1.07 version (i dont know if you can check this) this kind of energy/force calculations are not going to give the right thing.
Thank you for your time!!
Andy

Thanks Anton Gladky, that solved my question.