← Back to team overview

yade-users team mailing list archive

[Question #680139]: rigid (gts)Pfacets?

 

New question #680139 on Yade:
https://answers.launchpad.net/yade/+question/680139

Dear all,

recently I saw the pfacet element implemented in yade and had a look at [Effeindzourou2016].

The pfacets seem to be constructed to model deformable objects, like membranes. I was wondering if it is possible to use them for shape modelling of rigid particles (either a single pfacet object or a gtsPfacet consisting of several pfacets). 

Probably this was not done before, but I would appreciate your opinion if you expect minor or major changes in the code necessary. 

What I tried so far:

In  660289: how to simulate a rigid body when use gtsPfacet and CohFrictMat  it was advised to clumps several pfacets together with O.bodies.appendClumped(). 

My MWE is stated below. It’s a simplification from the trunk/examples/pfacet/gts-pfacet.py. At the beginning, I copied modified code  from trunk/py/gridpfacet.py  to create a clumped gtsPfacet. In my code, the gridNodes and gridConnections are added to the scene with O.bodies.append(). A list of pfacets is then added with O.bodies.appendClumped().

The clump of pfacets is created, but in the Viewer the pfacets disappear after the first time step. Then, the particle falls through the wall. I guess this is caused, because the pfacet elements are missing. Remaining are gridNodes and gridConnections and for them no interactions are defined with the wall. Any idea how to construct the clump correctly?

Thanks for your help and best regards,
Bettina

MWE: This reads trunk/examples/pfacet/box.gts, copy to working directory.
#---------------------------------------------------------------------------------------------------------------------
from yade import qt
from yade.gridpfacet import *
import gts, os.path, locale

locale.setlocale(locale.LC_ALL, 'en_US.UTF-8')   # Note: gts is locale-dependent. If, for example, german locale is used, gts.read()-function does not import floats normally

'''
if you get "Error: unsupported locale setting"
-> type as root: "dpkg-reconfigure locales"
-> choose "en_US.UTF-8" (press space to choose)
'''

#-------------------------------------------------------------------------------------
# adapted functions from trunk/py/gridpfacet.py to make clump of pfacets
#-------------------------------------------------------------------------------------

def ClumpedGtsPFacet(meshfile,shift=Vector3.Zero,scale=1.0,
							radius=1,wire=True,fixed=True,materialNodes=-1,material=-1,color=None):
	"""
	Imports mesh geometry from .gts file and automatically creates connected :yref:`PFacet3<PFacet>` elements. For an example see :ysrc:`examples/pfacet/gts-pfacet.py`.

	:param string filename: .gts file to read.
	:param [float,float,float] shift: [X,Y,Z] parameter shifts the mesh.
	:param float scale: factor scales the mesh.
	:param float radius: radius used to create the :yref:`PFacets<PFacet>`.
	:param materialNodes: specify :yref:`Body.material` of :yref:`GridNodes<GridNode>`. This material is used to make the internal connections.
	:param material: specify :yref:`Body.material` of :yref:`PFacets<PFacet>`. This material is used for interactions with external bodies.

	See documentation of :yref:`yade.utils.sphere` for meaning of other parameters.

	:returns: lists of :yref:`GridNode<GridNode>` ids `nodesIds`, :yref:`GridConnection<GridConnection>` ids `cylIds`, and :yref:`PFacet<PFacet>` ids `pfIds`
	"""
	import gts,yade.pack
	surf=gts.read(open(meshfile))
	surf.scale(scale,scale,scale)
	surf.translate(shift[0],shift[1],shift[2]) 
	nodesIds=[]; cylIds=[]; pfIds=[]
	pfacetList=[]
	for face in surf.faces():
		a=face.vertices()[0].coords()
		b=face.vertices()[1].coords()
		c=face.vertices()[2].coords()
		pfacetList.append(ClumpedPfacetCreator1([a,b,c],radius=radius,nodesIds=nodesIds,cylIds=cylIds,pfIds=pfIds,wire=wire,fixed=fixed,materialNodes=materialNodes,material=material,color=color))
		#print a,b,c
	pfIds=O.bodies.appendClumped(pfacetList)
	return nodesIds,cylIds,pfIds



def ClumpedPfacetCreator1(vertices,radius,nodesIds=[],cylIds=[],pfIds=[],wire=False, fixed=True,materialNodes=-1,material=-1,color=None):
	"""
	Create a :yref:`PFacet<PFacet>` element from 3 vertices and automatically append to simulation. The function uses the vertices to create :yref:`GridNodes<GridNode>` and automatically checks for existing nodes.
	
	:param [Vector3,Vector3,Vector3] vertices: coordinates of vertices in the global coordinate system.
	:param float radius: radius used to create the :yref:`PFacets<PFacet>`.
	:param list nodesIds: list with ids of already existing :yref:`GridNodes<GridNode>`. New ids will be added.
	:param list cylIds: list with ids of already existing :yref:`GridConnections<GridConnection>`. New ids will be added.
	:param list pfIds: list with ids of already existing :yref:`PFacets<PFacet>`. New ids will be added. 
	:param materialNodes: specify :yref:`Body.material` of :yref:`GridNodes<GridNode>`. This material is used to make the internal connections.
	:param material: specify :yref:`Body.material` of :yref:`PFacets<PFacet>`. This material is used for interactions with external bodies.
	
	See documentation of :yref:`yade.utils.sphere` for meaning of other parameters.
	"""
	n=len(nodesIds)
	k=[0,0,0]
	f=[0,0,0]
	u=0
	nod=0
	for i in vertices:
		u=0
		for j in nodesIds:
			if(i==O.bodies[j].state.pos):
				f[nod]=j
				k[nod]=1
				u+=1
		nod+=1
		test=True
		#if(u==0):
		for GN in nodesIds:
			if(i==O.bodies[GN].state.pos):
				  u=1
		if(u==0):
			nodesIds.append( O.bodies.append(gridNode(i,radius,wire=wire,fixed=fixed,material=materialNodes,color=color)) )

	if(k==[0,0,0]):
		pfObj=ClumpedPfacetCreator3(nodesIds[n],nodesIds[n+1],nodesIds[n+2],cylIds=cylIds,pfIds=pfIds,wire=wire,material=material,color=color,fixed=fixed )
	if(k==[1,0,0]):
		pfObj=ClumpedPfacetCreator3(f[0],nodesIds[n],nodesIds[n+1],cylIds=cylIds,pfIds=pfIds,wire=wire,material=material,color=color,fixed=fixed )
	if(k==[0,1,0]):
		pfObj=ClumpedPfacetCreator3(nodesIds[n],f[1],nodesIds[n+1],cylIds=cylIds,pfIds=pfIds,wire=wire,material=material,color=color,fixed=fixed )
	if(k==[0,0,1]):
		pfObj=ClumpedPfacetCreator3(nodesIds[n],nodesIds[n+1],f[2],cylIds=cylIds,pfIds=pfIds,wire=wire,material=material,color=color,fixed=fixed )
	if(k==[1,1,0]):
		pfObj=ClumpedPfacetCreator3(f[0],f[1],nodesIds[n],cylIds=cylIds,pfIds=pfIds,wire=wire,material=material,color=color,fixed=fixed )
	if(k==[0,1,1]):
		pfObj=ClumpedPfacetCreator3(nodesIds[n],f[1],f[2],cylIds=cylIds,pfIds=pfIds,wire=wire,material=material,color=color,fixed=fixed )
	if(k==[1,0,1]):
		pfObj=ClumpedPfacetCreator3(f[0],nodesIds[n],f[2],cylIds=cylIds,pfIds=pfIds,wire=wire,material=material,color=color,fixed=fixed )
	if(k==[1,1,1]):
		pfObj=ClumpedPfacetCreator3(f[0],f[1],f[2],cylIds=cylIds,pfIds=pfIds,wire=wire,material=material,color=color,fixed=fixed )
	return pfObj



def ClumpedPfacetCreator3(id1,id2,id3,cylIds=[],pfIds=[],wire=True,material=-1,color=None,fixed=True,mask=-1):
	"""
	Create a :yref:`PFacet` element from 3 already existing :yref:`GridNodes<GridNode>` which are not yet connected. The element is automatically appended to the simulation.
	
	:param int id1,id2,id3: id of the 3 :yref:`GridNodes<GridNode>` forming the :yref:`PFacet`.
	
	See documentation of :yref:`yade.gridpfacet.pfacetCreator1` for meaning of other parameters.
	"""  
	radius=O.bodies[id1].shape.radius
	try:
		cylIds.append(O.bodies.append(gridConnection(id1,id2,radius=radius,material=material,color=color,wire=wire,mask=mask)))
	except:
		pass
	try:
		cylIds.append(O.bodies.append(gridConnection(id2,id3,radius=radius,material=material,color=color,wire=wire,mask=mask)))
	except:
		pass	      
	try:
		cylIds.append(O.bodies.append(gridConnection(id3,id1,radius=radius,material=material,color=color,wire=wire,mask=mask)))
	except:
		pass	      
	#pfIds.append(O.bodies.appendClumped())
	return pfacet(id1,id2,id3,wire=wire,material=material,color=color,mask=mask)

#-------------------------------------------------------------------------------------

################
### ENGINES  ###
################

O.engines=[
	ForceResetter(),
	InsertionSortCollider([
		Bo1_Wall_Aabb(),
		Bo1_PFacet_Aabb(),
	],sortThenCollide=True),
	InteractionLoop(
	[
        Ig2_GridNode_GridNode_GridNodeGeom6D(),
		Ig2_Wall_PFacet_ScGeom(),Ig2_Wall_Sphere_ScGeom()
	],
	[
        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.5,label='ts'), 
	NewtonIntegrator(gravity=(0,-9.81,0),damping=0.1,label='newton')
]

#disabled GlobalStiffnessTimeStepper due to problems, adapt dt if necessary!
O.dt=1e-6

#################
### MATERIAL  ###
#################

O.materials.append(CohFrictMat(young=1e8,poisson=0.3,density=2650,frictionAngle=radians(20),normalCohesion=1e100,shearCohesion=1e100,momentRotationLaw=True,label='gridNodeMat'))
O.materials.append(FrictMat(young=1e8,poisson=0.3,density=2650,frictionAngle=radians(20),label='pFacetMat'))

###################
### IMPORT MESH ###
###################
radius=0.5
wire=False
fixed=False

z=-1.6
color=[0,0,1]

#--- original gtsPFacet construction --- WORKING
#nodesIds0,cylIds0,pfIds0 =  gtsPFacet('box.gts',shift=(0.,0.,0.),scale=2.,radius=radius,wire=wire,fixed=fixed,materialNodes='gridNodeMat',material='pFacetMat',color=color)

#--- try to make a clump out of it --- NOT WORKING
nodesIds0,cylIds0,pfIds0 =ClumpedGtsPFacet('box.gts',shift=(0,0,0),scale=2.,radius=radius,wire=wire,fixed=fixed,materialNodes='gridNodeMat',material='pFacetMat',color=color)

#--- mass and inertia should be updated?
O.bodies.updateClumpProperties(discretization=10)

#####################
#####   Wall      ###
#####################

O.bodies.append(utils.wall(position=z,sense=0, axis=1,color=Vector3(1,0,0),material='pFacetMat'))

##########
## VIEW ##
##########

qt.Controller()
qtv = qt.View()
qtr = qt.Renderer()
qtr.light2=True
qtr.lightPos=Vector3(1200,1500,500)
qtr.bgColor=[1,1,1]
qtv.ortho=True


O.saveTmp()



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