# i.geom.normal is not normalized (unit) vector from the branch vector (for the case of sphere-sphere contact)

I want to calculate the mechanical anisotropy, the branch vector is needed. However, I find the contact normal obtained by i.geom.normal is not normalized (unit) vector from the branch vector (for the case of sphere-sphere contact), The specific process is as follows:
(2) find all the interactions of particle 12 by "for i in O.bodies[12].intrs(): print (i.isReal,i.id1,i.id2)" and the result is as follows:
(True, 3, 12)
(True, 12, 97)
(True, 12, 211)
(True, 427, 12)
(True, 12, 1033)
(True, 1171, 12)
(True, 1196, 12)
(True, 12, 2736)
(True, 2777, 12)
(True, 12, 3222)
(True, 3453, 12)
(True, 4202, 12)
(True, 12, 7845)
(True, 12, 8119)
(True, 8496, 12)
(True, 10817, 12)
(True, 11570, 12)
(True, 12, 13485)
(True, 12, 15013)
(True, 12, 19529)
(True, 12, 21199)
(True, 12, 22291)
(True, 12, 22449)
(True, 12, 22803)
(True, 12, 23442)
(True, 12, 24412)
(True, 25196, 12)
(True, 12, 29244)
(True, 12, 30387)
(True, 12, 34784)
(True, 12, 35280)
(True, 36427, 12)
(True, 37206, 12)
(True, 12, 40297)
(True, 41898, 12)
(True, 12, 42413)
(True, 12, 44052)
(True, 12, 45495)
(True, 45948, 12)
(True, 12, 47578)
(True, 12, 48322)
(True, 12, 49112)
(True, 12, 49738)
(True, 49923, 12);
(3) then print the branch vector, contact normal vector and contact normal force vector of O.interactions[3,12] by the follow code:
i = O.interactions[3,12]
b1,b2 = [O.bodies[j] for j in (i.id1,i.id2)]
print("branch12",b2.state.pos - b1.state.pos)
print("branch21",b1.state.pos - b2.state.pos)
print("normal",i.geom.normal)
print("force",i.phys.normalForce + i.phys.shearForce)
print("force1",O.forces.f(i.id1))
print("force2",O.forces.f(i.id2))

the result is as follows:

('branch12', Vector3(0.005264805042519041,-0.0011961713201094834,-0.0024246928681093648))
('branch21', Vector3(-0.005264805042519041,0.0011961713201094834,0.0024246928681093648))
('normal', Vector3(-0,-1,-0))
('force', Vector3(0,-0.6724835758663493,-0))
('force1', Vector3(0,28.188489135826057,0))
('force2', Vector3(-0.0003050240398732389,0.00005915155744011741,0.00009299475981276988))

It can be seen that the branch vector, contact normal vector and contact normal force vector are not in the same direction. Why is this?
Thanks!

## Question information

Language:
English Edit question
Status:
Solved
For:
Assignee:
No assignee Edit question
Solved by:
Zhicheng Gao
Solved:
Last query:
 Revision history for this message Jan Stránský (honzik) said on 2021-11-22: #1

Hello,

there might be several reasons, difficult to tell without actual data [1]. For a general case:

A) contact normal vector vs. contact normal force:
> ... contact normal force vector ...
> print("force",i.phys.normalForce + i.phys.shearForce)

you are printing TOTAL force, not normal force.
If i.geom.shearForce is non-zero, then you naturally get non-parallel vector

O.forces.f(i.id1) result is "independent" of constituent interactions

B) branch vector vs. contact normal vector:

i.geom.normal is computed in InteractionLoop.
Body positions are THEN updated by NewtonIntegrator.
So b2.state.pos - b1.state.pos and i.geom.normal (in general) are not parallel.

Cheers
Jan

 Revision history for this message Zhicheng Gao (zhichenggao) said on 2021-11-22: #2

Thanks, Jan

>you are printing TOTAL force, not normal force. If i.geom.shearForce is non-zero, then you naturally get non-parallel vector
I have modified it to print("force",i.phys.normalForce）and I found out that it was my mistake, body 3 is a wall, not a particle. I checked the other contacts and they all met the parallel conditions. The question now is why this phenomenon occurs on the wall.

Best wishes
Gao

 Revision history for this message Jan Stránský (honzik) said on 2021-11-22: #3

> The question now is why this phenomenon occurs on the wall.

just print the wall's state.pos and it should answer you itself :-)
Basically, wall's state.pos has nothing to do with the contact point, which determines geom.normal

cheers
Jan

 Revision history for this message Zhicheng Gao (zhichenggao) said on 2021-11-22: #4

For i = O.interactions[A,B]，is i.phys.normalForce the force exerted on grain A by grain B or grain B by grain A?

 Revision history for this message Jan Stránský (honzik) said on 2021-11-22: #5

a) see documentation [2]
b) have a try yourself
Cheers
Jan

 Revision history for this message Zhicheng Gao (zhichenggao) said on 2021-11-22: #6