yade-users team mailing list archive

[Xiaolong ZHAO <question692626@xxxxxxxxxxxxxxxxxxxxx>]

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

Hello, dear professors,
I want to use an HCP packing in a sphere to simulate a single-particle crushing test. I chose the JCF model since it is more suitable for a rock particle. Recently, I check the particle size effect on the agglomerate. I kept all the model parameters the same include the ratio of particle size L to elementary particle size R.  It could be inferred that all the samples should have the same fabric. However, it is not. The coordination number and bonds number are different, which leads to different simulation results. 
An interesting thing was that when L and R were reduced or enlarged by a multiple of 2, the CN or bonds number can be constant, and the same simulation results can be obtained.  I think maybe the error of float number arouses this problem, but I have no idea to handle it.  
Below is the part of the script I used. You can change the "enlarge_factor", if enlarge_factor=1/2,2,4,8,16...., CN is constant, but not for other cases, such as enlarge_factor=3,5,6,7.......
Thank you in advance!
Xiaolong
###############################
from yade import utils,pack,plot,timing
from math import *
import numpy as np
readParamsFromTable(LL=1.1e-2,RR=4e-4)   #LL=3.0e-2,RR=1e-3   LL=1.1e-2,RR=4e-4
enlarge_factor=1
young=3.3e9
poisson=.3
frictionAngle=radians(30)
tensileStrength=22.5e6
cohesion=225e6
gamma_int=1.5
rate=0.02
from yade.params import table


rayon_s=float(table.RR*enlarge_factor)*1.0    #elementary sphere
rayon_g=float(table.LL/2.0*enlarge_factor)    #0.05  #large sphere
 
O.materials.append(JCFpmMat(type=1,young=young,poisson=poisson,frictionAngle=frictionAngle,density=2680,
                   tensileStrength=tensileStrength,cohesion=cohesion,
                   #jointNormalStiffness=1e7,jointShearStiffness=1e7,jointCohesion=1e6,jointFrictionAngle=radians(20),jointDilationAngle=0.0,
                   label='rock'))
sp=pack.regularHexa(yade._packPredicates.inSphere(center=(0,0,0),radius=rayon_g),gap=0,radius=rayon_s,material='rock')

O.bodies.append(sp)
count=0 
for i in O.bodies:
    if isinstance(O.bodies[i.id].shape,Sphere):  
        count=count+1

print('ball number=',count)

Rx=(max([O.bodies[b.id].state.pos[0] for b in O.bodies if isinstance(b.shape,Sphere)])
   -min([O.bodies[b.id].state.pos[0] for b in O.bodies if isinstance(b.shape,Sphere)]))*0.5+rayon_s
Ry=(max([O.bodies[b.id].state.pos[1] for b in O.bodies if isinstance(b.shape,Sphere)])
   -min([O.bodies[b.id].state.pos[1] for b in O.bodies if isinstance(b.shape,Sphere)]))*0.5+rayon_s
Rz=(max([O.bodies[b.id].state.pos[2] for b in O.bodies if isinstance(b.shape,Sphere)])
   -min([O.bodies[b.id].state.pos[2] for b in O.bodies if isinstance(b.shape,Sphere)]))*0.5+rayon_s
R=(Rx+Ry+Rz)/3
print('Rx=',Rx,'Ry=',Ry,'Rz=',Rz,'R=',R)

O.engines=[
           ForceResetter(),
           InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=gamma_int,label='is1aabb'),Bo1_Facet_Aabb(),Bo1_Wall_Aabb(),Bo1_Box_Aabb()]),
           InteractionLoop(
           [Ig2_Sphere_Sphere_ScGeom6D(interactionDetectionFactor=gamma_int,label='iss2sc'),Ig2_Sphere_Sphere_ScGeom(),
            Ig2_Wall_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom6D()],
           [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1),            
            Ip2_FrictMat_FrictMat_FrictPhys()],
           [Law2_ScGeom_FrictPhys_CundallStrack(),
            Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=False,smoothJoint=False)]),
            NewtonIntegrator(damping=0.8),
          ]

O.dt=0
O.step() #for CoheFrictM to release stress within agglomerate; for JCF model to create initial contacts
CN=avgNumInteractions(skipFree=True)
nbCohIntrs=0
for i in O.interactions:
    if isinstance(O.bodies[i.id1].shape,Sphere) and isinstance(O.bodies[i.id2].shape,Sphere):
        if i.phys.isCohesive:
            nbCohIntrs+=1
print ("nbCohIntrs=",nbCohIntrs)
print('CN=',CN)
#now reset the interaction radius and go ahead
is1aabb.aabbEnlargeFactor=1.
iss2sc.interactionDetectionFactor=1.
O.dt=utils.PWaveTimeStep()
#O.run()
utils.waitIfBatch()

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.