corrupted double-linked list

Hi all,

I am running my Yade model. when the Yade model is small (8000 particles), the code works well. But when the same model has 1000000 particles, I get the following error:

corrupted double-linked list
Aborted (core dumped)

Can anyone help me understand why I get this error? The model compiles well, but when the first iteration is solved, the error appears.

Thank you for your help.

MY CODE:

#MATERIAL PROPERTIES
intR = 0.20 #connectivity parameter
DENS = 2500.00 #Density (kg/mm^3)
YOUNG = 1800.00 #Young's modulus (N/mm^2)
FRICT = 18.00 #Friction angle (degrees)
ALPHA = 0.20 #Poisson's ratio
TENS = 100.00 #Tensile strength (N/mm^2)
COH = 1000.00 #Shear strength (N/mm^2)

mn,mx=Vector3(-0.500000,-0.500000,-0.500000),Vector3(100.500000,86.525190,86.236515) # corners of the initial packing

#CREATING NODES
with open('Nodes.txt') as f:
  Nodes = [[float(x) for x in line.split()] for line in f]

O.materials.append(JCFpmMat(type=1,density=DENS,young=YOUNG,poisson=ALPHA,frictionAngle=radians(FRICT),tensileStrength=TENS,cohesion=COH,label='mat1'))

walls=aabbWalls([mn,mx],thickness=0,material='mat1')
wallIds=O.bodies.append(walls)

# PARTICLES
for i in range(1000000):
  if Nodes[i][3] == 1:
    O.bodies.append([ sphere(center=(Nodes[i][0],Nodes[i][1],Nodes[i][2]),radius=0.50,color = (1,0,0),material='mat1'),])
  else:
    O.bodies.append([ sphere(center=(Nodes[i][0],Nodes[i][1],Nodes[i][2]),radius=0.50,color = (0,1,0),material='mat1'),])

nbSpheres = 1000000.00

# PRINTING
from yade import plot
from pprint import pprint

# SIMULATION LOOP
O.engines=[
    ForceResetter(),
    InsertionSortCollider([Bo1_Sphere_Aabb()]),
    InteractionLoop(
        [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=intR)],
        [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1)],
        [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM()]
    ),
    GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=10,timestepSafetyCoefficient=0.5, defaultDt=utils.PWaveTimeStep()),
    DFNFlowEngine(dead=1,label="flow"),
    NewtonIntegrator(damping= 0.80),
    PyRunner(command='BreakingBonds()',iterPeriod=1),
    PyRunner(command='Writting()',iterPeriod= 500),
]

#WRITTING RESULTS
number = 0
def Writting():
  global number
  number = number + 1
  f= open("Info%d.txt" %(number),"w+")
  f.write("%5.10f \n" % (O.time))
  f.write("%5.10f \n" % (flow.getPorePressure((50.000,43.013,42.868))))
  for i in range (20000):
    f.write("%5.10f \n" % (SaveDataNodes[i][0]))
  f.close()

def history():
  plot.addData(t=O.time,p=flow.getPorePressure((50.000,43.013,42.868)))

#Apply Qin (volume flow rate) in the center of the model.
flow.dead = 0
flow.defTolerance=0.3
flow.meshUpdateInterval=1
flow.useSolver=3
flow.permeabilityFactor=1
flow.viscosity=10
flow.bndCondIsPressure=[1,1,1,1,1,1]
flow.bndCondValue=[0,0,0,0,0,0]
flow.boundaryUseMaxMin=[0,0,0,0,0,0]
O.dt= 0.0010000000
O.dynDt=False

ti = O.time
O.engines=O.engines+[PyRunner(iterPeriod=2,command='history()')]
O.engines=O.engines+[PyRunner(iterPeriod=1,command='Pumping()')]
def Pumping():
  if O.time<(ti+0.2):
    flow.imposeFlux(Vector3(50.000,43.013,42.868),-1.0000000000)
  else:
    flow.imposeFlux(Vector3(50.000,43.013,42.868),0)

def BreakingBonds():
  global Bonds0
  CurrentBonds = 0
  for ii in O.interactions:
    CurrentBonds+=1
  if Bonds0 > CurrentBonds:
    Bonds0 = CurrentBonds
    print("BROKEN!")
    print(O.time)
    Writting()

#PLOTTING
plot.plots={'t':('p')}
plot.plot()

#CREATING CONNECTIVITY
with open('Connectivity.txt') as f:
  Connectivity = [[int(x) for x in line.split()] for line in f]

#CREATING INTERACTIONS
for i in range(5905449):
  if Connectivity[i][2] == 1:
    createInteraction(Connectivity[i][0],Connectivity[i][1])

Connectivity.clear()

#CREATING SAVEDATA
with open('SaveDataNodes.txt') as f:
  SaveDataNodes = [[int(x) for x in line.split()] for line in f]

#BOUND COUNT
Bonds0=0
for i in O.interactions:
  Bonds0+=1
print (Bonds0)