# Polyhedra splitting with Mohr-Coulomb-Weibull criterion

Asked by Karol Brzezinski on 2020-12-29

I am trying to understand and use an implementation of a combined Mohr-Coulomb–Weibull (MCW) failure criterion [1].

I have created an example based on the "splitter.py" example prepared by François Kneib. My example consists of a cube compressed between two walls. I created a polyhedraMat, which is splittable and has defined all the properties of Mohr-Coulomb–Weibull failure criterion. I set all the strengths to 1, probability P=0.63, and high Weibull modulus m=1000 (just in case). I expect to get strength close to 1 Pa. Instead of that, the strength is around 1.6 e5. It is independent of any parameter dedicated for MCW criterion. On the other hand, it depends on "strength" material property (it is not included in the example, but the default value is 100).

The major issue is:
- Do I have to activate the MCW criterion somehow (or deactivate the previous one [2])?

And there are two minor issues that I do not understand:
- Do I have to set "Wei_S0". From what I understood from [1] it is some general value that could be used if the tension and compression strengths were the same. But if I declare both strengths (sigmaCD, sigmaCZ), I expect that Wei_S0 is calculated based on them. Am I right?
- Since the default "strength" is 100, why do I get 1.6e5 in my "experiment"? The scale effect should even decrease strength. In the experiment, the contact force between the wall and the cube is divided by crosssection area (1m2) to calculate stress.

Cheers,
Karol

[1] Gladkyy, A., & Kuna, M. (2017). DEM simulation of polyhedral particle cracking using a combined Mohr-Coulomb–Weibull failure criterion. Granular Matter, 19(3), 41.
[2] Eliáš, J. (2014). Simulation of railway ballast using crushable polyhedral particles. Powder Technology, 264, 458-465.
---------------------------------------

# FUNCTION DENINITIONS
force = O.forces.f(1)[1]# cube wall interaction
sig = force/(d**2)
plot.addData(t = O.time, force = force, sig = sig)

d = 1 # cube length

#cube material
gravel1 = PolyhedraMat()
gravel1.IsSplitable = True
gravel1.Wei_P = 0.63
gravel1.Wei_m = 1000
gravel1.Wei_S0 = 1
gravel1.Wei_V0 = d**3
gravel1.sigmaCD = 1
gravel1.sigmaCZ = 1

#wall material
steel = PolyhedraMat()
steel.young = 1e10
steel.frictionAngle = 0

# create cube
p1 = polyhedra_utils.polyhedra(gravel1, v= [(0,0,0),(0,d,0),(d,d,0),(d,0,0),(0,0,d),(0,d,d),(d,d,d),(d,0,d)])

w1 = utils.wall(+d*1.001, axis=1, sense=-1, material=steel)
w2 = utils.wall(-0.001*d, axis=1, sense=+1, material=steel)

v = 5e-4
w1.state.vel = (0,-v,0)
O.bodies.append((p1,w1,w2))

O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Polyhedra_Aabb(),Bo1_Wall_Aabb()]),
InteractionLoop(
[Ig2_Wall_Polyhedra_PolyhedraGeom(), Ig2_Polyhedra_Polyhedra_PolyhedraGeom()],
[Ip2_PolyhedraMat_PolyhedraMat_PolyhedraPhys()],
[Law2_PolyhedraGeom_PolyhedraPhys_Volumetric()]
),
NewtonIntegrator(),
PolyhedraSplitter(iterPeriod=100),
]

O.dt= 0.1*polyhedra_utils.PWaveTimeStep()

plot.plots={'t':'sig'}
plot.plot(subPlots =False)

try:
qt.Controller()
v = qt.View()
v.ortho = True
except:
pass

O.run(50000),

## Question information

Language:
English Edit question
Status:
Solved
For:
Assignee:
No assignee Edit question
Solved by:
Karol Brzezinski
Solved:
2021-01-08
Last query:
2021-01-08