yade-users team mailing list archive

[xuelong14 <question671283@xxxxxxxxxxxxxxxxxxxxx>]

Question #671283 on Yade changed:
https://answers.launchpad.net/yade/+question/671283

Description changed to:
Hello everyone,

I did my triaxial compression tests under different confining pressures.
Sigma3 in the horizontal directions is constant all through the test and
the cubic sample is compressed in the vertical direction with a constant
strain rate.  Sigma3 is equal to 100kPa or 400kPa in my tests. Both
samples contract first and dilate later as the compression goes on.

As well all know, larger volume contraction and lower volume dilation
afterwards will occur through the test under higher sigma3. But my
results are just reverse to this. In the case when sigma3 is 100kPa, the
volume contraction is larger and the volume dilation is maller
afterwards. I really don't know why this happens.

My sample is trimmed from a big gravity deposited sample. My friend
obtained his sample by expanding spheres inside a cubic box and found
the same weird deformation.  We use the same contact law.  We have
adjusted the elastic modulus and the poisson ratio, but it doesn't work.

I did triaxial tests with PFC with linear contact model under different
confining pressures and everything is fine.

Does anyone else meet this kind of problem?

Below is my code:

from yade import ymport
from yade import utils

############################################
###   DEFINING VARIABLES AND MATERIALS   ###
############################################

compFricDegree = 26.565 # initial contact friction during the confining phase 
finalFricDegree = 26.565 # contact friction during the deviatoric loading
rate=-0.01 # loading rate (strain rate)
damp=0.7 # damping coefficient
stabilityThreshold=0.01 # we test unbalancedForce against this value in different loops (see below)
young=5e8 # contact stiffness
mn,mx=Vector3(0,0,0),Vector3(1,1,1) # corners of the initial packing


## create materials for spheres and plates
O.materials.append(FrictViscoMat(betan=0.2, young=young,poisson=0.3,frictionAngle=radians(compFricDegree),density=2600,label='spheres'))
O.materials.append(FrictViscoMat(betan=0.2, young=young,poisson=0.3,frictionAngle=0,density=0,label='walls'))

## create walls around the packing
walls=aabbWalls([mn,mx],thickness=0,material='walls')
wallIds=O.bodies.append(walls)

print "Load clumps"
deposited_bodies = ymport.textClumps('cubic_sample_clump_miu0.5_dep_damp0.1_35608_hor_layer_plane.txt', shift=Vector3(0,0,0), material='spheres')
O.bodies.updateClumpProperties()#get more accurate clump masses/volumes/inertia


############################
###   DEFINING ENGINES   ###
############################

triax=TriaxialStressController(
	## TriaxialStressController will be used to control stress and strain. It controls particles size and plates positions.
	## this control of boundary conditions was used for instance in http://dx.doi.org/10.1016/j.ijengsci.2008.07.002
	maxMultiplier=1.+2e4/young, # spheres growing factor (fast growth)
	finalMaxMultiplier=1.+2e3/young, # spheres growing factor (slow growth)
	thickness = 0,
	## switch stress/strain control using a bitmask. What is a bitmask, huh?!
	## Say x=1 if stess is controlled on x, else x=0. Same for for y and z, which are 1 or 0.
	## Then an integer uniquely defining the combination of all these tests is: mask = x*1 + y*2 + z*4
	## to put it differently, the mask is the integer whose binary representation is xyz, i.e.
	## "100" (1) means "x", "110" (3) means "x and y", "111" (7) means "x and y and z", etc.
	stressMask = 7,
	internalCompaction=False, # If true the confining pressure is generated by growing particles
)

contact = Law2_ScGeom_FrictViscoPhys_CundallStrackVisco(traceEnergy=True)
newton=NewtonIntegrator(damping=damp)

O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
		[Ip2_FrictViscoMat_FrictViscoMat_FrictViscoPhys()],
		[contact]
	),
	## We will use the global stiffness of each body to determine an optimal timestep (see https://yade-dem.org/w/images/1/1b/Chareyre&Villard2005_licensed.pdf)
	GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
	triax,
	newton,
	# save data for Paraview
	VTKRecorder(fileName='post/3d-vtk-',recorders=['all'],iterPeriod = 500000)
]

O.dt=.5*utils.PWaveTimeStep()

# for the samples of this size (diametre is around 1.8cm), the O.dt is 8.8e-5 s, which is about a hundred times of 
# the small sample, and this leads to much faster calculation speed.

#######################################
###   APPLYING CONFINING PRESSURE   ###
#######################################

## we define the lateral stresses during the test, here the same 100kPa as for the initial confinement.
p_const = 100000

#the value of (isotropic) confining stress defines the target stress to be applied in all three directions
triax.goal1=triax.goal2=triax.goal3=-p_const

while 1:
  O.run(1000, True)
  ##the global unbalanced force on dynamic bodies, thus excluding boundaries, which are not at equilibrium
  unb=unbalancedForce()
  print 'unbalanced force:',unb,' mean stress: ',triax.meanStress
  if unb<stabilityThreshold and abs(-p_const-triax.meanStress)/p_const<0.001:
    break

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


##############################
###   DEVIATORIC LOADING   ###
##############################

## Change contact friction (remember that decreasing it would generate instantaneous instabilities)
setContactFriction(radians(finalFricDegree))

##set stress control on x and y, we will impose strain rate on z
triax.stressMask = 3
##now goal2 is the target strain rate
triax.goal3=rate


O.run(4000000,True)

O.save('final.yade.gz')

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.