# Interparticle spacing

Hello,

I would like to compute the spacing between adjacent individual particles (spheres) in a given packing (i.e. interparticle spacing). To do that, I thought about getting the length of the sides from the Delaunay triangulation process in the FlowEngine. So I created a packing of spheres, then run the flow engine to do the triangulation, then I want to get the length of the side. My code is shown below and I am using Yadedaily in running this code. My questions are:

1- Is there a more efficient approach to get the interparticle spacing?

2- I found two functions in Yade documentation [1]: getVertices and printVertices, are those vertices belongs to the Delaunay triangulation (i.e. the vertices of the tetrahedra that forms the mesh for flow modeling)?

3- What are the values that I get from getVertices function? are they the Ids of the cells?

4- printVertices gives a txt file with six columns: id, x, y, z, alpha, fictitious. What does alpha and id mean?

Thank you,

Othman

[1] https:/

-----------

from yade import pack

P=1 #Pa

visc=1e-3 #Pa.sec - taken from Catalano, Chareyre, Bathelemy (2014)

density=1000 #kg/m3

g=9.81 #m/s2

sp_radius=2.36/1000

########## create walls ##########

mnx=0

mny=0

mnz=0

mxx=50.8/1000

mxy=50.8/1000

mxz=50.8/1000

mn,mx=Vector3(

walls=aabbWalls

wallIds=

yade.qt.View()

########## spheres ##########

sp=pack.

sp.makeCloud(

sp.toSimulation()

Height=

dP=P/Height #Pa/m

##check how to fix the particles

#for i in sp:

# body= O.bodies[i]

# body.state.

print ('porosity = ', utils.porosity())

########## flow engine ##########

flow=FlowEngine()

flow.useSolver=3

flow.permeabili

flow.viscosity=visc

flow.bndCondIsP

flow.bndCondVal

flow.boundaryUs

O.dt=1e-6

O.dynDt=False

flow.emulateAct

#get verticies positions

flow.printVerti

n=flow.nCells()

x=[]

for i in range (n):

l=flow.

x.append(l)

print (x)

## Question information

- Language:
- English Edit question

- Status:
- Solved

- For:
- Yade Edit question

- Assignee:
- No assignee Edit question

- Solved by:
- Jan Stránský

- Solved:
- 2021-01-27

- Last query:
- 2021-01-27

- Last reply:
- 2021-01-27

Robert Caulk (rcaulk) said : | #1 |

>1- Is there a more efficient approach to get the interparticle spacing?

Yes. Iterate on interactions and compute distance between particle centers.

>2- I found two functions in Yade documentation [1]: getVertices and printVertices, are those vertices belongs to the Delaunay triangulation (i.e. the vertices of the tetrahedra that forms the mesh for flow modeling)?

Probably yes, but you are linking to the entire doc. Better to link to individual functions like this [1] so we know what you are referencing.

>3- What are the values that I get from getVertices function? are they the Ids of the cells?

getVertices gives you the 4 vertex ids comprising the argument cell.

>4- printVertices gives a txt file with six columns: id, x, y, z, alpha, fictitious. What does alpha and id

id - vertex id

alpha - if the vertex is alpha boundary or not. Special boundary type. If you dont know what this means then you don't need to worry about it.

[1]http://

Othman Sh (othman-sh) said : | #2 |

Hi Robert,

Thanks for your answers. Would you please elaborate on #1 (Iterate on interactions and compute distance between particle centers)? how do I get the ids of adjacent/

Thanks,

Othman

Robert Caulk (rcaulk) said : | #3 |

You would iterate on the interaction container and query the constituent body ids:

separations =[]

for i in O.interactions:

id1,id2 = i.id1,i.id2

x = np.linalg.

separations.

print(separations)

Cheers,

Robert

Othman Sh (othman-sh) said : | #4 |

Thanks Robert. I did the code below but strangely I am not getting any outcome. The separation array is still empty!.

I am copying what I get in the terminal

Welcome to Yade 20201215-

Using python version: 3.8.5 (default, Jul 28 2020, 12:59:40)

[GCC 9.3.0]

TCP python prompt on localhost:9000, auth cookie `sycsuk'

XMLRPC info provider on http://

Running script inter.py

porosity = 0.7337488803402193

number of interactions = 3776

[]

[[ ^L clears screen, ^U kills line. F12 controller, F11 3D view (press "h" in 3D view for help), F10 both, F9 generator, F8 plot. ]]

My code is below.

Thanks

Othman

-------

from yade import pack

import numpy as np

sp_radius=2.36/1000

mnx=0

mny=0

mnz=0

mxx=50.8/1000

mxy=50.8/1000

mxz=50.8/1000

mn,mx=Vector3(

yade.qt.View()

########## spheres ##########

sp=pack.

sp.makeCloud(

sp.toSimulation()

print ('porosity = ', utils.porosity())

O.engines=[

ForceResetter(),

InsertionSortC

Bo1_Sphere_

Bo1_Facet_Aabb()

]),

InteractionLoop(

[

Ig2_

Ig2_

],

[

Ip2_

],

[

Law2_

],

),

NewtonIntegrat

]

O.step() #creating interactions

print('number of interactions = ',len(O.

separations =[]

for i in O.interactions:

id1,id2 = i.id1,i.id2

x = np.linalg.

separations.

print(separations)

Jérôme Duriez (jduriez) said : | #5 |

Maybe your interactions in O.interactions are not real. Non-real interactions are accounted for if you ask for len(O.interacti

Also note that for any Vector3 "vec" in yade terminal, like O.bodies[

Jan Stránský (honzik) said : | #6 |

Hello,

> I would like to compute the spacing between adjacent individual particles ... getting the length of the sides from the Delaunay triangulation process in the FlowEngine.

Isn't FlowEngine an overkill here?

Would "plain" Delaunay triangulation / Voronoi tessellation (using e.g. numpy) be an option?

Or Yade's TesselationWrapper [1]?

cheers

Jan

[1] https:/

Othman Sh (othman-sh) said : | #7 |

Hello all,

As an update about using Delaunay triangulation of the FlowEngine, printVerticies [1] gave me a list of the vertices but they are not grouped for each tetrahedron (i.e Delaunay tetrahedon). Instead, they are listed based on the bodies ids. So I cannot use them to calculate the distances between the adjacent particle centers because I don't know if body #1 is actually the neighbor of body #2.

To solve this issue, I used getCellVolume [2] to compute the volume of the cells which I think it is the volume of the tetrahedrons. Then from each tetrahedron volume, I can get the length of the sides of that tetrahedron (volume= side length^3 / 6*sqrt(2))) which should be equal to the distance between the centers of the particles forming the tetrahedron. My question here is:

- Does "cell" in getCellVolume means the tetrahedron element? I just want to confirm if my approach is correct

- Is there any error in this approach that I am not aware of?

Jerome,

That seems right because I checked O.interactions.

Jan,

That would be a good idea too. However, I couldn't find a function that can give the positions of the Delaunay tetrahedrons corners. There is computeVolumes [3] to compute the volume of Voronoi cells but this will not be useful because I think Voronoi cells are not the same as the Delaunay tetrahedron (right?), so I can't use that in calculating the interparticle spacing.

Thanks,

Othman

[1] https:/

[2] https:/

[3] https:/

Jan Stránský (honzik) said : | #8 |

> I couldn't find a function that can give the positions of the Delaunay tetrahedrons corners

position of Delaunay tetrahedron corners are centers of particles

> Voronoi cells are not the same as the Delaunay tetrahedron (right?)

right, they are dual to each other [4]

A MWE using scipy.spatial.

###

from scipy.spatial import Delaunay

from yade import pack

sp = pack.SpherePack()

sp.makeCloud(

sp.toSimulation()

pts = [b.state.pos for b in O.bodies]

#

t = Delaunay(pts)

for i1,i2,i3,i4 in t.simplices:

for j1,j2 in ((i1,i2)

p1,p2 = pts[j1],pts[j2]

# TODO: remove duplicites?

# TODO: consider radius of particles?

###

cheers

Jan

[4] https:/

[5] https:/

Robert Caulk (rcaulk) said : | #9 |

>- Is there any error in this approach that I am not aware of?

Yes. You are assuming all sides of the tetrahedron are equal.

Othman Sh (othman-sh) said : | #10 |

Thank you Jan and Robert.

Jan, your MWE is great. This solved my problem.

Othman Sh (othman-sh) said : | #11 |

Thanks Jan Stránský, that solved my question.