yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #12985
[Branch ~yade-pkg/yade/git-trunk] Rev 4001: move PBC example scripts
------------------------------------------------------------
revno: 4001
committer: bchareyre <bruno.chareyre@xxxxxxxxxxxxxxx>
timestamp: Tue 2017-02-07 17:46:13 +0100
message:
move PBC example scripts
added:
examples/PeriodicBoundaries/
examples/PeriodicBoundaries/cellFlipping.py
examples/PeriodicBoundaries/peri3dController_example1.py
examples/PeriodicBoundaries/peri3dController_shear.py
examples/PeriodicBoundaries/peri3dController_triaxialCompression.py
examples/PeriodicBoundaries/periodic-compress.py
examples/PeriodicBoundaries/periodic-grow.py
examples/PeriodicBoundaries/periodic-shear.py
examples/PeriodicBoundaries/periodic-simple-shear.py
examples/PeriodicBoundaries/periodic-simple.py
examples/PeriodicBoundaries/periodic-triax-settingHsize.py
examples/PeriodicBoundaries/periodic-triax.py
examples/PeriodicBoundaries/periodicSandPile.py
--
lp:yade
https://code.launchpad.net/~yade-pkg/yade/git-trunk
Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-pkg/yade/git-trunk/+edit-subscription
=== added directory 'examples/PeriodicBoundaries'
=== added file 'examples/PeriodicBoundaries/cellFlipping.py'
--- examples/PeriodicBoundaries/cellFlipping.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/cellFlipping.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,48 @@
+# coding: utf-8
+# 2017 Bruno Chareyre <bruno.chareyre~a~grenoble-inp.fr>
+"Demonstrate cell flipping in periodic boundary conditions"
+from yade import pack,qt,plot
+
+O.periodic=True
+O.cell.Hsize=Matrix3(0.1,0,0, 0,0.1,0, 0,0,0.1)
+sp=pack.SpherePack()
+radius=5e-3
+num=sp.makeCloud((0,0,0),(.1,.1,.1),radius,.2,100,periodic=True)
+O.bodies.append([sphere(s[0],s[1]) for s in sp])
+
+
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()],verletDist=.05*radius),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()]
+ ), PeriTriaxController(dynCell=True,mass=0.2,maxUnbalanced=0.01,relStressTol=0.02,goal=[-1e4,-1e4,-1e4],stressMask=7,globUpdate=5,maxStrainRate=[1.,1.,1.],doneHook='triaxDone()',label='triax'),
+ NewtonIntegrator(damping=.2),
+]
+O.dt=2.e-5
+phase=0
+def triaxDone():
+ global phase
+ if phase==0:
+ print 'Here we are: stress',triax.stress,'strain',triax.strain,'stiffness',triax.stiff
+ print 'Now shearing.'
+ O.cell.velGrad=Matrix3(0,1,0, 0,0,0, 0,0,0)
+ triax.stressMask=7
+ triax.goal=[-1e4,-1e4,-1e4]
+ phase+=1
+ O.saveTmp()
+ O.pause()
+
+O.run(-1,True)
+
+def addPlotData():
+ plot.addData(t=O.time,gamma=O.cell.trsf[0,1], unb = unbalancedForce(),s=getStress()[0,1])
+O.engines=O.engines+[PyRunner(command='addPlotData()',iterPeriod=20)]
+plot.plots={'gamma':'unb'}
+plot.plot()
+
+O.engines=O.engines+[PyRunner(command='if O.cell.hSize[0,1]>O.cell.hSize[0,0]: flipCell()',iterPeriod=20)]
+
+#now click play and watch the flips happening, the evolution of stress or unbalanced force should not show any discontinuity
\ No newline at end of file
=== added file 'examples/PeriodicBoundaries/peri3dController_example1.py'
--- examples/PeriodicBoundaries/peri3dController_example1.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/peri3dController_example1.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,67 @@
+# peri3dController_example1.py
+# script, that explains funcionality and input parameters of Peri3dController
+
+from yade import pack, plot
+
+# create some material
+O.materials.append(CpmMat(neverDamage=True,young=25e9,frictionAngle=.7,poisson=.2,sigmaT=3e6,epsCrackOnset=1e-4,relDuctility=30))
+
+# create periodic assembly of particles
+initSize=1.2
+sp=pack.randomPeriPack(radius=.05,initSize=Vector3(initSize,initSize,initSize),memoizeDb='/tmp/packDb.sqlite')
+sp.toSimulation()
+
+# plotting
+#plot.live=False
+plot.plots={'progress':('sx','sy','sz','syz','szx','sxy',),'progress_':('ex','ey','ez','eyz','ezx','exy',)}
+def plotAddData():
+ plot.addData(
+ progress=p3d.progress,progress_=p3d.progress,
+ sx=p3d.stress[0],sy=p3d.stress[1],sz=p3d.stress[2],
+ syz=p3d.stress[3],szx=p3d.stress[4],sxy=p3d.stress[5],
+ ex=p3d.strain[0],ey=p3d.strain[1],ez=p3d.strain[2],
+ eyz=p3d.strain[3],ezx=p3d.strain[4],exy=p3d.strain[5],
+ )
+
+# in how many time steps should be the goal state reached
+nSteps=4000
+
+O.dt=PWaveTimeStep()/2
+EnlargeFactor=1.5
+EnlargeFactor=1.0
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=EnlargeFactor,label='bo1s')]),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=EnlargeFactor,label='ig2ss')],
+ [Ip2_CpmMat_CpmMat_CpmPhys()],[Law2_ScGeom_CpmPhys_Cpm()]),
+ NewtonIntegrator(),
+ Peri3dController(
+ goal=(20e-4,-6e-4,0, -2e6,3e-4,2e6), # Vector6 of prescribed final values (xx,yy,zz, yz,zx,xy)
+ stressMask=0b101100, # prescribed ex,ey,sz,syz,ezx,sxy; e..strain; s..stress
+ nSteps=nSteps, # how many time steps the simulation will last
+ # after reaching nSteps do doneHook action
+ doneHook='print "Simulation with Peri3dController finished."; O.pause()',
+
+ # the prescribed path (step,value of stress/strain) can be defined in absolute values
+ xxPath=[(465,5e-4),(934,-5e-4),(1134,10e-4)],
+ # or in relative values
+ yyPath=[(2,4),(7,-2),(11,0),(14,4)],
+ # if the goal value is 0, the absolute stress/strain values are always considered (step values remain relative)
+ zzPath=[(5,-1e7),(10,0)],
+ # if ##Path is not explicitly defined, it is considered as linear function between (0,0) and (nSteps,goal)
+ # as in yzPath and xyPath
+ # the relative values are really relative (zxPath gives the same - except of the sign from goal value - result as yyPath)
+ zxPath=[(4,2),(14,-1),(22,0),(28,2)],
+ xyPath=[(1,1),(2,-1),(3,1),(4,-1),(5,1)],
+ # variables used in the first step
+ label='p3d'
+ ),
+ PyRunner(command='plotAddData()',iterPeriod=1),
+]
+
+O.step()
+bo1s.aabbEnlargeFactor=ig2ss.interactionDetectionFactor=1.
+
+O.run(); #O.wait()
+plot.plot(subPlots=False)
=== added file 'examples/PeriodicBoundaries/peri3dController_shear.py'
--- examples/PeriodicBoundaries/peri3dController_shear.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/peri3dController_shear.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,61 @@
+# peri3dController_shear.py
+# script, that explains funcionality and input parameters of Peri3dController on the example of
+# shear test with rotated periodic cell (that enables strain localization).
+# The simulation is run on rotated cell to enable localization and strain softening
+# (you can also try simulation with different angles of rotation to pbtain different results.
+
+
+from yade import pack,plot,qt
+
+# define material
+O.materials.append(CpmMat(young=25e9,poisson=.2,sigmaT=3e6,epsCrackOnset=1e-4,relDuctility=1e-6))
+
+# create periodic assembly of particles
+initSize=1.2
+sp=pack.randomPeriPack(radius=.05,initSize=Vector3(initSize,initSize,initSize),memoizeDb='/tmp/packDb.sqlite')
+angle=0.
+#angle=pi/4
+rot=Matrix3(cos(angle),-sin(angle),0, sin(angle),cos(angle),0, 0,0,1)
+sp.toSimulation(rot=rot)
+
+# plotting
+plot.live=False
+plot.plots={'iter':('sx','sy','sz','syz','szx','sxy',),'iter_':('ex','ey','ez','eyz','ezx','exy',),'exy':('sxy',)}
+def plotAddData():
+ plot.addData(
+ iter=O.iter,iter_=O.iter,
+ sx=p3d.stress[0],sy=p3d.stress[1],sz=p3d.stress[2],
+ syz=p3d.stress[3],szx=p3d.stress[4],sxy=p3d.stress[5],
+ ex=p3d.strain[0],ey=p3d.strain[1],ez=p3d.strain[2],
+ eyz=p3d.strain[3],ezx=p3d.strain[4],exy=p3d.strain[5],
+ )
+
+O.dt=PWaveTimeStep()/2
+
+# define the first part of simulation, hydrostatic compression
+enlargeFactor=1.5
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=enlargeFactor,label='bo1s')]),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=enlargeFactor,label='ig2ss')],
+ [Ip2_CpmMat_CpmMat_CpmPhys()],[Law2_ScGeom_CpmPhys_Cpm()]),
+ NewtonIntegrator(),
+ Peri3dController( goal=(0,0,0, 0,0,5e-3), # Vector6 of prescribed final values
+ stressMask=0b011111,
+ nSteps=2000,
+ doneHook='print "Simulation with Peri3dController finished."; O.pause()',
+ maxStrain=.5,
+ label='p3d'
+ ),
+ PyRunner(command='plotAddData()',iterPeriod=1),
+]
+O.step()
+bo1s.aabbEnlargeFactor=ig2ss.interactionDetectionFactor=1.0
+
+renderer=qt.Renderer()
+renderer.intrPhys,renderer.shape=True,False
+Gl1_CpmPhys.dmgLabel=False
+qt.View()
+O.run()
+#plot.plot()
=== added file 'examples/PeriodicBoundaries/peri3dController_triaxialCompression.py'
--- examples/PeriodicBoundaries/peri3dController_triaxialCompression.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/peri3dController_triaxialCompression.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,79 @@
+# peri3dController_triaxialCompression.py
+# script, that explains funcionality and input parameters of Peri3dController on the example of
+# triaxial compression.
+# Firstly, a hydrostatic preassure is applied, then a strain along z axis is increasing
+# while x- and y- stress is constant
+# The simulation is run on rotated cell to enable localization and strain softening
+# (you can also try simulation with command sp.toSimulation() with no rotation,
+# in this case there is almost no difference, but in script peri3dController_shear,
+# the cell rotation has significant effect)
+
+from yade import pack,plot,qt
+
+# define material
+O.materials.append(FrictMat())
+
+# create periodic assembly of particles
+initSize=1.2
+sp=pack.randomPeriPack(radius=.05,initSize=Vector3(initSize,initSize,initSize),memoizeDb='/tmp/packDb.sqlite')
+angle=0
+rot=Matrix3(cos(angle),0,-sin(angle), 0,1,0, sin(angle),0,cos(angle))
+sp.toSimulation(rot=rot)
+
+# plotting
+plot.live=False
+plot.plots={'iter':('sx','sy','sz','syz','szx','sxy',),'iter_':('ex','ey','ez','eyz','ezx','exy',),'ez':('sz',)}
+def plotAddData():
+ plot.addData(
+ iter=O.iter,iter_=O.iter,
+ sx=p3d.stress[0],sy=p3d.stress[1],sz=p3d.stress[2],
+ syz=p3d.stress[3],szx=p3d.stress[4],sxy=p3d.stress[5],
+ ex=p3d.strain[0],ey=p3d.strain[1],ez=p3d.strain[2],
+ eyz=p3d.strain[3],ezx=p3d.strain[4],exy=p3d.strain[5],
+ )
+
+O.dt=PWaveTimeStep()/2
+
+# define the first part of simulation, hydrostatic compression
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()]),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],[Law2_ScGeom_FrictPhys_CundallStrack()]),
+ NewtonIntegrator(),
+ Peri3dController( goal=(-1e7,-1e7,-1e7, 0,0,0), # Vector6 of prescribed final values
+ stressMask=0b000111,
+ nSteps=500,
+ doneHook='print "Hydrostatic load reached."; O.pause()',
+ youngEstimation=.5e9, # needed, when only nonzero prescribed values are stress
+ maxStrain=.5,
+ label='p3d'
+ ),
+ PyRunner(command='plotAddData()',iterPeriod=1),
+]
+O.run(); O.wait()
+
+# second part, z-axis straining and constant transversal stress
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()]),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],[Law2_ScGeom_FrictPhys_CundallStrack()]),
+ NewtonIntegrator(),
+ Peri3dController( goal=(-1e7,-1e7,-4e-1, 0,0,0), # Vector6 of prescribed final values
+ stressMask=0b000011,
+ nSteps=1000,
+ xxPath=[(0,1),(1,1)], # the first (time) zero defines the initial value of stress considered nonzero
+ yyPath=[(0,1),(1,1)],
+ doneHook='print "Simulation with Peri3dController finished."; O.pause()',
+ maxStrain=.5,
+ label='p3d',
+ strain=p3d.strain, # continue from value reached in previous part
+ stressIdeal=Vector6(-1e7,-1e7,0, 0,0,0), # continue from value reached in previous part
+ ),
+ PyRunner(command='plotAddData()',iterPeriod=1),
+]
+O.run();O.wait()
+plot.plot(subPlots=False)
=== added file 'examples/PeriodicBoundaries/periodic-compress.py'
--- examples/PeriodicBoundaries/periodic-compress.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/periodic-compress.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,48 @@
+'''
+This example shows compression of a packing with a periodic cell.
+'''
+
+O.periodic=True
+O.cell.setBox(20,20,10)
+
+from yade import pack,timing
+
+O.materials.append(FrictMat(young=30e9,density=2400))
+p=pack.SpherePack()
+p.makeCloud(Vector3(0,0,0),Vector3(20,20,10),1,.5,700,True)
+for sph in p:
+ O.bodies.append(sphere(sph[0],sph[1]))
+
+
+O.timingEnabled=True
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()],allowBiggerThanPeriod=True),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()],
+ ),
+ PeriIsoCompressor(charLen=.5,stresses=[-50e9,-1e8],doneHook="print 'FINISHED'; O.pause() ",keepProportions=True),
+ NewtonIntegrator(damping=.4)
+]
+O.dt=PWaveTimeStep()
+O.saveTmp()
+#print O.cell.refSize
+from yade import qt; qt.Controller(); qt.View()
+O.run()
+O.wait()
+timing.stats()
+#while True:
+# O.step()
+
+
+# now take that packing and pad some larger volume with it
+#sp=pack.SpherePack()
+#sp.fromSimulation() # take spheres from simulation; cellSize is set as well
+#O.reset()
+#print sp.cellSize
+#sp.cellFill((30,30,30))
+#print sp.cellSize
+#for s in sp:
+# O.bodies.append(sphere(s[0],s[1]))
=== added file 'examples/PeriodicBoundaries/periodic-grow.py'
--- examples/PeriodicBoundaries/periodic-grow.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/periodic-grow.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,39 @@
+"""Script that shrinks the periodic cell progressively.
+It prints strain and average stress (computed from total volume force)
+once in a while."""
+
+from yade import timing
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()],allowBiggerThanPeriod=True),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()],
+ ),
+ NewtonIntegrator(damping=.6)
+]
+import random
+for i in xrange(250):
+ O.bodies.append(sphere(Vector3(10*random.random(),10*random.random(),10*random.random()),.5+random.random()))
+cubeSize=20
+# absolute positioning of the cell is not important
+O.periodic=True
+O.cell.setBox(cubeSize,cubeSize,cubeSize)
+O.dt=PWaveTimeStep()
+O.saveTmp()
+from yade import qt
+qt.Controller(); qt.View()
+O.run(200,True)
+rate=-1e-3*cubeSize/(O.dt*200)*Matrix3.Identity
+O.cell.velGrad=rate
+
+print 'Please be patient...'
+
+for i in range(0,25):
+ O.run(2000,True)
+ F,stiff=totalForceInVolume()
+ dim=O.cell.refSize; A=Vector3(dim[1]*dim[2],dim[0]*dim[2],dim[0]*dim[1])
+ avgStress=sum([F[i]/A[i] for i in 0,1,2])/3.
+ print 'strain',(cubeSize-dim[0])/cubeSize,'avg. stress ',avgStress,'unbalanced ',unbalancedForce()
+#O.timingEnabled=True; timing.reset(); O.run(200000,True); timing.stats()
=== added file 'examples/PeriodicBoundaries/periodic-shear.py'
--- examples/PeriodicBoundaries/periodic-shear.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/periodic-shear.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,41 @@
+O.periodic=True
+O.cell.setBox(.55,.55,.55)
+O.bodies.append(facet([[.4,.0001,.3],[.2,.0001,.3],[.3,.2,.2]]))
+O.bodies.append(sphere([.3,.1,.4],.05,fixed=False))
+O.bodies.append(sphere([.200001,.2000001,.4],.05,fixed=True))
+O.bodies.append(sphere([.3,0,0],.1,fixed=True))
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()],label='isc'),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()]
+ ),
+ NewtonIntegrator(gravity=(0,0,-10)),
+ #PyRunner(command='doCellFlip()',realPeriod=5)
+]
+
+def doCellFlip():
+ flip=1 if O.cell.trsf[1,2]<0 else -1;
+ flipCell(Matrix3(0,0,0, 0,0,flip, 0,0,0))
+
+#g=0.
+#while False:
+# O.cellShear=Vector3(.2*sin(g),.2*cos(pi*g),.2*sin(2*g)+.2*cos(3*g))
+# time.sleep(0.001)
+# g+=1e-3
+O.cell.trsf=Matrix3(1,0,0, 0,1,.5, 0,0,1)
+O.dt=2e-2*PWaveTimeStep()
+O.step()
+O.saveTmp()
+rrr=yade.qt.Renderer()
+rrr.intrAllWire,rrr.bound=True,True
+#isc.watch1,isc.watch2=0,-1
+
+#import yade.qt,time
+#v=yade.qt.View()
+#v.axes=True
+#v.grid=(True,True,True)
+
+
=== added file 'examples/PeriodicBoundaries/periodic-simple-shear.py'
--- examples/PeriodicBoundaries/periodic-simple-shear.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/periodic-simple-shear.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,57 @@
+# coding: utf-8
+# 2009 © Václav Šmilauer <eudoxos@xxxxxxxx>
+# 2011 ©Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>
+"Test and demonstrate use of PeriTriaxController."
+from yade import pack,qt
+
+O.periodic=True
+O.cell.setBox(.1,.1,.1)
+#O.cell.Hsize=Matrix3(0.1,0,0, 0,0.1,0, 0,0,0.1)
+sp=pack.SpherePack()
+radius=5e-3
+num=sp.makeCloud((0,0,0),(.1,.1,.1),radius,.2,500,periodic=True) # min,max,radius,rRelFuzz,spheresInCell,periodic
+O.bodies.append([sphere(s[0],s[1]) for s in sp])
+
+
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()],verletDist=.05*radius),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()]
+ ),
+ PeriTriaxController(dynCell=True,mass=0.2,maxUnbalanced=0.01,relStressTol=0.02,goal=[-1e4,-1e4,0],stressMask=3,globUpdate=5,maxStrainRate=[1.,1.,1.],doneHook='triaxDone()',label='triax'),
+ NewtonIntegrator(damping=.2),
+]
+
+phase=0
+def triaxDone():
+ global phase
+ if phase==0:
+ print 'Here we are: stress',triax.stress,'strain',triax.strain,'stiffness',triax.stiff
+ print 'Now shearing.'
+ O.cell.velGrad=Matrix3(0,6,0, 0,0,0, 0,0,0)
+ triax.stressMask=7
+ triax.goal=[-1e4,-1e4,-1e4]
+ phase+=1
+ O.pause()
+ #elif phase==1:
+ #print 'Here we are: stress',triax.stress,'strain',triax.strain,'stiffness',triax.stiff
+ #phase+=1
+ ##print 'Done, pausing now.'
+ ##O.pause()
+
+O.dt=PWaveTimeStep()
+#O.run(7000);
+#qt.View()
+##r=qt.Renderer()
+##r.bgColor=1,1,1
+#O.wait()
+#O.saveTmp()
+#O.cell.velGrad=Matrix3(0,0,0, -6,0,0, 0,0,0)
+#O.run(5000);
+#O.wait()
+
+#O.cell.velGrad=Matrix3(0,6,0, 0,0,0, 0,0,0)
+#O.run(5000);
=== added file 'examples/PeriodicBoundaries/periodic-simple.py'
--- examples/PeriodicBoundaries/periodic-simple.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/periodic-simple.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,33 @@
+"""Simple test of periodic collider.
+A few spheres falling down in gravity field and one moving accross.
+Includes a clump.
+"""
+
+from yade import timing
+
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()],label='collider'),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()],
+ ),
+ TranslationEngine(translationAxis=(1,0,0),velocity=10,ids=[0]),
+ NewtonIntegrator(damping=.4,gravity=[0,0,-10])
+]
+O.bodies.append(sphere([-4,0,11],2,fixed=True))
+O.bodies.append(sphere([0,-2,5.5],2))
+O.bodies.append(sphere([0,2,5.5],2))
+O.bodies.appendClumped([sphere([0,4,8],.8),sphere([0,5,7],.6)])
+# sets up the periodic cell
+O.periodic=True
+O.cell.setBox(10,10,10)
+# normally handled in by the simulation... but we want to have the rendering right before start
+#O.cell.postProcessAttributes()
+O.dt=.1*PWaveTimeStep()
+O.saveTmp()
+from yade import qt
+qt.Controller()
+qt.View()
+#O.timingEnabled=True; timing.reset(); O.run(200000,True); timing.stats()
=== added file 'examples/PeriodicBoundaries/periodic-triax-settingHsize.py'
--- examples/PeriodicBoundaries/periodic-triax-settingHsize.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/periodic-triax-settingHsize.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,46 @@
+# coding: utf-8
+# 2011 ©Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>
+"Demonstrate the compression of a periodic cell with non-trivial initial geometry."
+from yade import pack,qt
+O.periodic=True
+
+O.cell.hSize=Matrix3(1.0, -0.15, -0.10,
+ -0.2 ,1.5, 0.3,
+ 0.3, -0.3, 1.0)
+sp=pack.SpherePack()
+num=sp.makeCloud(hSize=O.cell.hSize, rMean=-0.01,rRelFuzz=.2, num=500,periodic=True, porosity=0.52,distributeMass=False)
+O.bodies.append([sphere(s[0],s[1]) for s in sp])
+
+
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()]),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()]
+ ),
+ #We compress the packing isotropicaly first
+ PeriTriaxController( dynCell=True,mass=0.2,maxUnbalanced=0.01, relStressTol=0.01,goal=(-1e4,-1e4,-1e4), stressMask=7,globUpdate=5, maxStrainRate=(1.,1.,1.), doneHook='triaxDone()', label='triax'),
+ NewtonIntegrator(damping=.2),
+ #PyRunner(iterPeriod=500,command='print "strain: ",triax.strain," stress: ",triax.stress')
+]
+O.dt=0.5*PWaveTimeStep()
+qt.View()
+
+phase=0
+def triaxDone():
+ global phase
+ if phase==0:
+ print 'Here we are: stress',triax.stress,'strain',triax.strain
+ #Here we reset the transformation, the compressed packing corresponds to null strain
+ O.cell.trsf=Matrix3.Identity
+ print 'Now εxx will go from 0 to .4 while σzz and σyy will be kept the same.'
+ triax.stressMask=6
+ triax.goal=(-0.4,-1e4,-1e4)
+
+ phase+=1
+ elif phase==1:
+ print 'Here we are: stress',triax.stress,'strain',triax.strain
+ print 'Done, pausing now.'
+ O.pause()
=== added file 'examples/PeriodicBoundaries/periodic-triax.py'
--- examples/PeriodicBoundaries/periodic-triax.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/periodic-triax.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,48 @@
+# coding: utf-8
+# 2009 © Václav Šmilauer <eudoxos@xxxxxxxx>
+"Test and demonstrate use of PeriTriaxController."
+from yade import pack,qt
+O.periodic=True
+
+O.cell.hSize=Matrix3(0.1, 0, 0,
+ 0 ,0.1, 0,
+ 0, 0, 0.1)
+
+sp=pack.SpherePack()
+radius=5e-3
+num=sp.makeCloud(Vector3().Zero,O.cell.refSize,radius,.2,500,periodic=True) # min,max,radius,rRelFuzz,spheresInCell,periodic
+O.bodies.append([sphere(s[0],s[1]) for s in sp])
+
+
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()],verletDist=.05*radius),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()]
+ ),
+ #PeriTriaxController(maxUnbalanced=0.01,relStressTol=0.02,goal=[-1e4,-1e4,0],stressMask=3,globUpdate=5,maxStrainRate=[1.,1.,1.],doneHook='triaxDone()',label='triax'),
+ #using cell inertia
+ PeriTriaxController(dynCell=True,mass=0.2,maxUnbalanced=0.01,relStressTol=0.02,goal=(-1e4,-1e4,0),stressMask=3,globUpdate=5,maxStrainRate=(1.,1.,1.),doneHook='triaxDone()',label='triax'),
+ NewtonIntegrator(damping=.2),
+]
+O.dt=PWaveTimeStep()
+O.run();
+qt.View()
+
+phase=0
+def triaxDone():
+ global phase
+ if phase==0:
+ print 'Here we are: stress',triax.stress,'strain',triax.strain,'stiffness',triax.stiff
+ print 'Now εz will go from 0 to .2 while σx and σy will be kept the same.'
+ triax.goal=(-1e4,-1e4,-0.2)
+ phase+=1
+ elif phase==1:
+ print 'Here we are: stress',triax.stress,'strain',triax.strain,'stiffness',triax.stiff
+ print 'Done, pausing now.'
+ O.pause()
+
+
+
=== added file 'examples/PeriodicBoundaries/periodicSandPile.py'
--- examples/PeriodicBoundaries/periodicSandPile.py 1970-01-01 00:00:00 +0000
+++ examples/PeriodicBoundaries/periodicSandPile.py 2017-02-07 16:46:13 +0000
@@ -0,0 +1,54 @@
+#!/usr/bin/python
+# -*- coding: utf-8 -*-
+# © 2012 Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx>
+"""Script showing how large bodies can be combined with periodic boundary conditions using InsertioSortCollider::allowBiggerThanPeriod=True (1)."""
+
+from yade import pack
+from pylab import rand
+from yade import qt
+
+O.periodic=True
+length=0.4
+height=0.6
+width=0.2
+thickness=0.01
+
+O.cell.hSize=Matrix3(length, 0, 0,
+ 0 ,3.*height, 0,
+ 0, 0, width)
+
+O.materials.append(FrictMat(density=1,young=1e5,poisson=0.3,frictionAngle=radians(30),label='boxMat'))
+lowBox = box( center=(length/2.0,height-thickness/2.0,width/2.0), extents=(length*1000.0,thickness/2.0,width*1000.0) ,fixed=True,wire=False)
+O.bodies.append(lowBox)
+
+radius=0.01
+O.materials.append(FrictMat(density=1000,young=1e4,poisson=0.3,frictionAngle=radians(30),label='sphereMat'))
+sp=pack.SpherePack()
+sp.makeCloud((0.*length,height+1.2*radius,0.25*width),(0.5*length,2*height-1.2*radius,0.75*width),-1,.2,2000,periodic=True)
+O.bodies.append([sphere(s[0],s[1],color=(0.6+0.15*rand(),0.5+0.15*rand(),0.15+0.15*rand())) for s in sp])
+
+O.dt=0.2*PWaveTimeStep()
+O.usesTimeStepper=True
+newton=NewtonIntegrator(damping=0.6,gravity=(0,-10,0))
+
+O.engines=[
+ ForceResetter(),
+ #(1) This is where we all big bodies, else it would crash due to the very large bottom box:
+ InsertionSortCollider([Bo1_Box_Aabb(),Bo1_Sphere_Aabb()],allowBiggerThanPeriod=True),
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
+ [Ip2_FrictMat_FrictMat_FrictPhys()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()]
+ ),
+ GlobalStiffnessTimeStepper(timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8,defaultDt=O.dt),
+ newton
+]
+
+Gl1_Sphere.stripes=1
+
+from yade import qt
+qt.View()
+print('Press PLAY button')
+
+
+
