yade-dev team mailing list archive

[Anton Gladky <gladky.anton@xxxxxxxxx>]

Hi, Vaclav,

I have attached the script, which shows different results on 2091 and 2092
revisions. It is formal uniaxial-test script, I have made it smaller,
getting rid of unnecessary stuff.
It makes 20 steps and calculates the number of isCohesive contacts.

isCohesive contacts can be created in RockPM during first 10 steps, but in
99% cases they are creating only on a first step.

It seems, that other parameters did not change efficiently after r2091. Just
a little bit. It can be connected with fixing regularHexa packing recently.

I appreciate your help and ready to help you in fixing this problem.

Thank you.
______________________________
Anton Gladkyy


2010/4/1 Václav Šmilauer <eudoxos@xxxxxxxx>

> > There is an output from 2091:
> > #iter  Cohesive
> > 0 17716
> > 1 17716
>
> > r2092:
> > #iter  Cohesive
> > 0 17716
> > 1 18535
>
> This looks serious.
>
> Can you send me over the script (in private if you dont' want on the
> list) and I will investigate it. There is no reason why those number
> should differ.
>
> Are results (other than number of interactions) also different?
>
> Vaclav
>
>
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-dev<https://launchpad.net/%7Eyade-dev>
> Post to     : yade-dev@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~yade-dev<https://launchpad.net/%7Eyade-dev>
> More help   : https://help.launchpad.net/ListHelp
>

# -*- encoding=utf-8 -*-
from yade import pack,ymport
from yade import utils
from math import *


#Read parameters from table_____________________________________________
utils.readParamsFromTable(grainSizeInput=4,grainSizeMinMax=2,specimeSizeDiametrInput=50,specimeSizeLengthInput=100,strengthInput=100,correctionCoefCompress=3.00,correctionCoefStretch=0.04,correctionCoefShear=0.08,pressTypeInput=0,dampingInput=0.1,intRadiusInput=1.1,velocityPressInput=0.05,frictionAngleInput=0.52,dTimeInput=0,packingAlgorithmInput=0,generateAgainInput=0,noTableOk=True)
#_______________________________________________________________________


PI = pi;
a = specimeSizeDiametrInput #Size a, [mm] X
b = specimeSizeDiametrInput #Size b, [mm] Y
h = specimeSizeLengthInput #Size h, [mm] Z

grainSize = grainSizeInput #Grain size, [mm]
distBetweenGrains = 0 #Distance between grains, [mm]

rho=2500 #Grain density, [kg/m^3]
poisson = 0.3 #Poisson coefficient
young = 37e9 #Young modulus
frictionAngle = frictionAngleInput # Friction angle
stressCompressMax = strengthInput*1000000
intRadius = intRadiusInput

boxThick = 3 #Box thickness, [mm]
boxHeight = 10 #Box height, [mm]
boxCoef = 3 #In how many times box larger, than an example
pressCoef = 5 #In how many times press larger, than an example

#_______________________________________________________________________
a = 0.001*a
b = 0.001*b
h = 0.001*h
grainSize = 0.001*grainSize
radiusGrain = grainSize/2
distBetweenGrains = 0.001*distBetweenGrains
boxThick = 0.001*boxThick
boxHeight = 0.001*boxHeight


n2Mpa=PI*((a)**2)/4*1000000					#Compression
predictedForce = n2Mpa/1000000*stressCompressMax/correctionCoefCompress

videoOn=1
#Material Properties____________________________________________________
RockId=O.materials.append(RpmMat(young=young,frictionAngle=frictionAngle,poisson=poisson,density=rho,exampleNumber=1,stressCompressMax=100000000,Brittleness=0.15,initCohesive=True))
RockPress=O.materials.append(RpmMat(young=young,frictionAngle=frictionAngle,poisson=poisson,density=rho,exampleNumber=5))
kw_press={'color':[0.757,0.278,0.976],'dynamic':False,'wire':False,'material':RockPress}

ids_spheres=O.bodies.append(pack.regularHexa(pack.inCylinder((0,0,-boxHeight/2*boxCoef),(0,0,-boxHeight/2*boxCoef+h),a/2),radius=radiusGrain,gap=distBetweenGrains,color=(0.282,0.875,1),material=RockId))

print str(len(ids_spheres)) + ' particles generated'


#Find the hihghest and the lowest sphere________________________________
maxHight=O.bodies[ids_spheres[0]].state.pos[2]+O.bodies[ids_spheres[0]].shape.radius
minHight=O.bodies[ids_spheres[0]].state.pos[2]-O.bodies[ids_spheres[0]].shape.radius

maxLength=O.bodies[ids_spheres[0]].state.pos[1]+O.bodies[ids_spheres[0]].shape.radius
minLength=O.bodies[ids_spheres[0]].state.pos[1]-O.bodies[ids_spheres[0]].shape.radius

maxWidth=O.bodies[ids_spheres[0]].state.pos[0]+O.bodies[ids_spheres[0]].shape.radius
minWidth=O.bodies[ids_spheres[0]].state.pos[0]-O.bodies[ids_spheres[0]].shape.radius

for i in ids_spheres:
	if ((O.bodies[i].state.pos[2]+O.bodies[ids_spheres[i]].shape.radius)>maxHight):
		maxHight=O.bodies[i].state.pos[2]+O.bodies[ids_spheres[i]].shape.radius
	if ((O.bodies[i].state.pos[2]-O.bodies[ids_spheres[i]].shape.radius)<minHight):
		minHight=O.bodies[i].state.pos[2]-O.bodies[ids_spheres[i]].shape.radius
    
	if ((O.bodies[i].state.pos[1]+O.bodies[ids_spheres[i]].shape.radius)>maxLength):
		maxLength=O.bodies[i].state.pos[1]+O.bodies[ids_spheres[i]].shape.radius
	if ((O.bodies[i].state.pos[1]-O.bodies[ids_spheres[i]].shape.radius)<minLength):
		minLength=O.bodies[i].state.pos[1]-O.bodies[ids_spheres[i]].shape.radius
    
	if ((O.bodies[i].state.pos[0]+O.bodies[ids_spheres[i]].shape.radius)>maxWidth):
		maxWidth=O.bodies[i].state.pos[0]+O.bodies[ids_spheres[i]].shape.radius
	if ((O.bodies[i].state.pos[0]-O.bodies[ids_spheres[i]].shape.radius)<minWidth):
		minWidth=O.bodies[i].state.pos[0]-O.bodies[ids_spheres[i]].shape.radius

realHight = maxHight-minHight
realLength = maxLength-minLength
realWidth = maxWidth-minWidth

#Press preparation______________________________________________________
O.bodies.append(utils.facet([[(minWidth+realWidth/2-2*realWidth),0,minHight],[(minWidth+realWidth/2+realWidth),(minLength+realLength/2+realLength*2),(minHight)],[(minWidth+realWidth/2+realWidth),(minLength+realLength/2-realLength*2),(minHight)]],**kw_press))

#Engines________________________________________________________________
O.engines=[
	ForceResetter(),
	BoundDispatcher([Bo1_Sphere_Aabb(aabbEnlargeFactor=intRadius,label='is2aabb'),Bo1_Facet_Aabb()]),
	InsertionSortCollider(),
	InteractionDispatchers(
		[Ig2_Sphere_Sphere_Dem3DofGeom(distFactor=intRadius,label='ss2d3dg'),Ig2_Facet_Sphere_Dem3DofGeom()],
		[Ip2_RpmMat_RpmMat_RpmPhys()],
		[Law2_Dem3DofGeom_RockPMPhys_Rpm()]
	),
  GravityEngine(gravity=[0,0,-9.81]),
	NewtonIntegrator(damping=dampingInput),
	PeriodicPythonRunner(iterPeriod=1,command='cohesiveContacts()',label='plotDataCollector')
]

def cohesiveContacts():
	print len([1 for i in O.interactions if i.phys.isCohesive])
	

#Start Simulation_______________________________________________________
O.dt=utils.PWaveTimeStep()*0.2
O.saveTmp('init')

# to create initial contacts
O.step()
# now reset the interaction radius and go ahead
ss2d3dg.distFactor=1.
is2aabb.aabbEnlargeFactor=1.

# This is for multi-tasks
O.run(20)
O.wait()