one dimensional element that interacts with spheres
Goodmorning,
i need help with yade DEM simulation.
i need to know how to insert a one dimensional element (like a cylinder) that interacts with the spheres (that simulated the soil).
i have the script of direct shear test and i would like to know how and where to add this one dimensional element in the script.
i write the script below :
###
from yade import pack,plot,export
import math
sp=pack.
O.periodic=True
# dimensions of sample (fixed by particle size such as L/D~15)
RADIUS=0.001
length=0.06
height=length/3
width=length
thickness=RADIUS
# friction angles
compFRIC=40. # during compaction (controls porosity)
FRIC=40. # during shear
# boundary conditions
PI=25e2 # sample preparation: pressure applied for the isotropic compaction
SN=25e3 # normal stress
RATE=0.1 # shear velocity (top plate)
# simulation control
DAMPSHEAR=0.3
ITER=4e4
VTK=20
OUT='PTx'
####
O.cell.
O.materials.
O.materials.
upBox = utils.box(
lowBox = utils.box(
O.bodies.
sp.makeCloud(
O.bodies.
effCellVol=
volRatio=
#print 'volRatio=
O.engines=[
ForceResetter()
,InsertionSort
,InteractionLoop(
[Ig2_
[Ip2_
[Law2_
)
,GlobalStiffne
,PeriTriaxCont
,NewtonIntegra
,PyRunner(
,VTKRecorder(
]
def dataRecorder():
h=vol=
h=O.bodies[
vol=h*
contactStress=
for o in O.bodies:
if isinstance(
nb_s+=1
vol_s += 4.*pi/3.
n = 1-vol_s/vol
nbFrictCont=0.
for i in O.interactions:
if i.isReal and i.phys.
nbFrictCont+=1
plot.addData(
iter=O.iter
,contactStres
,xW=(
,height=h
,volume=vol
,porosity=n
)
def triaxDone():
global phase
volRatio=
h=O.bodies[
vol=h*
contactStress=
vol_s=
Rmin=1e6
for o in O.bodies:
if isinstance(
nbSph+=1
Rmean+
if o.shape.
if o.shape.
vol_s += 4.*pi/3.
Rmean=Rmean/nbSph
n = 1-vol_s/vol
print ('DONE! iter=',O.iter,'| sample generated: nb spheres',nbSph,', Rmean=',Rmean,', Rratio=
print ('Changing contact properties now')
tt=TriaxialCom
tt.setContactP
triax.dead=True
O.pause()
#### Initialization
print ('SAMPLE PREPARATION!')
recData.dead=True # uncomment to record what is happening during stress initialization
O.run(600000,1)
saveSolid.
#saveSolid.
O.step()
saveSolid.dead=True
O.save(
print ('Normal stress (platen) = ',O.forces.
print ('Normal stress (contacts) = ',getStress(
#### Applying normal stress
print ('NORMAL LOADING! iter=',O.iter)
recData.dead=True
stage=0
stiff=fnPlaten=
def servo():
global stage,stiff,
if stage==0:
currentSN=
unbF=
#print 'SN=',SN,'| current SN = ',currentSN,' | unbF=',unbF
boundaryVel=
O.bodies[
if ( (abs(currentSN-
stage+=1
fnPlaten=
print ('Normal stress =',currentSN,' | unbF=',unbF)
## the following computes the stiffness of the plate (used for stress control of the top plate)
for i in O.interactions.
stiff+
print ('DONE! iter=',O.iter)
O.pause()
if stage==1:
fnDesired=
#boundaryVel=
boundaryVel=
O.bodies[
O.engines = O.engines[
O.run(35000,1)
print ('STABILIZING! iter=',O.iter)
O.run(10000,1)
# coloring particles to see vertical stripes in material
dxi=4*(2*RADIUS) # can be a function of cell length dxi=O.cell.
n=int(length/dxi)
xmin=1e6
for i in range(0,n):
for o in O.bodies:
if o.id>1:
if o.state.
if (o.state.
o.shape.
for o in O.bodies:
if (o.state.
o.shape.
saveSolid.dead=True
#saveSolid.
O.step()
saveSolid.dead=True
O.save(
if recData.
print ('Normal stress (platen) = ',O.forces.
print ('Normal stress (contacts) = ',getStress(
#### preparing for shearing
Gl1_Sphere.
print ('Gluing spheres to boundary walls')
gluedSpheres=[]
## gluing particles in contact with the walls
for i in O.interactions:
if i.isReal:
if isinstance(
O.bodies[
gluedSpheres += [O.bodies[i.id2]]
if isinstance(
O.bodies[
gluedSpheres += [O.bodies[i.id1]]
for i in O.interactions:
if i.isReal and ( O.bodies[
O.bodies[
O.bodies[
O.bodies[
O.bodies[
i.phys.
print ('nb of glued spheres=
#### shearing
print ('SHEARING! iter=',O.iter)
for i in O.interactions:
if i.isReal:
if isinstance(
O.bodies[
O.bodies[
saveSolid.dead=True
O.step()
saveSolid.dead=True
O.save(
recData.dead=False
newton.
saveSolid.dead=True
saveSolid.
shearVel=0
iterShear=O.iter
while 1:
O.run(int(10),1)
if shearVel<RATE:
shearVel+
# only top wall moves
O.bodies[
## top and bottom walls move
#O.bodies[
#O.bodies[
if ( O.iter >= (int(iterShear+
print ('iter=',O.iter,' -> FINISHED!')
plot.
O.save(
sys.exit(0)
thanks for the answer and sorry for my english.
Alessandro
Question information
- Language:
- English Edit question
- Status:
- Answered
- For:
- Yade Edit question
- Assignee:
- No assignee Edit question
- Last query:
- Last reply:
Can you help with this problem?
Provide an answer of your own, or ask alessandro solari for more information if necessary.