← Back to team overview

yade-users team mailing list archive

Re: [Question #689434]: How to generate many polyhedra randomly without specific data in Potential Blocks code?

 

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

    Status: Answered => Open

weijie is still having a problem:
Hi Vasileios,
      I tried to make a random polyhedron with reference to the information you provided. However, I encountered a difficulty. A part of the converted polyhedron is displayed in black and it does not feel closed.How can I solve it?
    Following is my code:
#####################
from yade import polyhedra_utils,pack,plot,utils,export,qt,ymport
import numpy as np
import math
import random
#########################################################################Define physical parameters
powderDensity = 2500
O.materials.append(FrictMat(young=-1,poisson=-1,frictionAngle=radians(0.0),density=powderDensity,label='frictionless'))
young = 1e9
FricDegree1= atan(.6)
m = PolyhedraMat()
m.density = powderDensity 
m.young = young
m.poisson = 0.25
m.frictionAngle = FricDegree1
#########################################################################Forming polyhedra with polyhedra function
t = polyhedra_utils.polyhedra(m,size = (0.06,0.06,0.06),seed =5)
#t.state.pos = (0.,0.,0.5)
O.bodies.append(t)
########################################################################Polyhedron transformation
def dvalue(vecn1,pp1):
 dd1=-1*(vecn1[0]*pp1[0]+vecn1[1]*pp1[1]+vecn1[2]*pp1[2])
 return dd1  
chosenR=0.0001
for b in O.bodies:
 aa=[]
 bb=[]
 cc=[]
 dd=[]
 if not isinstance(b.shape,Polyhedra): # skip non-polyhedra bodies
  continue
 vs = [b.state.pos + b.state.ori*v for v in b.shape.v] # vertices in global coords
 face2=b.shape.GetSurfaces()

 id1=0
 while id1<len(face2):
  face11=face2[id1]
  if len(face11)>2:
   vec1=vs[face11[2]]-vs[face11[1]]; vec1.normalize() 
   vec2=vs[face11[0]]-vs[face11[1]]; vec2.normalize() #Normalize this object in-place.
   vects=vec1.cross(vec2) 

   dvalue2=dvalue(vects,vs[face11[0]]) # Some dvalue2 values return equal to 0. Check this part of your script once more.
   dv=dvalue2-chosenR

   if dv<0:
    vects=[-1*vects[0],-1*vects[1],-1*vects[2]]
    dv=-1*dv

   aa.append(vects[0])
   bb.append(vects[1])
   cc.append(vects[2])
   dd.append(dv)                     

  id1=id1+1
#######################################################################################
O.engines=[
	ForceResetter(),
	InsertionSortCollider([PotentialBlock2AABB()],verletDist=0.01, avoidSelfInteractionMask=2),
	InteractionLoop(
		[Ig2_PB_PB_ScGeom(twoDimension=False, unitWidth2D=1.0, calContactArea=True)],
		[Ip2_FrictMat_FrictMat_KnKsPBPhys(kn_i=1e8, ks_i=1e7, Knormal=1e8, Kshear=1e7, useFaceProperties=False, viscousDamping=0.2)],
	),
	#GlobalStiffnessTimeStepper(),
	NewtonIntegrator(damping=0.0,exactAsphericalRot=True,gravity=[0,-9.81,0]),
	#PotentialBlockVTKRecorder(fileName='./vtk/cubePBscaled',iterPeriod=10000,twoDimension=False,sampleX=50,sampleY=50,sampleZ=50,maxDimension=0.2,label='vtkRecorder')
]
#print(aa,bb,cc,dd)
bbb=Body()
wire=False
color=[125,2,1]
highlight=True
bbb.shape=PotentialBlock(k=0.0, r=chosenR, R=0.0, a=aa, b=bb, c=cc, d=dd, id=len(O.bodies)-1)
utils._commonBodySetup(bbb, bbb.shape.volume, bbb.shape.inertia,material='frictionless', pos=[0,0,0], fixed=False)
bbb.state.pos = [0,0.1,0]
# bbb.state.ori=Quaternion((0,1,0),pi)
lidID = O.bodies.append(bbb)

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