keeping constant confinement

Hello all,

I'm trying to create a cylinder full of spheres with a specific porosity and save the data in a text file to be used in another simulation. To do that:
1- I create a cloud of spheres in a box
2- apply triaxial stress until the desired porosity is reached
3- extract a cylinder by crop/filter the packed spheres

The problem is that when I extract the cylinder, the porosity is set back to the original cloud porosity, not the porosity that was reached after applying the triaxial stress. So I'm guessing that I need to apply constant confinement stress on the cylinder to keep it packed to the desired porosity. Can anyone help with that? Below is my code. (note: I'm new in using Yade and still learning how it works)

Thanks

############################ Material properties #############################

sigmaIso=-1e3
frictangel=.52
targetp=0.65 #this is the targeted porosity
O.materials.append(FrictMat(young = 5e10, poisson = 0.15,frictionAngle = radians(frictangel), density=1920))

############################ generate loose packing #############################

from yade import pack, qt, plot, export
sp=pack.SpherePack()
sp.makeCloud((0,0,0),(.15,.15,.3),rMean=.007,rRelFuzz=.2,periodic=True)

sp.toSimulation()
#yade.qt.View()

############################ Triaxial compaction #############################

O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom()],
      [Ip2_FrictMat_FrictMat_FrictPhys()],
      [Law2_ScGeom_FrictPhys_CundallStrack()]
   ),
   PeriTriaxController(label='triax',
      goal=(sigmaIso,sigmaIso,sigmaIso),stressMask=7,
      dynCell=True,maxStrainRate=(100,100,100),
      maxUnbalanced=.1,relStressTol=1e-3,
      doneHook='compactionFinished()'
   ),
   NewtonIntegrator(damping=.2),
   PyRunner(command='P()',iterPeriod=500)
]
O.dt=.5*PWaveTimeStep()
## stop conditions
O.run()
def P():

 if utils.porosity()>targetp:
  print ('porosity = ', utils.porosity())
 if utils.porosity()<targetp:
  O.pause()
  print 'Finished'
  #### cylinder extraction
  pred=pack.inCylinder((.05,.05,.05),(.05,.05,.2503),.0508)
  sp=SpherePack()
  sp.fromSimulation()
  O.switchScene()
  spFilter=filterSpherePack(pred,sp,Material=Material, returnSpherePack=True)
  spFilter.toSimulation()
  #### export data
  export.textExt('myspheres','x_y_z_r')
  print ('porosity = ', utils.porosity()) #check the porosity of the extracted cylinder

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,-.2)
   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()