Question #376480 on Yade changed: https://answers.launchpad.net/yade/+question/376480
Weimin Song gave more information on the question: Code above has bugs. Below is the modified code: readParamsFromTable(rMean=.05,rRelFuzz=.3,maxLoad=2e8,minLoad=10) from yade.params.table import * from yade import pack, plot aggregate=O.materials.append(FrictMat(density=3e3, young=30e9,poisson=.3,frictionAngle=.5,label="aggregate")) O.bodies.append(geom.facetBox((.5,.5,.5),(.5,.5,.5),wallMask=31)) sp=pack.SpherePack() sp.makeCloud((0,0,0),(1,1,1),rMean=rMean,rRelFuzz=rRelFuzz) sp.toSimulation(material=aggregate) #sp.toSimulation() wallMat=O.materials.append(FrictMat(density=4e3, young=500e9, poisson=0.1, frictionAngle=0.5, label="wallMat")) wall=O.bodies.append(wall((.5,.5,1),axis=2,sense=-1, color=(0,0,1), material=wallMat)) O.engines=[ ForceResetter(), InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]), InteractionLoop( [Ig2_Sphere_Sphere_L3Geom(),Ig2_Facet_Sphere_L3Geom(),Ig2_Wall_Sphere_L3Geom()], [Ip2_FrictMat_FrictMat_FrictPhys()], [Law2_L3Geom_FrictPhys_ElPerfPl()] ), NewtonIntegrator(gravity=(0,0,-9.81),damping=0.5), PyRunner(command='checkUnbalanced()',realPeriod=2,label='checker'), ] O.dt=.5*PWaveTimeStep() def checkUnbalanced(): if O.iter<10000: return if unbalancedForce()>0.2: return # O.bodies.append(wall(max([b.state.pos[2]+b.shape.radius for b in O.bodies if isinstance(b.shape,Sphere)]),axis=2,sense=-1, material=wallMat)) global plate # without this line, the plate variable would only exist inside this function plate=O.bodies[-1] # the last particles is the plate plate.state.vel=(0,0,-.4) O.engines=O.engines+[PyRunner(command='addPlotData()',iterPeriod=50)] checker.command='unloadPlate()' def unloadPlate(): if abs(O.forces.f(plate.id)[2])>maxLoad: plate.state.vel*=-1 checker.command='stopUnloading()' def stopUnloading(): if abs(O.forces.f(plate.id)[2])<minLoad: plot.saveDataTxt('1.txt') O.pause() ######remove the facet and append the box O.bodies.erase(wall) boxMat=O.materials.append(FrictMat(density=4e3, young=500e9, poisson=0.1, frictionAngle=0.5, label="boxMat")) O.bodies.append(box((.5,.5,1.05),(0.1,0.2,0.05),orientation=Quaternion((1, 0, 0), 0), dynamic=True, fixed=False, wire=False, color=(0.5,0.5,0.5), highlight=False, material=boxMat)) O.engines=[ ForceResetter(), # sphere, facet, box InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Box_Aabb()]), InteractionLoop( # the loading plate is a wall, we need to handle sphere+sphere, sphere+facet, sphere+wall [Ig2_Sphere_Sphere_L3Geom(),Ig2_Facet_Sphere_L3Geom(),Ig2_Box_Sphere_ScGeom()], [Ip2_FrictMat_FrictMat_FrictPhys()], [Law2_ScGeom_FrictPhys_CundallStrack()] ), NewtonIntegrator(gravity=(0,0,-9.81),damping=0.4), PyRunner(command='checkUnbalanced1()',realPeriod=2,label='checker1'), ] O.dt=.05*PWaveTimeStep() def checkUnbalanced1(): global box # without this line, the box variable would only exist inside this function box=O.bodies[-1] # the last particles is the box box.state.vel=(0,0,-1) O.engines=O.engines+[PyRunner(command='addPlotData()',iterPeriod=50)] checker.command='unloadBox()' def unloadBox(): if abs(O.forces.f(box.id)[2])>maxLoad: box.state.vel*=-1 checker.command='stopUnloading1()' def stopUnloading1(): if abs(O.forces.f(plate.id)[2])<minLoad: plot.saveDataTxt('2.txt') O.pause() def addPlotData(): if not isinstance(O.bodies[-1].shape,Wall): plot.addData(); return Fz=O.forces.f(plate.id)[2] plot.addData(Fz=Fz,w=plate.state.pos[2]-plate.state.refPos[2],unbalanced=unbalancedForce(),i=O.iter) # besides unbalanced force evolution, also plot the displacement-force diagram plot.plots={'i':('unbalanced',),'w':('Fz',)} plot.plot() O.run() waitIfBatch() -- You received this question notification because your team yade-users is an answer contact for Yade. _______________________________________________ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp