O.engines[].dead=True

We need an MWE to properly assess [1].

[1]https://www.yade-dem.org/wiki/Howtoask

Hi,
Please provide a working script.
Bruno

Le ven. 17 mai 2019 09:37, Francesco <email address hidden>
a écrit :

> New question #680892 on Yade:
> https://answers.launchpad.net/yade/+question/680892
>
> Hi yade users,
>
> I'm using the rotation and translation engines to move bodies and pack of
> particles.
> The simulation starts with the engines set to dead= false.
> After I use Pyrunners in order to stop the movements. For example these
> are the engines:
>
>
> RotationEngine(angularVelocity=-wt,dead=False,ids=nozzle,rotateAroundZero=True,rotationAxis=(1,0,0),zeroPoint=[0,y_rotation,0]),
>
>
> RotationEngine(angularVelocity=-wt,dead=False,ids=pack_rotation,rotateAroundZero=True,rotationAxis=(1,0,0),zeroPoint=[0,y_rotation,0]),
>
> This is the function of the pyrunner:
>
> def motion_5():
> O.engines[11].dead=True
> O.engines[12].dead=True
>
> And this is the pyrunner:
>
>
> PyRunner(command='motion_5()',iterPeriod=1,firstIterRun=iterrotation_2,nDo=2)
>
> I do not understand why the rotations don't stop.
> Does somebody find the same problem?
>
> Thank you so much,
>
> Francesco
>
>
> --
> You received this question notification because your team yade-users is
> an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>
>
>

Hi,

I attached the code. I commented some functions related to data storing that could create problems.

Thank you so much.

Francesco

Code:

from numpy import linspace
from yade import pack, qt
import numpy as np
import os
import glob
import re
import smtplib
from email.MIMEMultipart import MIMEMultipart
from email.MIMEText import MIMEText
from email.MIMEBase import MIMEBase
from email import encoders

######################################################
# Parameters
######################################################
mi=10**-3 # Millimetres transformation
mu=10**-6 # Micron transformation
x1=10*mi # Box half side
z1=10*mi # Box heigth
dist=100*mi # Nozzle target distance
h=z1+(dist) # Nozzle heigth from z=0
rIn=6*mi # Inner radius
rOut=12*mi # Cylinder radius
eps=0.1*mi # Tollerance
itereps=5
rPack=rOut-eps # Pack radius
rd=rOut-rIn # Difference between radius
l2=80*mi # Cylinder length
l_cloud=50*mi # Particles cloud length
c1=(4*l_cloud)/9 # External particles pack range
c2=(1*l_cloud)/9 # Middle particles pack range
l3=1*mi # Distance from the nozle indentation
v=57 # Particles initial velocity
vi_y=-v*cos(pi*angle/180)
vi_z=-v*sin(pi*angle/180)
vf1=57 # Fluid velocity 1
vf2=220 # Fluid velocity 2
rMean=40*mu # Mean particles radius
Rho=3.61 # Air density 1
Cd=0.47 # Sphere drag coefficient
wt=100000 # Box rotation
angle=22.5 # Nozzle angle
box_angle=45
vx=10000
vz=10000

######################################################
# Nozzle creation
######################################################
delta=((pi*0.5)/16)
poly=((rIn,h),
 (rOut-(rd*cos(delta)),h+(rIn*sin(delta))),
 (rOut-(rd*cos(2*delta)),h+(rIn*sin(2*delta))),
 (rOut-(rd*cos(3*delta)),h+(rIn*sin(3*delta))),
 (rOut-(rd*cos(4*delta)),h+(rIn*sin(4*delta))),
 (rOut-(rd*cos(5*delta)),h+(rIn*sin(5*delta))),
 (rOut-(rd*cos(6*delta)),h+(rIn*sin(6*delta))),
 (rOut-(rd*cos(7*delta)),h+(rIn*sin(7*delta))),
 (rOut-(rd*cos(8*delta)),h+(rIn*sin(8*delta))),
 (rOut-(rd*cos(9*delta)),h+(rIn*sin(9*delta))),
 (rOut-(rd*cos(10*delta)),h+(rIn*sin(10*delta))),
 (rOut-(rd*cos(11*delta)),h+(rIn*sin(11*delta))),
 (rOut-(rd*cos(12*delta)),h+(rIn*sin(12*delta))),
 (rOut-(rd*cos(13*delta)),h+(rIn*sin(13*delta))),
 (rOut-(rd*cos(14*delta)),h+(rIn*sin(14*delta))),
 (rOut-(rd*cos(15*delta)),h+(rIn*sin(15*delta))),
 (rOut,h+rIn),(rOut,h+rIn+(l2/8)),(rOut,h+rIn+(l2/4)),
 (rOut,h+rIn+((3*l2)/8)),(rOut,h+rIn+((4*l2)/8)),(rOut,h+rIn+((5*l2)/8)),
 (rOut,h+rIn+((6*l2)/8)),(rOut,h+rIn+((7*l2)/8)),(rOut,h+rIn+(l2))
)
thetas=linspace(0,2*pi,num=32,endpoint=True)
Nozzle=O.materials.append(FrictMat(density=7.8*10**3,young=210e9,poisson=0.3,label="Nozzle"))
pts=pack.revolutionSurfaceMeridians([[(pt[0],pt[1]) for pt in poly] for theta in thetas],thetas)
surf=pack.sweptPolylines2gtsSurface(pts)
nozzle=O.bodies.append(pack.gtsSurface2Facets(surf))

######################################################
# Creating Target
######################################################

Aluminum=O.materials.append(FrictMat(density=2.81*10**3,young=71e9,poisson=0.33,label="Aluminum"))
p=utils.box((x1,0,0),(x1,x1,z1),fixed=True,material="Aluminum")
box=O.bodies.append(p)
l=len(O.bodies)
t_id=l-1
p1=utils.box((-x1,0,0),(x1,x1,z1),fixed=True,material="Aluminum")
box1=O.bodies.append(p1)
l_1=len(O.bodies)
t_id1=l_1-1
p2=utils.box((2*x1,0,0),(2*x1,x1,z1),fixed=True,material="Aluminum")
box2=O.bodies.append(p2)
l_2=len(O.bodies)
t_id2=l_2-1
p3=utils.box((-2*x1,0,0),(2*x1,x1,z1),fixed=True,material="Aluminum")
box3=O.bodies.append(p3)
l_3=len(O.bodies)
t_id3=l_3-1
p4=utils.box((x1,-2*x1,0),(x1,x1,z1),fixed=True,material="Aluminum")
box4=O.bodies.append(p4)
l_4=len(O.bodies)
t_id4=l_4-1

######################################################
# Creating shots
######################################################
V_cloud=((2*rOut)*(2*rOut))*l_cloud*1.1
V_cyl=(pi*rOut*rOut*(l_cloud))
n_shot=1287
dens=(n_shot/V_cyl)
n_shot_cloud=int(dens*V_cloud)
corner=(rOut)
Ceramic=O.materials.append(FrictMat(density=2.3*10**3,young=300e9,poisson=0.30,label="Ceramic"))
sp=pack.SpherePack()
sp.makeCloud((-corner,-corner,h+rIn+l3),(corner,corner,h+rIn+l3+(1.1*l_cloud)),num=n_shot_cloud,psdSizes=[53*10**-6,63*10**-6,90*10**-6,106*10**-6,125*10**-6],psdCumm=[0,0.42,0.76,0.92,1])
pred_1=pack.inCylinder((0,0,h+rIn+l3),(0,0,h+rIn+l3+(c1)),radius=rPack)
pred_2=pack.inCylinder((0,0,h+rIn+l3+c1),(0,0,h+rIn+l3+(c1+c2)),radius=rPack)
pred_3=pack.inCylinder((0,0,h+rIn+l3+c1+c2),(0,0,h+rIn+l3+(c1+c2+c1)),radius=rPack)
sp2=pack.filterSpherePack(pred_1,sp,returnSpherePack=True)
sp3= pack.filterSpherePack(pred_2,sp,returnSpherePack=True)
sp4=pack.filterSpherePack(pred_3,sp,returnSpherePack=True)
sp2.toSimulation()
j1=len(O.bodies)
sp3.toSimulation()
j2=len(O.bodies)
sp4.toSimulation()
j=len(O.bodies)

pack_rotation=range(l_4+1,j)

O.dt=utils.PWaveTimeStep()

([utils.setBodyVelocity(i,(0,vi_y,vi_z),axis="xyz") for i in range(l_4,j)])

#periodrecord=20
#dzrecord=1*mi
#nrecord=int((l_cloud+rIn+l3)/dzrecord)
#startrecordtime=0
#endrecordtime=int((1*mi)/(O.dt*v))
#totalrecordtime=endrecordtime*nrecord
#deltarecordtime=(endrecordtime-startrecordtime)
#ndorecord=(int(deltarecordtime/periodrecord))*nrecord*2
#deltarecord=(endrecordtime-startrecordtime)*O.dt
#firstrecord=40
#lastrecord=45
#hin=h+rIn+l3
#shot=[]
#mass=[]
#Vnorm=[]
#def txt():
# for ii in range(l_4,j):
# if O.bodies[ii].state.pos[2]>h:
# if O.bodies[ii].state.pos[2]<hin:
# vb=O.bodies[ii].state.vel[2]
# radius=O.bodies[ii].shape.radius
# Fd=0.5*Rho*Cd*(pi*(radius)**2)*(vb-(-vf2))**2
# O.forces.setPermF(ii,(0,0,-Fd))
# else:
# if O.bodies[ii].state.vel[2]>0:
# vb=O.bodies[ii].state.vel[2]
# radius=O.bodies[ii].shape.radius
# Fd=0.5*Rho*Cd*(pi*(radius)**2)*(vb-(-vf2))**2
# O.forces.setPermF(ii,(0,0,-Fd))
# else:
# O.forces.setPermF(ii,(0,0,0))
# if O.bodies[ii].state.pos[2] < h:
# if O.bodies[ii].state.pos[2] > h-eps:
# name="/home/michelangelo/Scrivania/YADE_PIONA/NOZZLE_FLOW_EDGE_37/RESULTS_2/velocity"+str(ii)+".txt"
# with open(name,"a") as f:
# if os.stat(name).st_size == 0 :
# vout=O.bodies[ii].state.vel
# V=(((vout[0])**2)+((vout[1])**2)+((vout[2])**2))**(0.5)
# z=O.bodies[ii].state.pos[2]
# m=O.bodies[ii].state.mass
# f.write("{:.2f}\t{:.2f}\t{:.2f}\t{:.3f}\t{}\n".format(vout[0],vout[1],vout[2],z,ii))
# mass.append(m)
# Vnorm.append(V)
# if O.time>(firstrecord*deltarecord):
# if O.time<(lastrecord*deltarecord):
# shot.append(ii)
# else:
# pass
# else:
# pass
# else:
# pass

#start=int((dist)/(O.dt*v))
#
#rmax=((125*10**-6)/2)
#
#def txt_3():
# global yref0
# global vref0
# for ii in shot:
# name="/home/michelangelo/Scrivania/YADE_PIONA/NOZZLE_FLOW_EDGE_37/RESULTS_3/contact"+str(ii)+".txt"
# with open(name,"a+") as f:
# if O.bodies[ii].state.vel[2]<-1:
# if os.stat(name).st_size == 0 :
# if O.bodies[ii].state.pos[0]<0:
# if O.bodies[ii].state.pos[2]<z1+(O.bodies[ii].state.pos[0]+rmax)*tan((pi*(angle))/180):
# cp=O.bodies[ii].state.pos
# vc=O.bodies[ii].state.vel
# t=O.time
# f.write("{:.6f}\t{:.6f}\t{:.6f}\t{:.2f}\t{:.2f}\t{:.2f}\t{}\t{}\n".format(cp[0],cp[1],cp[2],vc[0],vc[1],vc[2],ii,t))
# else:
# pass
# else:
# if O.bodies[ii].state.pos[2]<(z1+rmax)-(O.bodies[ii].state.pos[0])*tan((pi*(90-angle))/180):
# cp=O.bodies[ii].state.pos
# vc=O.bodies[ii].state.vel
# t=O.time
# f.write("{:.6f}\t{:.6f}\t{:.6f}\t{:.2f}\t{:.2f}\t{:.2f}\t{}\t{}\n".format(cp[0],cp[1],cp[2],vc[0],vc[1],vc[2],ii,t))
# else:
# pass
# else:
# if O.bodies[ii].state.pos[0]<0:
# if O.bodies[ii].state.pos[2]<z1+(O.bodies[ii].state.pos[0]+rmax)*tan((pi*(angle))/180):
# lines=f.readlines()
# for x in lines:
# old=(re.split('\t|\n',x))
# yref0=float(old[1])
# vref0=float(old[4])
# if yref0-mi<=round(O.bodies[ii].state.pos[1],6)<=yref0+mi:
# pass
# else:
# if round(O.bodies[ii].state.vel[1],2)==vref0:
# pass
# else:
# cp=O.bodies[ii].state.pos
# vc=O.bodies[ii].state.vel
# t=O.time
# f.write("{:.6f}\t{:.6f}\t{:.6f}\t{:.2f}\t{:.2f}\t{:.2f}\t{}\t{}\n".format(cp[0],cp[1],cp[2],vc[0],vc[1],vc[2],ii,t))
# else:
# pass
# else:
# if O.bodies[ii].state.pos[2]<(z1+rmax)-(O.bodies[ii].state.pos[0])*tan((pi*(90-angle))/180):
# lines=f.readlines()
# for x in lines:
# old=(re.split('\t|\n',x))
# yref0=float(old[1])
# vref0=float(old[4])
# if yref0-mi<=round(O.bodies[ii].state.pos[1],6)<=yref0+mi:
# pass
# else:
# if round(O.bodies[ii].state.vel[1],2)==vref0:
# pass
# else:
# cp=O.bodies[ii].state.pos
# vc=O.bodies[ii].state.vel
# t=O.time
# f.write("{:.6f}\t{:.6f}\t{:.6f}\t{:.2f}\t{:.2f}\t{:.2f}\t{}\t{}\n".format(cp[0],cp[1],cp[2],vc[0],vc[1],vc[2],ii,t))
# else:
# pass
#
#
#
#periodflow=deltarecordtime
#nDoflow=nrecord
#massflow=[]
#stepmass=[]
#Vmean=[]
#NewVnorm=[]
#c=0
#def flowfunction():
# global c
# if len(massflow) == 0:
# M=sum(mass)
# Q=(M/deltarecord)
# massflow.append(Q)
# stepmass.append(M)
# Vm=np.mean(Vnorm)
# Vmean.append(Vm)
# c=len(Vnorm)
# else:
# M=sum(mass)
# ii=len(stepmass)-1
# Mb=stepmass[ii]
# Q=((M-Mb)/deltarecord)
# massflow.append(Q)
# stepmass.append(M)
# NewVnorm=Vnorm[c:]
# Vm=np.mean(NewVnorm)
# Vmean.append(Vm)
# c=len(Vnorm)
#
#startflowprint=(startrecordtime+(nrecord*(deltarecordtime)))
#def flowprint():
# with open("/home/michelangelo/Scrivania/YADE_PIONA/NOZZLE_FLOW_EDGE_37/RESULTS_1/massflow.txt","a") as f:
# step=len(massflow)
# for ii in range(step):
# mf=massflow[ii]
# vm=Vmean[ii]
# f.write("{}\t{}\n".format(mf,vm))
#
#
#startend=int((dist)/(O.dt*v))
#zp=[]
#impacted=0
#def end():
# global impacted
# for ii in shot:
# name="/home/michelangelo/Scrivania/YADE_PIONA/NOZZLE_FLOW_EDGE_37/RESULTS_3/contact"+str(ii)+".txt"
# with open(name,"a+") as f:
# if os.stat(name).st_size == 0 :
# zp.append(0)
# else:
# zp.append(1)
# zplength=len(zp)
# impacted=zp.count(1)
# percentage=100*(impacted)*(zplength)**-1
# if percentage>95:
# O.pause()
# else:
# pass

rotationtime=(((pi*box_angle)/180)/(wt))
iterrotation=int(rotationtime/O.dt)
def motion():
 O.engines[4].angularVelocity=0
 O.engines[5].angularVelocity=0
 O.engines[10].angularVelocity=0

translationtime_x=(2*x1*(sin(pi*box_angle/180)))/vx
itertranslation_x=int(translationtime_x/O.dt)
stoptranslation_x=iterrotation+itertranslation_x

translationtime_z=(2*x1*(1-sin(pi*(90-box_angle)/180)))/vz
itertranslation_z=int(translationtime_z/O.dt)
stoptranslation_z=stoptranslation_x+itertranslation_z

def motion_1():
 O.engines[6].dead=False
 O.engines[7].dead=False

def motion_2():
 O.engines[6].dead=True
 O.engines[7].dead=True

def motion_3():
 O.engines[8].dead=False
 O.engines[9].dead=False

def motion_4():
 O.engines[8].velocity=0
 O.engines[9].velocity=0

y_rotation=x1-x1*(1-sin(pi*(90-angle)/180))
rotationtime_2=(((pi*angle)/180)/(wt))
iterrotation_2=int(rotationtime/O.dt)

def motion_5():
 O.engines[11].dead=True
 O.engines[12].dead=True

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Box_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_MindlinPhys(en=MatchMaker(matches=[(Ceramic,Ceramic,0.69),(Ceramic,Aluminum,0.43),(Ceramic,Nozzle,1)]))],
  [Law2_ScGeom_MindlinPhys_Mindlin()]
 ),
 NewtonIntegrator(gravity=(0,0,-9.81),damping=0.13,label="newton"),
 RotationEngine(angularVelocity=-wt,dead=False,ids=[t_id],rotateAroundZero=True,rotationAxis=(0,0,1),zeroPoint=[0,-x1,0]),
 RotationEngine(angularVelocity=wt,dead=False,ids=[t_id1],rotateAroundZero=True,rotationAxis=(0,0,1),zeroPoint=[0,-x1,0]),
 TranslationEngine(velocity=vx,dead=True,ids=[t_id2],translationAxis=(1,0,0)),
 TranslationEngine(velocity=-vx,dead=True,ids=[t_id3],translationAxis=(1,0,0)),
 TranslationEngine(velocity=-vz,dead=True,ids=[t_id2],translationAxis=(0,1,0)),
 TranslationEngine(velocity=-vz,dead=True,ids=[t_id3],translationAxis=(0,1,0)),
 RotationEngine(angularVelocity=-wt,dead=False,ids=[t_id4],rotateAroundZero=True,rotationAxis=(0,0,1),zeroPoint=[0,-x1,0]),
 RotationEngine(angularVelocity=-wt,dead=False,ids=nozzle,rotateAroundZero=True,rotationAxis=(1,0,0),zeroPoint=[0,y_rotation,0]),
 RotationEngine(angularVelocity=-wt,dead=False,ids=pack_rotation,rotateAroundZero=True,rotationAxis=(1,0,0),zeroPoint=[0,y_rotation,0]),
   #PyRunner(command='txt()',iterPeriod=periodrecord,nDo=ndorecord),
   #PyRunner(command='txt_3()',firstIterRun=start,iterPeriod=periodrecord),
   #PyRunner(command='flowfunction()',iterPeriod=periodflow,firstIterRun=endrecordtime,nDo=nDoflow),
   #PyRunner(command='flowprint()',iterPeriod=100,firstIterRun=startflowprint,nDo=1),
   PyRunner(command='motion()',iterPeriod=1,firstIterRun=iterrotation,nDo=2),
   PyRunner(command='motion_1()',iterPeriod=1,firstIterRun=iterrotation,nDo=2),
   PyRunner(command='motion_2()',iterPeriod=1,firstIterRun=stoptranslation_x,nDo=2),
   PyRunner(command='motion_3()',iterPeriod=1,firstIterRun=stoptranslation_x,nDo=2),
   PyRunner(command='motion_4()',iterPeriod=1,firstIterRun=stoptranslation_z,nDo=2),
   PyRunner(command='motion_5()',iterPeriod=1,firstIterRun=iterrotation_2,nDo=2),
]

O.saveTmp()

The script is not Working for me (and a tough sell on "Minimal" too :-) )

I get: NameError: name 'angle' is not defined

Sifting through your script, I will take a guess at your problem:

RotationEngine allows you to add an angularVelocity to a body. The engine is not responsible for the dynamics of the system. Turning the engine to dead means the RotationEngine is no longer setting an angularVelocity to the body. But the body still exists and still has an angular velocity.

I guess you are looking for something more along the lines of:

def motion_5():
 O.engines[11].angularVelocity=0
 O.engines[11].angularVelocity=0

If you want to set the angularVelocity to 0 for a set of bodies.

Cheers,

Robert

https://www.yade-dem.org/wiki/Howtoask

Thank you Robert!

Your idea looks to be good.
I used this option for the translation of other bodies and it actually works.

I'll do what you suggested me.

Cheers,

Francesco