New question #701341 on Yade:
https://answers.launchpad.net/yade/+question/701341

Hello,
I have encountered a problem similar to that in [1], but my script is different 
from it. I have read the relevant answers in [1].

My question Description:
My script is about fluid-solid coupling triaxial compression,relevant research 
is carried out by changing confining pressure and osmotic pressure.I have 
successfully run some scripts (under the condition of confining pressure of 20, 
40 and 60MPa respectively, the osmotic pressure is 10, 20, 30, 40 and 50MPa 
respectively).
However, when the confining pressure is 80MPa, the script can run normally when 
the osmotic pressure is 10 and 20MPa. When it increases to 30 and 40MPa, an 
error prompt as described in the title appears.

Some problems related to locking some degrees of freedom and walls are 
mentioned in [1],In my script, I did turn off the movement of the wall in the X 
and Y directions (the script is shown below).I have tried to set the thickness 
of the wall to a positive value, the error does not appear, but the result is 
completely different from the law obtained from my previous work (the peak 
strength and peak strain of the sample become larger)

I'm confused and don't know how to solve this problem. In other words, if it's 
because of the wall thickness problem, why can the previous scripts run 
normally?

Thanks for kindly help!!!(Maybe my problem description is too long..sorry)

---------------------------------------------script-----------------------------------------
from yade import pack, ymport, plot, utils, export, timing
from builtins import range 

Sy=40e6
Wy=-80e6
rate=-0.2
damp=0.4
stabilityThreshold=0.001
key='_triax_base_'
young=3000e9
name='JCFPM_triax'
compFricDegree=30
poisson=0.4
OUT=str(Wy)+'_HM-coupling'

mn,mx=Vector3(0,0,0),Vector3(0.05,0.05,0.05)
O.materials.append(JCFpmMat(type=1,density=2640,young=young,poisson=poisson,tensileStrength=160e6,cohesion=1600e6,frictionAngle=radians(10),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)

O.bodies.append(ymport.text('packing-JCFPM.spheres',scale=1,shift=Vector3(0,0,0),material='sphere',color=(0.6,0.5,0.15)))

triax=TriaxialStressController(
        maxMultiplier=1.+4e8/young, 
        finalMaxMultiplier=1.+4e7/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],#0+边界id为right
        #s22=-triax.stress(triax.wall_top_id)[1],#1+边界id为top
        s33=-triax.stress(triax.wall_front_id)[2],#2+边界id为front
        numberTc=interactionLaw.nbTensCracks,
        numberSc=interactionLaw.nbShearCracks,
        unb=unbalancedForce(),
        flux1=flow.getBoundaryFlux(4),
        flux2=flow.getBoundaryFlux(5),
        timestep=O.dt
)
        plot.saveDataTxt(OUT)

def stop_condition():
        extremum=max(abs(s) for s in plot.data['s33'])
        s=abs(plot.data['s33'][-1])
        e=abs(plot.data['e33'][-1])
        if  e < 0.001:
                return
        if abs(s)/abs(extremum) < 0.9 and abs(s)>abs(Wy) :
                print('Max stress and strain:',extremum,e)
                presentcohesive_count = 0
                for i in O.interactions:
                        if hasattr(i.phys, 'isCohesive'):
                                if i.phys.isCohesive == True:
                                        presentcohesive_count+=1
                print('the number of cohesive bond now 
is:',presentcohesive_count)
                O.wait()



O.engines=[
        ForceResetter(),
        
InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=1.5,label='is2aabb'),Bo1_Box_Aabb()]),
        InteractionLoop(
                
[Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=1.5,label='ss2sc'),Ig2_Box_Sphere_ScGeom()],
                [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys()],
                
[Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=True,Key=OUT+'_Crack',label='interactionLaw')]
        ),
        
GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
        triax,
        FlowEngine(dead=1,label='flow'),
        
PyRunner(iterPeriod=int(1000),initRun=True,command='recorder()',label='data',dead=0),
        PyRunner(iterPeriod=1000,command='stop_condition()',dead=0),
        
VTKRecorder(iterPeriod=500,initRun=True,fileName='triax/JFFPM-',recorders=['spheres','jcfpm','cracks'],Key='triax',label='vtk',dead=1),
        newton,
]

O.step()
ss2sc.interactionDetectionFactor=-1
is2aabb.aabbEnlargeFactor=-1

cohesiveCount = 0
for i in O.interactions:
 if hasattr(i.phys, 'isCohesive'):
              if i.phys.isCohesive == True:
                     cohesiveCount+=1
print('the origin total number of cohesive bond is:',cohesiveCount)

triax.goal1=triax.goal2=triax.goal3=Wy
while 1:
        O.run(1000,1)
        unb=unbalancedForce()
        print('unbalanced force:',unb,'mean stres:',triax.meanStress)
        if unb<stabilityThreshold and abs(Wy-triax.meanStress)/abs(Wy)<0.001:
                break
cohesiveCount = 0
for i in O.interactions:
 if hasattr(i.phys, 'isCohesive'):
              if i.phys.isCohesive == True:
                     cohesiveCount+=1
print('the first total number of cohesive bond is:',cohesiveCount)

yade.timing.reset()
flow.dead=0
flow.debug=False
flow.fluidBulkModulus=2.2e9
flow.permeabilityFactor=10
flow.decoupleForces =  False
flow.meshUpdateInterval=1000
flow.useSolver=4
flow.viscosity=0.001
flow.bndCondIsPressure=[0,0,0,0,1,1]
flow.bndCondValue=[0,0,0,0,0,Sy]

triax.internalCompaction=False
flow.emulateAction()
triax.stressMask=3
triax.wall_left_activated=False
triax.wall_right_activated=False
triax.wall_bottom_activated=False
triax.wall_top_activated=False
triax.goal1=Wy
triax.goal2=Wy
triax.goal3=rate

plot.plots={'e33':(('s33','g'),None,('unb','b')),'i':(('numberTc','b'),('numberSc','r'),None,('s33','y'))}
plot.plot(subPlots=False)
O.run()
-----------------------------------------------------------------------------

[1]https://answers.launchpad.net/yade/+question/693071

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.

_______________________________________________
Mailing list: https://launchpad.net/~yade-users
Post to     : yade-users@lists.launchpad.net
Unsubscribe : https://launchpad.net/~yade-users
More help   : https://help.launchpad.net/ListHelp

Reply via email to