← Back to team overview

yade-users team mailing list archive

[Question #258651]: mistakes in calculating shearInc

 

New question #258651 on Yade:
https://answers.launchpad.net/yade/+question/258651

Hi all,

I have problems in well understanding the code inside 'ScGeom.cpp' when calculating shearInc,
the results given by my python code is quite different from what given by yade simulations.
Does anyone know what's wrong in my python code? Here it is:


import numpy as np
i=O.interactions[237,259] ## choose a random contact for instance
n=np.hstack(i.geom.normal) ## transform a vector into an array to simply the programming process
p=i.geom.penetrationDepth
ra=O.bodies[i.id1].shape.radius
rb=O.bodies[i.id2].shape.radius
vel1=np.hstack(O.bodies[i.id1].state.vel)
angVel1=np.hstack(O.bodies[i.id1].state.angVel)
vel2=np.hstack(O.bodies[i.id2].state.vel)
angVel2=np.hstack(O.bodies[i.id2].state.angVel)

## c++ code given by ScGeom.cpp
#alpha = (radius1+radius2)/(radius1+radius2-penetrationDepth);
#relativeVelocity = (rbp2->vel-rbp1->vel)*alpha + rbp2->angVel.cross(-radius2*normal) - rbp1->angVel.cross(radius1*normal);
#relativeVelocity = relativeVelocity-normal.dot(relativeVelocity)*normal;
#shearInc = relativeVelocity*scene->dt;

alpha=(ra+rb)/(ra+rb-p)
relVel=(vel2-vel1)*alpha+np.cross(rb*n,angVel2)-np.cross(ra*n,angVel1) ##getIncidentVel
relVel=relVel-np.dot(relVel,n)*n
shearInc=relVel*O.dt

shearInc     ### =array([ -4.46043319e-07,  -8.65634198e-07,  -1.11027027e-07])
i.geom.shearInc ## =Vector3(3.662833531906604e-07,-1.6118379736082308e-06,2.7062622606841303e-09)
O.step()
i.geom.shearInc ## =Vector3(3.0053927015226938e-07,-1.5987564318912085e-06,-8.2929936709194397e-09)

Thanks in advance!
Best regards.
Lingran 


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