← Back to team overview

yade-dev team mailing list archive

Re: computation of MomentBending

 

Luc Sibille said:     (by the date of Mon, 25 Oct 2010 19:51:42 +0200)

> Hi all,
> 
> I am still investigating the different moment-rotation contact laws in 
> Yade (version 0.50 and last update version).
> 
> I have two interrogations:
> 
> 1/ First in every case relative particle rotations is totaly recomputed 
> at each time step from the initial orientation with a code like that:
> 
> Quaternionr delta((rbp1.ori * (initialOrientation1.conjugate())) * 
> (initialOrientation2 * (rbp2.ori.conjugate())));
> 
> is it safe to make such a computation, where after a long time 
> simulation, relative rotation can become very large more than 2*pi 
> whereas there is after a correction for the computed relative rotation 
> if it is greater than 2*pi.

Remember that this is a relative rotation. Not rotation of both
spheres. I was doing tests, like rotating whole beam (a line of
spheres) around various axes, and the beam still behaved properly,
even though it was doing hundreds of 2*pi rotations. I mean - it was
working, because whole sample was rotating, so relative rotations
between spheres remained very small.

If you have a sample, where in a single contact, one of the spheres "rotates
away", then indeed troubles are possible. That needs checking.

> In addition:
> > Just to note that this formula will "fail" if the relative rotation of a
> > particle will be more than 180 degrees (the quaternion will take it the
> > other way around, as it goes always the shortest way). 
> 
> 
> I imagined a priori to compute this relative rotation incrementally (as 
> the shear force), to avoid to deal with great relative rotation. What is 
> your opinion? 

> why did you choose this formula and not an incremental 
> formula?

It is simpler. Just calculate the rotation difference.

If you want to track 2 pi rotations, you can store previously
calculated angle, and compare it with a new one. When you spot a 2*pi
jump in values, then increment some integer counter of 'total
rotations'. In fact I am doing that in lattice model, so it works if I
spin a rod around another rod several times, and then I can see it
unspinning back the exact same amount of full 2*pi rotations.


I recall that Vaclav has also implemented a 'total distance, no
increments' approach on shear force.

 
> 2/ Then I rediscovered that aa.angle() can return nan, as commented in 
> source file, it happens when the quaternion is close to identity. Please 
> apologise for my lack of culture, but what does represent physicaly a 
> quaternion close to identity?

hmm. All quaternions that we use are identity quaternions:
a^2+b^2+c^2+d^2=1. If that condition is not kept, then it is not a
rotation, but just some random quaternions.

To explain by analogy: think about complex numbers: a^2+b^2=1 means
that we are using only complex numbers from a circle of radius=1. If
that condition is not kept, then it is just some random complex
number.

So perhaps you mean, that if a quaternion is not a unit quaternion
then angle()==nan. 
This is correct, 
 - you cannot have an angle if quaternion does not represent rotation.
 - you cannot have an angle if a^2+b^2+c^2+d^2 is different from one.

-- 
Janek Kozicki                               http://janek.kozicki.pl/  |



Follow ups

References