yade-users team mailing list archive

[Ruidong LI <question707148@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Needs information => Open

Ruidong LI gave more information on the question:
Hi! Jan
Thanks for your reply.

> The loading is simply such that a constant velocity is prescribed (s.state.vel = ...)
But which velocity value should I define to achieve quasi-equilibrium compression? Can you show me the code? I am sorry for my stupidity.

>Therefore I ask for a MWE and better description of "are destroyed" (anyway useful even with the code).
Actually, what I mean is that when you run the 'triax' py without any change, the triangle facets will break. And my current problem is not with the 'triax' py.
My problem is that I just finished the part about generating flexible membranes using Pfacet and don't know what to do next. What I want to achieve is isotropic compression with a confining pressure equal to 50 kPa. And then, maintain the confining pressure in the x and y directions and apply a displacement-controlled movement of the top wall. 

my MWE is presented as below



################################
###            1             ###
###   LOAD INITIAL PACKING   ###
################################
dir_inipacking = '/home/kyle/Yade/install/bin/lentils/xml/iniPacking.xml'
O.load(dir_inipacking)

# Define engine
O.engines = [
 	ForceResetter(),
 	InsertionSortCollider([Bo1_Sphere_Aabb(),
 				Bo1_PFacet_Aabb(),], 
 				sortThenCollide=True),
 	InteractionLoop(
 		[
 			Ig2_GridNode_GridNode_GridNodeGeom6D(),
 			Ig2_GridConnection_GridConnection_GridCoGridCoGeom(),
 			Ig2_Sphere_Sphere_ScGeom6D(),
 			Ig2_Sphere_PFacet_ScGridCoGeom(),
 		],
 		[
        		Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True,setCohesionOnNewContacts=True),
  			Ip2_FrictMat_FrictMat_FrictPhys()
  		],
 		[
        		Law2_ScGeom6D_CohFrictPhys_CohesionMoment(),
  			Law2_ScGeom_FrictPhys_CundallStrack(),
  			Law2_ScGridCoGeom_FrictPhys_CundallStrack(),
  			Law2_GridCoGridCoGeom_FrictPhys_CundallStrack(),
 		]
 	),
    	GlobalStiffnessTimeStepper(timestepSafetyCoefficient=0.8,label='ts'),
 	NewtonIntegrator(gravity=(0,-9.81,0),damping=0.1,label='newton')
]

# encoding: utf-8
from yade import qt
from yade.gridpfacet import *
import math

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

# The following lines are used parameter definitions
readParamsFromTable(
	# directory	
	dir_vtk = '/home/kyle/Yade/install/bin/lentils/vtk/20', # directory of storing vtk files
	dir_txt = '/home/kyle/Yade/install/bin/lentils/txt/20/Num20_servo.txt', # directory of storing txt files
	dir_xml = '/home/kyle/Yade/install/bin/lentils/xml/20', # directory of storing xml files
        
        # material parameters
        young_w = 5e9, # Young's modulus of plates and walls
        young_g = 1.8e6, # Young's modulus of grains [1]
        den_w = 7874, # density of plates and walls
        den_g = 980, # density of grains 
        poi_w = 0.3, # possion's ratio of plates and walls
        poi_g = 0.25, # possion's ratio of grains
        friAngle_w = 0.5,#radians(15), # friction angle of plates and walls 
        friAngle_g = 0.5,#radians(29), # friction angle of grains

        # Parameters of cylindrical walls
        x_cyl = 0.0547, # x-coordinate of center of cylindrical walls
        y_cyl = 0.0535, # y-coordinate of center of cylindrical walls
        z_cyl = 0, # z-coordinate of center of cylindrical walls
        r_cyl = 0.0358, # radius of the cylinder
        h_cyl = 0.14, # height of the cylinder
        
        # confining pressure
        preStress = -50000,
        
        # axial strain rate
        strainRate=-100
              
)
from yade.params.table import *

########################################
###                 3                ###
###   CREATE TOP AND BOTTOM PLATES   ###
########################################

# erase rigid cylindrical facets
for b in O.bodies:
   if isinstance(b.shape,Facet): O.bodies.erase(b.id) 

# create top and bottom plates
h_cyl = max([b.state.pos[2] + b.shape.radius for b in O.bodies if isinstance(b.shape, Sphere)])
top_plate = O.bodies.append(yade.geom.facetCylinder((x_cyl,y_cyl,h_cyl),radius=r_cyl*1.2,height=0,segmentsNumber=40,color=(1,1,1)))
bottom_plate = O.bodies.append(yade.geom.facetCylinder((x_cyl,y_cyl,0),radius=r_cyl*1.2,height=0,segmentsNumber=40,color=(1,1,1)))

'''
vel = strainRate * (h - rParticle * 2 * bcCoeff)
for s in bot:
	s.shape.color = (1, 0, 0)
	s.state.blockedDOFs = 'xyzXYZ'
	s.state.vel = (0, 0, -vel)
for s in top:
	s.shape.color = Vector3(0, 1, 0)
	s.state.blockedDOFs = 'xyzXYZ'
	s.state.vel = (0, 0, vel)
'''

######################################
###                4               ###
###   GENERATE FLEXIBLE MEMBRANE   ###
######################################

# Materials
O.materials.append(CohFrictMat(young=1e7,poisson=1,density=2650,frictionAngle=radians(30),normalCohesion=3e7,shearCohesion=3e7,momentRotationLaw=True,label='gridNodeMat'))
O.materials.append( FrictMat( young=1e7,poisson=0.1,density=2650,frictionAngle=radians(30),label='pFacetMat' ))

# Generate flexible membrane
nodesIds=[]
cylIds=[]
pfacets=[]
width=2*r_cyl #diameter of cylinder
height=h_cyl #height of cylinder
nw=40 # no of nodes along perimeter
nh=25 # no of nodes along height
rNode=width/100 #radius of grid node
color1=[255./255.,102./255.,0./255.]
color2=[0,0,0]
color3=[248/255,248/255,168/255]
color4=[156/255,160/255,98/255]
#color3=[165/255,245/255,251/255]
#color4=[101/255,153/255,157/255]
rCyl2 = 0.5*width / cos(pi/float(nw))
vCenter = Vector3(x_cyl, y_cyl, 0)
for r in range(nw):
	for h in range(nh):
  		v1 = Vector3(rCyl2*cos(2*pi*(r+0)/float(nw)),rCyl2*sin(2*pi*(r+0)/float(nw)), height*(h+0)/float(nh))+vCenter
  		v2 = Vector3(rCyl2*cos(2*pi*(r+1)/float(nw)),rCyl2*sin(2*pi*(r+1)/float(nw)), height*(h+0)/float(nh))+vCenter
  		v3 = Vector3(rCyl2*cos(2*pi*(r+1)/float(nw)),rCyl2*sin(2*pi*(r+1)/float(nw)), height*(h+1)/float(nh))+vCenter
  		v4 = Vector3(rCyl2*cos(2*pi*(r+0)/float(nw)),rCyl2*sin(2*pi*(r+0)/float(nw)), height*(h+1)/float(nh))+vCenter
  		V1=(O.bodies.append(gridNode(v1,rNode,wire=False,fixed=False, material='gridNodeMat',color=color2)) )
  		V2=(O.bodies.append(gridNode(v2,rNode,wire=False,fixed=False, material='gridNodeMat',color=color2)) )
  		V3=(O.bodies.append(gridNode(v3,rNode,wire=False,fixed=False, material='gridNodeMat',color=color2)) )
  		V4=(O.bodies.append(gridNode(v4,rNode,wire=False,fixed=False, material='gridNodeMat',color=color2)) )
  		# append node ids to seperate matrix for later use
  		nodesIds.append(V1)
  		nodesIds.append(V2)
  		nodesIds.append(V3)
  		nodesIds.append(V4)
  		#create grid connection
  		cylIds.append(O.bodies.append(gridConnection(V1,V2,rNode,material='gridNodeMat',color=color3)))
  		cylIds.append(O.bodies.append(gridConnection(V2,V3,rNode,material='gridNodeMat',color=color3)))
  		cylIds.append(O.bodies.append(gridConnection(V1,V3,rNode,material='gridNodeMat',color=color3)))
  		cylIds.append(O.bodies.append(gridConnection(V3,V4,rNode,material='gridNodeMat',color=color3)))
  		cylIds.append(O.bodies.append(gridConnection(V4,V1,rNode,material='gridNodeMat',color=color3)))
  		#create Pfacets
  		pfacets.append(O.bodies.append(pfacet(V1,V2,V3,wire=True,material='pFacetMat',color=color3)))
  		pfacets.append(O.bodies.append(pfacet(V1,V3,V4,wire=True,material='pFacetMat',color=color4)))
  	
# calculate void ratio
def myClumpVoidRatio(VolTot, density):
	mass=sum([ b.state.mass for b in O.bodies if b.isClump])
	VolSolid = mass/density
	poro = (VolTot-VolSolid)/VolSolid
	return poro	
h_cyl = max([b.state.pos[2] + b.shape.radius for b in O.bodies if isinstance(b.shape, Sphere)])
Vol = math.pi * h_cyl * r_cyl * r_cyl
poro_ini = myClumpVoidRatio(Vol, den_g)
print('Initial void ratio:', str(poro_ini))

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