Question #675926 on Yade changed: https://answers.launchpad.net/yade/+question/675926
yangjunjie posted a new comment: Hi Robert: Thank you for your reply !! But I meet some problems: I added the code you provided to mine, there seems to be some conflict between plt.show() and plot.plot(), Or my understanding is wrong. Help me plz.. This is my whole code: readParamsFromTable(rMean=.05,rRelFuzz=.3,maxLoad=1e6,minLoad=1e4) from yade.params.table import * from yade import pack, plot 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) sp.toSimulation() O.engines=[ ForceResetter(), InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]), InteractionLoop( [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()], [Ip2_FrictMat_FrictMat_FrictPhys()], [Law2_ScGeom_FrictPhys_CundallStrack()] ), NewtonIntegrator(gravity=(0,0,-9.81),damping=0.5), PyRunner(command='checkUnbalanced()',realPeriod=2,label='checker'), PyRunner(command='replace()',iterPeriod=100,label ='rep'), #PyRunner(command='getRadiusArray()',iterPeriod=100,label ='psd'), ] O.dt=.5*PWaveTimeStep() ############################################################## def selectBodiesToCrush(): stresses = bodyStressTensors() # e.g. based on overall particle stress ret = [] for i,stress in enumerate(stresses): if stress[0,0] > 50: # Particle breakage criterion ret.append(i) return ret def createChildren(b): p = b.state.pos r = b.shape.radius # Particle replacement mode s1 = sphere(p+(.5*r,0,0),.5*r) s2 = sphere(p-(.5*r,0,0),.5*r) return s1,s2 def replace(): bodiesToDelete = selectBodiesToCrush() for i in bodiesToDelete: if isinstance(O.bodies[i].shape,Sphere): # Limitation of the smallest particles that can break b = O.bodies[i] children = createChildren(b) O.bodies.erase(i) O.bodies.append(children) ############################################################## #import matplotlib.pyplot as plt #import numpy as np #def getRadiusArray(): # radii=[] # for b in O.bodies: # if isinstance(b.shape,Sphere): # radii.append(b.shape.radius) # return radii # #radii=getRadiusArray() #values, base = np.histogram(radii, bins=40) #cumulative = np.array(np.cumsum(values),dtype=float) #perc = cumulative/max(cumulative)*100 #plt.semilogx(base[:-1], cumulative, c='blue') #plt.show() ############################################################## def checkUnbalanced(): if O.iter<5000: return if unbalancedForce()>.1: 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)) global plate plate=O.bodies[-1] plate.state.vel=(0,0,-.1) O.engines=O.engines+[PyRunner(command='addPlotData()',iterPeriod=200)] 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(O.tags['d.id']+'.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) plot.plots={'i':('unbalanced',),'w':('Fz',)} plot.plot() O.run() #from yade import utils #binsSizes, binsProc, binsSumCum= utils.psd(bins=5, mass=True) 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