Improved quaternion application

[infusion]

Description

I am the author of Quaternion.js ( https://github.com/infusion/Quaternion.js ) and thought it might be useful for a wider audience to benefit from an improvement I've derived here https://raw.org/proof/vector-rotation-using-quaternions/

In it's core it just makes use of basic algebra to improve vector rotation using quaternions, which is heavily used in Three.js.

[github-actions]

📦 Bundle size

Full ESM build, minified and gzipped.

Filesize dev	Filesize PR	Diff
644.2 kB (159.7 kB)	644.2 kB (159.7 kB)	-30 B

🌳 Bundle size after tree-shaking

Minimal build including a renderer, camera, empty scene, and dependencies.

Filesize dev	Filesize PR	Diff
437.3 kB (105.8 kB)	437.3 kB (105.8 kB)	-30 B

[Mugen87]

Why did you remove these unit tests?

[infusion]

These unit tests are silly. They don't actually proof the correctness of the actual rotation. A non-normalized quaternion that is applied to the rotation algorithm can only serve as a kind of fingerprint for the formulas if the arbitrary numbers that get were put in, generate a certain output. As the formulas have changed, they don't provide the same output for non-normalized quaternions. All other unit-tests that work on actual versors work perfectly fine.

[Mugen87]

Do you mind quickly explaining what kind of improvement the PR introduces?

[infusion]

Certainly! The rotation of a vector using a quaternion was implemented in a simple way, like most textbooks provide. Using some assumptions (like those of a normalized quaternion), these equations can be simplified, resulting in less multiplications and thus faster execution. The proof of correctness was provided with the link and all unit tests (except the fingerprint) of Three.js also pass.

[Mugen87]

@WestLangley What do you think about this PR?

[WestLangley]

This PR produces the wrong answer if the quaternion does not have unit length. The current implementation does not assume the quaternion is normalized.

I would only pursue such a change if there was a compelling reason to do so, and we are willing to accept the consequences of a behavior change.

I should also mention that the algorithm proposed in this PR is used in other engines when the quaternion is assumed to be a unit quaternion. For example, Unity.

[infusion]

Where do you use non-normalized quaternions and where would it be useful in a 3D engine except the made-up test case, I removed? A quaternion only represents a SO(3) group in case of length 1. What consequences do you see here?

I think the reason is similar to why Unity chose to use a more optimized algorithm to apply quaternions: performance. And what your original applyQuaternion() method is calculating is Q * A(v) * Q' with A(v) being the augmented vector of v and Q' the conjugate of Q. I never came across an application where this formula was applied for |Q| != 1, since you would introduce a strange scaling at the same time. So in my opinion it should be save as probably no one is using applyQuaternion() to pass something different than a versor.

[Mugen87]

I should also mention that the algorithm proposed in this PR is used in other engines when the quaternion is assumed to be a unit quaternion. For example, Unity.

Are you referring on this code?

https://github.com/Unity-Technologies/UnityCsReference/blob/7765d52c6cc13c363796ca00437f1a3209943991/Runtime/Export/Math/Quaternion.cs#L96-L117

BabylonJS assumes unit quaternions as well (see comment) although their implementation matches the current three.js code.

https://github.com/BabylonJS/Babylon.js/blob/24a61865ef7051b297be0bfe6bcaf6bb133fa090/packages/dev/core/src/Maths/math.vector.ts#L1240

Right now, I tend to towards a more optimized solution since assuming unit-quaternions seems to be a common thing in other 3D projects. Unity states in its documentation page:

Note that Unity expects Quaternions to be normalized.

We could add something similar to the three.js docs.

[infusion]

I did not check/verify the Unity code. But compared to my derivation, Unity has also quite a few more multiplications.

[WestLangley]

This PR produces the wrong answer if the quaternion does not have unit length. The current implementation does not assume the quaternion is normalized.

Where do you use non-normalized quaternions

three.js does not normalize user-set quaternions to the precision of the machine -- although the user-set length may be approximately one.

This PR would be correct if three.js called normalize() on all quaternions.

//

It is fine to improve the performance of a method, but you are also returning a different result.

[Mugen87]

This PR would be correct if three.js called normalize() on all quaternions.

I'm afraid doing this more or less equalizes the performance benefit from rewriting the method. The idea is to just assume unit quaternions.

I don't think the mentioned precision is necessary in context of 3D. Quaternions are normally set up with library methods like decompose() or setFromEuler(). These output unit quaternions and I think it's wasted resources to call normalized() everywhere.

[infusion]

@WestLangley The current implementation also expects the quaternion to be normalized and we also have to live with numerical instability here. As @Mugen87 says, users mostly set up quaternions/versors using setFromEuler() or from rotation matrices and thus are okay to work in the 3D space.

As I said I never came across an application of the formula used for vector rotation using quaternions that were not normalized. And since they are expected to be okayish normalized, the mathematical proof I've given shows the formulas are equivalent. So there is no "different" result.

Since all unit tests pass, the new calculation does not introduce an error large enough to be recognized by the current coverage (which is pretty good IMO).

[WestLangley]

The current implementation is exact.

With the current algorithm, the direction of the rotated vector is always correct, even if the quaternion is not perfectly normalized.

With this PR, the direction of the rotated vector is incorrect unless the quaternion length is exactly 1. (The direction varies as the scale changes.)

• If you want to assume every quaternion has length 1 -- which is what we are doing now -- it is better to retain the existing algorithm.

• If you want to modify the library and ensure every quaternion has length 1, then it is OK to use the algorithm in this PR. (I would advise against this option because it involves changing the users data.)

• If you want to make it clear in the docs that a quaternion set by the user must have unit length, then it is OK to use the algorithm in this PR. (Quaternions whose components are computed by three.js are not a problem.)

Edited for clarity.

[mrdoob]

It would also be good to know what the performance benefits are. Have you measured it?

[RemusMar]

If you want to modify the library and ensure every quaternion has length 1, then it is OK to use the algorithm in this PR.

+1

It would also be good to know what the performance benefits are.

I think this is a MUST for this PR

[Mugen87]

FYI: The folks from BabylonJS merged this implementation in their project: BabylonJS/Babylon.js#14075

BTW: When optimizing math code, the number of operations per type is often used as a performance metric. Simply because less basic math operations like multiplications or additions require less computation by the CPU.

[RemusMar]

The folks from BabylonJS merged this implementation in their project: BabylonJS/Babylon.js#14075

Did you notice the Popov72's comment?

I'm not seeing any measurable difference in the two PGs (javascript optimization can be frustrating...), but the new method is not slower and we say in the comment of the method that the quaternion must be a unit quaternion

[Mugen87]

@RaananW says one comment above:

I don't see a big difference in performance, with slight improvement in the new PR.

This is what I would expect in a performance measurement as well (simply because the method performs less computations).

[RemusMar]

Comparing the code changes, the result is obvious: there is a (minor) speed boost.
But is this enough to limit the entire logic to normalized quaternions?

[Mugen87]

I don't think this limitation matters. Like mentioned above you essentially always work with unit-quaternions. It doesn't matter if they are perfectly normalized since minor floating point inaccuracies do not visibly influence the result.

Let's change the discussion into a different direction. I think this PR can only be blocked if someone shares a live example that demonstrates how orienting a vector with the new routine fails. Then the limitation is visible for everyone (including myself) and can easily be understood. Right now, I have the feeling the objections against this change are too theoretical and from a pure mathematical point of view (which does not necessarily reflect real world applications in 3D). BTW: The quaternion for the test case has to be setup in a common way. Not manually but via euler angles or other math routines like decompose().

I've also seen no other engine so far that uses the same approach like three.js (well, BabylonJS but they have recently changed their code). Something like this is always raises doubts, tbh.

[infusion]

Just to give another argument: GLMatrix also uses this approach: https://glmatrix.net/docs/quat.js.html#line652

[WestLangley]

What other engines do is not relevant. Other engines have different limitations and use cases than three.js.

//

The issue has nothing to do with quaternions computed by three.js, the issue pertains to quaternions set by the user.

Object3D.quaternion represents a rotation and must be of unit length.

However, the math library is used not only by three.js internally, but also by users in their apps.

Users can instantiate a quaternion in their app, and the quaternion does not have to be of unit length.

Neither the quaternion constructor, nor Quaternion.set() require the quaternion to be of unit length.

//

As an example, consider an app that uses three.js to demonstrate the properties of quaternions.

The app could do something like this:

// arbitrary quaternion
const q = new THREE.Quaternion( 2, 0, 0, 0 );

// current - correct
new THREE.Vector3( 0, 0, 1 ).applyQuaternion( q ); // { x: 0, y: 0, z: - 4 }

// this PR - wrong
new THREE.Vector3( 0, 0, 1 ).applyQuaternion( q ); // { x: 0, y: 0, z: - 7 }

//

As I said in my first post, this PR involves a behavior change.

Which of the following two options is the most appropriate?

• make it clear in the docs that a quaternion set by the user must represent a pure rotation (have unit length).

• ensure that the Quaternion class works correctly for any quaternion (so the class is not limited).

I am open to either.

[Mugen87]

make it clear in the docs that a quaternion set by the user must represent a pure rotation (have unit length).

I vote for this option.

[infusion]

Any update on this or anything we can do to bring this into three.js?

[Mugen87]

make it clear in the docs that a quaternion set by the user must represent a pure rotation (have unit length).

Updating the documentation page of Quaternion like suggested should be done in this PR.

[zeux]

Just a quick note on the suggested application method: as far as I know it has been a de-facto standard in game industry for a couple decades. (the Unity3D link above is very different - it's just a quat->matrix conversion followed by matrix*vector multiplication, so I guess they didn't get the memo)

I think the original formulation is attributed to NVidia SDK (not sure who the author would be); you can find a copy of this on the internet in a few places and it runs like so (I unfortunately don't remember which NVidia SDK this is referring to and whether it's archived anywhere...):

Vector3D<T> Quaternion<T>::operator * (Vector3D<T> V)
{
	Vector3D<T> uv, uuv;
	uv = Axis.Cross(V);
	uuv = Axis.Cross(uv);
	uv *= (2.0f * Angle);
	uuv *= 2.0f;
	return V + uv + uuv;
}

You would note that this is slightly less efficient than the proposed implementation, because it misses the opportunity to move the 2* multiplication to the Cross. Also, Axis and Angle in this snippet are poorly named and refer to q.xyz and q.w. The implementation above can be found in Ogre (https://github.com/OGRECave/ogre/blob/master/OgreMain/src/OgreQuaternion.cpp#L359-L371), Irrlicht 3D (https://github.com/zaki/irrlicht/blob/master/include/quaternion.h#L690-L702), Godot (https://github.com/godotengine/godot/blob/master/core/math/quaternion.h#L91-L98), Urho3D (https://github.com/urho3d/urho3d/blob/master/Source/Urho3D/Math/Quaternion.h#L258-L288), CryEngine (https://github.com/MergHQ/CRYENGINE/blob/release/Code/CryEngine/CryCommon/CryMath/Cry_Quat.h#L857-L873), definitely others. Some engines like O3DE & Bullet Physics didn't get the memo as well I think, ah well. Some of the above links use a slightly more efficient formulation that is also suggested here that folds multiplications by 2 together.

Fabian Giesen also wrote about this a few years ago, if you're interested in a separate algebraic derivation:
https://fgiesen.wordpress.com/2019/02/09/rotating-a-single-vector-using-a-quaternion/. His end result is the same as the one proposed here.

Of course, the same assumption applies as it's the same derivation - quaternion must be unit length.

I personally think that quaternions with length that's substantially different from 1 aren't particularly important, however it might be interesting to look at numerical stability for quaternions that are slightly off-unit: given a quaternion with length 1+eps, what happens when you rotate a unit vector many times with either method?

I've created a test that does just that, and looks at the impact of introducing a small (1e-5) error in the original quaternion's length on the length and the direction of the result:
https://gist.github.com/zeux/21dc26af93139028ba1b278596415575

Running this, I get:

0.20428159582151195 -0.5773688016646978 0.8159196421859614 length 1.0202012380074728
0.20117753205335792 -0.5670648091976522 0.7988224679466441 length 1.000076113368286
0.20023656922873762 -0.5659361899936409 0.7997634307713477 length 1.0000000000000848
0.20023656922872274 -0.565936189993595 0.7997634307712778 deviation 1
0.20116222091915187 -0.5670216512698829 0.798761671505368 deviation 0.9999984807105694
0.20023656922872063 -0.565936189993593 0.7997634307712799 deviation 1

(these are, in order, repeated application of quaternion with an injected error using current code, repeated application using new code, repeated application of a unit quaternion for reference, and then the same 3 vectors normalized, with deviation as a dot-product compared to reference unit quaternion transformation)

From this we can see that the length gets distorted much more severely by the current method three.js is using; the "new" method is more forgiving. The old method might be a little bit more precise wrt preserving the vector's direction, but with the cosine error only reaching 2e-6 after 1000 transformations this probably is not a problem.

I would value slightly improved performance and improved length preservation due to floating-point numerical issues more highly than preserving the old behavior for quaternions with a large length; non-unit quaternions can be useful in context of algorithms like dual-quaternion skinning with scaling, but that requires custom math anyhow. I do agree that performance improvement is going to be very difficult to assess in context of JS though.

[aras-p]

the Unity3D link above is very different - it's just a quat->matrix conversion followed by matrix*vector multiplication, so I guess they didn't get the memo

FWIW the "new" Unity mathematics library (Unity.Mathematics package) is also doing the proposed method (link):

public static float3 mul(quaternion q, float3 v)
{
    float3 t = 2 * cross(q.value.xyz, v);
    return v + q.value.w * t + cross(q.value.xyz, t);
}

[Mugen87]

@WestLangley Are you okay if we merge the PR? I'll update the documentation with a note that quaternions set by users must have unit length. I'll also add a note to the migration guide.

[WestLangley]

@WestLangley Are you okay if we merge the PR? I'll update the documentation with a note that quaternions set by users must have unit length. I'll also add a note to the migration guide.

Yes.