Quaternion algebra and geometry

Revision en6, by adamant, 2018-06-01 01:47:45

Hi everyone! As you may already know (if you don't then I advice you to learn it), in 2D geometry it is convenient to use complex numbers to handling points operations and rotations. Now I would like to tell you about similar construction, which allows you to work efficiently with 3D geometry, in particular to use vectors and linear operations on them, calculate many popular operations (like dot or cross products) and maintain rotations in 3D space.

Let — are vectors in 3D space. Following notation from analytical geometry will be used:

  1. Dot product: , where — is the angle between and in their common plane, — length of vector .

  2. Cross product: , where — unit vectors, directed along , and axes. In 3D space cross product is a vector having length and directed orthogonally to their common plane in such way that if you'll watch on vectors from the endpoint of cross product, shortest rotation from to will be counter-clockwise.

So, quaternion — is hypercomplex number, which can be represented as , where — are real numbers and — are imaginary units. We can define multiplication operation on quaternions with following equality . The whole quaternionic multiplication table can be derived from this identity:

Assuming that quaternion multiplication is distributive over addition (), we can multiply the quaternions, "opening parenthesis". Note that quaternionic multiplication will be associative (), but in common case not commutative (there are such that ).

Quaternions can be also represented as , where — is vector of 3D space with unit vectors . Components and are called scalar and vector parts of quaternions. Assume we have quaternions and . Then their multiplication can be represented as . Consider multiplication of two quaternions with zero scalar parts.

.

Note that it can be shortly written as . Thus we have following formula for quaternionic multiplication:

Finally note that since , quaternion can be represented in form , where — are complex numbers. Also note that since , we have identity or, introducing notation : . Thus multiplication of quaternions and can be represented as . Number , introduced earlier, is called conjugated to complex number . So we have another one formula for quaternion multiplication:

Let's show that any non-zero quaternion is invertible, i.e. for each quaternion there is inverse such that . In this case we can divide quaternions by multiplicating with inverse element. Like in complex numers we can introduce for quaternion conjugated quaternion . From that formula we can see that . Thus we can introduce quaternionic norm (which is generalisation of length): and inverse element . Note that . Indeed from multiplication formula we get . From here we can see and .

Now we are ready to learn about the most useful property of quaternions which is handling rotations in 3D space. From now on we will consider vectors to be quaternions with zero scalar part. Let's introduce conjugation operation of vector a by quaternion g the result of which will be vector [1] . This identity is equivalent to the following: . Let , then we can get . Considering scalar and vector parts separately we will get:

Note that because ||ab|| = ||a||·||b|| we can get that norms of vectors and are equal. First equation of the system means that orthogonal projections on axis are equal for and . Set of vectors with same norm and projection on some direction form a circumference around that direction. Thus we see that can be obtained from by its rotation around on some angle. Let's find it. From now on we will say that rotation is clockwise or counter-clockwise depending on what it's like if we watch it from the endpoint of vector around which rotation is happening.

Let's assume that . If it's not true, we can divide by square root of its norm, it will not affect . Now we can assume that , where . Also using first equation and using notation and , we can transform second equation in following way: . (On the other hand we subtracted from both sides, on the other one since we can see that )

Note that and are lying in the plane orthogonal to , thus, this is the rotation plane. Now we should see that is vector rotated around vector clockwise. Thus left part of second equation is , rotated around clockwise on angle , and the right one is , rotated around counter-clockwise on . Thus we see that is rotated around counter-clockwise on .


To sum up: quaternion with unit norm has such property that for each expression specifies vector rotated around on counter-clockwise.

[1] . Scalar part of this product equals to .

And now the implementation. Ideone: #BBOx9H

// To begin, we introduce the notation for the class and its elements, which we will use.
typedef double ftype;
typedef complex<ftype> point;
typedef complex<point> quater;
#define qA real()
#define qB imag()
#define qs qA.real()
#define qx qA.imag()
#define qy qB.real()
#define qz qB.imag()

const ftype pi = acos(-1.);
const ftype eps = 1e-9;

// Multiplication via formulas in complex terms.
quater operator * (quater a, quater b)
{
    return {{a.qA * b.qA - a.qB * conj(b.qB)},
            {a.qA * b.qB + a.qB * conj(b.qA)}};
}

// Conjugated element
quater conj(quater a)
{
    return {conj(a.qA), - a.qB};
}

// Quaternionic norm.
ftype norm(quater a)
{
    return (a * conj(a)).qs;
}

// Length.
ftype abs(quater a)
{
    return sqrt(norm(a));
}

// Dividing via inverse element.
quater operator / (quater a, quater b)
{
    return a * conj(b) / point(norm(b));
}

// Quaternion from vector coordinates.
quater vec(ftype x, ftype y, ftype z)
{
    return {{0, x}, {y, z}};
}

// Basis vectors.
const quater ex = vec(1, 0, 0);
const quater ey = vec(0, 1, 0);
const quater ez = vec(0, 0, 1);

// Vector part of quaternion.
quater vec(quater a)
{
    return a -= a.qs;
}

// Dot product.
ftype dot(quater a, quater b)
{
    return -(a * b).qs;
}

// Cross product.
quater cross(quater a, quater b)
{
    return vec(a * b);
}

// Triple (mixed) product.
ftype mix(quater a, quater b, quater c)
{
    return dot(a, cross(b, c));
}

// Conjugation of vector a by quaternion g.
quater conj(quater a, quater g)
{
    return g * a / g;
}

// Quaternion representing rotation around i on phi counter-clockwise.
quater rotation(quater i, ftype phi)
{
    return point(cos(phi / 2)) + i * point(sin(phi / 2));
}

// Rotate a around i counter-clockwise.
quater rotate(quater a, quater i, ftype phi)
{
    return conj(a, rotation(i, phi));
}

// Check if two vectors are equal.
bool cmp(const quater &a, const quater &b)
{
    return abs(a - b) < eps;
}

// Any vector which is orthogonal to v.
quater ortho(quater v)
{
    if(abs(v.qy) > eps)
        return vec(v.qy, -v.qx, 0);
    else
        return vec(v.qz, 0, -v.qx);
}

// Quaternion representing rotation from a to b via minimal angle.
quater min_rotation(quater a, quater b)
{
    a /= abs(a);
    b /= abs(b);
    if(cmp(a, -b)) // Degenerate case :(
        return rotation(ortho(b), pi);
    else
        return conj(a * (a + b));
}

// Angle of rotation from [-pi; pi] defined by quaternion.
ftype get_angle(quater a)
{
    a /= abs(a);
    return remainder(2 * acos(a.qs), 2 * pi);
}

// Quaternion representing rotation which turns nx in x, ny in y and [nx, ny] in z.
quater basis_rotation(quater nx, quater ny)
{
    nx /= abs(nx);
    ny /= abs(ny);
    quater a = min_rotation(nx, ex);
    ny = conj(ny, a);
    quater b = min_rotation(ny, ey);
    if(cmp(ny, -ey))
        b = rotation(ex, pi);
    return b * a;
}
Tags geometry, quaternions

History

 
 
 
 
Revisions
 
 
  Rev. Lang. By When Δ Comment
en6 English adamant 2018-06-01 01:47:45 264
ru24 Russian adamant 2018-06-01 01:47:21 264
en5 English adamant 2016-09-19 18:17:09 1145
ru23 Russian adamant 2016-09-19 18:13:54 1163
ru22 Russian adamant 2016-09-19 16:52:44 524
en4 English adamant 2016-09-19 16:49:30 550
en3 English adamant 2016-09-19 16:20:39 379
ru21 Russian adamant 2016-09-19 16:12:01 379
ru20 Russian adamant 2016-09-19 01:18:09 29
en2 English adamant 2016-09-19 01:14:50 0 Initial revision for English translation
en1 English adamant 2016-09-19 01:14:40 10492 Initial revision for English translation
ru19 Russian adamant 2016-09-18 17:50:20 307
ru18 Russian adamant 2016-09-18 01:48:24 58
ru17 Russian adamant 2016-09-18 01:36:21 2 Мелкая правка: 'x ab \neq ab$).\n\nКва' -> 'x ab \neq ba$).\n\nКва'
ru16 Russian adamant 2016-09-18 01:31:26 273
ru15 Russian adamant 2016-09-18 01:23:03 0
ru14 Russian adamant 2016-09-18 01:22:59 2368 (опубликовано)
ru13 Russian adamant 2016-09-18 00:41:39 3474
ru12 Russian adamant 2016-08-24 09:09:49 189 Мелкая правка: ' = a^2 + (b, b)$. Таки' -
ru11 Russian adamant 2016-08-23 23:44:16 15 Мелкая правка: ' z_{2} k) =$ $=\relax' -
ru10 Russian adamant 2016-08-23 23:28:55 448 Мелкая правка: '\vec a$.\n\nПодводя ' -hr
ru9 Russian adamant 2016-08-23 23:21:10 875 Мелкая правка: 'ормируем $g$, подели' -
ru8 Russian adamant 2016-08-23 22:50:46 431 Мелкая правка: ' [\vec b, vec i] \si' -> ' [\vec b, \vec i] \si'
ru7 Russian adamant 2016-08-23 22:36:24 690
ru6 Russian adamant 2016-08-23 19:10:36 117
ru5 Russian adamant 2016-08-23 17:31:22 1025 Мелкая правка: '|q|| = q \bar q$ и обра' -
ru4 Russian adamant 2016-08-23 16:53:34 23
ru3 Russian adamant 2016-08-23 16:38:59 358 Мелкая правка: 'н в виде $(a+bi)+(c+' -
ru2 Russian adamant 2016-08-23 16:27:50 1592 Мелкая правка: '_1 a_2)$. Заметим, ч' -
ru1 Russian adamant 2016-08-23 15:22:49 1260 Первая редакция (сохранено в черновиках)