New question #677263 on Yade:
https://answers.launchpad.net/yade/+question/677263

Hi,

I am conducting a simulation of particle breakage in a triaxial test, but at 
present I'm facing a strange problem: Normally, after the consolidation is 
completed, the deviatoric loading is entered, but deviatoric pressure is 
negative at the initial stage,  which means that the axial pressure at this 
stage is less than the confining pressure. In fact, for example, the set 
confining pressure is 100 kPa, and after the consolidation is completed, the 
pressures in all three directions are both 100 kpa, then entering deviatoric 
loading, the axial stress should be increased from 100 kpa. But in my test, 
after consolidation,  confining pressure is 100 kpa in three directions, but in 
deviatoric loading , axial stress is changed to 80 kpa , as the axial strain 
increases, the axial stress will gradually increase. When the axial strain is 
about 5%, the deviatoric stress will be greater than 0.

Do you have any idea about why deciatoric stress  is negative in the initial 
phase of the deviatoric loading?

I urgently need your help . Thanks advance!

Regards!
Cloud

this is my MWE(MWE1.py and MWE2.py are two python file, run MWE1.py first ):
# MWE1.py
from yade import pack,export,ymport,bodiesHandling
import random
mn=(0,0,0)
mx = (40e-3,40e-3,80e-3)
sp=pack.SpherePack()
sp.makeCloud(mn,mx,-1,0.15,500,False,0.95,seed=1) # 500 particles
sp.toSimulation()
export.text('/tmp/123.txt')

import random
random.seed(1)
sp = ymport.text('/tmp/123.txt')
colors = [randomColor() for s in sp]

def Children(b):
    p = b.state.refPos
    r = b.shape.radius
    s1 = sphere(p,0.5*r,color=colors[si])
    s2 = sphere(p+(.75*r,0,0),.25*r,color=colors[si])
    s3 = sphere(p+(-.75*r,0,0),.25*r,color=colors[si])
    s4 = sphere(p+(0,.75*r,0),.25*r,color=colors[si])
    s5 = sphere(p+(0,-.75*r,0),.25*r,color=colors[si])
    return s1,s2,s3,s4,

for si,i in enumerate(sp):
        oriBody = Quaternion(Vector3(1,0,0),random.random()*(pi/2))
        oriBody *= Quaternion(Vector3(0,1,0),random.random()*(pi/2))
        oriBody *= Quaternion(Vector3(0,0,1),random.random()*(pi/2))
        b=O.bodies[si]
        children = Children(b)
        O.bodies.erase(si)
        a=bodiesHandling.spheresModify(children,orientation=oriBody,copy=True)
        ids=O.bodies.append(a)
        for j in ids:
                O.bodies[j].agglomerate=si

export.textExt('/tmp/456.txt','x_y_z_r_attrs',attrs=['b.agglomerate'])





#MWE2.py
from yade import ymport
key='haxi_'
mn = (0,0,0)
mx = (40e-3,40e-3,80e-3)
young = 300e6
compFricDegree = 30
damp = 0.6
stabilityThreshold = 0.01
rate = -0.02
finalFricDegree = 20
targetPorosity = 0.45
poisson = 0.22
confiningS = 100e3

O.materials.append(FrictMat(
    young=young,
    poisson=poisson,
    frictionAngle=0,
    density=0,
    label='walls')
    )

walls = aabbWalls([mn,mx],thickness=0,material='walls',wire=True)
wallIds = O.bodies.append(walls)

O.materials.append(CohFrictMat(
    young=young,
    poisson=poisson,
    isCohesive=True,
    frictionAngle=radians(compFricDegree),
    density=2650,
    momentRotationLaw=True,
    normalCohesion=1e20, # Cn
    shearCohesion=1e20, # Cs
    alphaKr=0.2, # rolling stiffness
    label='spheres')
    )

attrs = []
sp = ymport.textExt('/tmp/456.txt',format='x_y_z_r_attrs',attrs=attrs)
n = max(int(a[0]) for a in attrs)+1
colors = [randomColor() for _ in xrange(n)]
for s,a in zip(sp,attrs):
        aa = int(a[0])
        s.agglomerate = aa
        s.shape.color = colors[aa]
O.bodies.append(sp) # the last defined material is usedself.


triax=TriaxialStressController(
        maxMultiplier=1.+2e4/young, # spheres growing factor (fast growth)
        finalMaxMultiplier=1.+2e3/young, # spheres growing factor (slow growth)
        thickness = 0,
        stressMask = 7,
        internalCompaction=False,
)


factor=1.05
O.engines=[
        ForceResetter(),
        
InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=factor,label='bo1aabbs'),Bo1_Box_Aabb()]),
        InteractionLoop(
                
[Ig2_Sphere_Sphere_ScGeom6D(interactionDetectionFactor=factor,label='ig2sss'),Ig2_Box_Sphere_ScGeom()],
        [Ip2_FrictMat_FrictMat_FrictPhys() 
,Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True)],
        
[Law2_ScGeom_FrictPhys_CundallStrack(),Law2_ScGeom6D_CohFrictPhys_CohesionMoment(
            useIncrementalForm=True,
            label='conhesiveLaw')]
        ),
    
GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=1000,timestepSafetyCoefficient=0.5),
    triax,
    TriaxialStateRecorder(iterPeriod=100,file='WallStresses'+key,truncate=1),
    NewtonIntegrator(damping=damp),
]

O.run(1,True)

ig2sss.interactionDetectionFactor = bo1aabbs.aabbEnlargeFactor = 1
for i in O.interactions:
    if not isinstance(i.phys,CohFrictPhys):
        continue
    b1,b2 = [O.bodies[ii] for ii in (i.id1,i.id2)]
    if b1.agglomerate != b2.agglomerate:
        O.interactions.erase(i.id1,i.id2)
    else:
        i.phys.normalAdhesion = 2.9  # normalAdhesion
        i.phys.shearAdhesion = 6.5  # shearAdhesion
        i.phys.unp = i.geom.penetrationDepth

O.run(1,True)

yade.qt.Controller(), yade.qt.View()

#######################################
###   APPLYING CONFINING PRESSURE   ###
#######################################
triax.goal1=triax.goal2=triax.goal3=-confiningS

while 1:
 O.run(1000, True)
 unb=unbalancedForce()
 print 'unbalanced force:',unb,' mean stress: ',triax.meanStress
 if unb<stabilityThreshold and 
abs(-confiningS-triax.meanStress)/confiningS<0.001:
    break

O.save('confinedState'+key+'.yade.gz')
print "###      Isotropic state saved      ###"

triax.width0=triax.width
triax.depth0=triax.depth
triax.height0=triax.height
#############################
##   DEVIATORIC LOADING   ###
###############################
setContactFriction(radians(finalFricDegree))
triax.stressMask = 3
triax.goal3= rate
triax.goal1=-confiningS
triax.goal2=-confiningS


O.saveTmp()


from yade import plot

## a function saving variables
def history():
        plot.addData(
            e11=-triax.strain[0],
            e22=-triax.strain[1],
            e33=-triax.strain[2],
                    ev=triax.volumetricStrain,
                    s11=-triax.stress(triax.wall_right_id)[0],
                    s22=-triax.stress(triax.wall_top_id)[1],
                    s33=-triax.stress(triax.wall_front_id)[2],
            devi = -triax.stress(triax.wall_front_id)[2] + 
triax.stress(triax.wall_top_id)[1],
                    i=O.iter)

def stoploading():
    axial_strain = abs(triax.strain[2])
    if axial_strain > 0.2:
        print 'axial strain reach Limitation, stop loading!'
        O.pause()

if 1:
    
O.engines=O.engines+[PyRunner(iterPeriod=500,command='history()',label='recorder')]
    O.engines=O.engines+[PyRunner(command='stoploading()',iterPeriod=1000)]

##O.run(100,True)

plot.plots={'e33':('ev',),' e33 ':('devi')}

# display on the screen (doesn't work on VMware image it seems)
plot.plot()

# In that case we can still save the data to a text file at the the end of the 
simulation, with:
plot.saveDataTxt('results'+key)
#or even generate a script for gnuplot. Open another terminal and type  
"gnuplot plotScriptKEY.gnuplot:
plot.saveGnuplot('plotScript'+key)


-- 
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

Reply via email to