← Back to team overview

yade-users team mailing list archive

Re: [Question #663930]: pfacetCreator rapidly expanding problem?

 

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

Description changed to:
I have a box filled with "grass" using pfacetCreator. The problem is,
the pfacetCreator rapidly expands(does crazy things) after the DOF's are
opened up for the nodes. Why is this happening? I would like the grass
to be able to move, but stay the same length. And constrained at the
bottom nodes only.


Code: (I have a box to keep the balls in, but commented it out, the problem stills occurs )

from yade.gridpfacet import *
import numpy as np
from yade import utils
from yade import ymport
import sys,time,math,random

# Spheres information
sphereRadius = .004 # [m]
nu = .48 
#G = 300000 # [Pa]
den_rub = 24299.1 # [kg/m^3]
yng_rub = 3000000 # [Pa] 
fric_rub = radians(38)  # [degrees] 
O.materials.append(FrictMat(frictionAngle=fric_rub,young=yng_rub,density=den_rub,poisson=nu,label='rubber'))

### Steps
spheres_fall = 200
Trans = 105280+spheres_fall
All_steps = 252510+spheres_fall


blen = 0.21230  ## Length of Ball pit
bhei = .025     ## Height of Ball pit
pred = pack.inAlignedBox((-blen/2,(0+bhei+.001),-blen/2),(blen/2,(bhei+bhei+.001),blen/2))
sp = pack.randomDensePack(pred, spheresInCell=2732, radius=sphereRadius, memoizeDb='/home/justin/Desktop/real/tmp/ball_packing65_r4_h25.sqlite', returnSpherePack=True)
sp.toSimulation(color=(255,0,0),wire=False,material='rubber')


## Steel box informaion
#steel = O.materials.append(FrictMat(young=200e9,density=8050,poisson=.3)) # Pa, kg/m^3
#metal_base45 = O.bodies.append(ymport.gmsh("geo/fliped/based_70.mesh",scale=.001,shift=Vector3(0,-.00635,0), material=steel,color=(0,0,1)))


## Grass Information
O.materials.append(CohFrictMat(young=3e9,density=1150,poisson=.3,frictionAngle=10,normalCohesion=1e7,shearCohesion=1e7,momentRotationLaw=True,label='grass')) 
rCyl = 0.0006  ## Grass was about 1.2 [mm] wide
nL = 2        ## No exact Number here, just trial and error
L = .022       ## Height of spheres after leveling

## Time step set to 20% of Rayleigh Wave
#O.dt=  0.000407535525016 #.2*utils.RayleighWaveTimeStep() # 0.000407535525016
O.dt= .2*utils.RayleighWaveTimeStep()
print O.dt


O.engines=[
	###Reset all forces stored in Scene::forces (O.forces in python). Typically, this is the first engine to be run at every step. 		In addition, reset those energies that should be reset, if energy tracing is enabled.
	## Resets forces and moments that act on bodies
	ForceResetter(),

	## Using bounding boxes find possible body collisions.
	InsertionSortCollider([
		Bo1_Sphere_Aabb(), 
		Bo1_Facet_Aabb(), 
		Bo1_PFacet_Aabb(),
	]),
	InteractionLoop([
		Ig2_Facet_Sphere_ScGeom(),
		Ig2_Sphere_Sphere_ScGeom(),
		Ig2_GridNode_GridNode_GridNodeGeom6D(),
		Ig2_Sphere_PFacet_ScGridCoGeom(),
		Ig2_PFacet_PFacet_ScGeom(),
		],
		[
		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
	]),
	NewtonIntegrator(damping=.2,gravity=[0,-9.81,0],label='newtonInt'),
	PyRunner(command='switchTranslationEngine()',firstIterRun=Trans, nDo = 6, label='switchEng', dead = False),
	PyRunner(command='turnongrass()',firstIterRun=spheres_fall, nDo = 3, dead = False),
]

color=[255./255.,102./255.,0./255.]

### Grass Creation
### Create all nodes first :
nodesIds=[]
nodesIds2=[]
idxc = -1
x_gap = 0.009 ## Between lumps is roughly 9 [mm]
z_gap = 0.01905 ## Between lines of backing is .75 inch apart
range_x = int(math.floor(blen/x_gap)) ## finding the range for x
range_z = int(math.floor(blen/z_gap)) ## finding the range for z
cen_z = -(range_z/2)*z_gap  ## Allows the "box" of grass to be center in Z
#sys.exit()

for ii in range(0,range_z):
  cen_x = -(range_x/2)*x_gap  # Allows the "box" of grass to be center in X
  for jj in range(0,range_x):
    for i in np.linspace(0,L,nL):
      nodesIds.append(O.bodies.append(gridNode([cen_x,i,cen_z],rCyl,wire=False,fixed=False,material='grass')))
    for i in np.linspace(0,L,nL):
      nodesIds2.append(O.bodies.append(gridNode([cen_x+.001,i,cen_z],rCyl,wire=False,fixed=False,material='grass')))
    idxc += 1
    d = idxc*nL ## Start of grass fiber
    cen_x += x_gap
  cen_z += z_gap 


## Creates Nodes Connetions
move_up = 0   ## Moving up to node above 
for k in range(0,nL-1):
  for i in range(0,len(nodesIds),nL):
    pfacetCreator3(nodesIds[i+move_up],nodesIds2[i+move_up],nodesIds2[i+1+move_up],cylIds=[],pfIds=[],wire=False,material='grass',color=color) 
  for i in range(0,len(nodesIds2),nL):
    pfacetCreator3(nodesIds2[i+1+move_up],nodesIds[i+1+move_up],nodesIds[i+move_up],cylIds=[],pfIds=[],wire=False,material='grass',color=color) 
  move_up += 1

## Creates BC's for all Nodes 
for kk in range(0,len(nodesIds)):
  O.bodies[nodesIds[kk]].state.blockedDOFs='xyzXYZ' 
  O.bodies[nodesIds2[kk]].state.blockedDOFs='xyzXYZ' 

## Blocking DOF's for top and bottom Nodes
btm_nodes = range(0,len(nodesIds),nL)
top_nodes = range(nL-1,len(nodesIds),nL)
def turnongrass():
	print "Turning on Grass"
	for kk in range(0,len(nodesIds)):
#		if any(x == kk for x in top_nodes) == False:              ## Allows only the top nodes to have ZERO DOF's
		  if any(xx == kk for xx in btm_nodes) == False:          ## Allows only the bottom nodes to have ZERO DOF's
 			O.bodies[nodesIds[kk]].state.blockedDOFs=''       ## Unblocking all DOF's - Except Above
 			O.bodies[nodesIds2[kk]].state.blockedDOFs=''      ## Unblocking all DOF's - Except Above

def run(self):
	iterToRun=All_steps
	O.run(iterToRun, True)

def switchTranslationEngine():
	print "Switch from TranslationEngine engine to RotationEngine"
	transeng.dead = True
	rotengine.dead = False
	# Force Recorders
	trans_force.dead = True
	rotation_force.dead = False
	switchEng.dead = True

from yade import qt
qt.View()

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