yade-users team mailing list archive

[Yor1 <question269691@xxxxxxxxxxxxxxxxxxxxx>]

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

Hi everyone,

I'm trying to determine the uniaxial strain of the beam in 3-points bending test. 
When i execute my code i obtain this error message "Erreur de segmentation (core dumped)".
This is a part of my code:

from yade import ymport, utils, pack, export
from yade import plot
from pylab import *
import math

#---------------- SIMULATIONS DEFINED HERE (assembly, material, boundary conditions)

#### packing (previously constructed)
PACKING='poutre400_notch'
OUT=PACKING+'_flexion_3_points'

#### Simulation Control
DAMP=0.4 # numerical damping
saveData=100 # data record interval
iterMax=200000 # maximum number of iteration (to be adjusted)
saveVTK=10000 # Vtk files record interval

#### Boundary Conditions
strainRate=0.002

#### Material microproperties -> Lac du Bonnet granite (cf. A DEM model for soft and hard rock, Scholtes & Donze, JMPS 2013)
intR=1.5# allows near neighbour interaction and coordination number K=13 (determined with coordinationNumber.py -> to be adjusted for every packing)
DENS=4000 # could be adapted to match material density: dens_DEM=dens_rock*(V_rock/V_particles)=dens_rock*1/(1-poro_DEM) -> poro?
YOUNG=65e9 
FRICT=10
ALPHA=0.4
TENS=8e6 
COH=160e6
#### material definition
sphereMat = JCFpmMat(type=1,density=DENS,young=YOUNG,poisson = ALPHA,frictionAngle=radians(FRICT),tensileStrength=TENS,cohesion=COH)
wallMat = JCFpmMat(type=0,density=DENS,young=YOUNG,frictionAngle=radians(0))

for mat in (sphereMat,wallMat):
   O.materials.append(mat) # then wallMat will be used if material is not specified

#### preprocessing to get dimensions
O.bodies.append(ymport.text(PACKING+'.spheres',scale=1.,shift=Vector3(0,0,0)))

dim=utils.aabbExtrema()
xinf=dim[0][0]
xsup=dim[1][0]
X=xsup-xinf
yinf=dim[0][1]
ysup=dim[1][1]
Y=ysup-yinf
zinf=dim[0][2]
zsup=dim[1][2]
Z=zsup-zinf

r=1
#### preprocessing to get dimensions
R=0
Rmax=0
numSpheres=0.
for o in O.bodies:
 if isinstance(o.shape,Sphere):
   numSpheres+=1
   R+=o.shape.radius
   if o.shape.radius>Rmax:
     Rmax=o.shape.radius
Rmean=R/numSpheres

O.reset() # all previous lines were for getting dimensions of the packing to create walls at the right positions (below) because walls have to be genrated after spheres for FlowEngine


O.bodies.append(geom.facetCylinder(center=(xinf+X/5.,yinf-0.9*r,Z/2.),radius=r,height=Z,orientation=Quaternion((1, 0, 0), 0),segmentsNumber=20,wire=False,material=wallMat)) # bas gauche
O.bodies.append(geom.facetCylinder(center=(0.8*X,yinf-0.9*r,Z/2.),radius=r,height=Z,orientation=Quaternion((1, 0, 0), 0),segmentsNumber=20,wire=False,material=wallMat)) # bas droite
piston=O.bodies.append(geom.facetCylinder(center=(0.5*X,Y+r,Z/2.),radius=r,height=Z,dynamic=False,orientation=Quaternion((1, 0, 0), 0),segmentsNumber=20,wire=False,material=wallMat)) # haut

### packing
beam=O.bodies.append(ymport.text(PACKING+'.spheres',scale=1.,shift=Vector3(0,0,0),material=sphereMat))


### set a color to the spheres
for o in O.bodies:
 if isinstance(o.shape,Sphere):
   o.shape.color=(0.7,0.5,0.3)

crossSectionArea=X*Z

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

#### simulation is defined here (DEM loop, interaction law, servo control, recording, etc...)
##### simulation piston's motion

for i in range(0,len(piston)):
	O.bodies[piston[i]].state.vel[1]=-1

##### simulation beam's grains motion
O.engines=[
        ForceResetter(),
        InsertionSortCollider([Bo1_Box_Aabb(),Bo1_Sphere_Aabb(aabbEnlargeFactor=intR,label='Saabb'),Bo1_Facet_Aabb()]),
	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')]
	),
	UniaxialStrainer(strainRate=strainRate,axis=1,asymmetry=0,crossSectionArea=crossSectionArea,blockDisplacements=1,blockRotations=1,setSpeeds=0,stopStrain=0.1,label='strainer'),
        GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=10,timestepSafetyCoefficient=0.4, defaultDt=0.1*utils.PWaveTimeStep()),
        NewtonIntegrator(damping=DAMP,label="newton"),
        PyRunner(iterPeriod=int(saveData),initRun=True,command='recorder()',label='data'),
        VTKRecorder(iterPeriod=int(saveVTK),initRun=True,fileName=OUT+'-',recorders=['spheres','jcfpm','cracks'],Key=OUT,label='vtk')
]

Best regards

-- 
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.