VoxelPorosity

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

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

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)

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

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

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

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