calculation of stress tensor and fabric tensor with rigid boundary

Asked by ceguo

Hi yade-developers/users,

To calculate macro stress tensor from micro quantities, we usually the formula sigma_ij = 1/V sum(L_i x F_j), where L is the branch vector and F is the contact force and the summation is over all contacts. The branch vector is the connection between two contacting particles. My question is how about the boundary contact when rigid boundaries are used. I think these contacts should be considered (please correct me if I'm wrong). So these boundary branch vectors should be connections between particle center and the contacting points, right?
Then regarding to the calculation of fabric tensor, phi_ij = 1/N_c sum(n_i x n_j), where N_c is the total contact number, n is the unit normal of the contact plane. So with rigid boundary, each inter-particle contact should be counted twice while the boundary contact should only be considered once. Is my understanding correct? Or all the contacts should be treated in the same manner?
Also in the calculation of branch vector tensor, do we have to treat the inter-particle branch vector and the boundary branch vector differently?

I'm quite confused about these calculations. Hope anyone can help.

Best regards!
Ning

Question information

Language:
English Edit question
Status:
Expired
For:
Yade Edit question
Assignee:
No assignee Edit question
Last query:
Last reply:
Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#1

1/ Yes it works with branch = particle center -> contact point.
2/ If the fabric statistics imply a division by N_c, then why would some
contacts be counted twice? I don't get your idea.
3/ There is no difference between sphere-sphere and sphere-boundary
branches for computing stress.

HTH

Bruno

Revision history for this message
ceguo (hhh-guo) said :
#2

Hi Bruno,

Thanks for your reply! By saying count twice I mean for one physical contact N_c=2, and the counted normals are n and (-n). I have this consideration because I think in this way: for each particle-particle contact, it supports TWO particles; while for each particle-wall contact, it supports only ONE particle. For example, when we consider coordination number, Z=N_c/N_p, where N_c is the contact number (twice the number of physical particle-particle contacts) and N_p is the particle number.

I'm not sure my understanding is correct. So please correct me if I'm wrong.

Best regards!
Ning

Revision history for this message
Bruno Chareyre (bruno-chareyre) said :
#3

I see. Well, I would consider the wall as a particle, hence no difference in the counting, but it's a matter of convention after all.

Revision history for this message
ceguo (hhh-guo) said :
#4

So you mean we don't need to discriminate the two kinds of contacts, i.e. each physical contact, no matter particle-particle or particle-wall, should both be counted once (or both twice).

I want to make clear because the two different treatments (with or without discrimination) will generate different results for fabric tensor or coordination number. There is no problem for stress tensor calculation.

Revision history for this message
Chiara Modenese (chiara-modenese) said :
#5

On 26 September 2011 10:50, ceguo <email address hidden>wrote:

> Question #172345 on Yade changed:
> https://answers.launchpad.net/yade/+question/172345
>
> Status: Answered => Open
>
> ceguo is still having a problem:
> So you mean we don't need to discriminate the two kinds of contacts,
> i.e. each physical contact, no matter particle-particle or particle-
> wall, should both be counted once (or both twice).
>
It is not of a choice to count contacts twice or once. When you compute
coordination number you should count them twice. Eventually you can avoid to
consider rattles in the calculation (see Thornton 2000). Have a look at
utils.avgNumInteractions(). Note that from the total number of contacts only
1*N (not twice) of rattles (particles with one contact only) is subtracted.
I do not think there should be any distinction for wall-sphere contacts or
at least I do not see any reason for that.

Chiara

> I want to make clear because the two different treatments (with or
> without discrimination) will generate different results for fabric
> tensor or coordination number. There is no problem for stress tensor
> calculation.
>

> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Revision history for this message
Chiara Modenese (chiara-modenese) said :
#6

On 26 September 2011 10:50, ceguo <email address hidden>wrote:

> Question #172345 on Yade changed:
> https://answers.launchpad.net/yade/+question/172345
>
> Status: Answered => Open
>
> ceguo is still having a problem:
> So you mean we don't need to discriminate the two kinds of contacts,
> i.e. each physical contact, no matter particle-particle or particle-
> wall, should both be counted once (or both twice).
>
> I want to make clear because the two different treatments (with or
> without discrimination) will generate different results for fabric
> tensor or coordination number. There is no problem for stress tensor
> calculation.
>
As to the fabric tensor, there is only one possible definition in the code.
However you may want to implement other definitions. Again, I do not think
you should distinguish between the types of contact.

> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Revision history for this message
ceguo (hhh-guo) said :
#7

Hi Bruno and Chiara,

Here is a small script testing my concern:

from yade import pack

def numInteract():
   print O.interactions.countReal()

O.materials.append(FrictMat(young=150e6,poisson=.3,frictionAngle=.4,density=2600))
O.bodies.append(
   [utils.sphere(center=[0,0,0],radius=0.5),
   utils.sphere(center=[1,0,0],radius=0.5),
   utils.sphere(center=[0,1,0],radius=0.5),
   utils.sphere(center=[1,1,0],radius=0.5),
   utils.sphere(center=[0,0,1],radius=0.5),
   utils.sphere(center=[1,0,1],radius=0.5),
   utils.sphere(center=[0,1,1],radius=0.5),
   utils.sphere(center=[1,1,1],radius=0.5)]
)
walls = utils.aabbWalls(thickness=0.1)
O.bodies.append(walls)

O.dt=.5*utils.PWaveTimeStep()
O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()],sweepLength=.025),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
      [Ip2_FrictMat_FrictMat_FrictPhys()],
      [Law2_ScGeom_FrictPhys_CundallStrack()]
   ),
   NewtonIntegrator(),
   TriaxialCompressionEngine(
      internalCompaction=False,
      sigmaIsoCompaction=100e3,
      max_vel=10
   ),
   PyRunner(command='numInteract()',iterPeriod=10)
]
O.run(100); O.wait()

For this special lattice (regularly aligned eight particles), the number of particle-particle contacts is 12 and the particle-wall contacts is 8. So I suppose the coordination number is (12x2+8)/8=4 because we can see there are exactly 4 contacts for each particle. This means each p-p contact should be counted twice while each p-w contact once. I'm confused why the countReal() gives me a value of 9 and avgNumInteractions() gives me 1.2857 in the test?

I'm sorry I'm new to YADE and may code the script wrongly. Sorry for the trouble.

Best regards!
Ning

Revision history for this message
Chiara Modenese (chiara-modenese) said :
#8

By running your script I get:
O.interactions.countReal() = 29 # this is the total number of interactions
which are real
utils.avgNumInteractions() = 4.1428571428571432 # this is the coordination
number (so called apparent or average CN)

Did you find a problem while using the above functions?

Chiara

On 26 September 2011 13:20, ceguo <email address hidden>wrote:

> Question #172345 on Yade changed:
> https://answers.launchpad.net/yade/+question/172345
>
> Status: Answered => Open
>
> ceguo is still having a problem:
> Hi Bruno and Chiara,
>
> Here is a small script testing my concern:
>
> from yade import pack
>
> def numInteract():
> print O.interactions.countReal()
>
>
> O.materials.append(FrictMat(young=150e6,poisson=.3,frictionAngle=.4,density=2600))
> O.bodies.append(
> [utils.sphere(center=[0,0,0],radius=0.5),
> utils.sphere(center=[1,0,0],radius=0.5),
> utils.sphere(center=[0,1,0],radius=0.5),
> utils.sphere(center=[1,1,0],radius=0.5),
> utils.sphere(center=[0,0,1],radius=0.5),
> utils.sphere(center=[1,0,1],radius=0.5),
> utils.sphere(center=[0,1,1],radius=0.5),
> utils.sphere(center=[1,1,1],radius=0.5)]
> )
> walls = utils.aabbWalls(thickness=0.1)
> O.bodies.append(walls)
>
> O.dt=.5*utils.PWaveTimeStep()
> O.engines=[
> ForceResetter(),
>
> InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()],sweepLength=.025),
> InteractionLoop(
> [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
> [Ip2_FrictMat_FrictMat_FrictPhys()],
> [Law2_ScGeom_FrictPhys_CundallStrack()]
> ),
> NewtonIntegrator(),
> TriaxialCompressionEngine(
> internalCompaction=False,
> sigmaIsoCompaction=100e3,
> max_vel=10
> ),
> PyRunner(command='numInteract()',iterPeriod=10)
> ]
> O.run(100); O.wait()
>
> For this special lattice (regularly aligned eight particles), the number
> of particle-particle contacts is 12 and the particle-wall contacts is 8.
> So I suppose the coordination number is (12x2+8)/8=4 because we can see
> there are exactly 4 contacts for each particle. This means each p-p
> contact should be counted twice while each p-w contact once. I'm
> confused why the countReal() gives me a value of 9 and
> avgNumInteractions() gives me 1.2857 in the test?
>
> I'm sorry I'm new to YADE and may code the script wrongly. Sorry for the
> trouble.
>
> Best regards!
> Ning
>
> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Revision history for this message
ceguo (hhh-guo) said :
#9

Hi Chiara,

What I paste here is exactly what I used (only I add print utils.avgNumInteractions()). The output is:
Welcome to Yade 0.60.3
TCP python prompt on localhost:9000, auth cookie `suakcs'
XMLRPC info provider on http://localhost:21000
Running script test.py
9 1.28571428571
9 1.28571428571
9 1.28571428571
9 1.28571428571
8 1.14285714286
9 1.28571428571
9 1.28571428571
9 1.28571428571
9 1.28571428571
So what is the cause of the problem? Is there any conflict between yade and esys-particle? I have the latter installed.
Besides, as you got countReal=29 and particle number is 8, why CN=4.143. I think the formula should be CN=2xC/N.

Best regards!
Ning

Revision history for this message
Chiara Modenese (chiara-modenese) said :
#10

Let's avoid errors from numerical precision. Do like this and let me know
what you get:
1) Increase the radius of the particles to 0.51
2) print [(i.id1, i.id2) for i in O.interactions]

Let me know, thanks.
Chiara

[(i.id1, i.id2) for i in O.interactions]

On 27 September 2011 03:05, ceguo <email address hidden>wrote:

> Question #172345 on Yade changed:
> https://answers.launchpad.net/yade/+question/172345
>
> Status: Answered => Open
>
> ceguo is still having a problem:
> Hi Chiara,
>
> What I paste here is exactly what I used (only I add print
> utils.avgNumInteractions()). The output is:
> Welcome to Yade 0.60.3
> TCP python prompt on localhost:9000, auth cookie `suakcs'
> XMLRPC info provider on http://localhost:21000
> Running script test.py
> 9 1.28571428571
> 9 1.28571428571
> 9 1.28571428571
> 9 1.28571428571
> 8 1.14285714286
> 9 1.28571428571
> 9 1.28571428571
> 9 1.28571428571
> 9 1.28571428571
> So what is the cause of the problem? Is there any conflict between yade and
> esys-particle? I have the latter installed.
> Besides, as you got countReal=29 and particle number is 8, why CN=4.143. I
> think the formula should be CN=2xC/N.
>
> Best regards!
> Ning
>
> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Revision history for this message
ceguo (hhh-guo) said :
#11

Hi Chiara,

After doing (1) without (2), I got:
35 5.0
35 5.0
36 5.14285714286
36 5.14285714286
36 5.14285714286
36 5.14285714286
36 5.14285714286
36 5.14285714286
36 5.14285714286
And doing (2) without (1), I got:
[(13, 4), (12, 0), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6)]
[(13, 4), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (12, 0)]
[(13, 4), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (12, 0)]
[(13, 4), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (12, 0)]
[(2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4)]
[(12, 0), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4)]
[(2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4), (4, 0)]
[(2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4), (4, 0)]
[(12, 0), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4)]

Is this due to precision error?
Ning

Revision history for this message
Chiara Modenese (chiara-modenese) said :
#12

Your previous results were very different from mine, that was not due to
numerical precision but I could not tell what it was (I also did compare the
output of your original script with another Yade user and more or less
results were the same although not exactly). Let's say that this is not the
best example to test but I get your same result now. Also to answer your
question, walls are treated like spheres in the calculation of CN.

Chiara

On 27 September 2011 08:41, ceguo <email address hidden>wrote:

> Question #172345 on Yade changed:
> https://answers.launchpad.net/yade/+question/172345
>
> Status: Answered => Open
>
> ceguo is still having a problem:
> Hi Chiara,
>
> After doing (1) without (2), I got:
> 35 5.0
> 35 5.0
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
> And doing (2) without (1), I got:
> [(13, 4), (12, 0), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6)]
> [(13, 4), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (12, 0)]
> [(13, 4), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (12, 0)]
> [(13, 4), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (12, 0)]
> [(2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4)]
> [(12, 0), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4)]
> [(2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4), (4, 0)]
> [(2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4), (4, 0)]
> [(12, 0), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4)]
>
> Is this due to precision error?
> Ning
>
> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Revision history for this message
Chiara Modenese (chiara-modenese) said :
#13

On 27 September 2011 08:41, ceguo <email address hidden>wrote:

> Question #172345 on Yade changed:
> https://answers.launchpad.net/yade/+question/172345
>
> Status: Answered => Open
>
> ceguo is still having a problem:
> Hi Chiara,
>
> After doing (1) without (2), I got:
> 35 5.0
> 35 5.0
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
> 36 5.14285714286
>
Sorry can you also print here the output of (2) when you increase the
radius? Thanks.

> And doing (2) without (1), I got:
> [(13, 4), (12, 0), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6)]
> [(13, 4), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (12, 0)]
> [(13, 4), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (12, 0)]
> [(13, 4), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (12, 0)]
> [(2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4)]
> [(12, 0), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4)]
> [(2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4), (4, 0)]
> [(2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4), (4, 0)]
> [(12, 0), (2, 6), (2, 3), (5, 7), (5, 1), (7, 3), (3, 1), (7, 6), (13, 4)]
>
> Is this due to precision error?
> Ning
>
> --
> You received this question notification because you are a member of
> yade-users, which is an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>

Revision history for this message
ceguo (hhh-guo) said :
#14

Hi Chiara,

Here is the output with radii increased:
[(13, 4), (13, 6), (13, 5), (13, 7), (12, 2), (12, 3), (12, 1), (11, 2), (11, 6), (11, 7), (11, 3), (10, 0), (10, 4), (10, 5), (10, 1), (8, 0), (8, 2), (0, 2), (0, 4), (0, 1), (2, 6), (2, 3), (4, 6), (4, 5), (6, 7), (5, 7), (5, 1), (7, 3), (3, 1), (9, 1), (8, 4), (8, 6), (9, 3), (9, 7), (9, 5)]
[(13, 4), (13, 6), (13, 5), (13, 7), (12, 0), (12, 2), (12, 3), (12, 1), (11, 2), (11, 6), (11, 7), (11, 3), (10, 0), (10, 4), (10, 5), (10, 1), (8, 0), (8, 2), (0, 2), (0, 1), (2, 6), (2, 3), (4, 6), (4, 5), (6, 7), (5, 7), (5, 1), (7, 3), (3, 1), (9, 1), (8, 4), (8, 6), (9, 3), (9, 7), (9, 5)]
[(13, 4), (13, 6), (13, 5), (13, 7), (12, 0), (12, 2), (12, 3), (12, 1), (11, 2), (11, 6), (11, 7), (11, 3), (10, 0), (10, 4), (10, 5), (10, 1), (8, 0), (8, 2), (0, 2), (0, 4), (0, 1), (2, 6), (2, 3), (4, 6), (4, 5), (6, 7), (5, 7), (5, 1), (7, 3), (3, 1), (9, 1), (8, 4), (8, 6), (9, 3), (9, 7), (9, 5)]
[(13, 4), (13, 6), (13, 5), (13, 7), (12, 0), (12, 2), (12, 3), (12, 1), (11, 2), (11, 6), (11, 7), (11, 3), (10, 0), (10, 4), (10, 5), (10, 1), (8, 0), (8, 2), (0, 2), (0, 4), (0, 1), (2, 6), (2, 3), (4, 6), (4, 5), (6, 7), (5, 7), (5, 1), (7, 3), (3, 1), (9, 1), (8, 4), (8, 6), (9, 3), (9, 7), (9, 5)]
[(13, 4), (13, 6), (13, 5), (13, 7), (12, 0), (12, 2), (12, 3), (12, 1), (11, 2), (11, 6), (11, 7), (11, 3), (10, 0), (10, 4), (10, 5), (10, 1), (8, 0), (8, 2), (0, 2), (0, 4), (0, 1), (2, 6), (2, 3), (4, 6), (4, 5), (6, 7), (5, 7), (5, 1), (7, 3), (3, 1), (9, 1), (8, 4), (8, 6), (9, 3), (9, 7), (9, 5)]
[(13, 4), (13, 6), (13, 5), (13, 7), (12, 0), (12, 2), (12, 3), (12, 1), (11, 2), (11, 6), (11, 7), (11, 3), (10, 0), (10, 4), (10, 5), (10, 1), (8, 0), (8, 2), (0, 2), (0, 4), (0, 1), (2, 6), (2, 3), (4, 6), (4, 5), (6, 7), (5, 7), (5, 1), (7, 3), (3, 1), (9, 1), (8, 4), (8, 6), (9, 3), (9, 7), (9, 5)]
[(13, 4), (13, 6), (13, 5), (13, 7), (12, 0), (12, 2), (12, 3), (12, 1), (11, 2), (11, 6), (11, 7), (11, 3), (10, 0), (10, 4), (10, 5), (10, 1), (8, 0), (8, 2), (0, 2), (0, 4), (0, 1), (2, 6), (2, 3), (4, 6), (4, 5), (6, 7), (5, 7), (5, 1), (7, 3), (3, 1), (9, 1), (8, 4), (8, 6), (9, 3), (9, 7), (9, 5)]
[(13, 4), (13, 6), (13, 5), (13, 7), (12, 0), (12, 2), (12, 3), (12, 1), (11, 2), (11, 6), (11, 7), (11, 3), (10, 0), (10, 4), (10, 5), (10, 1), (8, 0), (8, 2), (0, 2), (0, 4), (0, 1), (2, 6), (2, 3), (4, 6), (4, 5), (6, 7), (5, 7), (5, 1), (7, 3), (3, 1), (9, 1), (8, 4), (8, 6), (9, 3), (9, 7), (9, 5)]
[(13, 4), (13, 6), (13, 5), (13, 7), (12, 0), (12, 2), (12, 3), (12, 1), (11, 2), (11, 6), (11, 7), (11, 3), (10, 0), (10, 4), (10, 5), (10, 1), (8, 0), (8, 2), (0, 2), (0, 4), (0, 1), (2, 6), (2, 3), (4, 6), (4, 5), (6, 7), (5, 7), (5, 1), (7, 3), (3, 1), (9, 1), (8, 4), (8, 6), (9, 3), (9, 7), (9, 5)]

Besides, I don't like the idea to treat walls as particles because the walls (MASSLESS) don't have equilibrium problems. The equilibrium is required only for particles. We can analyze the average contact number for each particle using regular lattice. In this example with 8 particles (each has 6 contacts, 3 for p-p and 3 for p-w), and we can calculate CN=(12x2+24x1)/8=6, because we have 12 p-p contacts which should be counted twice and 24 p-w contacts which should be counted only once. We can also use a much simpler case: only one particle circumferenced by four walls, so there are 4 p-w contacts, CN=4x1/1=4. The model can be easily extended.

What do you think this kind of idea?

Ning

Revision history for this message
Launchpad Janitor (janitor) said :
#15

This question was expired because it remained in the 'Open' state without activity for the last 15 days.