New question #694002 on Yade: https://answers.launchpad.net/yade/+question/694002
Hello, I wish you and your family perfect health. I am modeling a triaxial test. The problem that I have is that at the beginning of the simulation, some particles come out from the walls. Therefore, it destroys the shape of micro-strain contours. In addition, the initial stresses which are sensed at the walls are so high (more than 1 MPa). I am using YADE 2020.01a, and Ubuntu 20.04. My python version is 3.8.5. Attached, you can find my script. I know it should not be too long, but this is the problem that I have. I was wondering if you could give me a solution. Thank you very much. Best Regards, Alireza from yade import utils, plot from yade import pack, qt from datetime import datetime import math import pylab from numpy import * qtr=qt.Renderer() qtr.bgColor=(1,1,1) #============================================================== #================= define the materials ======================= #============================================================== O.materials.append(CohFrictMat(normalCohesion= 1e20, shearCohesion= 1e20, isCohesive= True, young=6.81e8, density=1377.5, poisson=0.3, frictionAngle= 0.31, fragile=False, label='aggregate-48')) O.materials.append(CohFrictMat(normalCohesion= 1e20, shearCohesion= 1e20, isCohesive= True, young=6.81e8, density=1532.2, poisson=0.3, frictionAngle= 0.31, fragile=False, label='aggregate-814')) O.materials.append(CohFrictMat(normalCohesion= 1e20, shearCohesion= 1e20, isCohesive= False, young=2e11, density=0, poisson=0.3, frictionAngle= 0, fragile=False, label='wall-frictionless')) #=============================================================== #=================== define walls ============================ #=============================================================== walls=aabbWalls([(0,0,0),(10e-2,10e-2,10e-2)],thickness=0.001,oversizeFactor=1.5,material='wall-frictionless') wallIds=O.bodies.append(walls) #=============================================================== #=================== define packing ============================ #=============================================================== nums=['t','t','t'] mats=['aggregate-48','aggregate-814','aggregate-814'] coke=((1.875e-3,3*1000),(0.9475e-3,3*2734),(0.50125e-3,3*24970)) color=((0,0,1),(0,1,0),(1,0,0)) tolerance=[(2e-4),(1e-4),(5e-5)] for i in range(len(nums)): nums[i]=pack.SpherePack() nums[i].makeCloud((5e-3,5e-3,5e-3),(9.95e-2,9.95e-2,9.95e-2),rMean=coke[i][0],rRelFuzz=tolerance[i],num=coke[i][1]) O.bodies.append([utils.sphere(c,r,material=mats[i],color=color[i]) for c,r in nums[i]]) #=============================================================== #=================== define Engine ============================= #=============================================================== O.engines=[ ForceResetter(), InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb(),Bo1_Facet_Aabb()]), InteractionLoop([Ig2_Sphere_Sphere_ScGeom6D(),Ig2_Box_Sphere_ScGeom6D(),Ig2_Facet_Sphere_ScGeom6D()], [Ip2_CohFrictMat_CohFrictMat_CohFrictPhys()], [Law2_ScGeom_FrictPhys_CundallStrack(),Law2_ScGeom6D_CohFrictPhys_CohesionMoment()] ), TriaxialStateRecorder(iterPeriod=10000,file='WallStresses'), # VTKRecorder(fileName='3d-vtk-',recorders=['all'],iterPeriod=1000), # qt.SnapshotEngine(fileBase='3d-',iterPeriod=200,label='snapshot'), NewtonIntegrator(damping=0.4,gravity=[0,0,-9.81]), PyRunner(iterPeriod=10000,command='addPlotData()'), ] #yade.qt.Controller(), yade.qt.View() O.dt=5e-8 #O.dt=.2*PWaveTimeStep() def addPlotData(): plot.addData(exx=-triax.strain[0], eyy=-triax.strain[1], ezz=-triax.strain[2], ev=(triax.strain[0]+triax.strain[1]+triax.strain[2]), sxx=-triax.stress(triax.wall_right_id)[0], syy=-triax.stress(triax.wall_top_id)[1], szz=-triax.stress(triax.wall_front_id)[2], mass=sum(b.state.mass for b in O.bodies if isinstance(b.shape,Sphere)), q=-(triax.stress(triax.wall_front_id)[2]-triax.stress(triax.wall_right_id)[0]), p=triax.meanStress, R=triax.stress(triax.wall_front_id)[2]/sigmaIso, por=triax.porosity, i=O.iter,**O.energy,total=O.energy.total()), #================================================================= #================= APPLYING CONFINING PRESSURE ================ #================================================================= sigmaIso=-5e4 triax=TriaxialStressController( maxMultiplier=1.000, finalMaxMultiplier=1.000, thickness = 0, stressMask = 7, internalCompaction=False, max_vel=0.05, ) O.engines=O.engines+[triax] triax.goal1=sigmaIso triax.goal2=sigmaIso triax.goal3=sigmaIso triax.wall_back_activated=True while 1: O.run(100, True) unb=unbalancedForce() meanS=(triax.stress(triax.wall_right_id)[0]+triax.stress(triax.wall_top_id)[1]+triax.stress(triax.wall_front_id)[2])/3 print('mean Stress:',triax.meanStress,'porosity:', triax.porosity,'Unbalancing:',unb) if triax.porosity<0.37 and abs(sigmaIso-triax.meanStress)/abs(sigmaIso)<0.02: print('Isotropic strain1:',-triax.strain[0], 'Isotropic strain 2:',-triax.strain[1], 'Isotropic strain 3:', -triax.strain[2]) break print("### Isotropic state saved ###") triax.depth0=triax.depth triax.height0=triax.height triax.width0=triax.width O.save('RVE3000-sizeDis3-solid-Iso5e4-Rate005-Isopart.yade') #==================================================================== #===================== DEVIATORIC LOADING ======================= #==================================================================== triax.stressMask = 3 triax.goal1 = sigmaIso triax.goal2 = sigmaIso triax.goal3 = -0.05 triax.max_vel=0.05 O.trackEnergy=True O.saveTmp() #==================================================================== #==================== Micro Strain =========================== #==================================================================== TW=TesselationWrapper() TW.triangulate() TW.computeVolumes() TW.volume(10) TW.setState(0) step=0 while 1: step +=1 O.run(1000000,True) print('Dev strain1:',-triax.strain[0], 'Dev strain 2:',-triax.strain[1], 'Dev strain 3:', -triax.strain[2]) TW.setState(1) TW.defToVtk("strain_"+str(step)+".vtk") if abs(triax.strain[2]) > 0.5: # makeVideo(snapshot.snapshots,'3d.mpeg',fps=10,bps=10000) print(O.iter) print('time is ', O.time) print('porosity:',triax.porosity) print('strain 1:',-triax.strain[0], 'strain 2:',-triax.strain[1], 'strain 3:',-triax.strain[2]) O.save('RVE3000-sizeDis3-solid-Iso5e4-Rate005-Devpart.yade') plot.saveDataTxt('RVE3000-sizeDis3-solid-Iso5e4-Rate005-Devpart.txt') break O.pause() -- 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