yade-users team mailing list archive

Thread
Date
Re: [Question #686318]: random dense pack with clumps

To: yade-users@xxxxxxxxxxxxxxxxxxx
From: Grace Mejico <question686318@xxxxxxxxxxxxxxxxxxxxx>
Date: Fri, 20 Dec 2019 17:57:28 -0000
Reply-to: question686318@xxxxxxxxxxxxxxxxxxxxx
Sender: bounces@xxxxxxxxxxxxx
Question #686318 on Yade changed:
https://answers.launchpad.net/yade/+question/686318

    Status: Answered => Open

Grace Mejico is still having a problem:
thank you so much for your answer, i try to incorporate what you suggest . Even though, it continues appearing " In [1]: FATAL /build/yade-fDuCoe/yade-2018.02b/core/ThreadRunner.cpp:30 run: Exception occured:
Body #573: Body::material type JCFpmMat doesn't correspond to Body::state type State (should be JCFpmState instead)." 

I would like to know if i have to change also the engine? or the problem
its related to the engine with JCFPm type? or is there any other thing
im doing wrong in the script? please, hope i can solve this error. Thank
you in advance, im very thankful for the responses.

im copying the script with the modifications 
from yade import pack, qt, plot
from numpy import arange
import numpy as np
import random

#DATOS DE ENTRADA
r=0.07
H=1 #altura del muro
B=.45 #base mayor
b=.25 #base menor
L=4 #longitud del muro
#EN SERVICIO
EmpujeTop=2400#Pa2400 
EmpujeBottom=9100#Pa9100
damp=0.7
aceloutofplane=0
ii=0

# GEOMETRIA ********************************************************************
a=(B-b)*.5
Vol1=pack.inParallelepiped(o=(0,0,0), a=(L,0,0), b=(0,B,0), c=(0,0,H))
Vol2=pack.inParallelepiped(o=(L,B+a,0), a=(0,B+a,0), b=(0,a+b,H), c=(L,a+b,H))
Vol3=pack.inParallelepiped(o=(L,-a,0), a=(0,-a,0), b=(0,a,H), c=(L,a,H))
Vol4=Vol1-Vol2-Vol3-Vol3
Esferas=pack.randomDensePack(Vol4,radius=r,rRelFuzz=0)
O.bodies.append(Esferas)

relRadList1 = [.04,.03]
relPosList1 = [[0.035,0,0],[.083,0,0]]
#peanut:
relRadList2 = [.05,.05]
relPosList2 = [[0,0,0],[.073,0,0]]

relRadList4 = [0.025,0.025,0.025,0.025]
relPosList4 = [[0,0,0],[.043,0,0],[0.02,0.045,0],[0.02,0,045]]
#stick:
relRadList3 = [0.05,0.025,0.05,0.025]
relPosList3 = [[.04,0,0],[.123,0,0.022],[0.14,0,0],[0.123,0,-0.022]]
templates= []
templates.append(clumpTemplate(relRadii=relRadList1,relPositions=relPosList1))
templates.append(clumpTemplate(relRadii=relRadList2,relPositions=relPosList2))
templates.append(clumpTemplate(relRadii=relRadList3,relPositions=relPosList3))
templates.append(clumpTemplate(relRadii=relRadList4,relPositions=relPosList4))
O.bodies.replaceByClumps(templates,[.5,.3,.1,.1])

# DEFINICION DEL MATERIAL ******************************************************
mat1 = JCFpmMat()
mat1.cohesion =0 #Pa32.37e3
mat1.density = 2821 #kg/m3
mat1.frictionAngle = radians(58.74) #rad
mat1.poisson = 0.18
mat1.young = 37.04e9 #Pa

# DEFINICION DE LA JUNTA  ******************************************************
mat1.jointCohesion = 262000 #Pa
mat1.jointFrictionAngle = radians(40)
mat1.jointNormalStiffness= 1.75e8 #Pa/m
mat1.jointShearStiffness=0.7e8 #Pa/m
#mat1.jointDilationAngle = 0
#mat1.jointTensileStrength = 0
O.materials.append(mat1)


# APLICACION DE FUERZAS ********************************************************
listSphereA=[]
listSphereB=[]
listSphereC=[]
listSphereD=[]
for x1 in O.bodies:
	if isinstance(x1.shape,Sphere) and x1.state.pos[1]<3*r and x1.state.pos[2]<=H and x1.state.pos[2]>3*H/4.000:
		listSphereA.append(x1.id)
for x2 in O.bodies:
	if isinstance(x2.shape,Sphere) and x2.state.pos[1]<3*r and x2.state.pos[2]<=3*H/4.000 and x2.state.pos[2]>H/2.000:
		listSphereB.append(x2.id)
for x3 in O.bodies:
	if isinstance(x3.shape,Sphere) and x3.state.pos[1]<3.5*r and x3.state.pos[2]<=H/2.000 and x3.state.pos[2]>H/4.000:
		listSphereC.append(x3.id)
for x4 in O.bodies:
	if isinstance(x4.shape,Sphere) and x4.state.pos[1]<3*r and x4.state.pos[2]<=H/4.000 and x4.state.pos[2]>0.000:
		listSphereD.append(x4.id)
ha=1*H/8.000
hb=3*H/8.000
hc=5*H/8.000
hd=7*H/8.000
m=EmpujeBottom-EmpujeTop
Ea=m*ha+EmpujeTop
Eb=m*hb+EmpujeTop
Ec=m*hc+EmpujeTop
Ed=m*hd+EmpujeTop

Fa=Ea*(H/4.000)*L/len(listSphereA)
Fb=Eb*(H/4.000)*L/len(listSphereB)
Fc=Ec*(H/4.000)*L/len(listSphereC)
Fd=Ed*(H/4.000)*L/len(listSphereD)

def aplicarFuerzaA():
	for i in listSphereA:
		O.forces.setPermF(i,(0,Fa,0))
def aplicarFuerzaB():
	for i in listSphereB:
		O.forces.setPermF(i,(0,Fb,0))
def aplicarFuerzaC():
	for i in listSphereC:
		O.forces.setPermF(i,(0,Fc,0))
def aplicarFuerzaD():
	for i in listSphereD:
		O.forces.setPermF(i,(0,Fd,0.98))


# CONDICIONES DE BORDE *********************************************************

Borde1=utils.wall(0.43,axis=0, sense=0)
Borde2=utils.wall(L-0.43,axis=0, sense=0)
Borde3=utils.wall(0, axis=1, sense=0)
Borde4=utils.wall(0.05, axis=2, sense=0)

O.bodies.append(Borde1)
O.bodies.append(Borde2)
O.bodies.append(Borde3)
O.bodies.append(Borde4)

listSphereBase=[]
#listSphereBorde1=[]
#listSphereBorde2=[]
#restringir movimiento en la base
#for y1 in O.bodies:
#	if isinstance(y1.shape,Sphere) and y1.state.pos[2]<=2.5*r:
#		listSphereBase.append(y1.id)
#for y2 in O.bodies:
#	if isinstance(y2.shape,Sphere) and y2.state.pos[0]<=2.5*r:
#		listSphereBorde1.append(y2.id)
#for y3 in O.bodies:
#	if isinstance(y3.shape,Sphere) and y3.state.pos[0]>L-2.5*r:
#		listSphereBorde2.append(y3.id)
#for i in listSphereBase:
#		O.bodies[i].state.blockedDOFs='xyz'
#for i in listSphereBorde1:
#		O.bodies[i].state.blockedDOFs='xyzXYZ'
#for i in listSphereBorde2:
#		O.bodies[i].state.blockedDOFs='xyzXYZ'


# DEFINICION DE FUNCIONES ******************************************************
listSphereTodo=[]
for m in O.bodies:
	if isinstance(m.shape,Sphere) and m.state.pos[2]<2*H:
		listSphereBase.append(m.id)

for i in listSphereTodo:
	O.bodies[i].state.blockedDOFs='XYZ'

def checkUnbalanced():
	print ('iter %d, unbalanced forces = %f'%(O.iter, utils.unbalancedForce()))
	print 'Fzas en c/ esfera (N)', Fa, Fb, Fc, Fd
#	if unbalancedForce()<.00001:
	print ( 'Desplazamiento superior',O.bodies[w].state.pos[1],O.bodies[w].state.pos[2])

listTopBodies=[]
for z1 in O.bodies:
	if isinstance(z1.shape,Sphere) and z1.state.pos[2]<=H and z1.state.pos[2]>H-3.000*r and z1.state.pos[0]>(L/2.000)-1.1*r and z1.state.pos[0]<(L/2.000)+1.1*r and z1.state.pos[1]>B/2.000:
		listTopBodies.append(z1.id)
w=max(listTopBodies)

def addPlotData():
	plot.addData(Ypos=O.bodies[w].state.pos[1],Zpos=O.bodies[w].state.pos[2],iteracion=O.iter,unbal=unbalancedForce())
plot.plots={'Ypos':('Zpos')}


# ENGINE ***********************************************************************
O.engines=[
		ForceResetter(),
		InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Wall_Aabb()]),
		InteractionLoop(
			[Ig2_Sphere_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
			[Ip2_JCFpmMat_JCFpmMat_JCFpmPhys()],
			[Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM()]
		),
		NewtonIntegrator(gravity=(0,aceloutofplane,-9.81),damping=damp),
		PyRunner(firstIterRun=10000,command='aplicarFuerzaA()'),
		PyRunner(firstIterRun=10000,command='aplicarFuerzaB()'),
		PyRunner(firstIterRun=10000,command='aplicarFuerzaC()'),
		PyRunner(firstIterRun=10000,command='aplicarFuerzaD()'),
		PyRunner(iterPeriod=10000,command='checkUnbalanced()'),
		PyRunner(iterPeriod=1,command='addPlotData()'),
		VTKRecorder(iterPeriod=100,recorders=['spheres','jcfpm','facets','colors'],fileName='/tmp/p1-')
]

for b in O.bodies:
	b.mat=mat1
# DETALLES FINALES *************************************************************

O.dt=0.5*PWaveTimeStep() # establece el delta de tiempo como una fraccion del tiempo critico
plot.plot(subPlots=True) # llama a los sub plots
plot.live=True # ploteo en tiempo real
qt.Controller() # abre la ventana del controlador
V=qt.View() # abre la ventada de la vista
R=yade.qt.Renderer() # llama a la renderizacion
R.bgColor=(1.,1.,1.) # definel color blan
V.screenSize = (900,900) # tamano de pantalla

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