Triaxial total stress-strain curve using JCFpmMat

My code is as follows:

mn,mx=Vector3(0,0,0),Vector3(0.05,0.05,0.05)
O.materials.append(JCFpmMat(type=1,density=2000,young=9e9,poisson=0.3,tensileStrength=3e7,cohesion=3e7,frictionAngle=radians(20),label='sphere'))
O.materials.append(JCFpmMat(type=1,frictionAngle=0,density=0,label='wall'))

walls=aabbWalls([mn,mx],thickness=0,material='wall')
wallIds=O.bodies.append(walls)
sp=pack.SpherePack()
pred=pack.inAlignedBox(mn,mx)
sp=pack.randomDensePack(pred,radius=0.0015,rRelFuzz=0.333,returnSpherePack=True,memoizeDb='/tmp/triax.sqlite')
spheres=sp.toSimulation(material='sphere')

triax=TriaxialStressController(
 maxMultiplier=1.+2e4/young,
 finalMaxMultiplier=1.+2e3/young,
 thickness = 0,
 stressMask = 7,
 internalCompaction=True,
)
newton=NewtonIntegrator(damping=damp)
O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
  [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys()],
  [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=True,Key='triax_JCFPM',label='interactionLaw')]
 ),
 GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=10,timestepSafetyCoefficient=0.8),
 triax,
 TriaxialStateRecorder(iterPeriod=100,file='WallStresses'+key),
 newton,
]

triax.goal1=triax.goal2=triax.goal3=-1e6
while 1:
 O.run(500,1)
 unb=unbalancedForce()
 print('unbalanced force:',unb,'mean stres:',triax.meanStress)
 if unb<stabilityThreshold :
  break

triax.internalCompaction=False
triax.stressMask=3
triax.goal1=-1e6
triax.goal2=-1e6
triax.goal3=-0.02

O.run()

(I have omitted the post-processing parts such as plot,make it a MWE)

Hello,

Can you please limit your questions to 1 question per thread? This is our forum guideline, which is more precisely outlined here [1]. The reason we ask you to do this is because it helps keep the forum organized, easier to search, easier to read. Hence it is a strong knowledge base. We are happy to help you, but we also what to be as effective as possible with our volunteered time here. Thank you for helping keep the knowledge base organized and adhering to our guidelines [1]. Please read that.

>>1.The purpose of Isotropic compression stage is to make the particles dense,right?If so,in my script I use the randomDensePack,Is it necessary to retain this stage?

You still need to reach the stable unbalanced force, whether you seek a desired porosity or not.

2.In the deviatoric stage,the example use triax.goal2=rate to control the compression.I also try to use it in my script(triaxial with JCFpm),but after 4500000 iter,I can not get the ideal total stress-strain curve.(Is it related to my material parameters and strain rate?)

What is the "ideal total stress-strain curve?" How is it different from what you are computing?

3.With the JCFpmMat,I can record the number of ShearCracks and TensCracks.Can you give me a suggestion on how to control the end of the test (specimen compression failure)?

You can use PyRunner to insert logic into the engine list. Based on whatever logic you please, you can stop the simulation. Please review our documentation on PyRunners [2].

Cheers,

Robert

[1]https://www.yade-dem.org/wiki/Howtoask
[2]https://yade-dem.org/doc/tutorial-data-mining.html#pyrunner

Hello,Robert,
Sorry for violating the guideline.Sometimes I can't distinguish whether there is correlation between several problems...I will ask only a question per thread next time.Sorry again.

>>What is the "ideal total stress-strain curve?" How is it different from what you are computing?
In my understanding,in the traditional triaxial compression test of rock,the total stress-strain curve should include the following(at least):
1.Elastic stage(Linear rising curve ) 2.Plastic deformation and micro cracks appear (Nonlinear rising curve) 3.Strain softening stage (Descending curve) There should be a peak stress on behalf of the strength between 2. and 3.

As for my curve, I've run the script 1e8 iter so far,the stress-strain curve fluctuates all the time,and the stress is between 4e6 to 5.8e6... It seems that there is no damage to the sample?

>>You can use PyRunner to insert logic into the engine list. Based on whatever logic you please, you can stop the simulation.
Yes,I know it,but I have no idea about how to know the stop condition.Could you please give me some advice or reference.

Thanks again.

Hello,

>>Yes,I know it,but I have no idea about how to know the stop condition.Could you please give me some advice or reference.

I guess this depends on what you want to study...if you are trying to study the strain-softening phase then you'd need different stop logic than if you are only interested in the yield strength. Yield strength, look for when deviatoric stress is less than XX% of max historical stress. Strain softening, look for a max strain value to decide when to stop.

>>As for my curve, I've run the script 1e8 iter so far,the stress-strain curve fluctuates all the time,and the stress is between 4e6 to 5.8e6... It seems that there is no damage to the sample?

So there is no linear increase in stress before stabilizing? Sounds like the simulation is not stable. Plot unbalanced force to confirm.

Cheers,

Robert

Hello,

In fact,I want to study the complete stress-strain process,both yield strength and strain-softening phase.
In other words,It seems that the curve should have a rapid decline process, representing the failure of the sample?I have run the script with the triax.strain[2]=0.4(big enough?),but there is still no rapid decline of sigma,the stress still fluctuates up and down,unbalancedForce is very small (4e-5).
Another guess,I recorded nbTensCracks and nbShearCracks.For my whole process,they have been increasing.(Tc>>Sc).Cracks stands for the fracture of bond,right?If so,can the failure be judged by controlling the proportion of broken bonds to total bonded bonds?

Thanks for your patient help.

Hello,

>>I have run the script with the triax.strain[2]=0.4(big enough?)

Are you trying to set that value? It is just an output quantity, it is not a user input. [1]

I tried to run your code, but it does not work:

"NameError: name 'young' is not defined"

After I fix that, it then throws:

"NameError: name 'damp' is not defined"

Then I quit.

If you really want help, please review [2].

-rc

[1]https://yade-dem.org/doc/yade.wrapper.html#yade.wrapper.TriaxialStressController.strain
[2]https://www.yade-dem.org/wiki/Howtoask

Hello,
Sorry for post a MWE that can not run..(My original intention is to worry about the trouble caused by too long code),I will post my script that can be run.
>>Are you trying to set that value? It is just an output quantity, it is not a user input.
In fact, I didn't clearly express my meaning. I mean that I keep running, and the strain value has been large (=0.4), but there is still no stress-softening stage(the stress did not decrease rapidly,which represent failure,I guess)

Following is my script:
###########
from yade import pack,plot, utils, export, timing
import numpy as np
damp=0.2
stabilityThreshold=0.01
key='_triax_base_'
young=3.18e9

mn,mx=Vector3(0,0,0),Vector3(0.05,0.05,0.05)
O.materials.append(JCFpmMat(type=1,density=2000,young=9e9,poisson=0.3,tensileStrength=2e7,cohesion=2e7,frictionAngle=radians(20),label='sphere'))
O.materials.append(JCFpmMat(type=1,frictionAngle=0,density=0,label='wall'))

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

sp=pack.SpherePack()
pred=pack.inAlignedBox(mn,mx)
sp=pack.randomDensePack(pred,radius=0.0015,rRelFuzz=0.333,returnSpherePack=True,memoizeDb='/tmp/triax.sqlite')
spheres=sp.toSimulation(material='sphere')

triax=TriaxialStressController(
 maxMultiplier=1.+2e4/young,
 finalMaxMultiplier=1.+2e3/young,
 thickness = 0,
 stressMask = 7,
 internalCompaction=True,
)
newton=NewtonIntegrator(damping=damp)

def recorder():
 yade.plot.addData(
 i=O.iter,
 e11=-triax.strain[0],e22=-triax.strain[1],e33=-triax.strain[2],
 s11=-triax.stress(triax.wall_right_id)[0],
 s22=-triax.stress(triax.wall_top_id)[1],
 s33=-triax.stress(triax.wall_front_id)[2],
 numberTc=interactionLaw.nbTensCracks,
 numberSc=interactionLaw.nbShearCracks,
 unb=unbalancedForce()
)
 plot.saveDataTxt('triax_JCFPM_Wy=1.5Mpa,rate=-0.03')

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
  [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys()],
  [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=True,Key='triax_JCFPM_Wy=1.5Mpa,rate=-0.03',label='interactionLaw')]
 ),
 GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=10,timestepSafetyCoefficient=0.8),
 triax,
 TriaxialStateRecorder(iterPeriod=100,file='WallStresses'+key),
 newton,
 PyRunner(iterPeriod=int(1000),initRun=True,command='recorder()',label='data',dead=0),
]

triax.goal1=triax.goal2=triax.goal3=-1.5e6
while 1:
 O.run(500,1)
 unb=unbalancedForce()
 print('unbalanced force:',unb,'mean stres:',triax.meanStress)
 if unb<stabilityThreshold :
  break

triax.internalCompaction=False
triax.stressMask=3
triax.goal1=-1.5e6
triax.goal2=-1.5e6
triax.goal3=-0.03

plot.plots={'e33':('s33',None,'unb'),'i':('numberTc','numberSc',None,'s33')}
plot.plot()
O.run()
################
Thanks!

Hello,

I want to more clearly describe my current confusion about this simulation:
My current main problem is that my triaxial compression simulation does not have a strain softening process (that is, the stress does not drop significantly), and the stress is still fluctuating up and down when the strain is already large (0.5 or more).
I read other triaxial compression related papers that use bonded particle model (BPM), and there is such a strain-softening stage.
What should I do to make the curve become a full stress strain curve?

Thanks.

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