General strategies for testing for equality of dou

Question

So this seems to be a reoccurring problem for me- I'm trying to implement the line segment intersection and doubly connected edge list overlay algorithms in Computational Geometry by de Berg et al. Right now, I'm using the following function to test for equality of two values:

private static final double absTol = 1e-14;
private static final double relTol = 1e-10;        
public static final boolean areClose(double i, double j) {
    return (Math.abs(i-j) <= absTol) || 
           (Math.abs(i-j) <= relTol*Math.max(Math.abs(i), Math.abs(j)));
}

The main problem I'm having is that I'll occasionally experience a bug. For example, today I realized that one of my functions was failing because I had originally set absTol to 1e-16 and the above function failed when one of my functions should have identified that 0.0 and 1.535e-15 were equal. I fixed the problem by adjusting absTol up to 1e-14. So my question is, is there a better way to attack this problem then honing the tolerances? Or, should you modify the algorithm so that they simply output the wrong values instead of crashing? Not sure exactly how to proceed from here.

One thing I note is that in the line segment intersection algorithm, outlined here, two separate functions are necessary for calculating intersection points. Basically, the first time an intersection related event point is calculated, you calculate the intersection of two segments; I use the following algorithm:

public static Point2D findIntersection(Line2D line1, Line2D line2) {
    double s1X = line1.getX2()-line1.getX1();     
    double s1Y = line1.getY2()-line1.getY1();
    double s2X = line2.getX2()-line2.getX1();     
    double s2Y = line2.getY2()-line2.getY1();

    double s = (-s1Y*(line1.getX1()-line2.getX1())+s1X*(line1.getY1()-line2.getY1()))/(-s2X*s1Y+s1X*s2Y);
    double t = ( s2X*(line1.getY1()-line2.getY1())-s2Y*(line1.getX1()-line2.getX1()))/(-s2X*s1Y+s1X*s2Y);

    if (s >= 0 && s <= 1 && t >= 0 && t <= 1) {
        return new Point2D.Double(line1.getX1()+(t*s1X), line1.getY1()+(t*s1Y));
    } else {
        return null; 
    }
}

But, when determining the equality of lines in the status, I've been using:

private double getIntercept(Line2D line) {
    if (MathUtilities.areClose(line.getY1(), line.getY2())) {
        // line is horizontal. Set intersection to x value
        return this.x;
    } else {
        return line.getX1() + (line.getX2()-line.getX1())*((this.y-line.getY1())/(line.getY2()-line.getY1()));
    }   
}

to test where lines intersect the event line (and when two or more lines have the same intercept, they intersect the event line at the same position). So essentially, two different algorithms are used to find intersection points, which means the values I get differ slightly. Lastly, I've also come to realize that the test for equality used in areClose is not necessarily transitive, so that will also cause some problems.

Answer 1

Dealing with "exactness" in computational geometry is an issue that comes up more often than you might think. Other than fixed point arithmetic (which has already been suggested), another approach is to use adaptive precision arithmetic -- evaluating operations in better than double precision, but only when necessary to preserve correctness.

Implementing adaptive precision operations is very much a non-trivial task, but there are a few libraries that are available, i.e. Shewchuck's robust predicates and/or the MPFR library that is used by the CGAL geometry library.

Robustly detecting the intersection between two lines could be done via a few calls to Shewchuk's orientation routine.

Answer 2

You have an essential problem in the fact that you store coordinates to (potentially) higher precision than is meaningful to you. That is, if two values are always to be considered equal when they differ by less than absTol, then it is meaningless to store or use values that are not integer multiples of absTol, and performing computations with such values can produce anomalous results.

You also have an inherent inconsistency associated with using a relative tolerance in addition to an absolute tolerance. For example, if you scale a model either up or down then some of its areClose() relationships can change.

Have you considered using fixed-point arithmetic? That would avoid arithmetic anomalies and some scaling anomalies, plus it would give you some relief on the topic of evaluating equality. It wouldn't solve all your problems, but it might be worth a look.