yade-users team mailing list archive

[xujin <question701059@xxxxxxxxxxxxxxxxxxxxx>]

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

Description changed to:
Hi, I am studying ‘soil structure interaction’ and meeting some
questions again.

I have two scripts and they should run in sequence.

The first script is used to generate particles. The particles expand
until the wallstress reach 90kPa.

The second script will load the result of the first script. I delete all
the walls of the model and four walls were recreated. The front and back
walls aren’t recreated and I don’t know if it is correct. All walls are
fixed except the top one. I want to apply ‘-100000 Pa’ on the top wall
and I don’t know how to apply stress directly so I adopt this method
‘O.forces.setPermF(top wall id, (0,-100000*area(top wall),0))’. Is there
any method to apply stress directly on the wall? The top wall may slope
slightly(such as right higher than left) when add the stress. How do I
control it to stay level the whole time when the program is running (My
most pressing problem)?


The version of ubuntu is 21.10

The first script:

#########################################
### Defining parameters and variables ###
#########################################

#Material constants 
Density = 3000
FrictionAngle = 1.5
PoissonRatio = 0.5
Young = 300e6
Damp = 0.5
AvgRadius = 0.1
N_particles = 4000

#Wall constants
WDensity = 0
WFrictionAngle = 0.0
WPoissonRatio = 0.5
WYoung = 50e9

#Packing variables
mn = Vector3(1,0,0)
mx = Vector3(101,20,1)

#Confining variables
ConfPress1 = -90000  #pre-compression
ConfPress = -1.0e5

#time calculation
startT = O.time
endT = O.time
timeSpent = endT - startT

########################################
#import necessary packages
from yade import pack,plot,os,timing
import matplotlib; matplotlib.rc('axes',grid=True)
import pylab

########################################
### Sample Preparing ###################
########################################

#Create materials for spheres and plates
SphereMat = O.materials.append(FrictMat(young = Young, poisson = PoissonRatio, frictionAngle = radians(FrictionAngle), density = Density))
WallMat = O.materials.append(FrictMat(young = WYoung, poisson = WPoissonRatio, frictionAngle = radians(WFrictionAngle), density = WDensity))

#Create walls for packing
wallIds = O.bodies.append(aabbWalls([mn,mx],thickness=0.001,material=WallMat))

#Use SpherePack object to generate a random loose particles packing
#O.periodic = True
#O.cell.setBox(8,3,8)

sp = pack.SpherePack()

#psdSizes,psdCumm=[0.00142,0.00151,0.00160,0.00169,0.00179,0.00188,0.00197,0.00207,0.00216,0.00225,0.00235,0.00244,0.00253,0.00263,0.00272,0.00281],[0,0.067,0.133,0.2,0.267,0.333,0.4,0.467,0.533,0.6,0.667,0.733,0.8,0.867,0.933,1]

sp.makeCloud(Vector3(1,0,0.5),Vector3(101,20,0.5),num=N_particles,rMean=AvgRadius,rRelFuzz=.33,distributeMass=False)

sp.toSimulation(material = SphereMat)

#for b in O.bodies:
#   if isinstance(b.shape,Sphere):
#       b.state.blockedDOFs='zXY'
#       b.shape.color=(0.,0.,1.)


O.usesTimeStepper=True
O.trackEnerty=True
#################################
#####Defining triaxil engines####
#################################

###first step: compression#######
triax1=TriaxialStressController(
    wall_back_activated = False, #for 2d simulation
    wall_front_activated = False,
    thickness = 0.001,
    maxMultiplier=1.+1.5e5/Young, # spheres growing factor (fast growth)
    finalMaxMultiplier=1.+4e3/Young,
    #maxMultiplier = 1.002,
    internalCompaction = True, # If true the confining pressure is generated by growing particles
    #max_vel = 1.5,
    stressMask = 7,
    computeStressStrainInterval = 10,
    
    goal1 = ConfPress1,
    goal2 = ConfPress1,
    #goal3 = ConfPress,  
)


#O.dt=0.5*PWaveTimeStep()

newton=NewtonIntegrator(damping=Damp)


###engine
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, defaultDt=4*utils.PWaveTimeStep()),
    triax1,
    newton,
    PyRunner(realPeriod=10,command='checkUnbalanced()',label='check'),
]
# Simulation stop conditions defination
def checkUnbalanced():
    unb=unbalancedForce()
    if unb<0.001:
       O.pause()
       O.save('initial.yade.bz2')





the second script:

#Wall constants
WDensity = 0
WFrictionAngle = 0.0
WPoissonRatio = 0.5
WYoung = 50e9
Mat = O.materials.append(FrictMat(young = WYoung, poisson = WPoissonRatio, frictionAngle = radians(WFrictionAngle), density = WDensity))
FrictionAngle = 35

Damp = 0.25

#Confining variables
force = -10.0e7 #force=-100000*area(top wall)

#time calculation
startT = O.time
endT = O.time
timeSpent = endT - startT

#restart
O.load('initial.yade.bz2')

########################################
#import necessary packages
from yade import pack,plot,os,timing


for b in range(6):
    O.bodies.erase(b)

id_box1 = O.bodies.append(box((0.5,12.5,0.5),(0.5,12.5,0.5),fixed=True,material=Mat))#    left wall
id_box2 = O.bodies.append(box((101.5,12.5,0.5),(0.5,12.5,0.5),fixed=True,material=Mat))#  right wall
id_box3 = O.bodies.append(box((51,20.5,0.5),(50,0.5,0.5),fixed=False,material=Mat))#      top wall
id_box4 = O.bodies.append(box((51,-0.5,0.5),(70,0.5,0.5),fixed=True,material=Mat))#       bottom wall

print('bottom',O.bodies[-1].state.pos)
print('top',O.bodies[-2].state.pos)
print('right',O.bodies[-3].state.pos)
print('left',O.bodies[-4].state.pos)

#O.bodies[-2].state.blockedDOFs='yXZ'
#O.bodies[-1].state.blockedDOFs='xYZ'

for b in O.bodies:
   if isinstance(b.shape,Sphere):
#       b.state.blockedDOFs='zXY'
       b.shape.color=(0.,0.,1.)

O.forces.setPermF(O.bodies[-2].id,(0,force,0))

#O.bodies[-1].state.vel = (1,0,0)

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, defaultDt=4*utils.PWaveTimeStep()),
    NewtonIntegrator(damping=Damp),
    PyRunner(realPeriod=2,command='endCheck()',label='check')
#    PyRunner(realPeriod=10,command='checkUnbalanced()',label='check'),
]

# Simulation stop conditions defination
def endCheck():
    unb=unbalancedForce()
    if unb<0.0001:
        if O.iter<60000:
            return
        else:
            O.pause()
            print('bottom force:',O.forces.f(O.bodies[-1].id))
            print('top force:',O.forces.f(O.bodies[-2].id))
            print('right force',O.forces.f(O.bodies[-3].id))
            print('left force',O.forces.f(O.bodies[-4].id))
            O.save('Compression.yade.bz2')

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