yade-users team mailing list archive

[Xin Wu <question210746@xxxxxxxxxxxxxxxxxxxxx>]

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

Hallo All,
I ahve the Problem when i use the  [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
      [Law2_ScGeom_ViscElPhys_Basic()], come the Warning. And the simulationm breaks immediately.

Running script Rahmenscher82a.py
WARN  /home/wuxin/yade-trunk-2954021/pkg/dem/Shop.cpp:429 PWaveTimeStep: PWaveTimeStep has not found any suitable spherical body to calculate dt. dt is set to 1.0
ERROR /home/wuxin/yade-trunk-2954021/pkg/common/GravityEngines.cpp:19 action: GravityEngine is deprecated, consider using Newton::gravity instead (unless gravitational energy has to be tracked - not implemented with the newton attribute).
[[ ^L clears screen, ^U kills line. F12 controller, F11 3d view, F10 both, F9 generator, F8 plot. ]]

The Scripe ist :

# gravity deposition, continuing with oedometric test after stabilization
# shows also how to run parametric studies with yade-batch
############################################################################################
from yade import utils
from numpy import linspace
from numpy import arange
import gts
import itertools
from yade import pack
###########################################################################################
# The components of the batch are:
# 1. table with parameters, one set of parameters per line (ccc.table)
# 2. utils.readParamsFromTable which reads respective line from the parameter file
# 3. the simulation muse be run using yade-batch, not yade
#
# $ yade-batch --job-threads=1 03-oedometric-test.table 03-oedometric-test.py
#
rMean=.0096
rRelFuzz=.0016
maxLoad=4500
#minLoad=500
frictionAngleSt=radians(35)
frictionAngleBo=radians(42.97183463)
tc = 0.001
en = 0.3
es = 0.3
a=0.05
#id_SphereMat=O.materials.append(ViscElMat(kn=1e6,ks=1e6,cn=0.0,cs=0.0,density=1000,frictionAngle=angle))
#id_WallMat=O.materials.append(ViscElMat(kn=1e8,ks=1e8,cn=0.0,cs=0.0,density=1000,frictionAngle=angle))
id_SphereMat=O.materials.append(ViscElMat(kn=1e7,ks=1e7,cn=0.0,cs=0.0,density=2650,frictionAngle=frictionAngleBo))
id_WallMat=O.materials.append(ViscElMat(kn=1e8,ks=1e8,cn=0.0,cs=0.0,density=8000,frictionAngle=frictionAngleSt))


steel=O.materials[id_WallMat]
spheremat=O.materials[id_SphereMat]
# load parameters from file if run in batch
# default values are used if not run from batch
#utils.readParamsFromTable(rMean=.0096,rRelFuzz=.016,maxLoad=1e6,minLoad=1e4)
# make rMean, rRelFuzz, maxLoad accessible directly as variables later
#from yade.params.table import *

# create box with free top, and ceate loose packing inside the box
from yade import pack, plot,qt
fIDSI=O.bodies.append(utils.geom.facetBox((.15,.15,.135),(.15,.15,.045),wallMask=15,material=steel))
fIDSII=O.bodies.append(utils.geom.facetBox((.15,.15,.045),(.15,.15,.045),wallMask=31,material=steel))
sp=pack.SpherePack()
sp.makeCloud((0,0,0),(0.3,0.3,0.3250),rMean=rMean,rRelFuzz=rRelFuzz)
sp.toSimulation(material=spheremat)

O.engines=[
   ForceResetter(),
   # sphere, facet, wall
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
   InteractionLoop(
      # the loading plate is a wall, we need to handle sphere+sphere, sphere+facet, sphere+wall
      [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
      [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
      [Law2_ScGeom_ViscElPhys_Basic()]
   ),
   GravityEngine(gravity=(0,0,-9.81)),
   NewtonIntegrator(damping=0.7),
   # the label creates an automatic variable referring to this engine
   # we use it below to change its attributes from the functions called
   qt.SnapshotEngine(iterPeriod=60000,fileBase='/tmp/bulldozer-',label='snapshooter'),
   PyRunner(command='checkUnbalanced()',realPeriod=1,label='checker'),
   TranslationEngine(translationAxis=[1,0,0],velocity=0,ids=fIDSI,label='Transl'),
]
O.dt=.25*utils.PWaveTimeStep()
#Transl.velocity=0.01
#qt.Controller()
#qt.View()

# the following checkUnbalanced, unloadPlate and stopUnloading functions are all called by the 'checker'
# (the last engine) one after another; this sequence defines progression of different stages of the
# simulation, as each of the functions, when the condition is satisfied, updates 'checker' to call
# the next function when it is run from within the simulation next time

def AngVel():
    for b in O.bodies:
        if isinstance(b.shape,Sphere):
           b.state.blockedDOFs='Y'
           b.state.angVel[1]=0
           #b.state.angVel[1]=0


# check whether the gravity deposition has already finished
# if so, add wall on the top of the packing and start the oedometric test

def checkUnbalanced():
    # at the very start, unbalanced force can be low as there is only few contacts, but it does not mean the packing is stable
    if O.iter<240000: return 
    # the rest will be run only if unbalanced is < .1 (stabilized packing)
    if utils.unbalancedForce()>.05: return 
    # add plate at the position on the top of the packing
    # the maximum finds the z-coordinate of the top of the topmost particle
    #Transl.velocity=0.01
    fIDSIII=O.bodies.append(utils.wall(max([b.state.pos[2]+b.shape.radius for b in O.bodies if isinstance(b.shape,Sphere)]),axis=2,sense=-1,material=spheremat))
    global plate        # without this line, the plate variable would only exist inside this function
    plate=O.bodies[fIDSIII]  # the last particles is the plate
    # Wall objects are "fixed" by default, i.e. not subject to forces
    # prescribing a velocity will therefore make it move at constant velocity (downwards)
    plate.state.vel=(0,0,-.025)
    # start plotting the data now, it was not interesting before
    O.engines=O.engines+[PyRunner(command='addPlotData()',iterPeriod=1000)]
    # next time, do not call this function anymore, but the next one (unloadPlate) instead
    checker.command='unloadPlate()'
  

def unloadPlate():
   # if the force on plate exceeds maximum load, start unloading
   if abs(O.forces.f(plate.id)[2])>=maxLoad:
      O.engines=O.engines+[PyRunner(command='ServoContorl()',iterPeriod=1)]
      #O.engines=O.engines+[PyRunner(command='AngVel()',iterPeriod=1)]
      # next time, do not call this function anymore, but the next one (stopUnloading) instead
      checker.command='stopUnloading()'
   
      

def stopUnloading():
    if abs(O.forces.f(plate.id)[2])>(maxLoad-1):
       if abs(O.forces.f(plate.id)[2])<(maxLoad+1):
          Transl.velocity=0.0010


def ServoContorl():
    Fzs=O.forces.f(plate.id)[2]
    if Fzs==0 :
       Fzs=0.000000000000001
    global G
    global KN
    KN=0.000000000000000001
    for i in plate.intrs() :
        KN=KN+i.phys.kn
    G=a/(KN*(O.dt))
    plate.state.vel[2]=(G*(Fzs-maxLoad))
    






def addPlotData():
    Fz=O.forces.f(plate.id)[2]
    F = 0
    global i
    for i in fIDSI:
        F += O.forces.f(i)[0]
    plot.addData(t=O.time,Fx=-F,PX=(O.bodies[i].state.pos[0]-O.bodies[i].state.refPos[0]),Fz=Fz,w=plate.state.pos[2]-plate.state.refPos[2],i=O.iter)
    if (O.bodies[i].state.pos[0]-O.bodies[i].state.refPos[0])>0.050:
       plot.saveDataTxt('KnKs1e65f075D07.txt',vars=('i','PX','Fx','w','Fz'))
       O.pause()

# besides unbalanced force evolution, also plot the displacement-force diagram
plot.plots={'i':('w','Fz',),'PX':('Fx',)}
plot.plot()

############################################################################################################################
#O.saveTmp()
qt.Controller()
qt.View()
r=qt.Renderer()
#r.lightPos=Vector3(0,0,50)

##############################################################################################################################

O.run()
# when running with yade-batch, the script must not finish until the simulation is done fully
# this command will wait for that (has no influence in the non-batch mode)
utils.waitIfBatch()

-- 
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.