yade-dev team mailing list archive

[Václav S(milauer <eudoxos@xxxxxxxx>]

oh well.... (it is a complete simple simulation you can run)

# encoding: utf-8

from yade.wrapper import *
def facet(vertices,relative=False,young=30e9,poisson=0.2,color=[1,1,1]):
	"""First vertex is the reference point for this Facet (always in absolute coordinates). If relative is True, other vertices are given relatively to the reference point."""
	if not relative:
		for i in range(1,len(vertices)):
			vertices[i]=[vertices[i][j]-vertices[0][j] for j in range(3)]
	print vertices
	vStr='['+' '.join(['{%g %g %g}'%(v[0],v[1],v[2]) for v in vertices[1:]])+']'
	b=Body()
	b.shape=GeometricalModel("Facet",{'wire':True,'visible':True,'diffuseColor':color})
	b.shape.setRaw('vertices',vStr)
	b.mold=InteractingGeometry("InteractingFacet",{'diffuseColor':color})
	b.mold.setRaw('vertices',vStr)
	b.phys=PhysicalParameters("BodyMacroParameters",{'se3':[vertices[0][0],vertices[0][1],vertices[0][2],1,0,0,0],'mass':0,'inertia':[0,0,0],'young':young,'poisson':poisson})
	b.bound=BoundingVolume('AABB',{'diffuseColor':[0,1,0]})
	b['isDynamic']=False
	return b

o=Omega()
from yade import utils
from math import *

if 1:
	o.bodies.append(facet(vertices=[[1,0,0],[0,1,0],[-1,0,0],[0,-1,0]],young=7e4))
elif 0:
	nEdges=20
	fAngle=2*pi/nEdges
	for i in range(0,nEdges):
		o.bodies.append(facet([0,0,0],[[sin(fAngle*i),cos(fAngle*i),-1],[sin(fAngle*(i+1)),cos(fAngle*(i+1)),-1]]))
		o.bodies.append(facet([i,0,0],[[0,0,0],[0,1,0],[1,1,0],[1,0,0]],relative=True))
else:
	stl=STLImporter()
	stl.import_vertices = True
	stl.import_edges = True
	stl.import_facets = True
	stl.facets_wire = True
	stl.open('cone.stl')
	for i in xrange(stl.number_of_all_imported):
		b=Body()
		b['isDynamic']=False
		b.phys=PhysicalParameters('BodyMacroParameters',{'se3':[0,0,0,1,0,0,0],'mass':1,'inertia':[1,1,1],'young':30e9})
		if(i>=stl.number_of_imported_vertices): b.bound=BoundingVolume('AABB',{'diffuseColor':[0,1,0]})
		o.bodies.append(b)
	stl.import_geometry(o.bodies,0)

o.bodies.append(utils.sphere([.1,0,-10],0.2,wire=False,young=7e4))
o.bodies.append(utils.sphere([0,-.1,-1.5],0.2,young=7e4))
o.bodies.append(utils.sphere([0,.1,-2],0.2,young=7e4))

for b in o.bodies: b.mold.postProcessAttributes()
o.dt=1e-5
#o.dt=1e2*utils.PWaveTimeStep()

o.initializers=[
	StandAloneEngine('PhysicalActionContainerInitializer'),
	MetaEngine('BoundingVolumeMetaEngine',[EngineUnit('InteractingSphere2AABB'),EngineUnit('InteractingEdge2AABB'),EngineUnit('InteractingFacet2AABB'),EngineUnit('MetaInteractingGeometry2AABB'),EngineUnit('InteractingBox2AABB')])
]

o.engines=[
	StandAloneEngine('PhysicalActionContainerReseter'),
	MetaEngine('BoundingVolumeMetaEngine',[
		EngineUnit('InteractingSphere2AABB'),
		EngineUnit('InteractingBox2AABB'),
		EngineUnit('InteractingEdge2AABB'),
		EngineUnit('InteractingFacet2AABB'),
		EngineUnit('InteractingBox2AABB'),
		EngineUnit('MetaInteractingGeometry2AABB')
	]),
	StandAloneEngine('PersistentSAPCollider'),
	MetaEngine('InteractionGeometryMetaEngine',[
		EngineUnit('InteractingSphere2InteractingSphere4SpheresContactGeometry'),
		EngineUnit('InteractingBox2InteractingSphere4SpheresContactGeometry'),
		EngineUnit('InteractingVertex2InteractingSphere4SpheresContactGeometry'),
		EngineUnit('InteractingEdge2InteractingSphere4SpheresContactGeometry'),
		EngineUnit('InteractingFacet2InteractingSphere4SpheresContactGeometry'),
	]),
	MetaEngine('InteractionPhysicsMetaEngine',[EngineUnit('SimpleElasticRelationships')]),
	StandAloneEngine('ElasticContactLaw'),
	DeusExMachina('GravityEngine',{'gravity':[0,0,10]}),
	DeusExMachina('NewtonsDampedLaw',{'damping':0.2}),
	#StandAloneEngine('PositionOrientationRecorder',{'interval':60,'outputFile':'/tmp/_facet-'})
]

o.run(5)
o.save('/tmp/a.xml'); #o.load('/tmp/a.xml')
#o.run()