# yade-users team mailing list archive

## [Question #700411]: Aborted (core dumped) occured when I use the Hertz-Minlin contact law with non-zero krot and eta

```New question #700411 on Yade:

Aborted (core dumped)
I am running a periodic triaxial test using hertz-mindlin contact law. However, the program abort suddenly as soon as I add the rolling stiffness coefficient and rolling plastic limit (as shown above). If these values are set to 0, the simulation can continue.
My code is as follows:
# encoding: utf-8

num=2500,
)

import glob, os
import numpy as np

isBatch=runningInBatch() #check if run in batch mode
if isBatch:
print('Running: '+O.tags['description']) #mark the current simulation

num = table.num # number of soil particles
dScaling = 1e3 # density scaling to increase crticial time step
sigmac= -100000 # confining pressure
rate = 0.1 # strain rate
damp = 0.5 # damping coefficient
stabilityRatio = 1.e-3 # threshold for quasi-static condition
stressTolRatio = 1.e-3 # tolerance for stress goal

# corners to define specimen size
m_1=0.004*pow(num/2500,1/3)
m_2=0.008*pow(num/2500,1/3)
mn,mx=Vector3.Zero,Vector3(m_1,m_1,m_2)

# Soil sphere parameters
E=7e9 # micro Young's modulus
v=0.3 # micro Poisson's ratio
kr=0.5 # rolling/bending stiffness coefficient
eta=0.5 # rolling/bending plastic limit
mu = 26 # contact friction angle
rho = 2650*dScaling # soil density

# create materials
spMat = O.materials.append(

# create a cloud of ramdomly positioned spheres
O.periodic = True
O.cell.hSize = Matrix3(mx[0],0,0, 0,mx[1],0, 0,0,mx[2])
sp=pack.SpherePack()
sizes=[.00015,.000225,.0003,.0003625,.00045,.0006,0.00118]
cumm=[0,.215,.3807312826,.601,.8865351132,.9965177017,1]
sp.makeCloud(mn,mx,psdSizes=sizes,psdCumm=cumm,distributeMass=True,seed=1,num=num)

spIds=sp.toSimulation(material=spMat)

#define engines
O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom()],
[Ip2_FrictMat_FrictMat_MindlinPhys(
krot=kr,
ktwist=kr,
eta=eta
)],
[Law2_ScGeom_MindlinPhys_Mindlin(includeMoment=True)]
),
PeriTriaxController(label='triax',
goal=(sigmac,sigmac,sigmac),
# whether they are strains or stresses
# type of servo-control
dynCell=True,
maxStrainRate=(10*rate,10*rate,10*rate),
# wait until the unbalanced force goes below this value
maxUnbalanced=stabilityRatio,
relStressTol=stressTolRatio,
doneHook='compactionFinished()'
),
NewtonIntegrator(damping=damp,label='newton'),
]

O.dt=.3*PWaveTimeStep()

def compactionFinished():
if unbalancedForce()<stabilityRatio:
print(porosity())
triax.doneHook='startShearing()'

def startShearing():
# set the current cell configuration to be the reference one
O.cell.trsf = Matrix3.Identity
setRefSe3()
# set loading type: constant pressure in x,y, 30% compression in z
triax.globUpdate=1
# allow faster deformation along x,y to better maintain stresses
triax.maxStrainRate=(10*rate,10*rate,rate)
triax.goal=(sigmac,sigmac,-.3)
triax.maxUnbalanced=stabilityRatio
triax.relStressTol=stressTolRatio
triax.doneHook='Finished()'

def Finished():
print(str(num)+'finished')
O.pause()

s22=-triax.stress[1],s33=-triax.stress[2],e11=triax.strain[0],\
e22=triax.strain[1],e33=triax.strain[2],n=porosity(),\
V=O.cell.volume,M=triax.mass,Unb=unbalancedForece())
if isBatch:
dataName = 'triax2_'+ O.tags['description'] + '.txt'
else:
dataName = 'TriaxREV2_'+str(num)+'.txt'
plot.saveDataTxt(dataName,vars=('s11','s22','s33','e11','e22','e33','n','V','M'))

O.run()
waitIfBatch()

--