New question #676609 on Yade:
https://answers.launchpad.net/yade/+question/676609
Hi all,
When I used createInteraction() command to manually create intrs between 2
spheres, I got an error:
RuntimeError: Interaction #39+#18 already exists.
I have then check these 2 spheres #39 and #18:
Sphere #39 has no intrs, len(O.bodies[39].intrs()) = 0, while sphere #18 has 7
intrs but #39 is not in the neighboring list of sphere #18.
Please see the MWE (Minimal Working Example™) generating this error down below.
Please let me know why I have this RuntimeError error?
Best regards,
Thai-Son
PS: I have also tried aabbEnlargeFactor and interactionDetectionFactor but I
could not obtain what I want, that's why I have to try doing this
createInteraction manually. To be more specified, what I want is that all
spheres in the packing has at least 1 intrs.
#-------------- MWE
------------------------------#
# encoding: utf-8
#
------------------------------------------------------------------------------ #
#######################
### INPUT PARAMETERs ###
#######################
seed0=1 # = 1 -> same packing everytime
distributeMass0=False # True -> Mass fraction, False -> PARTICLE COUNT fraction
num_spheres0=200
mx=Vector3(.4e-2,.6e-2,.4e-2) # size of aggregate
mn=Vector3(0,0,0)
porosity0=0.4
rRelFuzz0=0.7 # for the case of polydisperse
compFricDegree = 15.0 # initial contact friction during the confining phase
(will be decreased during the REFD compaction process)
confiningS=-1e5
surfaceTension=0.07274 # [N/m]
young0= 50e9 # 50 GPa
poisson0= 0.2
density0 = 2530 # 2.53 g/cm3 = 2530 kg/m3
tensileStrength = 3.3e9 # 3.3 GPa
shearStrength = 3.3e9 # 3.3 GPa
####################################
## GENERATING PARTICLE PACKING ###
####################################
from yade import pack
import pylab
O.materials.append(CohFrictMat(young=young0,poisson=poisson0,density=density0,frictionAngle=radians(compFricDegree),
normalCohesion=tensileStrength, shearCohesion=shearStrength,
momentRotationLaw=True,etaRoll=0.1,label='spheres'))
O.materials.append(FrictMat(young=young0,poisson=poisson0,frictionAngle=0,density=0,label='frictionlessWalls'))
walls=aabbWalls((mn,mx),thickness=0,material='frictionlessWalls')
wallIds=O.bodies.append(walls)
### Aggregate generation ###
num_spheres1D=int(round(num_spheres0**(1./3))) # for calculation of rMean0
sp=pack.SpherePack()
rMean0=(1-rRelFuzz0) * ((mx[0]/num_spheres1D)/2) # Irregular case, in order to
obtain: num_spheres_real not less than num_spheres
sp.makeCloud(minCorner=mn,maxCorner=mx,rMean=rMean0,rRelFuzz=rRelFuzz0,porosity=porosity0,num=num_spheres0,distributeMass=distributeMass0,seed=seed0)
O.bodies.append([utils.sphere(center,rad,material='spheres') for center,rad in
sp])
### mechanical loading #### To confer a solid-like nature to the packing
triax=TriaxialStressController(
internalCompaction=True,
goal1=confiningS,
goal2=confiningS,
goal3=confiningS,
max_vel=10,
label="triax"
)
newton=NewtonIntegrator(damping=0.4)
unsat=UnsaturatedEngine(dead=1,label="unsat")
# ---------------------------------------------------------------------------
O.engines=[
ForceResetter(),
### Collsion detector
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
### Interaction handling
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom6D(),Ig2_Box_Sphere_ScGeom()],
[Ip2_FrictMat_FrictMat_FrictPhys(),Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(label="cohesiveIp")],
[Law2_ScGeom_FrictPhys_CundallStrack(),Law2_ScGeom6D_CohFrictPhys_CohesionMoment(
useIncrementalForm=True, #useIncrementalForm is turned
on as we want plasticity on the contact moments
always_use_moment_law=False, #if we want "rolling"
friction even if the contact is not cohesive (or cohesion is broken), we will
have to turn this true somewhere
label='cohesiveLaw')]
),
GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=1,timestepSafetyCoefficient=0.5),
triax,
newton
]
# ---------------------------------------------------------------------------
numDEMiter1stDEMloop=3000
numDEMiter2ndDEMloop=2000
numDEMiter3rdDEMloop=500
numDEMiterPerDryIter=250
while 1:
O.run(numDEMiter1stDEMloop,True)
unb=unbalancedForce()
if unb<0.01 and
abs(triax.goal1-triax.meanStress)/abs(triax.goal1)<0.001:
break
## REACH NEW EQU. STATE ###
finalFricDegree = 30 # contact friction during the deviatoric loading
# We move to deviatoric loading, let us turn internal compaction off to keep
particles sizes constant
triax.internalCompaction=False
# Change contact friction (remember that decreasing it would generate
instantaneous instabilities)
setContactFriction(radians(finalFricDegree))
while 1:
O.run(numDEMiter2ndDEMloop,True)
unb=unbalancedForce()
if unb<0.001 and
abs(triax.goal1-triax.meanStress)/abs(triax.goal1)<0.001:
break
triax.depth0=triax.depth
triax.height0=triax.height
triax.width0=triax.width
O.run(numDEMiter3rdDEMloop,True)
ei0=-triax.strain[0];ei1=-triax.strain[1];ei2=-triax.strain[2]
si0=-triax.stress(0)[0];si1=-triax.stress(2)[1];si2=-triax.stress(4)[2]
### counting the real number of spheres generated by YADE
num_spheres_real=0
for b in O.bodies:
if not isinstance(b.shape,Box):
num_spheres_real =num_spheres_real+1
num_spheres=num_spheres_real
### --------------------------------------------------------------------###
### ---- GLUEING THE PARTICLE manually ---- ###
### --------------------------------------------------------------------###
### calculating number of no-intrs bodies to see the effect of manuallly
creating bond between particles
NoIntrsBodies=[]
for i in range(6,len(O.bodies)): # exclude walls, id<6
if len(O.bodies[i].intrs()) == 0 and not i in NoIntrsBodies :
NoIntrsBodies+=[i]
print 'NoIntrsBodies before creating additional bonds: ', len(NoIntrsBodies)
#-----------------------------------------------------------#
# CREATE BONDS BETWEEN PARTICLES MANUALLY
#-----------------------------------------------------------#
# List of no-bond particles
EnlargeFactorManual=150./100.
ParticleWithManualBondCreatedList=[]
for i in NoIntrsBodies:
if not i in ParticleWithManualBondCreatedList:
for j in range(6,len(O.bodies)): # exclude walls, id<6, only
manually create bonds between spheres
if not i==j and not i in
ParticleWithManualBondCreatedList: # skip the sphere i itself and skip i if i
has been already added in this "for" loop
if
(float((O.bodies[i].state.pos-O.bodies[j].state.pos).norm())-O.bodies[j].shape.radius)
/ O.bodies[i].shape.radius < EnlargeFactorManual:
# --- MANUALLY CREATE INTERATION (BOND) BETWEEN
2 SPHERES --- #
createInteraction(i,j)
# -- this
ParticleWithManualBondCreatedList is used for avoiding checking a sphere which
already turned from no-bond to bonded sphere
if not i in
ParticleWithManualBondCreatedList:
ParticleWithManualBondCreatedList+=[i]
if not j in
ParticleWithManualBondCreatedList:
ParticleWithManualBondCreatedList+=[j]
NoIntrsBodiesAfter=[]
for i in range(6,len(O.bodies)): # exclude walls, id<6
if len(O.bodies[i].intrs()) == 0 and not i in NoIntrsBodiesAfter :
NoIntrsBodiesAfter+=[i]
print 'EnlargeFactorManual: ', EnlargeFactorManual
print 'NoIntrsBodies after creating additional bonds: ', len(NoIntrsBodiesAfter)
### --- DEBUGGING section --- # Please use this pack of codes after getting the
error "RuntimeError: Interaction #39+#18 already exists."
#%cpaste
print ''
print "i in NoIntrsBodies >>>", i in NoIntrsBodies
print "j in NoIntrsBodies >>>", j in NoIntrsBodies
print "i in ParticleWithManualBondCreatedList >>>",i in
ParticleWithManualBondCreatedList
print "j in ParticleWithManualBondCreatedList >>>",j in
ParticleWithManualBondCreatedList
print "EnlargeFactorManual: ",EnlargeFactorManual
print "len(O.bodies[i].intrs()):", len(O.bodies[i].intrs())
print "len(O.bodies[j].intrs()):", len(O.bodies[j].intrs())
if len(O.bodies[i].intrs()) >0:
for k in range(len(O.bodies[i].intrs())):
print "O.bodies[i].intrs()[",k,"].dict(): id1=",
O.bodies[i].intrs()[k].id1, " id2=", O.bodies[i].intrs()[k].id2
print " ", O.bodies[i].intrs()[k].dict()
if len(O.bodies[j].intrs()) >0:
for k in range(len(O.bodies[j].intrs())):
print "O.bodies[j].intrs()[",k,"].dict(): id1=",
O.bodies[j].intrs()[k].id1, " id2=", O.bodies[j].intrs()[k].id2
print " ", O.bodies[j].intrs()[k].dict()
--
You received this question notification because your team yade-users is
an answer contact for Yade.
_______________________________________________
Mailing list: https://launchpad.net/~yade-users
Post to : [email protected]
Unsubscribe : https://launchpad.net/~yade-users
More help : https://help.launchpad.net/ListHelp
