yade-users team mailing list archive

[Hanying Zhang <question694472@xxxxxxxxxxxxxxxxxxxxx>]

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

hi,

I have these errors "double free or corruption (out)/segment error/bus error" a few seconds after I click the run button of my script. It is basicly about generating a cylindrical specimen consisted of spheres and cylinders. My yade version is 2018.02b and I use ubuntu18.04. Below is my code.

#====================================================================
from yade.gridpfacet import *
from yade import pack, plot
#from random import random as r
import numpy as np
from numpy import *
import math

#parameters
# target size after equilibrium (before shearing): cylinder with hight=80e-3, radius=20e-3
# target size before equilibrium (initial size): cylinder with hight=800e-3, radius=20e-3
height = 500e-3
width = 40e-3
#rParticle = 0.11*1e-3/2
#rRelFuzz=.0005
mi,ma = (-20e-3,-20e-3,0),(20e-3,20e-3,500e-3)
nCyls,nSphs = 80,140
frictionAngleSph=30
frictionAngleCyl=30


#=============================meterials========================================
O.materials.append(FrictMat(young=4.0e6,poisson=.3,frictionAngle=30,density=2630e0,label='sphereMat'))#for sphere

O.materials.append(FrictMat(young=3e9,poisson=.15,frictionAngle=20,density=910e+0,label='extcylMat'))#for sphere-cylinder
O.materials.append(CohFrictMat(young=3e9,poisson=.15,density=910e0,frictionAngle=20,normalCohesion=1e40,shearCohesion=1e40,momentRotationLaw=True,label='intcylMat'))#for cylinder-cylinder

O.materials.append( CohFrictMat( young=3e8,poisson=0.15,density=910e0,frictionAngle=20,normalCohesion=3e100,shearCohesion=3e100,momentRotationLaw=True,label='gridNodeMat' ) )#for gridNodes
#O.materials.append(CohFrictMat(young=3e9,poisson=.15,density=910e6,frictionAngle=20,normalCohesion=1e40,shearCohesion=1e40,momentRotationLaw=True,label='gridNodeMat'))#for gridNodes
O.materials.append(FrictMat(young=4e6,poisson=0.3,density=1000,frictionAngle=20,label='pFacetMat')) #for pfacet


#==============================Engines=========================================
O.engines=[
	ForceResetter(),
	InsertionSortCollider([
		Bo1_Sphere_Aabb(),
		Bo1_GridConnection_Aabb(),
		Bo1_Facet_Aabb()
	]),
	InteractionLoop([
		Ig2_Sphere_Sphere_ScGeom(),
		Ig2_Box_Sphere_ScGeom(),
		Ig2_GridNode_GridNode_GridNodeGeom6D(),
		Ig2_Sphere_GridConnection_ScGridCoGeom(),
		Ig2_GridConnection_GridConnection_GridCoGridCoGeom(),
		Ig2_Facet_Sphere_ScGeom(),
	],
	[
		Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True,setCohesionOnNewContacts=False),	# internal cylinder physics
		Ip2_FrictMat_FrictMat_FrictPhys(),	# physics for external interactions, i.e., cylinder-cylinder, sphere-sphere, cylinder-sphere

		Ip2_FrictMat_CpmMat_FrictPhys(),
	],
	[
		Law2_ScGeom_FrictPhys_CundallStrack(),	# contact law for sphere-sphere
		Law2_ScGridCoGeom_FrictPhys_CundallStrack(),	# contact law for cylinder-sphere
		Law2_ScGeom6D_CohFrictPhys_CohesionMoment(),	# contact law for "internal" cylinder forces
		Law2_GridCoGridCoGeom_FrictPhys_CundallStrack(),	# contact law for cylinder-cylinder interaction
		#Law2_GridCoGridCoGeom_CohFrictPhys_CundallStrack(),
		Law2_ScGeom_CpmPhys_Cpm(),
	]
	),
	GlobalStiffnessTimeStepper(timestepSafetyCoefficient=0.5,label='ts'), 
	NewtonIntegrator(gravity=(0,0,-9.81),damping=0.5,label='newton'),
]

#=======================================sphere=================================================
sp = yade.pack.SpherePack()
sp.makeCloud(mi,ma,psdSizes=[0.5e-3,0.75e-3,0.8e-3,0.9e-3,1.2e-3,2e-3],psdCumm=[0.,0.2,0.4,0.6,0.8,1.0],num=1080)
#sp.makeCloud(mi,ma,porosity=0.66/1.66,psdSizes=[0.5e-3,0.75e-3,0.8e-3,0.9e-3,1.2e-3,2e-3],psdCumm=[0.,0.2,0.4,0.6,0.8,1.0],num=108000)
# e min =0.56, and e max =0.89, e=0.66
spheres=sp.toSimulation(color=(0,0.5,0.7),material='sphereMat')

num=0
vol=0
mass=0
for b in O.bodies:
    r=b.shape.radius
    x=b.state.pos[0]
    y=b.state.pos[1]
    z=b.state.pos[2]
    ID=b.id
    if x**2+y**2>=(20e-3-r)**2:
        O.bodies.erase(b.id)
    else:
        num+=1
        vol+=4./3.*np.pi*r**3
        mass+=2.69*4./3.*np.pi*r**3
print 'num, e, voidRatio =',num, 0.66, ((80e-3*(40e-3*0.5)**2*np.pi)-vol)/vol


#=======================================fiber=================================================
length=25.4e-3
diameter=0.3e-3
Xw=0.25e-2 # fiber content
fiber_vol=(80e-3*(40e-3*0.5)**2*np.pi)*Xw
Vfiber=length*np.pi*(0.5*diameter)**2
numFiber=fiber_vol/Vfiber
print "the numer of fibres =",int(numFiber)
#numFiber=156
Ne=int(length/(0.9e-3))
print 'Ne=',Ne
nodesIds=[]
cylIds=[]
numFibre=0
target_num=int(numFiber)
target_num=int(8)
####fiberfilewrite
fiberwr=open('fibers.txt','w')
for n in range(500):
    if numFibre<target_num:
		random.seed(n*5)
		z0=random.uniform(0, 80e-3 )
		random.seed( n+2)
		sita0=random.uniform(-np.pi,np.pi)
		random.seed( n+9)
		l_y=random.uniform(0,20e-3)
		x0=l_y*sin(sita0)
		y0=l_y*cos(sita0)
		
		random.seed(n+8)
		sita=random.uniform(-0.5*np.pi, 0.5*np.pi)
		random.seed(n+9)
		phi=random.uniform(0, 2*np.pi)  
		
		hz=length*sin(sita)
		hx=length*cos(sita)*cos(phi)
		hy=length*cos(sita)*sin(phi)
		
		x1=x0+hx
		y1=y0+hy
		z1=z0+hz

		if (x1**2+y1**2)<(40e-3*0.5)**2 and 0<(z1)<80e-3 and (x0**2+y0**2)<(40e-3*0.5)**2 and 0<z0<80e-3 :
		    fiberwr.write(str(numFibre)+'\t'+str(x0)+'\t'+str(y0)+'\t'+str(z0)+'\t'+str(x1)+'\t'+str(y1)+'\t'+str(z1)+'\n')
		    numFibre+=1
		    vertices=[]
		    for i in range(0, Ne+1):
		        px=float(i)*hx/float(Ne)+x0; py=float(i)*hy/float(Ne)+y0; pz=float(i)*hz/float(Ne)+z0;
			vertices.append([px,py,pz])        
		    cylinderConnection(vertices,0.15e-3,nodesIds,cylIds,color=[0.5,0.5,0],fixed=True,highlight=False,intMaterial='intcylMat',extMaterial='extcylMat')

#=======================================facet=================================================
nw = 24
nh = 60
R=0.00002
color=[255./255.,102./255.,0./255.]
rCyl2 = .5*width / cos(pi/float(nw))
# facets
facets = []
rCyl2 = .5*width / cos(pi/float(nw))
for r in xrange(nw):
	for h in xrange(nh):
		v1 = Vector3( rCyl2*cos(2*pi*(r+0)/float(nw)), rCyl2*sin(2*pi*(r+0)/float(nw)), height*(h+0)/float(nh) )
		v2 = Vector3( rCyl2*cos(2*pi*(r+1)/float(nw)), rCyl2*sin(2*pi*(r+1)/float(nw)), height*(h+0)/float(nh) )
		v3 = Vector3( rCyl2*cos(2*pi*(r+1)/float(nw)), rCyl2*sin(2*pi*(r+1)/float(nw)), height*(h+1)/float(nh) )
		v4 = Vector3( rCyl2*cos(2*pi*(r+0)/float(nw)), rCyl2*sin(2*pi*(r+0)/float(nw)), height*(h+1)/float(nh) )
		f1 = facet((v1,v2,v3),color=[255./255.,102./255.,0./255.],material='gridNodeMat')
		f2 = facet((v1,v3,v4),color=[255./255.,102./255.,0./255.],material='gridNodeMat')
		facets.extend((f1,f2))
O.bodies.append(facets)
#===============================================================================

Thanks for any help!!!!!!!



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