yade-users team mailing list archive

[nie jiayan <question679311@xxxxxxxxxxxxxxxxxxxxx>]

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

Hello, everyone. Now I want to do 2D DEM simulation with clumps in periodic boundary condition. But some confusions exits. Why when I run the simulation, the interactions between the constituent balls in clumps are also active. And the PeriTriaxController() failed. Anyone can help me. Many thanks in advance.
The following are my codes.
from yade import pack, qt, plot
import time
import numpy as np

############################################
###   Read external file information     ###
############################################

filename1='pebblex.txt'
filename2='pebbley.txt'
filename3='pebbler.txt'
filename4='clumpn.txt'

pbx=np.loadtxt(filename1)
pby=np.loadtxt(filename2)
pbr=np.loadtxt(filename3)
cmn=np.loadtxt(filename4)

cn=cmn.size-1

## create materials for spheres for initial consolidation
O.materials.append(FrictMat(density=2650,young=6.0e8,poisson=0.8,frictionAngle=0.01))

## process the external file information to get the clump configuration list

for ii in range(0,cn):
    aa=int(np.sum(cmn[:ii+1]))
    bb=int(np.sum(cmn[:ii+2]))
    bodyList=[]
    for jj in range(aa,bb):
        bodyList.append(O.bodies.append(sphere([pbx[jj],pby[jj],0.05],pbr[jj])))
    O.bodies.clump(bodyList,discretization=10)
############################################
###   DEFINING VARIABLES AND MATERIALS   ###
############################################

ts = time.time()
pressure = -1e5
size=0.30

## setup the periodic boundary
O.periodic=True
O.cell.hSize = Matrix3(size,0,0, 0,size,0, 0,0,0.1)

for p in O.bodies:
    if isinstance(p.shape, Clump):
        p.shape.color=(0.25,0.5,0.75)
        p.state.blockedDOFs = 'zXY'

        

print(len(O.bodies))
############################
###   DEFINING ENGINES   ###
############################

triax=PeriTriaxController(
		# specify target values and whether they are strains or stresses
		goal=(pressure,pressure,0),stressMask=3,
		# type of servo-control
		dynCell=True,maxStrainRate=(0.5,0.5,0),
		# wait until the unbalanced force goes below this value
		maxUnbalanced=0.001,relStressTol=1e-3,
                # call this function when goal is reached and the packing is stable
		doneHook='consolidationFinished()'
	)

newton=NewtonIntegrator(damping=.1)

O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb()]),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom(),],
		[Ip2_FrictMat_FrictMat_FrictPhys(),],
		[Law2_ScGeom_FrictPhys_CundallStrack(),]
	),
	triax,
	# TriaxialStateRecorder(iterPeriod=100,file='WallStresses'+table.key),
	newton,
        PyRunner(command='addPlotData()',iterPeriod=100),
]

O.dt=.5*PWaveTimeStep()

def addPlotData():
   
   pp = utils.porosity()  #this is the porosity of the cell. 
   ee = pp / (1-pp)  #this is the void ratio of the 3D cell.  
   a = [i for i in O.bodies if isinstance(i.shape, Clump) and i.intrs()==[]]
   plot.addData(unbalanced=unbalancedForce(),i=O.iter,   
                e11=-triax.strain[0], e22=-triax.strain[1], e33=-triax.strain[2],
		s11=-triax.stress[0],
		s22=-triax.stress[1],
		s33=-triax.stress[2],
                void = ee, floatPt=len(a),
   )

# define what to plot//since they have the same x (i), the latter i should have a space in it 'i '
plot.plots={'i':('unbalanced',),'i ':('s11','s22','s33'),' i':('e11','e22','e33'),
   ' i ':('void',),'i  ':('floatPt',),
   #each 'i' should have different number or location of space. Or the latter will overwrite the former one.
}
# show the plot
plot.plot()

yade.qt.Controller(), yade.qt.View()


def consolidationFinished():
        # set the current cell configuration to be the reference one
	O.cell.trsf=Matrix3.Identity
	print 'consolidationFinished'
        print('time: '+str((time.time()-ts)/60.)+'min')
        O.save('dense_packing_100kPa.yade.gz')
	O.pause()

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