New question #673124 on Yade: https://answers.launchpad.net/yade/+question/673124
I tried to make the displacement field of the second stage of triaxial compression. I have this answer in the following question, but I have not implemented this method. Here is my code, please tell me where the error is. https://answers.launchpad.net/yade/+question/478212. from yade import pack,plot,qt,export nRead=readParamsFromTable( num_spheres=1000, compFricDegree = 30, key='_triax_base_', unknownOk=True ) from yade.params import table num_spheres=table.num_spheres# number of spheres key=table.key targetPorosity = 0.42 compFricDegree = table.compFricDegree finalFricDegree = 30 rate=-0.02 damp=0.2 stabilityThreshold=0.01 young=5e8 mn,mx=Vector3(0,0,0),Vector3(1,2,1) # corners of the initial packing weiya=100000 poisson=0.1 O.materials.append(FrictMat(young=young,poisson=poisson,frictionAngle=radians(compFricDegree),density=2650,label='spheres')) O.materials.append(FrictMat(young=young,poisson=poisson,frictionAngle=0,density=0,label='walls')) walls=aabbWalls([mn,mx],thickness=0,material='walls') wallIds=O.bodies.append(walls) psdSizes,psdCumm =[0.19,0.191,0.195,0.196,0.2],[0.,0.1,0.5,0.6,1.] #[0.01,0.03,0.11,0.13,0.2],[0.,0.1,0.5,0.6,1.] # [0.05,0.07,0.14,0.155,0.2],[0.,0.1,0.5,0.6,1.] #[0.1,0.11,0.16,0.17,0.2],[0.,0.1,0.5,0.6,1.] #[0.15,0.16,0.184,0.19,0.2],[0.,0.1,0.5,0.6,1.] #[0.19,0.191,0.195,0.196,0.2],[0.,0.1,0.5,0.6,1.] sp=pack.SpherePack(); sp.makeCloud(mn,mx,num=num_spheres, psdSizes=psdSizes, psdCumm=psdCumm,distributeMass=True,porosity=0.42) O.bodies.append([sphere(center, rad, material='spheres') for center, rad in sp]) triax=TriaxialStressController( maxMultiplier=1.+2e4/young, finalMaxMultiplier=1.+2e3/young, thickness = 0, stressMask = 7, internalCompaction=False, ) newton=NewtonIntegrator(damping=damp) O.engines=[ ForceResetter(), InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]), InteractionLoop( [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()], [Ip2_FrictMat_FrictMat_FrictPhys()], [Law2_ScGeom_FrictPhys_CundallStrack()] ), GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8), triax, TriaxialStateRecorder(iterPeriod=100,file='WallStresses'+table.key), newton, ] yade.qt.Controller(), yade.qt.View() triax.goal1=triax.goal2=triax.goal3=-weiya while 1: O.run(1000, True) unb=unbalancedForce() print 'unbalanced force:',unb,' mean stress: ',triax.meanStress if unb<stabilityThreshold and abs(-weiya-triax.meanStress)/weiya<0.001: break print "### Isotropic state saved ###" import sys while triax.porosity>targetPorosity: compFricDegree = 0.95*compFricDegree setContactFriction(radians(compFricDegree)) print "\r Friction: ",compFricDegree," porosity:",triax.porosity, sys.stdout.flush() O.run(500,1) print 'Box Volume: ', triax.boxVolume print "### Compacted state saved ###" triax.internalCompaction=False setContactFriction(radians(finalFricDegree)) triax.stressMask = 5 triax.goal1=-weiya triax.goal2=rate triax.goal3=-weiya newton.damping=0.1 O.saveTmp() vtk = export.VTKExporter('/home/ma/dspl/jan-') vtk.exportSpheres(useRef=True,what=[('dspl','b.state.displ()')]) -- 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 : [email protected] Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
