# stl geometry problem

Hello everyone,
I am trying to simulate the rock cutting process with a disk cutter.
During the cut the disk translates and rotates. I have done that by using
Kinematic Engine.
The problem is with the geometry.
When disk rotates, deforms. The facets forming the disk seem to detach.
I have created the geometry (Solid), with 3 different ways (Blender, 3d builder and Autocad) without solving the problem.
Any suggestions on how can I fix it;
I am wondering it might affect forces exerted on the disk.
I am using Ubuntu 18.04 and Yade 2018.02b.
Thank you
I am sending the following: A file of the disk geometry in .stl format, and my code.

https://1drv.ms/u/s!AiKcYaY3wMDqgV9IE7I2RFUWVTxp

from __future__ import print_function
from __future__ import division
from yade import pack, plot, geom, pack,utils,ymport,export,bodiesHandling
import math

################# SIMULATIONS DEFINED HERE
#### Simulation
OUT='LinearCutTest_JCFPM'

#### Simulation Control
saveVTK=10000 # saving output files for paraview

#### Material microproperties
intR=1.25, # allows near neighbour interaction (GammaInt)
DENSITY=2700, # Density of the material (kg/m^3),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=68e9, # Eeq(GPa),elastic modulus
FRICTION_ANGLE=12, #(Degrees),Contact friction angle
POISSON=1,
TENS=22e6, #(Tensile Strenght T (MPa)),Defines the maximum admissible normal force in traction in the matrix (FnMax = tensileStrength * crossSection). (Pa)
COH=220e6, #(Local Cohesion C (MPa)),Defines the maximum admissible tangential force in shear, for Fn=0, in the matrix (FsMax = cohesion * crossSection). (Pa)

### Cutting_Specimen_dimensions
X = 54.85e-3,
Y=54.85e-3,
Z=180e-3,

)

#### material definition
Sample = O.materials.append(JCFpmMat(
young=YOUNG,
poisson=POISSON,
cohesion=COH,
tensileStrength=TENS,
density = DENSITY,
label='spheres'
))

#### create the specimen
specimen=pack.randomDensePack(
pack.inAlignedBox((0,0,0),(X,Y,Z)),
spheresInCell = 1000,
material=Sample,
rRelFuzz=0.34,
memoizeDb = '/tmp/Linear-Cutting_Test.db',
returnSpherePack = True,
)
specimen.toSimulation()

#### Disk Creation
Disk=ymport.stl('./p(5mm).stl')
O.bodies.append(Disk)
ids = [b.id for b in Disk]

################# ENGINES DEFINED HERE

O.engines=[
ForceResetter(),
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.5,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'),
CombinedKinematicEngine(ids=ids,comb=[
TranslationEngine(translationAxis=[0,0,1],velocity=-580),
RotationEngine(angularVelocity=6700,rotationAxis=[1,0,0],rotateAroundZero=False)],label='Eng')
]

################# 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 Disk))
Fx = abs(sum(O.forces.f(facet.id)[0] for facet in Disk))
Fy = abs(sum(O.forces.f(facet.id)[1] for facet in Disk))
'FS':Fx,
'FN':Fy,
'FR':Fz,
'dz': abs(Disk[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 Disk[2].state.pos[2] == 0:
O.pause()

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

#### 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.

# run simulation
plot.plot()
#O.run()

## Question information

Language:
English Edit question
Status:
Solved
For:
Assignee:
No assignee Edit question
Solved by:
Jan Stránský
Solved:
2019-11-06
Last query:
2019-11-06
2019-11-04
 Jan Stránský (honzik) said on 2019-11-04: #1

Hello,

- try to make the code MWE

M in MWE stands for minimal. In this case of moving a disk, specimen, material, recorder... are irrelevant.
In this light, you can create a few facets manually, still preserving the core of the problem and breaking the need of external files.
Creating a MWE and localizing the problem, often you can solve the problem yourself.

Solution:
- RotationEngine(...,rotateAroundZero=True) # [2]
- specify zeroPoint of RotationEngine [3]
- (probably) you will need to update zeroPoint each time step according to the translation of the cutter

a MWE:
###
O.bodies.append((
facet(((0,0,0),(0,-1,0),(0,0,-1))),
facet(((0,0,0),(0,-1,0),(0,0,+1))),
facet(((0,0,0),(0,+1,0),(0,0,-1))),
facet(((0,0,0),(0,+1,0),(0,0,+1))),
))
ids = (0,1,2,3)
O.engines = [
ForceResetter(),
CombinedKinematicEngine(ids=ids,comb=[
TranslationEngine(translationAxis=[0,0,1],velocity=-580,label="trans"),
RotationEngine(angularVelocity=6700,rotationAxis=[1,0,0],rotateAroundZero=True,zeroPoint=(0,0,0),label="rot")
],label='Eng'),
NewtonIntegrator(),
PyRunner(iterPeriod=1,command="rot.zeroPoint+=(0,0,trans.velocity*O.dt)")
]
O.dt = 1e-8
###

cheers
Jan

 Panos (panoschristo) said on 2019-11-06: #2

Thank you Jan Stránský for your response and sorry for replying late to it, i was trying to
make it work.
I did what you point me and it seems that the problem is fixed.
I only changed O.dt=PWaveTimeStep()

I am having a different problem a mismatch between Disk’s position and
the data plotted on a graph. But I wil open a new thread if I don’t figure it.
Thank you.

 Panos (panoschristo) said on 2019-11-06: #3

Thanks Jan Stránský, that solved my question.