yade-users team mailing list archive

[kalogeropoulos <question678136@xxxxxxxxxxxxxxxxxxxxx>]

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

Hi to all, 

I have written a code for the calculation of the cutting force, and the two others during a cutting test,

My question is : I am not sure i have put the friction angle between rock-cutter right.

As you will see i defined the two materials, and i put there the friction angle, one for sphere-sphere, and the other for facet-sphere. Is it right or i have to define it inside Ip2 or anywhere else??

Help!!!!!

### JCFPM script

# -*- coding: utf-8 -*-
# -*- encoding=utf-8 -*-
from __future__ import division
from yade import pack, plot, geom, pack,utils
from yade import _packPredicates


################# SIMULATIONS DEFINED HERE

#### packing (previously constructed)
OUT='LinearCutTest_JCFPM'

#### Simulation Control
#rate=0.01#-0.01 #deformation rate
#iterMax=800000 # maximum number of iterations Doesnt Need anymore
saveVTK=4000 # saving output files for paraview

#### Material microproperties
intR=1.25 # allows near neighbour interaction (can be adjusted for every packing)
DENS=5290.26 #2674.3 toy petrwmatos## could be adapted to match material density: dens_DEM=dens_rock*(V_rock/V_particles)=dens_rock*1/(1-poro_DEM) -> packing porosity as to be computed?
YOUNG=42e9
FRICT=43.94
ALPHA=2.5
TENS=11e6
COH=1100e6

#### Cutter Properties
x_cutter = 6.35e-3     # Width along x-axis
y_cutter = 1e-3     # Height along y-axis
z_cutter = 5e-3  # Length along z-axis
DOC = 5e-3          # Depth of Cut


#### material definition
def sphereMat(): return JCFpmMat(type=1,density=DENS,young=YOUNG,poisson=ALPHA,frictionAngle=radians(FRICT),tensileStrength=TENS,cohesion=COH)
sample = O.materials.append(sphereMat())
#### create the specimen
X = 54.85e-3
Y=54.85e-3
Z=146.5e-3
pred= pack.inAlignedBox((0,0,0),(X,Y,Z))
sps=SpherePack()
sp=pack.randomDensePack(pred,spheresInCell=10000,radius=1.3725e-3,rRelFuzz=.34,memoizeDb='/tmp/07012019.sqlite',returnSpherePack=True)
sp.toSimulation(material=sample) # up to you to define another sample here, e.g., with randomDensePack or anything else.


#### create the cutter
O.materials.append(JCFpmMat(young=600e9,poisson=0.225,frictionAngle=radians(30),density=1470000,label='bx'))
bx = geom.facetBox(center=(27.425e-3,Y+5e-3,-5e-3),extents=(x_cutter,-(5e-3+DOC),z_cutter),orientation=Quaternion((1,0,0),pi/36),wallMask=8+32)
O.bodies.append(bx)

for facet in bx:
    facet.state.blockedDOFs='xyzXYZ' 
    facet.state.vel=(0.0,0.,5.0)        # velocity along x-axis
      #block e.w na metakineitai mono kata

R=0
Rmax=0
nbSpheres=0.
for o in O.bodies:
 if isinstance(o.shape,Sphere):
   nbSpheres+=1
   R+=o.shape.radius
   if o.shape.radius>Rmax:
     Rmax=o.shape.radius
Rmean=R/nbSpheres


################# ENGINES DEFINED HERE
## damping = 0.4

O.engines=[
 ForceResetter(),
        InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=intR,label='Saabb'),Bo1_Facet_Aabb()],verletDist=.05*1.3725e-3),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=intR,label='SSgeom'),Ig2_Facet_Sphere_ScGeom()],
  [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1,label='interactionPhys')],
  [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=True,Key=OUT,label='interactionLaw')]
 ),
 GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=10,timestepSafetyCoefficient=0.5,defaultDt=utils.PWaveTimeStep()),
 NewtonIntegrator(damping=0.6,label='newton'),
 PyRunner(iterPeriod=int(100),initRun=True,command='recorder()',label='data'),
     #   VTKRecorder(iterPeriod=int(saveVTK),initRun=True,fileName=OUT+'-',recorders=['spheres','jcfpm','cracks','facets'],Key=OUT,label='vtk'),
 PyRunner(iterPeriod=int(100), command='Stop()',label='data')
]


################# RECORDER DEFINED HERE
def recorder():
   global Fz, Fx, Fy
   global dz
   Fz = 0.0
   Fz = abs(sum(O.forces.f(facet.id)[2] for facet  in bx))
   Fx = abs(sum(O.forces.f(facet.id)[0] for facet  in bx))
   Fy = abs(sum(O.forces.f(facet.id)[1] for facet  in bx))
   yade.plot.addData({'i':O.iter,
        'Fx':Fx,
        'Fy':Fy,
        'Fz':Fz,
        'dz': abs(bx[2].state.pos[2]),
        'tc':interactionLaw.nbTensCracks,
        'sc':interactionLaw.nbShearCracks,
        'te':interactionLaw.totalTensCracksE,
        'se':interactionLaw.totalShearCracksE,
        'unbF':utils.unbalancedForce()})
   plot.saveDataTxt(OUT)
    
def Stop():
    if bx[2].state.pos[2] >= Z+1e-3:
        O.pause()

# if you want to plot during simulation
plot.plots={'dz':('Fz')}
plot.plot()

################# PREPROCESSING

#### manage interaction detection factor during the first timestep and then set default interaction range ((cf. A DEM model for soft and hard rock, Scholtes & Donze, JMPS 2013))
O.step();
### initializes the interaction detection factor
SSgeom.interactionDetectionFactor=-1.
Saabb.aabbEnlargeFactor=-1.

#### coordination number verification
numSSlinks=0
numCohesivelinks=0
for i in O.interactions:
    if isinstance(O.bodies[i.id1].shape,Sphere) and isinstance(O.bodies[i.id2].shape,Sphere):
      numSSlinks+=1
    if i.phys.isCohesive :
      numCohesivelinks+=1

print "nblinks=", numSSlinks, " | nbCohesivelinks=", numCohesivelinks, " || Kcohesive=", 2.0*numCohesivelinks/nbSpheres
print 40*'-'
################# SIMULATION REALLY STARTS HERE
#O.run()

Thanks a Lot!!

A.K

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