New question #696648 on Yade: https://answers.launchpad.net/yade/+question/696648
Hi, How can I erase "walls" and the "plate" once my simulation is done (after o.pause). I also would like take a image of my simulation after o.pause. Here is my code: #!/usr/bin/env python #encoding: ascii # Testing of the Deformation Enginge with Luding Contact Law # Modified Oedometric Test # The reference paper [Haustein2017] from __future__ import print_function from yade import utils, plot, timing from yade import pack import pandas as pd o = Omega() # Physical parameters fr = 0.3 rho = 2000 Diameter = 16.5e-3 r1 = Diameter r2 = Diameter k1 = 1005.0 kp = 10.0*k1 kc = k1 * 0.0 ks = k1 * 0.1 DeltaPMax = Diameter/3.0 Chi1 = 0.34 o.dt = 1.0e-5 particleMass = 4.0/3.0*math.pi*r1*r1*r1*rho Vi1 = math.sqrt(k1/particleMass)*DeltaPMax*Chi1 PhiF1 = DeltaPMax*(kp-k1)*(r1+r2)/(kp*2*r1*r2) #************************************* # Add material mat1 = O.materials.append(LudingMat(frictionAngle=fr, density=rho, k1=k1, kp=kp, ks=ks, kc=kc, PhiF=PhiF1, G0 = 0.0)) # Spheres for compression sp=pack.SpherePack() sp.makeCloud((-11.0*Diameter,-11.0*Diameter,-25*Diameter),(11*Diameter,11*Diameter,30.0*Diameter), rMean=Diameter/2.0,rRelFuzz=.18) cyl = pack.inCylinder((0,0,0),(0,0,30*Diameter),11*Diameter-0.006) sp = pack.filterSpherePack(cyl,sp,True,material=mat1) sp.toSimulation(material=mat1) ##No of particles count = 0 for b in O.bodies: if isinstance(b.shape, Sphere): count +=1 ###################################################################### walls=O.bodies.append(yade.geom.facetCylinder((0,0,0),radius=11*Diameter,height=11*Diameter,segmentsNumber=20,wallMask=6,material=mat1)) df = pd.DataFrame(columns=['Cut_xy','Cut_z','Volume','Porosity']) df.to_csv('PH101_rev.csv') from csv import writer def append_list_as_row(file_name, list_of_elem): # Open file in append mode with open(file_name, 'a+', newline='') as write_obj: # Create a writer object from csv module csv_writer = writer(write_obj) # Add contents of list as last row in the csv file csv_writer.writerow(list_of_elem) # Add engines o.engines = [ ForceResetter(), InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=1.05), Bo1_Wall_Aabb(), Bo1_Facet_Aabb() ]), InteractionLoop( [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=1.05), Ig2_Facet_Sphere_ScGeom(), Ig2_Wall_Sphere_ScGeom()], [Ip2_LudingMat_LudingMat_LudingPhys()], [Law2_ScGeom_LudingPhys_Basic()] ), NewtonIntegrator(damping=0.1, gravity=[0, 0, -9.81]), #VTKRecorder(fileName='vtk-',recorders=['all'],iterPeriod=10000), PyRunner(command='checkForce()', realPeriod=1, label="fCheck"), DeformControl(label="DefControl") ] def checkForce(): # at the very start, unbalanced force can be low as there is only few # contacts, but it does not mean the packing is stable if O.iter < 20000: return # the rest will be run only if unbalanced is < .1 (stabilized packing) timing.reset() if unbalancedForce() > 0.3: return # add plate at upper box side highSphere = 0.0 for b in O.bodies: if highSphere < b.state.pos[2] and isinstance(b.shape, Sphere): highSphere = b.state.pos[2] else: pass O.bodies.append(wall(highSphere+0.5*Diameter, axis=2, sense=-1, material=mat1)) # without this line, the plate variable would only exist inside this # function global plate plate = O.bodies[-1] # the last particles is the plate # Wall objects are "fixed" by default, i.e. not subject to forces # prescribing a velocity will therefore make it move at constant velocity # (downwards) plate.state.vel = (0, 0, -.1) # start plotting the data now, it was not interesting before O.engines = O.engines + [PyRunner(command='addPlotData()', iterPeriod=1000)] # next time, do not call this function anymore, but the next one # (unloadPlate) instead fCheck.command = 'unloadPlate()' def unloadPlate(): # if the force on plate exceeds maximum load, start unloading # if abs(O.forces.f(plate.id)[2]) > 5e-2: if abs(O.forces.f(plate.id)[2]) > 4.3e3: plate.state.vel *= -1 # next time, do not call this function anymore, but the next one # (stopUnloading) instead fCheck.command = 'stopUnloading()' def stopUnloading(): if abs(O.forces.f(plate.id)[2]) == 0: # O.tags can be used to retrieve unique identifiers of the simulation # if running in batch, subsequent simulation would overwrite each other's output files otherwise # d (or description) is simulation description (composed of parameter values) # while the id is composed of time and process number # plot.saveDataTxt(O.tags['d.id'] + '.txt') plot.saveDataTxt('data'+ O.tags['id'] +'.txt') #print(timing.stats()) O.pause() #Want Erase wall and plate #Take snap of my simualation Diameter_cut=[0.2,0.5,1,1.5,2,2.5,3,3.5,4,4.5,5] height_tab=utils.aabbExtrema()[1][2]-utils.aabbExtrema()[0][2]/(2*len(Diameter_cut)) D=Diameter for i in range (0,len(Diameter_cut)): Cut_xy=Diameter_cut[i] Start_point=utils.aabbExtrema()[0]+(Cut_xy*D,Cut_xy*D,(i+1)*height_tab) End_point=utils.aabbExtrema()[1]-(Cut_xy*D,Cut_xy*D,(i+1)*height_tab) vol=(End_point[0]-Start_point[0])*(End_point[1]-Start_point[1])*(End_point[2]-Start_point[2]) Porosity_rev= [voxelPorosity(resolution=1600,start=Start_point,end=End_point)] row_contents = [i,Start_point,End_point,vol,Porosity_rev] append_list_as_row('PH101_rev.csv', row_contents) 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] - (-4*Diameter), unbalanced=unbalancedForce(), i=O.iter ) def defVisualizer(): with open("data.dat","a") as f: for b in O.bodies: if isinstance(b.shape, Sphere): rData = "{x},{y},{z},{r},{w}\t".format(x = b.state.pos[0], y = b.state.pos[1], z = b.state.pos[2], r = b.shape.radius + b.state.dR, w = plate.state.pos[2] ) f.write(rData) f.write("\n") O.timingEnabled=True O.run(1, True) plot.plots={'w':('Fz', None)} plot.plot() -- 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