yade-users team mailing list archive

[Jan Stránský <question234129@xxxxxxxxxxxxxxxxxxxxx>]

Question #234129 on Yade changed:
https://answers.launchpad.net/yade/+question/234129

    Status: Open => Answered

Jan Stránský proposed the following answer:
Hello Lingran,

actually it is very lucky that they are different :-) In c++ code, the
values are transfered from previous time step to the current one, while in
your python script, your values are always the current one. For example, Ft
= i.phys.shearForce is NOT the same as the "equivalent" in c++ ( Vector3r&
shearForce = geom->rotate(phys->shearForce); ). even neglecting the
rotation. shearForce in c++ is one time step earlier than in the python.
So, for one interaction, the correct approach should be:

Ft=i.phys.shearForce
O.step() # now you use the same value of Ft as in c++
du=i.geom.shearInc
kt=i.phys.ks
Ftt=Ft-kt*du
maxFs=Fn.norm()*i.phys.tangensOfFrictionAngle
if(Ftt.norm() > maxFs):
  ratio = maxFs / Ftt.norm()
  trialForce=Ftt
  Ftt *= ratio
  e=((1/kt)*(trialForce-Ftt)).dot(Ftt)
print Ftt,i.phys.shearForce # now both should be the same

For more interactions, you can store shearForces of the interactions in a
list and after O.step() use it in for loop over interaction.

A few commets to your script (not very important):
- you can access O.energy['plastDissip'] directly (without need of
O.energy.items() )
- during the simulation, you can store data to plot.data (with plot.addData
function, e.g.
   plot.addData(i=O.iter,t=O.time,e=E,p=O.energy['plastDissip'])
and then save it with plot.saveDataTxt function. For non-alphabetical
ordered of variables, use O.saveDataTxt(vars=('i','t','e','p'))
- I would use O.step() insteat of O.run(1,True), but it  really does not
matter :-)

cheers
Jan



2013/8/16 lingran <question234129@xxxxxxxxxxxxxxxxxxxxx>

> Question #234129 on Yade changed:
> https://answers.launchpad.net/yade/+question/234129
>
> lingran gave more information on the question:
> sorry for the mistakes. Here is the new version of my code:
>
> aa=dict(O.energy.items())
> E=aa['plastDissip']
> f=open('plastic energy dissipation.txt','w')
>
> for k in range (1,4000):
>         O.run(1,True)
>         aa=dict(O.energy.items())
>         for i in O.interactions:
>                 Ft=i.phys.shearForce
>                 du=i.geom.shearInc
>                 kt=i.phys.ks
>                 Fn=i.phys.normalForce
>                 Ftt=Ft-kt*du
>                 maxFs=Fn.norm()*i.phys.tangensOfFrictionAngle
>
>                 if(Ftt.norm() > maxFs):
>                         ratio = maxFs / Ftt.norm()
>                         trialForce=Ftt
>                         Ftt *= ratio
>                         e=((1/kt)*(trialForce-Ftt)).dot(Ftt)
>                         E=E+e
>
>         f.write('%d %f %f %f\n' %(O.iter,O.time,E,aa['plastDissip']))
>
> f.close()
>
> --
> 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     : yade-users@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~yade-users
> More help   : https://help.launchpad.net/ListHelp
>

-- 
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.