About periodic simple shear

Asked by Lei Hang on 2020-06-28

Hello, I am a new user for yade. When I study the yade documentation about periodic simple shear, I don't understand the meaning of "if 0:..., else:...". The detailed code is following:

from yade import pack,plot
# the "if 0:" block will be never executed, therefore the "else:" block will be
# to use cloud instead of regular packing, change to "if 1:" or something similar
if 0:
# create cloud of spheres and insert them into the simulation
# we give corners, mean radius, radius variation
sp=pack.SpherePack()
sp.makeCloud((0,0,0),(2,2,2),rMean=.1,rRelFuzz=.6,periodic=True)
# insert the packing into the simulation
sp.toSimulation(color=(0,0,1)) # pure blue
else:
# in this case, add dense packing
O.bodies.append(
pack.regularHexa(pack.inAlignedBox((0,0,0),(2,2,2)),radius=.1,gap=0,color=(0,0,
,!1))
)

Please help me understand the grammar of "if 0:..., else:...". Thanks a lot!

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Jan Stránský
Solved:
2020-07-02
Last query:
2020-07-02
Last reply:
2020-07-02

This question was reopened

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

Hello,

> I am a new user for yade

welcome :-)

next time please provide also a link to the script

This pattern is to easy switch between two different scenarios
Suppose:
###
if 0:
   createLoosePacking() # possibly a long code block
else:
   createDensePacking() # possibly a long code block
###

As it is, with "if 0:", is equivalent to (as the "if 0:" block is never executed and therefore "else" is always executed):
createDensePacking() # possibly a long code block

If you change it to
###
if 1:
   createLoosePacking() # possibly a long code block
else:
   createDensePacking() # possibly a long code block
###

then it is equivalent to (as "if 1" is always executed and therefore "else" is never executed):
createLoosePacking() # possibly a long code block

This way you can choose one scenario or the other just by altering 0/1 after the "if" keyword.

HTH
cheers
Jan

Lei Hang (h-stone) said : #2

Thanks for your quick help Jan! But I still don't understand (as the "if 0:" block is never executed and therefore "else" is always executed). The "if 0:" block is never executed, why we add the "if 0:"block? why we not just write the content of "else:" block? The same in "if 1:..., else:..."block, why we not just write the content of "if 1:"?

Thank you for your help!

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

Because you can easily choose one scenario or the other just by replacing 0 with 1 or vice versa.
The value "0" or "1" is not hardcoded, it is there to enable quick and easy selection.

If you know your final script needs just the "else" scenario, then just write that one.

But in the case of the chosen script, you can choose loose packing or you can choose densePacking.
If you let "if 0:", then you are using loose packing
If you change it to "if 1:", then you are using dense packing.
Then if you change it back to "if 0:", you use loose packing.
Then ...
everything by changing just one symbol. That is the meaning of it.

You could modify it e.g. in a more readable (but more verbose at the same time) code:
###
useLoosePacking = True
useDensePacking = False
# check that exactly one value is set to True
assert useLoosePacking or useDensePacking is True # at least one is True
assert useLoosePacking and useDensePacking is False # not both are True
#
if useLoosePacking:
   createLoosePacking() # possibly a long code block
elif useDensePacking:
   createDensePacking() # possibly a long code block
###

cheers
Jan

Lei Hang (h-stone) said : #4

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

Lei Hang (h-stone) said : #5

Hello,

When I study the yade script about periodic simple shear test, I don't understand some places. The periodic simple shear script link is following:

https://gitlab.com/yade-dev/trunk/blob/master/doc/sphinx/tutorial/04-periodic-simple-shear.py

I want to ask several questions:
1.What is the meaning of "O.cell.hSize=Matrix3(2,0,0, 0,2,0, 0,0,2)"? What is the meaning of the numbers in the Matrix?

2."pack.regularHexa(pack.inAlignedBox((0,0,0),(2,2,2)),radius=.1,gap=0,color=(0,0,1))". What is the meaning of the "color=(0,0,1)"? In some other scripts I also see "color=(1,1,1)"

3.What is the meaning of "O.cell.velGrad=Matrix3(-.1,0,0, 0,-.1,0, 0,0,-.1)"? What is the meaning of the numbers in the Matrix? In the "Checkstress()" part, it changes to "O.cell.velGrad=Matrix3(0,0,.1, 0,0,0, 0,0,0)"

4.What is the meaning of "O.cell.trsf[0,2]" in the "checkDistorsion" part? What do the numbers in the "trsf" mean?

5.In "plot.addData(exz=O.cell.trsf[0,2],szz=stress[2,2],sxz=stress[0,2],tanPhi=(stress[0,2]/stress[2,2]) if stress[2,2]!=0 else 0,i=O.iter)" part,What is the meaning of the number in the "streess[]"? For example, "stress[2,2]" or "stess[0,2]"

Thank you very much for your kindly help!

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

Hello,

first of all, I am not sure if periodic simulations are a good very first starting point to learn Yade..
The previous scripts are OK?

> 1.What is the meaning of "O.cell.hSize=Matrix3(2,0,0, 0,2,0, 0,0,2)"?

it assign Matrix3(2,0,0, 0,2,0, 0,0,2) object to O.cell [1] as its hSize [2] attribute

> What is the meaning of the numbers in the Matrix?

"Base cell vectors (columns of the matrix)" [2].
In this case, diagonal terms are "normal size" (i.e. 2x2x2), zero non-diagonal terms means no skew and no rotation
(non-diagonal terms control skew and/or rotation)

> 2."pack.regularHexa(pack.inAlignedBox((0,0,0),(2,2,2)),radius=.1,gap=0,color=(0,0,1))". What is the meaning of the "color=(0,0,1)"? In some other scripts I also see "color=(1,1,1)"

color of resulting particles

(in ideal case, the documentation should tell that **kw are keywords passed to utils.sphere function. You are welcome to improve the documentation! :-)

> 3.What is the meaning of "O.cell.velGrad=Matrix3(-.1,0,0, 0,-.1,0, 0,0,-.1)"?

you assign Matrix3(-.1,0,0, 0,-.1,0, 0,0,-.1) object to O.cell [1] as its hSize [3] attribute.
The components are components of velocity gradient (in matrix form).
In this case, only diagonal negative terms it is isotropic compaction.

> What is the meaning of the numbers in the Matrix?

components of velocity gradient [4]
The indexing of the components is very similar to those of deformation gradient (O.cell.trsf), see below.

> In the "Checkstress()" part, it changes to "O.cell.velGrad=Matrix3(0,0,.1, 0,0,0, 0,0,0)"

now it is changed to shear:
diagonal terms are zero = no change of size
non-diagonal terms are non-zero = some shear or rotation (in this case only one, xz, non-diagonal term is non-zero = shear in xz plane).

> 4.What is the meaning of "O.cell.trsf[0,2]" in the "checkDistorsion" part? What do the numbers in the "trsf" mean?

trsf is deformation gradient [5] of the cell
trsf[0,2] means its xz shear part (see below)

> stress=sum(normalShearStressTensors(),Matrix3.Zero)
> 5.In "plot.addData(exz=O.cell.trsf[0,2],szz=stress[2,2],sxz=stress[0,2],tanPhi=(stress[0,2]/stress[2,2]) if stress[2,2]!=0 else 0,i=O.iter)" part,What is the meaning of the number in the "streess[]"? For example, "stress[2,2]" or "stess[0,2]"

stress here is matrix form of Cauchy stress tensor [6]. Indices are 0=x, 1=y, 2=z
stress[2,2] means stress_zz, i.e. normal z component
stress[0,2] means stress_xz, i.e. shear xz component

The indices corresponds to those of hSize, velGrad and trsf (where it is related to strain/strain-rate counterpart)

cheers
Jan

[1] https://yade-dem.org/doc/yade.wrapper.html#cell
[2] https://yade-dem.org/doc/yade.wrapper.html#yade.wrapper.Cell.hSize
[3] https://yade-dem.org/doc/yade.wrapper.html#yade.wrapper.Cell.velGrad
[4] https://www.continuummechanics.org/velocitygradient.html
[5] https://en.wikipedia.org/wiki/Finite_strain_theory#Deformation_gradient_tensor
[6] https://en.wikipedia.org/wiki/Cauchy_stress_tensor

Lei Hang (h-stone) said : #7

Thank you for your answer! I have some other questions.

1. What is the meaning of the numbers(2,3,4,5) in "hsize[2]","hsize[3]","velocity gradient [4]","deformation gradient [5]"? Are they the order of attributes in this script? Does "hsize" mean the size of cell? Does "diagonal terms" mean the size of cell?

2.
"> In the "Checkstress()" part, it changes to "O.cell.velGrad=Matrix3(0,0,.1, 0,0,0, 0,0,0)"
now it is changed to shear:
diagonal terms are zero = no change of size
non-diagonal terms are non-zero = some shear or rotation (in this case only one, xz, non-diagonal term is non-zero = shear in xz plane)."
Here, how to determine it shears in xz? Why is it shearing rather than rotating? Why is it in the xz plane rather than in other planes?

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

Hello,

> O.cell [1]
> hSize [2]
> "Base cell vectors (columns of the matrix)" [2].
> O.cell [1]
> hSize [3]
> components of velocity gradient [4]
> trsf is deformation gradient [5]
> stress here is matrix form of Cauchy stress tensor [6].

This is the standard approach on this forum to "cite" references summarized at the end of the message (or at the end of a previous message in the thread)
Full links in the text looks messy and also this way you can easily reference one link at multiple places (like O.cell here)

> Does "hsize" mean the size of cell?

no, see [2] :-), where it states that hSize is "Base cell vectors (columns of the matrix)". I.e. vectors of the edges of the cell.

> Does "diagonal terms" mean the size of cell?

In case of diagonal matrix, the diagonal terms represent size. In general case not.

> O.cell.velGrad=Matrix3(0,0,.1, 0,0,0, 0,0,0)
> Here, how to determine it shears in xz?
> Why is it in the xz plane rather than in other planes?

because the only non-zero term is at index (0,2), corresponding to xz (0=x,2=z) shear component, which is shear in xz plane.

> Why is it shearing rather than rotating?

just because it is..
see e.g. [7], section "2D Transformations" (or google "2x2 matrix mapping scale skew rotation"). It is for 2D, which easier to understand and does not differ from 3D much.
velGrad is "roughly" (in "small strain sense") time derivative of deformation gradient (see [4] for full definition and relation).

cheers
Jan

[7] https://www.cs.brandeis.edu/~cs155/Lecture_06.pdf

Lei Hang (h-stone) said : #9

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

Lei Hang (h-stone) said : #10

Hello

Excuse me. I find I don't understand the non-zero term is at index(0,2) in the answer. How to associate the non-zero term with the index? Thanks so much!

"> O.cell.velGrad=Matrix3(0,0,.1, 0,0,0, 0,0,0)
> Here, how to determine it shears in xz?
> Why is it in the xz plane rather than in other planes?

because the only non-zero term is at index (0,2), corresponding to xz (0=x,2=z) shear component, which is shear in xz plane."

Lei Hang (h-stone) said : #11

If the O.cell.velGrad=Matrix3(.1,0,0, 0,0,0, 0,0,0), Which is the shear plane?

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

> How to associate the non-zero term with the index?

The Matrix3 is (xx, xy, xz, yx, yy, yz, zx, zy, zz), "one row" of the form:
xx xy xz
yx yy yz
zx zy zz

> If the O.cell.velGrad=Matrix3(.1,0,0, 0,0,0, 0,0,0), Which is the shear plane?

you have no shear here.
the first (the only non-zero) element is xx = normal x direction, so this correspond to stretching (enlargement) along x axis.

cheers
Jan

Lei Hang (h-stone) said : #13

Thank you very much Jan!

Lei Hang (h-stone) said : #14

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