yade-users team mailing list archive

[Huihuang Xia <question611087@xxxxxxxxxxxxxxxxxxxxx>]

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

Hi,

https://answers.launchpad.net/yade/+question/432617, in this question, damage ratio was defined and implemented in YADE. 
According to this definition, I built a impact crushing simulation using damage ratio in YADE, but my code got a wrong result, as damage ratio in this simulation was greater than 1. 

Here is my code:

#!/usr/bin/python
# -*- coding: utf-8 -*-

#Created by Huihuang Xia from Huaqiao University,Xiamen,P.R.C
#My E-Mail:huihuangxia@xxxxxxxx

###########################
# IMPORT MODULES
###########################

from yade import pack
from yade import plot
from yade import qt

###########################
# DEFINE MATERIALS
###########################

rock=CohFrictMat(young=5.98e7,poisson=0.3,alphaKr=3000,alphaKtw=3000,density=2678,frictionAngle=0.5,isCohesive=True,normalCohesion=7.9e6,shearCohesion=7.9e6,momentRotationLaw=True)
O.materials.append(rock)
steel=CohFrictMat(young=3.06e11,poisson=0.29,density=7861,frictionAngle=0.545,normalCohesion=0,shearCohesion=0)
O.materials.append(steel)
    
#################################
# CREATE SAMPLE & RIGID_WALL
#################################
wall=O.bodies.append(geom.facetBox(center=(0,0,0),extents=(0.015,0.015,0.0005),color=(1,1,0),material=steel,fixed=True))
pred=pack.inSphere(center=(0,0,0.006),radius=0.005)
assembly=pack.randomDensePack(pred,radius=0.0002,rRelFuzz=0.5,spheresInCell=2500,material=rock)
O.bodies.append(assembly)
for b in assembly:
	b.state.vel=(0,0,-20)

###########################
# ENGINES
###########################
O.engines=[
    ForceResetter(),
    InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=1.0),Bo1_Facet_Aabb()]),
    InteractionLoop(
        [Ig2_Sphere_Sphere_ScGeom6D(interactionDetectionFactor=1.0),Ig2_Facet_Sphere_ScGeom6D()],
        [Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True)],
        [Law2_ScGeom6D_CohFrictPhys_CohesionMoment()]
    ),
    VTKRecorder(fileName='post/impact-',recorders=['all'],iterPeriod=250),
    GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=5,timestepSafetyCoefficient=0.8,defaultDt=PWaveTimeStep()),
    NewtonIntegrator(gravity=(0,0,-9.81)),
    PyRunner(command='damageRatio()',realPeriod=5),
    PyRunner(command='addPlotData()',realPeriod=5),
]
O.trackEnergy=True
O.step()

###########################
# DEFINE FUNCATIONS
###########################

#define a function to save damage ratio
global sumCohBonds
sumCohBonds=0
for b in O.interactions: 
    if b.phys.cohesionBroken==True:
        continue
    sumCohBonds+=1              # calculate total intact bonds

numBroCohBonds=0
def damageRatio():
    global numBroCohBonds
    for br in O.interactions: 
        if br.phys.cohesionBroken==False:
        	continue
        numBroCohBonds+=1      # calculate broken bonds
    damageRatio=numBroCohBonds/float(sumCohBonds)
    print sumCohBonds,numBroCohBonds
    plot.saveDataTxt('data/damageratio.txt')
def addPlotData():
    damageRatio=numBroCohBonds/float(sumCohBonds)
    plot.addData(t=O.time,dr=damageRatio)
plot.plots={'t':('dr')}
plot.plot()

#set an optimal timestep
O.dt=utils.PWaveTimeStep()
O.usesTimeStepper=True

#3D view and controller
qt.View()
qt.Controller()

and here is the damage ratio vs time:

# dr		t
0.0	            9.49679622677e-09
0.0044209618287 1.60495856089e-06
0.0155272805693	4.9193404254e-06
0.0266335993099	8.39516779415e-06
0.0362303213284	1.10352770982e-05
0.0474444684063	1.44826140458e-05
0.0585507871469	1.73886336283e-05
0.0679318524908	2.02186788283e-05
0.0788225145568	2.3419099093e-05
0.091654086694	2.71228495603e-05
0.104485658831	3.03992441997e-05
0.117856372655	3.34097285591e-05
0.130687944792	3.48247511657e-05
0.144058658615	3.669561998e-05
0.158291999137	3.75218412351e-05
0.172633167997	3.86994439409e-05
0.190748328661	4.02569184773e-05
0.210696571059	4.18618769955e-05
0.236144058659	4.52902203506e-05
0.278089281863	4.98581792485e-05
0.339120120768	5.47350982335e-05
0.42009920207	       5.83200136999e-05
0.520379555747	6.13408059573e-05
0.640931636834	6.41705226598e-05
0.776148371792	6.70709729927e-05
0.946409316368	7.14288492846e-05
1.15074401553	       7.49846615977e-05
1.39076989433	      7.87721800155e-05
1.67392710804	     8.25338810098e-05
1.99611817986	    8.60126620602e-05
2.32769031702	    8.71703956744e-05
2.67360362303	  8.81834109127e-05
3.04539573	         9.07206891265e-05
3.46236791029	9.43729049081e-05
3.92160879879	9.80223033349e-05
4.41718783696	0.000101386989899
4.94899719646	0.000104757207048
5.50323485012	0.000107345725132
6.08561569981	0.000110695198656
6.67565236144	0.000111623153105
7.272697865	0.000112554955106
7.87524261376	0.000113345795924
8.48781539789	0.000116016886332
9.10739702394	0.000118131382236
9.72881173172	0.000120064040351
10.3572352814	0.000121999979727
10.989756308	0.000123849996843
11.6188268277	0.000126003034127


PS: this simulation got no error and my Yade (V.2016.06a)

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