Minimizing total floating point arithmetic error

I have a 2D convex polygon in 3D space and a function to measure the area of ​​a polygon.

public double area() {
    if (vertices.size() >= 3) {
        double area = 0;
        Vector3 origin = vertices.get(0);
        Vector3 prev = vertices.get(1).clone();
        prev.sub(origin);
        for (int i = 2; i < vertices.size(); i++) {
            Vector3 current = vertices.get(i).clone();
            current.sub(origin);
            Vector3 cross = prev.cross(current);
            area += cross.magnitude();
            prev = current;
        }
        area /= 2;
        return area;
    } else {
        return 0;
    }
}

      

To test that this method works in all orientations of the polygon, my program rotated it several times per iteration and calculated the area. So ...

Face f = poly.getFaces().get(0);
        for (int i = 0; i < f.size(); i++) {
            Vector3 v = f.getVertex(i);
            v.rotate(0.1f, 0.2f, 0.3f);
        }
        if (blah % 1000 == 0)
            System.out.println(blah + ":\t" + f.area());

      

My method seems to be correct when tested with a 20x20 square. However, the rotate method (a method in the Vector3 class) seems to introduce some error in the position of each vertex in the polygon, which affects the area calculation. Here is the Vector3.rotate () method

public void rotate(double xAngle, double yAngle, double zAngle) {
    double oldY = y;
    double oldZ = z;
    y = oldY * Math.cos(xAngle) - oldZ * Math.sin(xAngle);
    z = oldY * Math.sin(xAngle) + oldZ * Math.cos(xAngle);

    oldZ = z;
    double oldX = x;
    z = oldZ * Math.cos(yAngle) - oldX * Math.sin(yAngle);
    x = oldZ * Math.sin(yAngle) + oldX * Math.cos(yAngle);

    oldX = x;
    oldY = y;
    x = oldX * Math.cos(zAngle) - oldY * Math.sin(zAngle);
    y = oldX * Math.sin(zAngle) + oldY * Math.cos(zAngle);
}

      

Here is the output for my program in Iteration: Scope format:

0:  400.0
1000:   399.9999999999981
2000:   399.99999999999744
3000:   399.9999999999959
4000:   399.9999999999924
5000:   399.9999999999912
6000:   399.99999999999187
7000:   399.9999999999892
8000:   399.9999999999868
9000:   399.99999999998664
10000:  399.99999999998386
11000:  399.99999999998283
12000:  399.99999999998215
13000:  399.9999999999805
14000:  399.99999999998016
15000:  399.99999999997897
16000:  399.9999999999782
17000:  399.99999999997715
18000:  399.99999999997726
19000:  399.9999999999769
20000:  399.99999999997584

      

Since this is ultimately intended for the physics engine, I would like to know how I can minimize the cumulative error, since the Vector3.rotate () method will be used on a regular basis.

Thank!

A few odd notes:

  • The error is proportional to the amount of rotation. i.e. more rotation per iteration → more error per iteration.

  • There are more errors when passing doubles to the rotation function than when passing them floating point.

+3


source to share


2 answers


You will always have a cumulative error with repetitive floating point triggers - what exactly do they work. To deal with this, you basically have two options:



  • Just ignore it. Note that in your example, after 20,000 iterations (!), The Scope remains accurate to 13 decimal places. That's not bad, considering that doubles can only store about 16 decimal places to start with .

    Indeed, by plotting the graph, the area of ​​your square appears to be more or less linear:
    PLOT
    This makes sense if we assume that the effective determinant of your approximate rotation matrix is about 1 & minus; 3.417825 x 10 -18 which is good for the normal double precision floating point error range. In this case, the area of ​​your square will continue a very slow exponential decay towards zero, so you needtwo billion (2 x 10 9 )7.3 10 14 to get the area up to 399. Assuming 100 iterations per second, about seven and a half months 230 thousand years.

    Edit: When I first calculated how long it would take for the area to reach 399, it seems I made a mistake and somehow managed to overestimate the decay rate at around 400,000 (!). I have fixed the error above.

  • If you still feel like you don't need the cumulative error, the answer is simple: don't iterate over floating point rotations. Instead, store your current object orientation in a member variable, and use this information to always rotate the object from its original orientation to its current orientation.

    It's easy in 2D, since you just need to keep the angle. In 3D, I suggest storing either a quaternion or a matrix, and sometimes rescaling it so that its norm / determinant remains roughly the same (and if you use a matrix to represent the orientation of a rigid body, it remains roughly orthogonal ).

    Of course, this approach will not eliminate the cumulative error in object orientation, but rescaling ensures that the volume, region, and / or shape of the object is not affected.

+8


source


You say there is a cumulative error, but I don't believe there is (note how your output is not always omitted) and the rest of the error is due to rounding and loss of precision in the float.

I was working on a 2d physics engine at university (in java) and found double to be more accurate (of course it sees the dimensions of the oracles data

In short, you will never get rid of this behavior, you just have to accept the limitations of precision

EDIT:



Now I am looking at your .area function, possibly cumulative due to

+= cross.magnitude

      

but I have to say that the whole function looks a little strange. Why do you need to know the previous vertices to calculate the current area?

+1


source







All Articles