shape mismatch

Hi,

I'm conducting a triaxial test to verify the mechanic properties of my sample. However, after running for 20 minute (iter=3e6), some error occurred.

The simplified script for the triaxial test is shown below:

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

# periodic triaxial test simulation

from yade import pack, qt, plot
sigmaIso=-5e4

O.periodic=True

sp=pack.randomPeriPack(radius=0.2e-3, initSize=(3e-3,3e-3,3e-3))
sp.toSimulation()

idSand=O.materials.append(CohFrictMat(alphaKr=1.8,etaRoll=1.8,density=2600,young=1e8,
                                      isCohesive=0,momentRotationLaw=1,
                                      poisson=.2,frictionAngle=0.2679492,label="sand"))

for b in O.bodies:
    b.material=O.materials[1]

#for c in O.

O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom6D()],
      [Ip2_CohFrictMat_CohFrictMat_CohFrictPhys()],
      [Law2_ScGeom6D_CohFrictPhys_CohesionMoment()]
   ),
   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='compressionFinished()'
   ),
   NewtonIntegrator(damping=.2),
# PyRunner(command='checkStress()',iterPeriod=100,label='checker'),
   PyRunner(command='addPlotData()',iterPeriod=100),
]
O.dt=.5*PWaveTimeStep()

def addPlotData():
   plot.addData(unbalanced=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],
      n=yade.utils.porosity()
   )

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

def compressionFinished():
    print
    'Compression finished'
    # 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()

---------------------------------------------------------------------------

---------------------------------------------------------------------------
In [1]:
Unhandled exception in thread started by <function liveUpdate at 0x7f69ce36e140>

ValueError Traceback (most recent call last)
/usr/lib/x86_64-linux-gnu/yade/py/yade/plot.py in liveUpdate(timestamp)
    506 for ax in axes:
    507 try:
--> 508 ax.relim() # recompute axes limits
    509 ax.autoscale_view()
    510 except RuntimeError: pass # happens if data are being updated and have not the same dimension at the very moment

/usr/lib/python2.7/dist-packages/matplotlib/axes/_base.pyc in relim(self, visible_only)
   1936 for line in self.lines:
   1937 if not visible_only or line.get_visible():
-> 1938 self._update_line_limits(line)
   1939
   1940 for p in self.patches:

/usr/lib/python2.7/dist-packages/matplotlib/axes/_base.pyc in _update_line_limits(self, line)
   1799 Figures out the data limit of the given line, updating self.dataLim.
   1800 """
-> 1801 path = line.get_path()
   1802 if path.vertices.size == 0:
   1803 return

/usr/lib/python2.7/dist-packages/matplotlib/lines.pyc in get_path(self)
    955 """
    956 if self._invalidy or self._invalidx:
--> 957 self.recache()
    958 return self._path
    959

/usr/lib/python2.7/dist-packages/matplotlib/lines.pyc in recache(self, always)
    665 y = self._y
    666
--> 667 self._xy = np.column_stack(np.broadcast_arrays(x, y)).astype(float)
    668 self._x, self._y = self._xy.T # views
    669

/usr/lib/python2.7/dist-packages/numpy/lib/stride_tricks.pyc in broadcast_arrays(*args, **kwargs)
    247 args = [np.array(_m, copy=False, subok=subok) for _m in args]
    248
--> 249 shape = _broadcast_shape(*args)
    250
    251 if all(array.shape == shape for array in args):

/usr/lib/python2.7/dist-packages/numpy/lib/stride_tricks.pyc in _broadcast_shape(*args)
    182 # use the old-iterator because np.nditer does not handle size 0 arrays
    183 # consistently
--> 184 b = np.broadcast(*args[:32])
    185 # unfortunately, it cannot handle 32 or more arguments directly
    186 for pos in range(32, len(args), 31):

ValueError: shape mismatch: objects cannot be broadcast to a single shape
---------------------------------------------------------------------------

The plot nolonger refresh to show the real time stress and strain state of the test.
I think that there may be some limits in the source code when the calculation amount go large.
Actually, I tried the FrictMat and the error didn't occur. So I guess it has something to do with the calculation burden?
Can anyone give me some suggestion? Thanks a lot.

Zhang Wenyue