About periodic simple shear

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

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!

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

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

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!

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

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?

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

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

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."

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

Thank you very much Jan!

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

Hi Jan Stránský!

I'm hoping to hi-jack this question thread and expand upon the amazing discussion you and Lei Hang have had above. Seriously, the above conversation was an amazing review of the opensource code.

To start off my questions, why does the code first simulate a constant isotropic compression of O.cell.velGrad=Matrix3(-.1,0,0, 0,-.1,0, 0,0,-.1)?

Does this step ultimately cause the periodic boundary condition to be flush with the particle regularHexa pack? Because without it, the periodic boundary condition seems to simulate a regularHexa pack within the representative volume element but it is suspended inside of it at iteration 0.

I am asking because I am interested in a periodic simple shear simulation as well ... I just cannot connect the dots on why this step is logical. Is isotropic compression necessary prior to the simple shear of O.cell.velGrad=Matrix3(0,0,.1, 0,0,0, 0,0,0)?

Another question I have is about the section:
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)
 # color particles based on rotation amount

The outputted data takes on the form (via columns):
# exz i sxz szz tanPhi

Based on the conversation above, I recognize that sxz = shear stress in the xz plane and szz = normal stress in the zz direction.
I am confused about exz, i, and tanPhi:

To my understanding, exz represents the current transformation matrix of the cell at the given iteration step. Would this definition be analogous to the cell strain or deformation at a given time step?

I am fairly lost on what the data column 'i' could represent - so no guesses there.

Finally, I am curious about tanPhi. For my simulations I have the particle rotations blocked via:
b.state.blockedDOFs='XYZ'
Why is a change in angle being outputted when the rotational DOFs are being blocked? Does tanPhi represent something else here?

Thanks for all of the help Jan. Your assistance has and always is, greatly appreciated.

Cheers,
Kieran

Hello,

> why does the code first simulate a constant isotropic compression of O.cell.velGrad=Matrix3(-.1,0,0, 0,-.1,0, 0,0,-.1)?
> Is isotropic compression necessary prior to the simple shear of O.cell.velGrad=Matrix3(0,0,.1, 0,0,0, 0,0,0)?

Because the author decided it :-)
Probably because it corresponds to simulated real situation.
Definitely it is not necessary.

> Does this step ultimately cause the periodic boundary condition to be flush with the particle regularHexa pack? Because without it, the periodic boundary condition seems to simulate a regularHexa pack within the representative volume element but it is suspended inside of it at iteration 0.

Sorry, here I am lost...
Could you please reword the question and idea behind?

> exz ... Would this definition be analogous to the cell strain or deformation at a given time step?

In this very specific case yes, exz=O.cell.trsf[0,2] is shear strain (in small strain sense)
BUT, in general, be careful with mapping O.cell.trsf to strain (!)

> I am fairly lost on what the data column 'i' could represent

do you mean i=O.iter?
If yes, then "i" is iteration number

> Finally, I am curious about tanPhi.
> Why is a change in angle being outputted when the rotational DOFs are being blocked? Does tanPhi represent something else here?

Even you have rotational DOFs of particle blocked, the cell as a whole deforms (after all you apply some shear).

Anyway, here,
tanPhi=(stress[0,2]/stress[2,2]) = shear stress / normal stress
so it is not a geometric angle, but some stress measure.
To me it seems like reformulated cohesion-less (c=0) Mohr-Coulomb condition
shear stress = normal stress * tan(phi)
tan(phi) = shear stress / normal stress
so you basically output tangent of friction angle.

cheers
Jan

[8] https://en.wikipedia.org/wiki/Lode_coordinates

Hi Jan,

Thank you for the explanations! Let's see how I can reword my question above.

> Does this step ultimately cause the periodic boundary condition to be flush with the particle regularHexa pack? Because without it, the periodic boundary condition seems to simulate a regularHexa pack within the representative volume element but it is suspended inside of it at iteration 0.

During the compression portion of the simulation, I observe how the representative volume element (RVE) shrinks inward to the particle clump / regularHexa pack. Once this is complete, or the limitMeanStress is reached, the xz-plane shear begins. So back to when I asked above if this step is necessary for the PBC to be flush with the clump, I really was asking if we need to isotopically compress the RVE in order to shear the cube.

It sounds like it doesn't need to occur - which is promising to hear. I was also able to start the RVE shear simulation by setting the limitMeanStress=0 so it starts right away.

The compression stage is not necessary, of course.
It depends on what you want to simulate. The original script probably simulated some experiment with prestress and one outcome is tanPhi = shear stress / normal stress (which is meaningless for zero normal stress).

If you need to model only shear, then model only shear.

For some cohesive materials, also pre-tension + shear could be applied :-)

It is really up to you.

cheers
Jan