VoxelPorosity

Asked by Luis Barbosa on 2016-09-07

Hi,

I'm using utils.voxelPorosity to calc the porosity of a specific volume of a pack (regularHexa(pack.inSphere()))
But for any value of V (specified volume), the value of porosity is always 1.0

Does anyone know if this function is working well?

Thanks
Luis

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Luis Barbosa
Solved:
2016-09-24
Last query:
2016-09-24
Last reply:
2016-09-13
Christian Jakob (jakob-ifgt) said : #1

Hi,

For voxelPorosity you need to specify a start point/vector and an end point/vector, not a volume.
It will divide the cuboid range into res*res*res cells and counts cells inside and outside particles, where res is the resolution (default: 200). The function returns the ratio, which is porosity.

https://yade-dem.org/doc/yade.utils.html?highlight=voxelporo#yade._utils.voxelPorosity

Hope it helps,

Christian

Luis Barbosa (luis-pires-b) said : #2

Hi Jakob,

That is exactly what I am doing, specifying resolution, start and end point, for instance 200,(0,0,0),(1,1,1).

Thanks,
Luis

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

Hi,

Just to be sure: before calling voxelPorosity(), did you actually include in your YADE simulation the list of spheres (SpherePack) returned by regularHexa() ? Using something like sp.toSimulation() ?

(A minimal working example may help here)

Luis Barbosa (luis-pires-b) said : #4

Hi, well I'm using O.bodies.append(), as follows:

###################################Aggregate Formation#######################################
rad,gap=.0002,0 #particle radius and gap
r=random.uniform(rmin, rmax) #aggregate radius random according to the class
ag = pack.regularHexa(pack.inSphere((0,0,r),r),radius=rad,gap=gap)
a = len(ag)
c = 0.4*a #percentage of maximum particle removal
n=random.randint(1, int(c)) #number of spheres to exclude randomically

#randomly remove spheres from pack
for i in range(0,len(ag)-n):
 b = ag[random.randint(0,len(ag)-1)]
 O.bodies.append(b)
 ag.remove(b)

########################################Interactions#########################################
O.engines=[
  ForceResetter(),

 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
  [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1)],
  [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(smoothJoint=False)]
 ),

  TranslationEngine(ids=p2,translationAxis=[0,0,-1],velocity=0.03),
  GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.5),
  NewtonIntegrator(damping=0.1,gravity=[0,0,-9.81]),

]

##########################MODULE FOR POROSITY MEASUREMENT############################

g=r/2
h=3*r/2
VP=utils.voxelPorosity(200,(-g,g,g),(g,-g,h))

Did you make basic checks like: how many spheres remain after removal? where are they located? what is their radius?
Most likely the problem is in your code, not in regularHexa() nor voxelPorosity().
Bruno

p.s. Note that you exemple is not self-contained, hence nearly useless.
rmin, rmax are not defined. See https://yade-dem.org/wiki/Howtoask

Luis Barbosa (luis-pires-b) said : #6

Hi Bruno, thanks.

Yes, I've made the checks, here follow my complete script:

#!/usr/bin/python
# -*- coding: utf-8 -*-

#REPRODUCAO COM DADOS EXPERIMENTAIS

##################################Diameter by Class##########################################
#Class1 = 2.38 to 4.76
#Class2 = 4.76 to 6.65
#Class3 = 6.65 to 9.52
#Class4 = 9.52 to 19.1
##################################Radius by Class##########################################
aggregateclass = 1

if aggregateclass==1:
 rmin=0.00119
 rmax=0.00238
 print "class 1"

if aggregateclass==2:
 rmin=0.00238
 rmax=0.00325
 print "class 2"

if aggregateclass==3:
 rmin=0.00325
 rmax=0.00476
 print "class 3"

if aggregateclass==4:
 rmin=0.00476
 rmax=0.00955
 print "class 4"
##################################Import Modulus#############################################
from yade import plot
from yade import pack
from yade import utils
from yade import pack
import math
import random

####################################Material#################################################
O.materials.append(JCFpmMat(type=1,young=70e9,poisson=0.3,frictionAngle=radians(30),density=2500,tensileStrength=5e6,cohesion=5e6,jointNormalStiffness=1e7,jointShearStiffness=1e7,jointCohesion=5e6,jointFrictionAngle=radians(20),jointDilationAngle=0.0,label='spheres'))

###################################Aggregate Formation#######################################
rad,gap=.0002,0 #particle radius and gap
r=random.uniform(rmin, rmax) #aggregate radius random according to the class
ag = pack.regularHexa(pack.inSphere((0,0,r),r),radius=rad,gap=gap)
a = len(ag)
c = 0.25*a #percentage of maximum particle removal
d = 0.135*a #percentage of minimum particle removal
n=random.randint(int(d), int(c)) #number of spheres to exclude randomically

#randomly append spheres from packing
for i in range(0,len(ag)-n):
 b = ag[random.randint(0,len(ag)-1)]
 O.bodies.append(b)
 ag.remove(b)

print "total de particulas", a
print "numero de esferas removidas", n
print "raio do agregado m", r
######################################Planes#################################################
p1=O.bodies.append(utils.geom.facetBox((0,0,0),(0.01,-0.01,0),wallMask=1))
p2=O.bodies.append(utils.geom.facetBox((0,0,r),(0.01,-0.01,r),wallMask=32))# os r's se somam
########################################Interactions#########################################
O.engines=[
  ForceResetter(),

 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
  [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1)],
  [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(smoothJoint=False)]
 ),

  TranslationEngine(ids=p2,translationAxis=[0,0,-1],velocity=0.03),
  GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.5),
  NewtonIntegrator(damping=0.1,gravity=[0,0,-9.81]),
  ForceRecorder(ids=p2,file='rep1.txt',iterPeriod=500),

]

##########################Aggregate Mass############################
MA=utils.getSpheresMass()
print "Mass Aggregate kg", MA
##########################MODULE FOR POROSITY MEASUREMENT############################
VE=utils.getSpheresVolume()
VT=(4*math.pi*r*r*r)/3
P=(VT-VE)/VT
print "Porosidade", P

g=r/2
h=3*r/2
VP=utils.voxelPorosity(200,(-g,g,g),(g,-g,h))
print "VP", VP

Christian Jakob (jakob-ifgt) said : #7

I investigated your script a little bit and tried to use another range in voxelPorosity() function:

VP=utils.voxelPorosity(200,(-g,-g,-g),(g,g,g))

So now I get

VP 0.757457375

which means the function is working.

I am not sure what g and h means in your script, but it seems you just did a wrong input. Also I dont know if voxelPorosity() is a good choice in this case since you have a spherical shaped packing ...

Christian

Luis Barbosa (luis-pires-b) said : #8

Thanks Christian,

I'm using this function because I extract the porosity of a cube inside the sphere, but I'm not sure if the function works when used in this way.

Thanks a lot,

Luis