Periodic triaxial test example in chapter 2 does not work

Asked by liucheng83

Hello,
Periodic triaxial test examplein chapter 2 (Tutorial) does not work.

------------when I use yade 0.70 it gives the message as fllow:
yade.plot: creating new line for kinetic
yade.plot: creating new line for elastPotential
yade.plot: creating new line for nonviscDamp
yade.plot: creating new line for velGradWork
yade.plot: creating new line for plastDissip
FATAL /build/buildd/yade-0.70.0/core/ThreadRunner.cpp:31 run: Exception occured:
Body #270 has velocity==NaN!

------------when I use yade 0.80 it gives the message as fllow:(Yade 0.80.0+1+9~precise1, from 2012-05-04
)

Yade [1]: FATAL /build/buildd/yade-stable-0.80.0+1+9~precise1/core/ThreadRunner.cpp:31 run: Exception occured:
PeriTriaxController run on aperiodic simulation.
yade.plot: creating new line for kinetic
yade.plot: creating new line for nonviscDamp

the scripts is in 2.6.6 Periodic triaxial test in yade document , for easy, list as fllowing:

----------------------------------begin
# encoding: utf-8

# periodic triaxial test simulation
#
# The initial packing is either
#
# 1. random cloud with uniform distribution, or
# 2. cloud with specified granulometry (radii and percentages), or
# 3. cloud of clumps, i.e. rigid aggregates of several particles
#
# The triaxial consists of 2 stages:
#
# 1. isotropic compaction, until sigmaIso is reached in all directions;
# this stage is ended by calling compactionFinished()
# 2. constant-strain deformation along the z-axis, while maintaining
# constant stress (sigmaIso) laterally; this stage is ended by calling
# triaxFinished()
#
# Controlling of strain and stresses is performed via PeriTriaxController,
# of which parameters determine type of control and also stability
# condition (maxUnbalanced) so that the packing is considered stabilized
# and the stage is done.
#

sigmaIso=-1e5

#import matplotlib
#matplotlib.use('Agg')

# generate loose packing
from yade import pack, qt, plot
sp=pack.SpherePack()
if 0:
   ## uniform distribution
   sp.makeCloud((0,0,0),(2,2,2),rMean=.1,rRelFuzz=.3,periodic=True)
elif 0:
   ## per-fraction distribution
   ## passing: cummulative percentage
   sp.particleSD2(radii=[.09,.1,.2],passing=[40,80,100],periodic=True,numSph=1000)
else:
   ## create packing from clumps
   # configuration of one clump
   c1=pack.SpherePack([((0,0,0),.1),((.15,0,0),.05),((0,.1,0),.05)])
   # make cloud using the configuration c1 (there could c2, c3, ...; selection between them would be random)
   sp.makeClumpCloud((0,0,0),(2,2,2),[c1],periodic=True)

# setup periodic boundary, insert the packing
sp.toSimulation()

O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_L3Geom()],
      [Ip2_FrictMat_FrictMat_FrictPhys()],
      [Law2_L3Geom_FrictPhys_ElPerfPl()]
   ),
   NewtonIntegrator(damping=.6),
   PeriTriaxController(label='triax',
      # specify target values and whether they are strains or stresses
      goal=(sigmaIso,sigmaIso,sigmaIso),stressMask=7,
      # type of servo-control
      dynCell=True,maxStrainRate=(.1,.1,.1),
      # wait until the unbalanced force goes below this value
      maxUnbalanced=.1,relStressTol=1e-3,
      # call this function when goal is reached and the packing is stable
      doneHook='compactionFinished()'
   ),
   PyRunner(command='addPlotData()',iterPeriod=100),
]
O.dt=.5*utils.PWaveTimeStep()

def addPlotData():
   plot.addData(unbalanced=utils.unbalancedForce(),i=O.iter,
      sxx=triax.stress[0],syy=triax.stress[1],szz=triax.stress[2],
      exx=triax.strain[0],eyy=triax.strain[1],ezz=triax.strain[2],
      # save all available energy data
      Etot=O.energy.total(),**O.energy
   )

# enable energy tracking in the code
O.trackEnergy=True

# define what to plot
plot.plots={'i':('unbalanced',),'i ':('sxx','syy','szz'),' i':('exx','eyy','ezz'),
   # energy plot
   ' i ':(O.energy.keys,None,'Etot'),
}
# show the plot
plot.plot()

def compactionFinished():
   # set the current cell configuration to be the reference one
   O.cell.trsf=Matrix3.Identity
   # change control type: keep constant confinement in x,y, 20% compression in z
   triax.goal=(sigmaIso,sigmaIso,-.3)
   triax.stressMask=3
   # allow faster deformation along x,y to better maintain stresses
   triax.maxStrainRate=(1.,1.,.1)
   # next time, call triaxFinished instead of compactionFinished
   triax.doneHook='triaxFinished()'
   # do not wait for stabilization before calling triaxFinished
   triax.maxUnbalanced=10

def triaxFinished():
   print 'Finished'
   O.pause()

-------------------------------------------------------end

Question information

Language:
English Edit question
Status:
Answered
For:
Yade Edit question
Assignee:
No assignee Edit question
Last query:
Last reply:
Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#1

Hi,
Thanks for reporting.
I can confirm the problem (found that yesterday while showing yade to
somebody).
You can find a relatively similar script in examples, called
periodic-triax.py. The main difference is that it does not plot curves.

So, if nobody can update the examples of chapter 2 (which were written
for an older version) soon, I would suggest you to check periodic-triax
in the meantime.
You may even find out what's wrong in the example of chapter 2 and let
us know.

Bruno

Revision history for this message
liucheng83 (lcheng83) said :
#2

Hello, Bruno
I try the fllowing scripts, but found the result of O.periodic is False, can anybody tell me why? Thank you !
---------------
sp=pack.SpherePack()
c1=pack.SpherePack([((0,0,0),.1),((.15,0,0),.05),((0,.1,0),.05)])
# make cloud using the configuration c1 (there could c2, c3, ...; selection between them would be random)
sp.makeClumpCloud((0,0,0),(2,2,2),[c1],periodic=True)
sp.toSimulation()
O.periodic
---------------
I think it maybe something wrong with the function makeClumpCloud(...,,periodic=True)

Revision history for this message
liucheng83 (lcheng83) said :
#3

Sorry , another scripts is fllowing, the result is true.
------------------------------
sp=pack.SpherePack()
c1=pack.SpherePack([((0,0,0),.1),((.15,0,0),.05),((0,.1,0),.05)])
sp.makeCloud((0,0,0),(2,2,2),rMean=.1,rRelFuzz=.3,periodic=True)
sp.toSimulation()
O.periodic
--------------------------------

Revision history for this message
liucheng83 (lcheng83) said :
#4

Is there another bug? when I run the Periodic triaxial test examplein chapter 2 (Tutorial) with
-----

 InteractionLoop(
  #[Ig2_Sphere_Sphere_L3Geom()],
  #[Ip2_FrictMat_FrictMat_FrictPhys()],
  #[Law2_L3Geom_FrictPhys_ElPerfPl()]
  [Ig2_Sphere_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_FrictPhys()],
  [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
-------------
the stress will always greater than 0, so the "goal=(-1e5,-1e5,-1e5),stressMask=7" can never be reached.
and I modified it with
----

 InteractionLoop(
  #[Ig2_Sphere_Sphere_L3Geom()],
  #[Ip2_FrictMat_FrictMat_FrictPhys()],
  #[Law2_L3Geom_FrictPhys_ElPerfPl()]
  [Ig2_Sphere_Sphere_ScGeom()],
  [Ip2_FrictMat_FrictMat_FrictPhys()],
  [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
-------------------
It give the stress little than 0, and It can be finished at last, Is there any bug for the "Law2_L3Geom_FrictPhys_ElPerfPl()"?
Anybody know can give me some suggestion, thank you !

Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#5

>Is there any bug for the "Law2_L3Geom_FrictPhys_ElPerfPl()"?

Yes.
Better don't use this law, it will be removed soon.
If you can make the script work properly, then add your name at the top of the script and send it to us please. It would replace the current example, and you would become author a yade's demo script. ;)

Can you help with this problem?

Provide an answer of your own, or ask liucheng83 for more information if necessary.

To post a message you must log in.