does iterperiodaffect the results of simulation?

To answer your question as ambiguously as you posed it:

Maybe.

(consider reading out posting guidelines to reduce the collective ambiguity) [1]

Cheers,

Robert

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

Hello,

please provide a MWE [1], without it it is just guessing.

> The simulation runs perfectly for the first. But after that, it shows the following error

what does "for the first" and "after that" mean?

> So I want to ask whether iterperiod affects the result of simulation?

in general, of course iterPeriod does affect results. If yes and how much depends on the code (please provide a MWE [1])

cheers
Jan

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

Hi Jan and Robert,

Sorry for not providing MWE before. Below is my MWE. I hope the information is enough. I am simulating a single grain compression experiment. Then get the load-displacement image. When using strain rate 1e-2 and iterperiod 1000, there was no problem with simulation. But when I use strain rate 1e-3 and iterperiod 1000 appeared shape mismatch: objects cannot be broadcast to a single shape. According to the previous question, I increased iterperiod to 10000. But the load-displacement image is not the same as iterperiod 1000. I think iterPeriod determines how often to save simulation data. It should not affect my load-displacement image. So I am very confused.

"""
A simple grain test
A sphere with radius 2.5e-3
Load in z direction by z-perpendicslar walls.

"""

    iterper=1000,

# engines
O.engines=[
                       ForceResetter(),
                       InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=intRadius,label='is2aabb'),Bo1_Wall_Aabb()]),
                       InteractionLoop(
                            [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=intRadius,label='ss2sc'),Ig2_Wall_Sphere_ScGeom()],
                             [Ip2_FrictMat_FrictMat_FrictPhys(),Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(
                             cohesiveTresholdIteration=10, label='jcf', xSectionWeibullShapeParameter=xSectionShape,weibullCutOffMin=weibullCutOffMin, weibullCutOffMax=weibullCutOffMax
                                            )],
                             [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(smoothJoint=True,Key=identifier, label='interactionLaw', neverErase=True,recordCracks=True,recordMoments=True,momentRadiusFactor=momentRadiusFactor),
                             Law2_ScGeom_FrictPhys_CundallStrack()],
    ),

           GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.4, defaultDt=0.1*utils.PWaveTimeStep()),
           VTKRecorder(dead=0,iterPeriod=iterper*2,initRun=True,fileName=(output+'-'),recorders=['jcfpm','cracks','facets','moments','intr'],Key=identifier,label='vtk'),
                       NewtonIntegrator(damping=0.7),
                       PyRunner(iterPeriod=100,command='addPlotData()',initRun=True),
                       PyRunner(iterPeriod=100,command='stopIfDamaged()'),
]

# plot stuff
def addPlotData():
 # forces of walls. f1 is "down", f2 is "up" (f1 needs to be negated for evlauation)
 f1,f2 = [O.forces.f(i)[2] for i in wallIDs]
 f1 *= -1
 # average force
 f = .5*(f1+f2)
 # displacement
 wall = O.bodies[wallIDs[1]]
 dspl = -1* wall.state.displ()[2]

 # store values
 plot.addData(
  t = O.time,
  i = O.iter,
  dspl = dspl,
  f1 = f1,
  f2 = f2,
  f = f,
 )
 plot.saveDataTxt('txt/data'+identifier+'.txt',vars= ('i', 'dspl','f'))

# plot dspl on x axis, loadon y1 axis and f,f1,f2 in y2 axis
plot.plots={'dspl':('f',)}

O.dt = 0
O.step() # to create initial contacts
# now reset the interaction radius and go ahead

# time step
O.dt = 0.004 * utils.PWaveTimeStep()

# run simulation
plot.plot()
O.run()

regards,
Li

Thanks for a ME, but please make it MWE, W=working, here
- iterper=1000, is syntax error
- output, intRadius, xSectionShape variables is not defined

cheers
Jan

Hi Jan,

below is new one.

from __future__ import print_function
from yade import plot,pack

"""
A simple grain test
A sphere with radius 2.5e-3
Load in z direction by z-perpendicslar walls.

"""

# load parameters from file if run in batch
# default values are used if not run from batch
readParamsFromTable(noTableOk=True, # unknownOk=True,
 weibullCutOffMax=10,
 weibullCutOffMin=0.1,
 xSectionShape = 2,
 xSectionScale = 0,

    momentRadiusFactor=4,

 young = 4e10,
    cohesion = 3.5e7,
    poisson = .4,
    density = 2706,
    frictionAngle=0.7,

    iterper=1000,
    intRadius =1.329,

 velocity = 1e-2,
 specimenRadius = .024,
 sphereRadius = 1.65e-3,
    rRelFuzz=0.33333,
)
from yade.params.table import *

identifier='shape-velocity '+str(velocity)+'-weibullshape'+str(xSectionShape)+'-momentRadiusFactor'+str(momentRadiusFactor)
outputDir='out_'+identifier
if not os.path.exists(outputDir):
 os.mkdir(outputDir)
if not os.path.exists('txt'):
 os.mkdir('txt')

output = './'+outputDir+'/'+identifier

# material
JCFmat = O.materials.append(JCFpmMat(young=young, cohesion = cohesion, density=density, frictionAngle = frictionAngle, poisson=poisson,
                 jointTensileStrength=3.5e7, jointNormalStiffness = 4e10, jointShearStiffness = 1.6e10, jointFrictionAngle= radians(30), label='JCFpmMat'))

# spheres
sp=pack.randomDensePack(
 pack.inSphere((0,0,0),specimenRadius),
 radius = sphereRadius, rRelFuzz = rRelFuzz, returnSpherePack = True
)
sp.toSimulation()

# walls
wallMat = O.materials.append(FrictMat(young=100e10,poisson=.05,frictionAngle=.25,density=70000, label='frictionlessWalls'))
zMin,zMax=[pt[2] for pt in aabbExtrema()]
wallIDs = O.bodies.append([wall((0,0,z),2) for z in (zMin,zMax)])
walls = wallMin,wallMax = [O.bodies[i] for i in wallIDs]
v = velocity
wallMin.state.vel = (0,0,0)
wallMax.state.vel = (0,0,-v)

# engines
O.engines=[
                       ForceResetter(),
                       InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=intRadius,label='is2aabb'),Bo1_Wall_Aabb()]),
                       InteractionLoop(
                            [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=intRadius,label='ss2sc'),Ig2_Wall_Sphere_ScGeom()],
                             [Ip2_FrictMat_FrictMat_FrictPhys(),Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(
                             cohesiveTresholdIteration=10, label='jcf', xSectionWeibullShapeParameter=xSectionShape,weibullCutOffMin=weibullCutOffMin, weibullCutOffMax=weibullCutOffMax
                                            )],
                             [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(smoothJoint=True,Key=identifier, label='interactionLaw', neverErase=True,recordCracks=True,recordMoments=True,momentRadiusFactor=momentRadiusFactor),
                             Law2_ScGeom_FrictPhys_CundallStrack()],
    ),

           GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.4, defaultDt=0.1*utils.PWaveTimeStep()),
           VTKRecorder(dead=0,iterPeriod=iterper*2,initRun=True,fileName=(output+'-'),recorders=['jcfpm','cracks','facets','moments','intr'],Key=identifier,label='vtk'),
                       NewtonIntegrator(damping=0.7),
                       PyRunner(iterPeriod=100,command='addPlotData()',initRun=True),
                       PyRunner(iterPeriod=100,command='stopIfDamaged()'),
]
# stop condition
def stopIfDamaged():
 if O.iter < 5000: # do nothing at the beginning
  return
 fMax = max(plot.data["f"])
 f = plot.data["f"][-1]
 if f/fMax < 0.6:
  print("Damaged, stopping.")
  print("ft = ",max(plot.data["f"]))
  O.pause()

# plot stuff
def addPlotData():
 # forces of walls. f1 is "down", f2 is "up" (f1 needs to be negated for evlauation)
 f1,f2 = [O.forces.f(i)[2] for i in wallIDs]
 f1 *= -1
 # average force
 f = .5*(f1+f2)
 # displacement
 wall = O.bodies[wallIDs[1]]
 dspl = -1* wall.state.displ()[2]

 # store values
 plot.addData(
  t = O.time,
  i = O.iter,
  dspl = dspl,
  f1 = f1,
  f2 = f2,
  f = f,
 )
 plot.saveDataTxt('txt/data'+identifier+'.txt',vars= ('i', 'dspl','f'))

# plot dspl on x axis, loadon y1 axis and f,f1,f2 in y2 axis
plot.plots={'dspl':('f',)}

O.dt = 0
O.step() # to create initial contacts
# now reset the interaction radius and go ahead

# time step
O.dt = 0.004 * utils.PWaveTimeStep()

# run simulation
plot.plot()
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)
waitIfBatch()

regards,
Li

Thanks. One more detail:
what should we change to get the error and at what time/iteration the error should occur?
cheers
Jan

Hi Jan,

when i change velocity = 1e-2 to 1e-3, I will get the error. This simulation in my computer take long time to run. About at what time/iteration I can't tell the time precisely. But in my computer after like 1 or 2 hours the error occur.

regards,
Li

Please copy and paste the full error from your terminal. Typing your interpretation of the error does not help us.

If it is the error I think it is (no way to know until you copy and paste it) I think this is a matplotlib library bug that we don't have control over.

Hi Robert,

below is the full error. I think it is a matplotlib library bug. But in the previous problem, we changed the iterperiod problem and solved it. But my problem is that when I change the iterperiod, my load-displacement image also changes accordingly. However, iterperiod should not change the load-displacement image? （my opinion）

Traceback (most recent call last):
File "/usr/lib/x86_64-linux-gnu/yade/py/yade/plot.py", line 514, in liveUpdate
ax.relimQ # recompute axes limits
File "/usr/lib/python3/dist-packages/matplotlib/axes/_base.py", line 2051, in
relim
self._update_line_limits(line)
File "/usr/lib/python3/dist-packages/matplotlib/axes/_base.py", line 1924, in
_update_line_limits
path = line.get_path()
File "/usr/lib/python3/dist-packages/matplotlib/lines.py", line 1027, in get_path
self.recache()
File "/usr/lib/python3/dist-packages/matplotlib/lines.py", line 679, in recache
self._xy = np.column_stack(np.broadcast_arrays(x, y)).astype(float)
File "<___ array_function___ internals?-", line 5, in broadcast_arrays
File "/usr/lib/python3/dist-packages/numpy/lib/stride_tricks.py", line 264, in
broadcast_arrays
shape = _broadcast_shape(*args)
File "/usr/lib/python3/dist-packages/numpy/lib/stride_tricks.py", line 191, in
_broadcast_shape
b = np.broadcast(*args[:32])
ValueError: shape mismatch: objects cannot be broadcast to a single shape

regards,
Li

It is a known bug [1], we cannot do much about it from what I can see.

> But in the previous problem, we changed the iterperiod problem and solved it

? what does that mean? Are you referring to a different question on launchpad?

> But my problem is that when I change the iterperiod, my load-displacement image also changes accordingly.

The number points shown in yade's plot depends on the iterperiod. So you are just plotting fewer points spaced further apart - the actual data points are not moving. You can collect data as often as you like using saveDataTxt (same way you are doing it now), then create your plot later using the saved txt file.

Cheers,

Robert

[1]https://gitlab.com/yade-dev/trunk/-/issues/110

> ? what does that mean? Are you referring to a different question on launchpad?

yes, this question [1].
And in [2] means this problem appear "randomly". So should I try more times？ Maybe it does not happen some times.

> The number points shown in yade's plot depends on the iterperiod. So you are just plotting fewer points spaced further apart - the actual data points are not moving. You can collect data as often as you like using saveDataTxt (same way you are doing it now), then create your plot later using the saved txt file.

I have the same view. I ran the program again. Still waiting for the result.

[1] https://answers.launchpad.net/yade/+question/670492
[2]https://gitlab.com/yade-dev/trunk/-/issues/110

This question was expired because it remained in the 'Open' state without activity for the last 15 days.