Request for suggestion regarding particle packing for CD triaxial test

Hi,

Please avoid in the future to post several questions in one Launchpad question. This being said,

1. See the notion of "Representative Elementary Volume (REV)" for DEM triaxial simulations in the scientific literature. There will be answer(s) there. ~ 8000 - 10000 in my opinion

3. This should be inconsequential. Note the triaxial test examples at https://gitlab.com/yade-dev/trunk/-/tree/master/examples/triax-tutorial (or somewhere on your computer)

4. YADE can do both.

Hi,
Thanks for the reply.

> Please avoid in the future to post several questions in one Launchpad question. - Sure

I am a new user to YADE. Based on your suggestions, I tried the following code for simple CD triaxial test stage but I think I have missed something in my code as the engine is not working after particle generation. Can you please let me know if I missed some thing.

from yade import pack
from numpy import arange
import itertools
import random
import yade.plot
## Extra things are also imported for future use
## corners of the initial packing
mn,mx=Vector3(0,0,0),Vector3(0.04,0.04,0.04)
finalFricDegree = 30
rate=0.02
stabilityThreshold=0.01
damp = 0.2
## create material #0, which will be used as default

id_Mat=O.materials.append(FrictMat(young=7e7,poisson=.4,frictionAngle=radians(5.8),density=2650,label='spheres'))
O.materials.append(FrictMat(young=7e6,poisson=.4,frictionAngle=0,density=3000,label='walls'))
Mat=O.materials[id_Mat]

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

#create assembly of spheres:
sp=pack.SpherePack()
sp.makeCloud(mn,mx,num=15000,rMean=.0005,rRelFuzz=0.0,periodic=False),

O.bodies.append([sphere(center,rad,material='spheres') for center,rad in sp])
print(len(sp),' particles generated.')
#sp.toSimulation()

from yade import qt
qt.Controller()
qt.View()

## hope that we got the ids right?!

triax = TriaxialStressController(
#maxMultiplier=1.+2e4/young, # spheres growing factor (fast growth)
#finalMaxMultiplier=1.+2e3/young, # spheres growing factor (slow growth)
 wall_bottom_id=wallIds[2],
 wall_top_id=wallIds[3],
 wall_left_id=wallIds[0],
 wall_right_id=wallIds[1],
 wall_back_id=wallIds[4],
 wall_front_id=wallIds[5],
 internalCompaction=False,

 thickness = 0,
 stressMask =7,
 max_vel=1,
)

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(timestepSafetyCoefficient=0.2),
 triax,
 #TriaxialStateRecorder(iterPeriod=100,file='WallStresses'+table.key),
 # you can add TriaxialStateRecorder and such here…
 newton,

]
#Gl1_Sphere.stripes=0
#if nRead==0: yade.qt.Controller(), yade.qt.View()

##############################
### Confining Stage ###
##############################

triax.goal1=triax.goal2=triax.goal3=-100000
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(-100000-triax.meanStress)/100000<0.001:
     break

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

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

##We move to deviatoric loading, let us turn internal compaction off to keep particles sizes constant
triax.internalCompaction=False

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

##set stress control on x and z, we will impose strain rate on y
triax.stressMask = 5
##now goal2 is the target strain rate
triax.goal2=-rate
## we define the lateral stresses during the test, here the same 10kPa as for the initial confinement.
triax.goal1=-100000
triax.goal3=-100000

##we can change damping here. What is the effect in your opinion?
newton.damping=0.1

##Save temporary state in live memory. This state will be reloaded from the interface with the "reload" button.
O.saveTmp()
O.run(100000,True)

Thanking you in advance,
Sarath

This question was expired because it remained in the 'Open' state without activity for the last 15 days.