Normal stiffness for sphere-polyhedra interaction

Asked by Cong Xu on 2019-02-05

Hello,

I have a question on how the normal stiffness between polyhedra and sphere is calculated. I'm using Yade 2018.02b.

I checked the code files:

In Ig2_Sphere_Polyhedra_ScGeom,

In Ip2_FrictMat_PolyhedraMat_FrictPhys,
void Ip2_FrictMat_PolyhedraMat_FrictPhys::go(const shared_ptr<Material>& pp1, const shared_ptr<Material>& pp2, const shared_ptr<Interaction>& interaction){
Ip2_FrictMat_FrictMat_FrictPhys().go(mat1,mat2,interaction);
}

And in Ip2_FrictMat_FrictMat_FrictPhys(),
Real Ea = mat1->young;
Real Eb = mat2->young;
Real Va = mat1->poisson;
Real Vb = mat2->poisson;

//half the harmonic average of the two stiffnesses, when (2*Ri*Ei) is the stiffness of a contact point on sphere "i"
Real Kn = 2*Ea*Ra*Eb*Rb/(Ea*Ra+Eb*Rb);
//same for shear stiffness
Real Ks = 2*Ea*Ra*Va*Eb*Rb*Vb/(Ea*Ra*Va+Eb*Rb*Vb);

I ran an example script, the sphere has radius 0.5, young 1e9, and polyhedra has young 1e10.
The kn between them in the simulation is 5e8, which is not the result by Kn = 2*Ea*Ra*Eb*Rb/(Ea*Ra+Eb*Rb)

Wondering if there is other equation to calculate the kn for the sphere-polyhedra interaction? Thank you!

Best regards,
Cong

The following is the example script.

import random

polyMat = PolyhedraMat(young=1e10,poisson=.05)
frictMat = FrictMat(young=1e9,poisson=.2)

O.materials.append((polyMat,frictMat))

poly = polyhedra_utils.polyhedra(polyMat,(1,2,3)); poly.wire=False
sph1 = sphere((2,2,2),.5,material=frictMat)
sph2 = sphere((-2,0,0),.5,material=frictMat)
sph3 = sphere((0,-2,0),.5,material=frictMat)
O.bodies.append((
poly,
sph1,
sph2,
sph3
))

O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Polyhedra_Aabb(),Bo1_Sphere_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Polyhedra_ScGeom()],
[Ip2_FrictMat_PolyhedraMat_FrictPhys()],
[Law2_ScGeom_FrictPhys_CundallStrack()],
),
NewtonIntegrator(),
]

O.dt = 1e-7
O.step()

poly.state.blockedDOFs = 'xyzXYZ'
for s in (sph1,sph2,sph3):
r=random.random
s.state.vel = -10*(s.state.pos+Vector3(r(),r(),r()))

qt.View()

Question information

Language:
English Edit question
Status:
Expired
For:
Assignee:
No assignee Edit question
Last query:
2019-02-23
2019-03-11
 Jan Stránský (honzik) said on 2019-02-07: #1

Hello,
strange, the result is as if both young values were 1e9..
what version of Yade do you use?
cheers
Jan

 Cong Xu (congxu) said on 2019-02-07: #2

Hello Jan,

Cheers
Cong

 Launchpad Janitor (janitor) said on 2019-02-23: #3

This question was expired because it remained in the 'Needs information' state without activity for the last 15 days.

 Cong Xu (congxu) said on 2019-02-23: #4

Hello,

Any result on this? Thank you!

Cheers,
Cong

 Jérôme Duriez (jduriez) said on 2019-03-01: #5

Hi,

Your script has a random behavior (in terms of an existence of a sphere-polyhedra interaction for instance) because of your use of polyhedra_utils.polyhedra().

I'd advice you show us a really minimal (you also do not need 4 bodies, 2 would be enough) working (always an interaction) script to increase your chances getting help from the (limited in size ?) polyhedra users community.

Jérôme

 Jan Stránský (honzik) said on 2019-03-01: #6

> Any result on this?

no, sorry, currently I do not have time for it.. But evidently it is a bug, so please open a bug report
cheers
Jan

 Launchpad Janitor (janitor) said on 2019-03-11: #7

This question was expired because it remained in the 'Open' state without activity for the last 15 days.