About transverse strain

Asked by AnnaMaria on 2019-11-14

Hello,
I am a new YADE user, trying to calibrate my DEM model
by performing uniaxial compression tests with UniaxialStrainer.
I need to compute transverse strain in a cylindrical specimen,
in order to match it with macroscopic poisson’s ratio.
Any ideas on how i could do it?
I think its not “part” of uniaxial strainer.
Thanks

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Jan Stránský
Solved:
2019-11-28
Last query:
2019-11-28
Last reply:
2019-11-28
Jérôme Duriez (jduriez) said : #1

Hi,

No, it's not part of UniaxialStrainer which considers only axial quantities. You need to measure yourself the transversal dimensions of your specimen. First possibility could be to use aabbDim()

Jan Stránský (honzik) said : #2

Hello,

> I think its not “part” of uniaxial strainer.

exactly. Because computation of it is not very straightforward. E.g. you can average "particle transverse strain" as newDistanceFromAxis/originalDistanceFromAxis-1. Then particles near the axis would be problematic because their original length is low (or possibly 0) and would corrupt the results.
Also it may happen that boundary particles "fly away", also corrupting the result.
So you can e.g. somehow choose a "central ring" and do computation on these particles.

Apart from computing transverse strain from uniaxial stress test, there is another option how to compute Poisson's ratio.
Poisson's ratio can also be computed by comparing Young's modulus (uniaxial stress modulus) to oedometric (uniaxial strain) modulus and/or bulk ("triaxial") modulus. See related discussion [1].

You can use both methods (transverse strain and another loading) and compare the results of Poisson's ratio :-)

cheers
Jan

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

AnnaMaria (annm5) said : #3

Thanks for your response Jérôme Duriez (jduriez).
I read in [1], as Jan Stránský (honzik) suggested that you have
to to define which particles are "suitable" for the computation.
That very near the axis had negative influence, as well as those "crushed away" from the surface and flying away. I choose those roughly half-way from the axis to surface.
Can I do this with aabbDim()? Choose those half-away?
Access a region within the cylinder?
utils.aabbDim(cutoff=0.0, centers=False)

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

Thanks

AnnaMaria (annm5) said : #4

Thanks Jan Stránský (honzik), for your response
it was while i was typing a comment on what Jérôme Duriez (jduriez) said
and i did not see it in time.

Im not sure how i can access this ring. With aabbDim?

Thanks

Jan Stránský (honzik) said : #5

> Can I do this with aabbDim()?

no, aabbDim computes aabb (axis aligned bounding box) dimensions of the sample (roughly the "maximum dimensions along axes").

To choose particles within a "ring" defined by radii r1 and r2, do:
### (not tested)
halfwayParticles = []
for b in O.bodies:
   x,y,z = b.state.pos # maybe b.state.refPos
   d = sqrt(pow(x,2)+pow(y,2)) # assuming z axis
   if d > r1 and d < r2:
      halfwayParticles.append(b)
###

cheers
Jan

AnnaMaria (annm5) said : #6

Thanks Jan Stránský (honzik)
i will try what you said ;)

An alternative approach, which may fit your needs or not, would be to use TriaxialStressController by prescribing a very small confinement on the other two directions. That would prevent particles from flying away and it would keep the boundaries in contact with the specimen, which makes defining transverse strain straightforward.
If the confinement is small enough it should approach the uniaxial response.
Cheers
Bruno

Jan Stránský (honzik) said : #8

Bruno's answer would work for prismatic shape. If you have cylinder, it would not work.

@Anna: what is the shape fo your specimen?

cheers
Jan

AnnaMaria (annm5) said : #9

Thanks a lot for your comments Jan and Bruno.
My specimen is cylindrical.
I am still trying to make it work. (the way that Jan Stránský (honzik) suggested)

AnnaMaria (annm5) said : #10

Hello everyone again:
I have created the hollow cylinder as Jan Stransky suggested, (MWE):
for p in O.bodies:
  x=p.state.refPos[0]
  z=p.state.refPos[2]
  d=sqrt(pow(x,2)+pow(z,2))
  if d>r1 and d<r2:
   p.shape.color=(1,1,1)
   #Mesure state.refPos of the particles, within this domain.
   xnew=p.state.refPos[0]
   znew=p.state.refPos[2]
   dnew=sqrt(pow(xnew,2)+pow(znew,2))
Then I am trying to record the lateral strain:
def recorder():
 global d
 global dnew
 elat=((d)-(dnew))/(d)
but it does not work. I am not getting any values!
I can’t understand what i have been doing wrong all this time.
Could you help me. Thanks in advance.

Jan Stránský (honzik) said : #11

Please provide a real MWE (W=working, a complete script which we can test ourselves).
Thanks
Jan

AnnaMaria (annm5) said : #12

Thanks Jan,
This is a working example.
I want to plot in the positive domain of the Cartesian coordinate system the stress vs strain, and in the negative the lateral strain. But right now it does not seem to work, even if only plot lateral strain. Stress seems to plot ok.

from yade import plot
#### Simulation Control
rate=-1 #deformation rate

#### Material
# default parameters or from table
readParamsFromTable(noTableOk=True, # unknownOk=True,
 Specimen_Radius = 0.027,
 Specimen_Height= 0.13,
 Sphere_Radius = 2e-3,
 intR=1.5,
 DENSITY=3000,
 YOUNG=60e9,
 FRICTION_ANGLE=50,
 POISSON=1,
 TENS=5e6,
 COHESIO=40e6,
 r1=0.010,
 r2=0.020,
 h1=0.061,
 h2=0.061,
)

from yade.params.table import *

#### material definition
Sample = O.materials.append(JCFpmMat(
young=YOUNG,
poisson=POISSON,
frictionAngle=radians(FRICTION_ANGLE),
cohesion=COHESIO,
tensileStrength=TENS,
density = DENSITY,
label='spheres'
))
#### Specimen Creation
sp=pack.randomDensePack(
 pack.inCylinder((0,0,0),(0,Specimen_Height,0),Specimen_Radius),
 radius=Sphere_Radius,
 spheresInCell=1000,
 material=Sample,
 rRelFuzz=0.12,
 memoizeDb='/tmp/Cpack.db',
 returnSpherePack=True,
)
sp.toSimulation()
for p in O.bodies:
  x=p.state.refPos[0]
  z=p.state.refPos[2]
  d=sqrt(pow(x,2)+pow(z,2))
  if d>r1 and d<r2:
   p.shape.color=(1,1,1)
   xnew=p.state.refPos[0]
   znew=p.state.refPos[2]
   dnew=sqrt(pow(xnew,2)+pow(znew,2))

#### boundary condition (utils.uniaxialTestFeatures)
bb=utils.uniaxialTestFeatures()
negIds,posIds,longerAxis,crossSectionArea=bb['negIds'],bb['posIds'],bb['axis'],bb['area']

################# ENGINES DEFINED HERE

O.engines=[
 ForceResetter(),
        InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=intR,label='Saabb')]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=intR,label='SSgeom')],
  [Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1,label='interactionPhys')],
  [Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(recordCracks=True,label='interactionLaw')]
 ),
 UniaxialStrainer(strainRate=rate,axis=longerAxis,asymmetry=0,posIds=posIds,negIds=negIds,crossSectionArea=crossSectionArea,blockDisplacements=1,blockRotations=1,setSpeeds=0,stopStrain=0.1,dead=1,label='strainer'),
 GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=10,timestepSafetyCoefficient=0.5, defaultDt=utils.PWaveTimeStep()),
 NewtonIntegrator(damping=0.2,label='newton'),
 PyRunner(iterPeriod=int(100),initRun=True,command='recorder()',label='data'),
]

################# RECORDER DEFINED HERE

def recorder():
  global d
  global dnew
  elat=(d-dnew)/d
  yade.plot.addData({'i':O.iter,
         'eps':abs(strainer.strain),
  'elat':elat,
         'sigma':abs(strainer.avgStress),})

#### Plot During Simulation
plot.plots={'elat':('sigma')}
plot.plot()

O.step();
### initializes the interaction detection factor
SSgeom.interactionDetectionFactor=-1.
Saabb.aabbEnlargeFactor=-1.
strainer.dead=0

Jan Stránský (honzik) said : #13

1) I have no problem running the script, plot.data shows the saved data.

2) Python/Yade reads the script and executes it line by line. In your case, it loops over the bodies and computes the d and dnew ONLY once at the beginning, then the value is unchanged during the whole simulation.
Put the code inside a function and call that function inside recorder() before plot.addData()

3) furthermore, d and dnew are just values form the last single body.
You want probably some average (maybe wighted?) of the values.

4) expressions for d and dnew does not differ, probably dnew should be computed using state.pos instead of state.refPos?

cheers
Jan

AnnaMaria (annm5) said : #14

Hi Jan,
I think I am getting close
2) I Changed the script as follows:
After creating the packing:
I added the following:
#sp.toSimulation()
def lateral():
   global elat
   global d
   global dnew
   for p in O.bodies:
    x=p.state.refPos[0]
    z=p.state.refPos[2]
    d=sqrt(pow(x,2)+pow(z,2))
    if d>r1 and d<r2:
     p.shape.color=(1,1,1)
     xnew=p.state.pos[0]
     znew=p.state.pos[2]
     dnew=sqrt(pow(xnew,2)+pow(znew,2))
     elat=(d-dnew)/d
     return elat

and then:
def recorder():
 lateral()
yade.plot.addData({'i':O.iter,'eps':abs(strainer.strain),'elat':(elat),'sigma':abs(strainer.avgStress)})

as you said.
The thing is that the function above does not seem to work correct. The spheres have the same color, while before they were white (color=(1,1,1)). Maybe the distance measurement is not working as well.
Also what do you mean:
3) furthermore, d and dnew are just values form the last single body.
You want probably some average (maybe wighted?) of the values.
Like this: elat=(avg(d)-avg(dnew))/avg(d)

Best Jan Stránský (honzik) said : #15

5)
Use more than one space for indentation. Here your recorder would be
###
def recorder():
    lateral()
yade.plot.addData({'i':O.iter,'eps':abs(strainer.strain),'elat':(elat),'sigma':abs(strainer.avgStress)})
###
meaning plot.addData outside recorder(), which is probably not what you want

6)
> The spheres have the same color
no, there is one sphere with different color.
In laterla(), you return directly when the very first suitable body is found.

6)+3) solution
###
def lateral():
    elatTot = 0.0 # sum of elat of all suitable particles
    nTot = 0 # number of suitable particles
    for p in O.bodies:
        x=...; z=...; d=...
        if d>r1 and d<r2:
            xnew=...; znew=...; dnew=...
            elat=(d-dnew)/d
            elatTot += elat
            nTot += 1
    return elatTot / nTot # average elat
###

cheers
Jan

AnnaMaria (annm5) said : #16

Jan Thank you very much for your effort!!!!
I really appreciate it.
I have only changed elatTot / nTot to elat_avg=elatTot / nTot
And then return elat_avg # average elat so that I can plot.

def lateral():
    elatTot = 0.0 # sum of elat of all suitable particles
    nTot = 0 # number of suitable particles
    for p in O.bodies:
 x=p.state.refPos[0]
 z=p.state.refPos[2]
 d=sqrt(pow(x,2)+pow(z,2))
        if d>r1 and d<r2:
           xnew=p.state.pos[0]
            znew=p.state.pos[2]
            dnew=sqrt(pow(xnew,2)+pow(znew,2))
            elat=(d-dnew)/d
            elatTot += elat
            nTot += 1
     elat_avg=elatTot / nTot
    return elat_avg # average elat

def recorder():
 lateral()
 yade.plot.addData({'i':O.iter,'eps':abs(strainer.strain),'elateral':(elat_avg),'sigma':abs(strainer.avgStress)})
 plot.saveDataTxt(OUT)

#### Plot During Simulation
plot.plots={'elateral':('sigma')}
plot.plot()

I also changed the recorder part. But when I run the code, I get the following Error:
File "<string>", line 1, in <module>
File "UCS.py", line 87, in recorder
yade.plot.addData({'i':O.iter,'eps':abs(strainer.strain),'elateral':(elat_avg),'sigma':abs(strainer.avgStress)})
NameError: global name 'elat_avg' is not defined

AnnaMaria (annm5) said : #17

I am using more than one spaces
But when i copy paste my code in lanchpad it just removes them.
def recorder():
 lateral()
 yade.plot.addData({'i':O.iter,'eps':abs(strainer.strain),'elateral':(elat_avg),'sigma':abs(strainer.avgStress)})

Jan Stránský (honzik) said : #18

def recorder():
 elat_avg = lateral() # here
 yade.plot.addData(...)

AnnaMaria (annm5) said : #19

Thanks Jan Stránský, that solved my question.