yade-users team mailing list archive

[liucheng83 <question233349@xxxxxxxxxxxxxxxxxxxxx>]

Question #233349 on Yade changed:
https://answers.launchpad.net/yade/+question/233349

liucheng83 gave more information on the question:
Hello, anybody can give some advice to my question? 
some strange things comes out, when I use the kinemCNLEngine,
(1) The fx of upperBox go from positive to negative in shear stage, why? Can fx negative?
(2) The gamma is not equal to shearSpeed*nCycShear, it is multiplied by 2, that is 2*shearSpeed*nCycShear
(3)The calculation is too slow, maybe a step per second, how to fasten the calculation ?
I change the m,kg,Pa to mm,kg, kPa as following.
Please test the script, thank you!
----

# -*- coding: utf-8 -*-

# J. Duriez duriez@xxxxxxxxxxxxxxx
# "PreProcessor" script to define a numerical model for simple shear box
# In order to illustrate various Kinem...Engine which allow to perform different loadings on the box
# There is here one § for each Kinem...Engine, comment/uncomment them to observe what you want

#NB : the run of the script is paused at each plot window (so that there
is time to observe it). Type "Return" in the Yade terminal to resume


from yade import plot
from yade.pack import *


O=Omega() 

#Def of the material which will be used
#O.materials.append(NormalInelasticMat(density=2600,young=4.0e9,poisson=.04,frictionAngle=.6,coeff_dech=3.0,label='Materiau1'))
idSoil=O.materials.append(FrictMat(density=2.7e-6,young=4e6,poisson=.3,frictionAngle=radians(30),label="matSoil"))

#Def of dimensions of the box
#length=5
#height=1
#width=2
length=5
height=1
width=2
thickness=0.05

porosity = 0.55
nSpheres=2000
rMean = pow( ((1-porosity) * length * height * width) / (nSpheres * 4.0/3.0 * pi) , 1.0/3.0 )

#Definition of bodies constituing the numerical model : six boxes corresponding to sids of the simple shear box, containing a particle samples
pred=inAlignedBox((0,0,-0.02),(0.1,0.02,0.02))

leftBox = box( center=(-thickness/2.0,(height)/2.0,0),
extents=(thickness/2.0,5*(height/2.0+thickness),width/2.0)
,fixed=True,wire=True)

lowBox = box( center=(length/2.0,-thickness/2.0,0),
extents=(length/2.0,thickness/2.0,width/2.0) ,fixed=True,wire=True)

rightBox = box( center=(length+thickness/2.0,height/2.0,0),
extents=(thickness/2.0,5*(height/2.0+thickness),width/2.0)
,fixed=True,wire=True)

upBox = box( center=(length/2.0,height+thickness/2.0,0),
extents=(length/2.0,thickness/2.0,width/2.0) ,fixed=True,wire=True)

behindBox = box(
center=(length/2.0,height/2.0,-width/2.0-thickness/2.0),
extents=(2.5*length/2.0,height/2.0+thickness,thickness/2.0),
fixed=True,wire=True)

inFrontBox = box(
center=(length/2.0,height/2.0,width/2.0+thickness/2.0),
extents=(2.5*length/2.0,height/2.0+thickness,thickness/2.0),
fixed=True,wire=True)

O.bodies.append([leftBox,lowBox,rightBox,upBox,behindBox,inFrontBox])

#memoizeDb='/tmp/simpleshear-triax-packings.sqlite'
#ListSph=randomDensePack(pred,radius=0.002,rRelFuzz=0.15,memoDbg=True,memoizeDb=memoizeDb)
#ListSph=randomDensePack(pred,radius=0.002,rRelFuzz=0.15,memoDbg=True,memoizeDb=memoizeDb,spheresInCell=100)
#O.bodies.append(ListSph)

sp=yade._packSpheres.SpherePack()
sp.makeCloud(Vector3(0,0.0,-width/2.0),Vector3(length,height,width/2.0),rMean,.15)
O.bodies.append([sphere(s[0],s[1]) for s in sp])


#---- Def of the engines ----#

O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
		[Ip2_FrictMat_FrictMat_FrictPhys()],
		[Law2_ScGeom_FrictPhys_CundallStrack()]
	),
	NewtonIntegrator(damping=.1)
	,PyRunner(iterPeriod=50,command='defData()')
	]

def defData():
	plot.addData(fy=O.forces.f(3)[1], # vertical component of the force sustained by the upper side of the shear box
		    fx=O.forces.f(3)[0], # horizontal component of the force sustained by the upper side of the shear box
		    step=O.iter,
		    gamma=O.bodies[3].state.pos[0] - length/2.0,
		    u=O.bodies[3].state.pos[1] - (height+thickness/2.0)
		    )

fy=O.forces.f(3)[1] # The fy above exists only as a label in plot.data


from yade import qt
qt.View()
#O.save("SimpleShearReady.xml")
#O.saveTmp("InitialState")


#---- Compression ----
O.engines = O.engines+[KinemCTDEngine(compSpeed=500,sigma_save=(),temoin_save=(),targetSigma=0.1,LOG=False)]

O.dt=.4*PWaveTimeStep()
print ''
print 'Be patient, running in progress (the shear box is being compressed)'
O.run(14000)
plot.plots={'step':('fy',)}
plot.plot(subPlots=False)
print 'Plotting curve (fy = normal force). Type Return to go ahead'
print ''
raw_input()
plot.plots={'u':('fy',)}
plot.plot(subPlots=False)
print 'Plotting curve (u = normal displacement). Type Return to go ahead'
print ''
raw_input()
#O.save('FinComp.xml')


##---- Shear at constant normal displacement ----
#nCycShear = 20000

#O.engines=O.engines[:5]+[KinemCNDEngine(shearSpeed=(length/7.0)/(nCycShear*O.dt),gamma_save=(),temoin_save=(),gammalim=length/7.0,LOG=False)]

#print 'Be patient, running in progress (the shear box is being sheared)'
#O.run(int(1.15*nCycShear),True)
#plot.plots={'step':('gamma',)}
#plot.plot(subPlots=False)
#print 'Plotting curve (gamma = tangential displacement). Type Return to go ahead'
#print ''
#raw_input()
#plot.plots={'gamma':('fx','fy',)}
#plot.plot(subPlots=False)
#print 'Plotting curve (fx = tangential force). Type Return to go ahead'
#print ''
#raw_input()


##---- A re-compression, from this initial sheared state ----
#O.engines=O.engines[:5]+[KinemCTDEngine(compSpeed=0.5,sigma_save=(),temoin_save=(),targetSigma=80000.0,LOG=False)]
#print 'Be patient, running in progress (the sample is being again compressed, from this sheared state)'
#O.run(10000,True)
#plot.plots={'u':('fx','fy',)}
#plot.plot(subPlots=False)
#print 'Plotting curve. Type Return to go ahead'
#print ''
#raw_input()
#print 'End of script'

#---- Shear at constant normal load/stress ----
nCycShear = 20000
O.engines=O.engines[:6]+[KinemCNLEngine(shearSpeed=(length/10.0)/(nCycShear*O.dt),gamma_save=(),temoin_save=(),gammalim=length/10.0,LOG=False)]
O.run(int(1.15*nCycShear))
plot.plots={'step':('gamma','u',)}
plot.plot()
raw_input()
plot.plots={'gamma':('fx','fy',)}
plot.plot()
raw_input()


#---- Shear at constant normal stifness ----
#nCycShear = 20000
#O.engines=O.engines[:6]+[KinemCNSEngine(shearSpeed=(length/10.0)/(nCycShear*O.dt),gammalim=length/10.0,LOG=False,KnC=1)]
#O.run(int(1.15*nCycShear),True)
#plot.plots={'step':('gamma','u',)}
#plot.plot()
#raw_input()
#plot.plots={'gamma':('fx','fy',)}
#plot.plot()
#raw_input()
#plot.plots={'u':('fy',)}
#plot.plot()
#raw_input()

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