yade-dev team mailing list archive

[Anton Gladky <gladky.anton@xxxxxxxxx>]

I did not know about that inspector :)

Ok, after the first step the acceleration is getting very high, after next
step velocity => positions etc.

I have attached the script. Please, check, do you have the same problem?

Thank you,

Anton



On Tue, Sep 7, 2010 at 1:51 PM, Vaclav Smilauer <eudoxos@xxxxxxxx> wrote:

> > I have installed Maverick 10.10 on virtual machine, I have also
> > catched the problem. So there is no problem with hardware.
> Do you mean you  still have the problem? In that case say what are
> positions of the particles. If you run the Inspector (from the main
> window), that should give you an idea. Are they insanely high values
> (for positions, for instance) or nans? Is it positions or orientations?
> Or velocities?
>
> Cheers, v.
>
>
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-dev<https://launchpad.net/%7Eyade-dev>
> Post to     : yade-dev@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~yade-dev<https://launchpad.net/%7Eyade-dev>
> More help   : https://help.launchpad.net/ListHelp
>

#!/usr/local/bin/yade-trunk -x
# -*- encoding=utf-8 -*-

import itertools, os, system
import math, numpy
from yade import pack,ymport,utils, plot

o=Omega()
X = 0
Y = 0
fiInput=0.1
radInput=2.5
NvInput=1
stepsNumber=1000000

RSpiral = 0.315/2.0
rSpiral = 0.159/2.0
hSpiral = 0.25

d = rSpiral*2.0
D = RSpiral*2.0
H = hSpiral
N=150

rad=radInput*1.0e-3
v_angular = N*2*math.pi/60

frictionAngle = 0.52
rho=2500


fi = fiInput
Nv = NvInput


#******************************************************************************************************
tc = 0.001
en = 0.008
es = 0.008
params=utils.getViscoelasticFromSpheresInteraction(10e3,tc,en,es)
facetMat=O.materials.append(ViscElMat(frictionAngle=frictionAngle,**params)) # **params sets kn, cn, ks, cs
sphereMat=O.materials.append(ViscElMat(density=rho,frictionAngle=frictionAngle,**params))
#******************************************************************************************************

kw_facets={'wire':True,'dynamic':False,'highlight':False}

#_______________________________________________
Vmax = math.pi*(D*D - d*d)*N*H/240
V1vitok = math.pi*(D*D - d*d)*H/4	
Vc = Nv*fi*V1vitok								
Voutput = Vmax*fi									
Vsph = 4.0/3.0*math.pi*rad*rad*rad
Ns = Vc/Vsph											
#_______________________________________________

HcylFill = 4*Vc/(math.pi*(D*D - d*d)*0.6)

ids_spheres=O.bodies.append(pack.regularHexa((pack.inCylinder((X,Y,0.0),(X,Y,-HcylFill),D/2.0)-pack.inCylinder((X,Y,0.0),(X,Y,-HcylFill),d/2.0)),radius=rad,gap=0.0,color=(0.929,0.412,0.412),mask=3,material=sphereMat))

s=O.bodies[ids_spheres[0]]
p=utils.getViscoelasticFromSpheresInteraction(s.state.mass,tc,en,es)
s.mat.kn,s.mat.cn,s.mat.ks,s.mat.cs=p['kn'],p['cn'],p['ks'],p['cs']

#Find the hihghest and the lowest sphere________________________________
maxHight=O.bodies[ids_spheres[0]].state.pos[2]+O.bodies[ids_spheres[0]].shape.radius
minHight=O.bodies[ids_spheres[0]].state.pos[2]-O.bodies[ids_spheres[0]].shape.radius

for i in ids_spheres:
	if ((O.bodies[i].state.pos[2]+O.bodies[ids_spheres[i]].shape.radius)>maxHight):
		maxHight=O.bodies[i].state.pos[2]+O.bodies[ids_spheres[i]].shape.radius
	if ((O.bodies[i].state.pos[2]-O.bodies[ids_spheres[i]].shape.radius)<minHight):
		minHight=O.bodies[i].state.pos[2]-O.bodies[ids_spheres[i]].shape.radius

realHight = 0-minHight
centerHight = minHight + realHight/2.0


o.dt=0.02*tc


#________________________________Sphere Transfer________________________

pressId=O.bodies.append(utils.facet(((X-2.0*D,Y-1.0*D,minHight),(X+2.0*D,Y-1.0*D,minHight),(X,Y+2.0*D,minHight)),mask=33,material=facetMat,**kw_facets))
O.bodies.append(utils.facetCylinder((X,Y,centerHight),D/2.0,realHight,wallMask=4,segmentsNumber=30,mask=33,material=facetMat,**kw_facets))
O.bodies.append(utils.facetCylinder((X,Y,centerHight),d/2.0,realHight,wallMask=4,segmentsNumber=30,mask=17,material=facetMat,**kw_facets))

VideoOn=1

o.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
	InteractionDispatchers(
		[Ig2_Sphere_Sphere_ScGeom(), Ig2_Facet_Sphere_ScGeom()],
		[Ip2_ViscElMat_ViscElMat_ViscElPhys()],
		[Law2_ScGeom_ViscElPhys_Basic()],
	),
	GravityEngine(gravity=[0.0,0.0,-9.810]),
	NewtonIntegrator(damping=.0),
	TranslationEngine(translationAxis=[0,0,1],velocity=0.1,subscribedBodies=[pressId])
]

O.saveTmp('init')

if (VideoOn):
	from yade import qt
	qt.Controller()
	r=qt.Renderer()
	r.bgColor=1,1,1
	qt.View()


#O.run(stepsNumber)
O.wait()