← Back to team overview

yade-users team mailing list archive

Re: [Question #194918]: change blockedDOFs of utils.box

 

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

onyourself posted a new comment:
hi

I encounter this problem when running my script about "2D triax"…i get
that "There should be a warning if mass = 0.0 and fixed = False ..."
just dont know how to correct it, below is my code.

#############################################################################

from yade.gridpfacet import *
from yade import pack, plot
from random import random
import numpy as np
from numpy import *
import math

#parameters
rParticle = 0.14e-3#diameter0.28
rRelFuzz=.0005
mi,ma = (0.,0.,0.),(1000e-3,0,250e-3)
#nCyls,nSphs = 8000,28216
frictionAngleSph=30
frictionAngleCyl=30

#meterials########
O.materials.append(FrictMat(young=4.0e6,poisson=.3,frictionAngle=20,density=910e+6,label='sphcylMat'))
O.materials.append(FrictMat(young=4.0e6,poisson=.3,frictionAngle=26,density=910e+6,label='cylinderMat'))
O.materials.append(FrictMat(young=4.0e6,poisson=.15,frictionAngle=30,density=2630e+6,label='sphereMat'))
O.materials.append(FrictMat(young=4.0e6,poisson=.15,frictionAngle=0,density=0,label='walls'))

## create walls around the packing
walls=aabbWalls([(0,-0.5,0),(1000e-3,0.5,250e-3)],thickness=0,material='walls')
wallIds=O.bodies.append(walls)

##sphere
sp = yade.pack.SpherePack()
sp.makeCloud(mi,ma,porosity=0.7/1.7,psdSizes=[.00014,0.00016,0.00022,.0003,.00035],psdCumm=[0.,0.1,0.3,0.6,1.],num=35000)
#psdSizes=[0.01,0.0125,0.015,0.018,0.02,0.023,0.025,0.027,0.029,0.03],psdCumm=[0,0.05,0.15,0.38,0.54,0.77,0.9,0.97,0.99,1]) # makeCloud for spheres
spheres=sp.toSimulation(color=(0,0.5,0.7),material='sphereMat')
print "cheers"

for i in O.bodies:
	i.state.blockedDOFs="yXZ"
	i.shape.color=(0.,0.5,0.7)

num=0
vol=0
mass=0
area=0
for b in O.bodies:
	if isinstance(b.shape,Sphere):
  	  r=b.shape.radius
  	  x=b.state.pos[0]
  	  y=b.state.pos[1]
  	  z=b.state.pos[2]
  	  ID=b.id
  	  num+=1
  	  #vol+=4./3.*np.pi*r**3
  	  mass+=2.63*4./3.*np.pi*r**3
	  area+=np.pi*r**2
print 'num,mass,voidRatio =',num,mass,((1000e-3*250e-3)-area)/area

length=6e-3
rfiber=0.023e-3/2
fibre_mass=np.pi*(rfiber)**2*length*0.91
Xw=0.5*1e-2 # fiber content
numFiber=Xw*mass/fibre_mass
print "the numer of fibres =",int(numFiber)

Ne=int(length/(0.28e-3/2))
print 'Ne=',Ne
nodesIds=[]
cylIds=[]
numFibre=0
#target_num=int(numFiber)
target_num=1200
deleted_sphere=0

####fiberfilewrite
fiberwr=open('fibers.txt','w')
cylNodes=[]

for n in range(2000):
	if numFibre<target_num:
		
		random.seed()
		l_x=random.uniform(0,1000e-3)
		l_z=random.uniform(0,250e-3)
		x0=l_x
		z0=l_z
		
		u1 = random.random()
		u2 = random.random()*pi
		u3 = random.random()*pi
		a = sqrt(1-u1)
		b = sqrt(u1)
		t=((a*cos(u2))**2+(b*sin(u3))**2)**0.5
		cosphi=a*cos(u2)/t
		sinphi=b*sin(u3)/t

		hx=length*cosphi
		hz=length*sinphi
		
		x1=x0+hx
		y1=0
		z1=z0+hz

		if 0<x1<1000e-3 and 0<z1<250e-3 and 0<x0<1000e-3 and 0<z0<250e-3:
			if len(cylNodes)<=1:
				fiberwr.write(str(numFibre)+'\t'+str(x0)+'\t'+str(0)+'\t'+str(z0)+'\t'+str(x1)+'\t'+str(0)+'\t'+str(z1)+'\n')
				cylNodes.append([x0,z0,x1,z1])
				numFibre+=1
				vertices=[]
				for i in range(0, Ne+1):
					px=float(i)*hx/float(Ne)+x0; py=0; pz=float(i)*hz/float(Ne)+z0;
					vertices.append([px,py,pz])        
				cylinderConnection(vertices,0.023e-3/2,nodesIds,cylIds,color=[0.5,0.5,0],highlight=False,intMaterial='cylinderMat',extMaterial='fMat')
				print "+1"
			else:
				valid=True
				#print len(cylNodes)
				for c in cylNodes:
					#if(max(c[0],c[2])>=min(x0,x1) 
					#and max(x0,x1)>=min(c[0],c[2]) 
					#and max(c[1],c[3])>=min(z0,z1)
					#and max(z0,z1)>=min(c[1],c[3])):
					x12=(x1-x0)
					z12=(z1-z0)
					x13=(c[0]-x0)
					z13=(c[1]-z0)
					x14=(c[2]-x0)
					z14=(c[3]-z0)
					cross1=x12*z13-x13*z12
					cross2=x12*z14-x14*z12
					x34=c[2]-c[0]
					z34=c[3]-c[1]
					x23=c[0]-x1
					z23=c[1]-z1
					x24=c[2]-x1
					z24=c[3]-z1
					cross3=x34*z23-x23*z34
					cross4=x34*z24-x24*z34
					if(cross1*cross2<=0
					and cross3*cross4<=0):
						valid=False
						print "next"
						break
				if valid:
					fiberwr.write(str(numFibre)+'\t'+str(x0)+'\t'+str(0)+'\t'+str(z0)+'\t'+str(x1)+'\t'+str(0)+'\t'+str(z1)+'\n')
					cylNodes.append([x0,z0,x1,z1])
					numFibre+=1
					#print "++1"
					vertices=[]
					for i in range(0, Ne+1):
						px=float(i)*hx/float(Ne)+x0; py=0; pz=float(i)*hz/float(Ne)+z0;
						vertices.append([px,py,pz])        
					cylinderConnection(vertices,0.023e-3/2,nodesIds,cylIds,color=[0.5,0.5,0],highlight=False,intMaterial='cylinderMat',extMaterial='fMat')

			###spherefileWrite
			##################################################################
			sphwr=open('spheres.txt','w')
			for b in O.bodies:
				if isinstance(b.shape,Sphere):
		    			r=b.shape.radius 
					if r>0.012e-3:                     
		    				x=b.state.pos[0]
			   			y=b.state.pos[1]
			    			z=b.state.pos[2]
						ID=b.id
						k=(z1-z0)/(x1-x0)
						#kx-z+z0-kx0=0
						A=k
						B=-1
						C=z0-k*x0
						#if (A*x+B*z+C)/(A^2+B^2)^0.5<=(0.14e-3+0.023e-3/2):
						if abs(k*x-1*z+C)/(k**2+1)**0.5<=(r+0.0115e-3):
							O.bodies.erase(b.id)
							deleted_sphere+=1
						#else:
						sphwr.write(str(ID)+'\t'+str(x)+'\t'+str(y)+'\t'+str(z)+'\t'+str(r)+'\n')
			sphwr.close()
			##################################################################
fiberwr.close()
print 'numFibre , deleted_sphere =',numFibre, deleted_sphere

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, mask = x*1 + y*2 + z*4
	## "100" (1) means "x", "110" (3) means "x and y", "111" (7) means "x and y and z", etc.
	stressMask = 2,
	goal1=-0.9,
	goal2=0,
	goal3=-0.9,
	internalCompaction=False # If true the confining pressure is generated by growing particles
)

O.engines=[
	ForceResetter(),
	InsertionSortCollider([
		Bo1_Box_Aabb(),
		Bo1_Sphere_Aabb(),
		Bo1_GridConnection_Aabb(),
	]),
	InteractionLoop([
		Ig2_Sphere_Sphere_ScGeom(),
		Ig2_Box_Sphere_ScGeom(),
		Ig2_GridNode_GridNode_GridNodeGeom6D(),
		Ig2_Sphere_GridConnection_ScGridCoGeom(),
		Ig2_GridConnection_GridConnection_GridCoGridCoGeom()
	],
	[
		Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True,setCohesionOnNewContacts=False),	# internal cylinder physics
		Ip2_FrictMat_FrictMat_FrictPhys()	# physics for external interactions, i.e., cylinder-cylinder, sphere-sphere, cylinder-sphere
	],
	[
		Law2_ScGeom_FrictPhys_CundallStrack(),	# contact law for sphere-sphere
		Law2_ScGridCoGeom_FrictPhys_CundallStrack(),	# contact law for cylinder-sphere
		Law2_ScGeom6D_CohFrictPhys_CohesionMoment(),	# contact law for "internal" cylinder forces
		Law2_GridCoGridCoGeom_FrictPhys_CundallStrack()	# contact law for cylinder-cylinder interaction
	]
	),
	GlobalStiffnessTimeStepper(timestepSafetyCoefficient=0.1,label='ts'), 
	#triax,
	#TriaxialStateRecorder(iterPeriod=100,file='WallStresses'+table.key),
	NewtonIntegrator(gravity=(0,0,10),damping=0.1,label='newton'),
]

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


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


Gl1_Sphere.stripes=1
yade.qt.Controller(), yade.qt.View()

print triax.porosity

###########################################################################

it is mainly about generating spheres and cylinders(fibers) then delete
objects if they interact each other."ValueError: cannot convert float
NaN to integer" occurs when i click the run button. thanks so much if
anyone can help!!!!!!!!

xxxe

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