Rotating vectors in space and high precision in C++

Question

This is my function to compute 3D rotation in C++ defined by an angle in radiant around axis.

 Vector rotate(const Vector& axis, const Vector& input, const double angle)    {
    double norm = 1/axis.norm();

    if(norm != 1)
        axis *= norm;

    double cos = std::cos(angle);
    double mcos = 1 - cos;
    double sin = std::sin(angle);

    double r1[3];
    double r2[3];
    double r3[3];

    double t_x, t_ym t_z;


    r1[0] = cos + std::pow(axis.x, 2)*mcos;
    r1[1] = axis.x*axis.y*mcos - axis.z * sin;
    r1[2] = axis.x*axis.z*mcos - axis.y * sin;

    r2[0] = axis.x*axis.y*mcos + axis.z*sin;
    r2[1] = cos + std::pow(axis.y, 2)*mcos;
    r2[2] = axis.x*axis.z*mcos - axis.x * sin;

    r3[0] = axis.x*axis.z*mcos - axis.y * sin;
    r3[1] = axis.z*axis.y*mcos - axis.x * sin;
    r3[2] = cos - std::pow(axis.z, 2) * mcos;

    return Vector(t_x, t_y, t_z);
}

The thing is if you try yourself to rotate the vector a n times pi/4 where n multiple of 4 (so you do complete revolution around the axis by doing four quarter of rotations) the error will propagate pretty fast.

Example (where err = input-output):

input: (1.265, 3.398, 3.333)
rotation axis: (2.33, 0.568, 2.689)

n: 8 (so two completes revolutions)
output: (1.301967, 1.533389, 4.138940)
error: (0.038697, -0.864611, 0.805940)


n: 400 (so 100 completes revolutions)
error: (472..., 166..., 673...)

What can I do ?

Constraints:

• Rotations are not predictable so not possible to do something like angle = pi/4 *n % 2*pi like @molbdnilo suggests. Because I have to chain translations and rotations to test if there's a collision.

Answer 1

This is one of the typical problems that you can run into with floats.

Floating point numbers are pretty exact on singular operations. In fact, for many operations you are guaranteed to get the most exact result that can be represented in the format, so any rounding errors that you get are solely due to fitting it into the representation.

However, as soon as you start chaining floating point operations, even those errors can accumulate. If you are lucky, you can make your algorithm numerically stable, so that the rounding errors cancel each other out in the end and you always stay in the ballpark of the correct result. Getting this right can be quite a challenge though, especially for complex computations. For instance, in your particular implementation, there is lots of potential for catastrophic cancellation introducing large rounding errors into the computation chain.

The easier solution is: Avoid chaining the float operations in the first place! Or to be more precise: Only chain those parts which you can keep numerically stable. Since you mentioned this is for a computer game: In a game you transform the geometry according to the camera matrix each frame. You never touch the geometry in memory, instead you simply adjust the camera matrix. That way, your source geometry is always fresh and the rounding error in each frame is simply the error from that single transformation.

Similarly, you usually don't update the camera matrix incrementally. Instead, you read the player's position and view and build the complete matrix from scratch from those vectors. Now the only challenge that you have left is make sure that you don't accumulate errors into the player position and view, but this is much easier than ensuring stability at the other end of the transformation pipeline.

Answer 2

Just store the original base vector and the rotation angle together, and do the calculation every time you need the current rotated value. You can cache this and invalidate every time the angle changes, but always work from the original base vector and total aggregate rotation.

Presto! No cumulative errors, because no chained calculations.

Also, if you're concerned about cumulative errors in the angle itself, store that in degrees and convert to radians when required. Again, pi is touched once in the degree->radian conversion, and you don't have a chain of approximate pi/n values contributing more errors.