Segment Intersection — Failing on Unknown Edge Cases

Question

I was recently at a programming competition. One of the problems was to determine whether two line segments are parallel, intersect, or do not intersect.

For all of the test cases I ran, my solution worked. However, when I submitted it, it was rejected (but no other information was given). I suspect there is either some edge case I am completely missing... or perhaps the judges messed up.

Input is entered as follows:

ax ay bx by cx cy dx dy

where the first segment is given by $$(a_x, a_y), (b_x, b_y)$$ and the second by $$(c_x, c_y), (d_x, d_y)$$

Some other rules/assumptions:

• We may assume that no segment has length zero; that is, each segment has distinct endpoints.

• Segments may share endpoint(s) with the other segment. Endpoints are included in the segment and thus are able to intersect (this also means that a segment whose endpoint lies within the other segment will intersect)

• Lines which both intersect and are parallel (either partial or complete overlap) should be considered parallel.

• All coordinates are integers on the open interval (-1000, 1000).

Here's a screenshot of the problem. Note that I'm ignoring the multiple tests option because I highly doubt that was the issue. enter image description here

import java.util.Scanner;
public class QP1505
{
 public static void main(String[] args)
 {
 Scanner in = new Scanner(System.in);
 int ax = in.nextInt();
 int ay = in.nextInt();
 int bx = in.nextInt();
 int by = in.nextInt();
 int cx = in.nextInt();
 int cy = in.nextInt();
 int dx = in.nextInt();
 int dy = in.nextInt();
 in.close();
 /*
 * Cross multiply to check equal slopes This takes care of
 * +/- infinity as well as zero
 */
 if ((by - ay) * (dx - cx) == (dy - cy) * (bx - ax))
 {
 System.out.println("PARALLEL");
 }
 else
 {
 double m1 = (double) (by - ay) / (bx - ax);
 double m2 = (double) (dy - cy) / (dx - cx);
 // Find point of intersection
 double x, y;
 if (Double.isInfinite(m1))
 {
 x = ax;
 y = m2 * (x - cx) + cy;
 }
 else if (Double.isInfinite(m2))
 {
 x = cx;
 y = m1 * (x - ax) + ay;
 }
 else
 {
 x = (m1 * ax - ay - m2 * cx + cy) / (m1 - m2);
 y = m1 * (x - ax) + ay;
 }
 // Check bounds
 if (((x >= ax && x <= bx) || (x <= ax && x >= bx))
 && ((y >= ay && y <= by) || (y <= ay && y >= by))
 && ((x >= cx && x <= dx) || (x <= cx && x >= dx))
 && ((y >= cy && y <= dy) || (y <= cy && y >= dy)))
 {
 System.out.println("INTERSECT");
 }
 else
 {
 System.out.println("DO NOT INTERSECT");
 }
 }
 }
}

And some test cases:

• Vertical Horizontal, Intersect

  0 0 0 5 -1 2 1 2
  INTERSECT
  

• Vertical Horizontal, No Intersect

  0 0 0 5 1 2 4 2
  DO NOT INTERSECT
  

• Horizontal Horizontal, Parallel

  0 0 5 0 1 1 5 1
  PARALLEL
  

• Vertical Regular, Intersect

  0 0 0 5 1 0 -1 4
  INTERSECT
  

• Vertical Regular, No Intersect

  0 0 0 5 1 0 4 6
  DO NOT INTERSECT
  

• Horizontal Regular, Intersect

  4 4 7 4 4 3 6 5
  INTERSECT
  

• Horizontal Regular, No Intersect

  2 2 8 2 3 3 9 6
  do not intersect
  

• Regular Regular, Intersect

  2 1 4 3 1 2 5 1
  INTERSECT
  

• Regular Regular, Parallel

  1 2 5 6 3 2 5 4
  PARALLEL
  

• Regular Regular, No Intersect

  1 2 5 6 3 2 8 5
  DO NOT INTERSECT
  

Let me know if anything else needs to be clarified (or any other test cases).

• baum: It does work, for every single test case I have run. But should it be migrated to stack overflow, perhaps?
• baum: That works. I'll post a bunch of test cases. Also, the problem is not available online as it was just released this morning.
• Gareth Rees: Don't use slope–intercept representation, use vectors! See this answer.
• Snowbody: Did the contest organizers say that the points were restricted in any way? Are they guaranteed to be integers? Are they in a specific range?
• baum: Problem was just released online; see edit.

Answer 1

Failure cases

I ran your program and was able to find an endless number of failures. My technique was to have one line touch the origin (0,0) and then generate all lines that pass directly through the origin by simply choosing mirrored coordinates. I was hoping to find roundoff errors and I did. Here are just a small sampling of failure cases, all of which print "Do not intersect" instead of "intersect":

0 0 0 5 7 29 -7 -29
0 0 0 5 11 -30 -11 30
0 0 0 5 11 -15 -11 15
0 0 0 5 11 25 -11 -25
0 0 0 5 13 -30 -13 30
0 0 0 5 13 -15 -13 15
0 0 0 5 13 27 -13 -27
0 0 0 5 14 29 -14 -29
0 0 0 5 19 21 -19 -21
0 0 0 5 21 23 -21 -23
0 0 0 5 21 27 -21 -27
0 0 0 5 22 -30 -22 30
0 0 0 5 22 -15 -22 15
0 0 0 5 22 25 -22 -25
0 0 0 5 23 -26 -23 26
0 0 0 5 23 -13 -23 13
0 0 0 5 23 27 -23 -27
0 0 0 5 25 -29 -25 29
0 0 0 5 25 7 -25 -7
0 0 0 5 25 14 -25 -14
0 0 0 5 25 28 -25 -28
0 0 0 5 26 -30 -26 30
0 0 0 5 26 -15 -26 15
0 0 0 5 26 27 -26 -27
0 0 0 5 27 29 -27 -29
0 0 0 5 28 29 -28 -29
0 0 0 5 29 15 -29 -15

Floating point?

My personal opinion about solving this problem correctly is that if your coordinates are constrained to be within -1000..1000, then you can do the whole problem using integers and not use floating point at all.

What you would need to do is create a fraction class that holds a numerator and denominator. The fraction class needs to be able to: multiply, divide, add, subtract, compare (i.e. everything you are doing with your doubles). Since the range of coordinates is so small, you shouldn't have a problem with overflow. By doing this you will be safe from roundoff errors.

• baum: ...and that must have been it. Amazing. One of the judges suggested we should have used fractions, but I didn't even think of that sort of test.
• JS1: @baum It's hard to get floating point numbers to not have roundoff errors. I think the parallel slopes test was ok (and it could have been done with integers). But finding the intersection point involved division and that's where you start to not be able to represent things exactly.
• baum: Out of curiosity, if the magnitude of the coordinates was not restricted, would there be any foolproof way of solving this problem? Perhaps BigInteger fractions? I'm thinking segments whose endpoints lie incredibly close to each other but are extremely far from the origin.
• JS1: @baum Depending on the situation, a long might suffice. If not, then BigInteger would work.
• Snowbody: The problem is that m can't be computed exactly, so x is not computed exactly either (it would have to be exactly 0 in order to work). Good job, enjoy your bounty.

Answer 2

As it was already pointed out, computing the slope of the line segments (as a double) is problematic because the slope can be "infinite" for vertical segments, and because of rounding errors.

My suggestion is to use a different algorithm which does not compute the slope. If the coordinates are integers then all intermediate values are integers as well and no rounding errors can occur. (I would choose the same algorithm for floating point coordinates because it does not need any special cases for vertical or nearly vertical segments.)

The idea is to describe the line segment from \$(a_x, a_y) \$ to \$ (b_x, b_y) \$ in its parametric form $$ (x, y) = (a_x, a_y) + u (b_x - a_x, b_y - a_y) ,円 \quad \text{where } 0 \le u \le 1 ,円. $$ Then the intersection of two line segments is a solution of the linear equation system $$ (b_x - a_x) u - (d_x - c_x) v = c_x - a_x \\ (b_y - a_y) u - (d_y - c_y) v = c_y - a_y $$ with \$ 0 \le u \le 1 \$ and \$ 0 \le v \le 1 \$. In the general case, the solution is given by $$ u = \frac {\Delta_u} \Delta ,円 , \quad v = \frac {\Delta_v} \Delta $$ with the determinants $$ \Delta = \begin{vmatrix} b_x - a_x & d_x - c_x \\ b_y - a_y & d_y - c_y \end{vmatrix} ,円 , \quad \Delta_u = \begin{vmatrix} c_x - a_x & d_x - c_x \\ c_y - a_y & d_y - c_y \end{vmatrix} ,円 , \quad \Delta_v = \begin{vmatrix} c_x - a_x & b_x - a_x \\ c_y - a_y & b_y - a_y \end{vmatrix} $$ If \$ \Delta = 0\$ then the line segments are parallel. Otherwise there is a unique solution \$u, v \$ to the linear equation system, and it easy to check if the solutions fall in the range \$ [0, 1] \$ without actually performing the division.

This leads to the following simple function:

static void checkIntersection(int ax, int ay, int bx, int by, int cx, int cy, int dx, int dy) {
 int det = (bx - ax)*(dy - cy) - (by - ay)*(dx - cx);
 if (det != 0) {
 /*
 * Lines intersect. Check if intersection point is on both segments:
 */
 int detu = (cx - ax)*(dy - cy) - (cy - ay)*(dx - cx);
 int detv = (cx - ax)*(by - ay) - (cy - ay)*(bx - ax);
 if (det < 0) {
 // Normalise to det>0 to simplify the following check.
 det = -det;
 detu = -detu;
 detv = -detv;
 }
 if (detu >= 0 && detu <= det && detv >= 0 && detv <= det) {
 System.out.println("INTERSECT");
 } else {
 System.out.println("NO NOT INTERSECT");
 }
 } else {
 /*
 * Lines are parallel (or identical):
 */
 System.out.println("PARALLEL");
 }
}

For integer coordinates in the range \$ -1000 \ldots 1000 \,ドル all the calculations can be done without overflow. Generally, for coordinates in a range \$ -M \ldots M \,ドル the computed values are bounded by \$ 8M^2\$ in absolute value, so you can choose the required type accordingly.

• Martin R: I just noticed that this is exactly what Gareth suggested in a comment to the question, slightly simpler here because we don't have to check if parallel segments overlap or not. I apologize for not noticing that earlier.
• Snowbody: Also note that in the event you're not restricted to (-1000,1000) for each coordinate, you do need to handle wraparound on addition/subtraction, and overflow on multiplication.
• Martin R: @Snowbody: Sure, that's what I tried to express (perhaps badly) in the last paragraph.

Answer 3

Abstraction

Code that's well abstracted is easier to maintain, ...and to debug. I'd start by defining what a segment is - it's an immutable data structure with 2 points, where each point is defined by X and Y values (I'm sure java has something like a Point class built-in somewhere).

So I would take this chunk:

 Scanner in = new Scanner(System.in);
 int ax = in.nextInt();
 int ay = in.nextInt();
 int bx = in.nextInt();
 int by = in.nextInt();
 int cx = in.nextInt();
 int cy = in.nextInt();
 int dx = in.nextInt();
 int dy = in.nextInt();
 in.close();

And extract it into its own getSegments method, so that this:

if ((by - ay) * (dx - cx) == (dy - cy) * (bx - ax))

Can look something like this:

if (isParallel(segment1, segment2))

And then you would have isParallel on its own, and you could write a dozen unit tests just for that method, covering every possible edge case.

---

Then you can extract another findIntersectionPoint method here, which would return a Point:

 double m1 = (double) (by - ay) / (bx - ax);
 double m2 = (double) (dy - cy) / (dx - cx);
 // Find point of intersection
 double x, y;
 if (Double.isInfinite(m1))
 {
 x = ax;
 y = m2 * (x - cx) + cy;
 }
 else if (Double.isInfinite(m2))
 {
 x = cx;
 y = m1 * (x - ax) + ay;
 }
 else
 {
 x = (m1 * ax - ay - m2 * cx + cy) / (m1 - m2);
 y = m1 * (x - ax) + ay;
 }

And then you can cover that method with another dozen unit tests, and - assuming the "check bounds" part is also extracted into its own method, simplify the calling code to this:

Point intersection = findIntersectionPoint(segment1, segment2);
if (intersects(intersection, segment1, segment2))
{
 // intersect
}
else
{
 // no intersect
}

And then you can bombard findIntersectionPoint and intersects with another dozen unit tests to cover all edge cases.

---

The key point in this answer is granularity: you have one single method that seems to work, but it's doing many things that can't be tested individually - so even if you do end up figuring out a failing test, you wouldn't know exactly where the problem is and where to start looking.

Defining a Segment object might be over-the-top, but I find it simplifies the code, at least how it reads: seeing segment2.point1.x instead of cx makes it much easier to mentally map the values.. at least in my poor little fried brain.

---

One thing that strikes me, as a c# programmer, is that your curly braces are... c#-style. Typically java code would look something like this:

if (Double.isInfinite(m1)) {
 x = ax;
 y = m2 * (x - cx) + cy;
} else if (Double.isInfinite(m2)) {
 x = cx;
 y = m1 * (x - ax) + ay;
} else {
 x = (m1 * ax - ay - m2 * cx + cy) / (m1 - m2);
 y = m1 * (x - ax) + ay;
}

• baum: That's useful... I'll try doing that. As for braces style, that's how my professor likes them (though he's not originally a Java programmer either)

Answer 4

Okay I found some actual examples. with integer coordinates.

EXAMPLE #1 - Wraparound

Let's try the most extreme case, where one of the segments is: (INTEGER.MIN_VALUE,0) (Integer.MAX_VALUE,1). The slope should be positive and nearly zero, but calculation fails spectacularly. In an effort to avoid integer division, the code casts the numerator to a double, but the problem is that the subtractions are still being performed using int arithmetic:

 double m1 = (double) (by - ay) / (bx - ax);
 m1 = (double)(1-0)/(Integer.MAX_VALUE-Integer.MIN_VALUE);

Here, Java wraps around:

 m1 = (1.0)/(-1);
 m1 = -1.0;

I don't think we need to continue here; the calculation of x will obviously be wrong. The wraparound could be avoided by making sure that the individual elements are each cast to double or java.Math.BigInteger, so that all arithmetic is performed with doubles or BigIntegers, but there are still problems with that approach -- see below.

EXAMPLE #2 - Roundoff error

In this example we make sure that wraparound doesn't happen, either by choosing an interval that doesn't overflow, or by casting to double. There's still a problem. Say we have (0,0) (Integer.MAX_VALUE,1) and (0,0) (Integer.MAX_VALUE-1,1).

 double m1 = (double) (by - ay) / (bx - ax);
 m1 = (double)(1 - 0)/(Integer.MAX_VALUE-0);
 m1 = (double)(1)/Integer.MAX_VALUE;
 m1 = 1.0/Integer.MAX_VALUE;
 m1 = 0b0.0000000000000000000000000000001000000000000000000000000000000100000000000000000000000;
 double m2 = (double) (dy - cy) / (dx - cx);
 m2 = (double)(1 - 0)/(Integer.MAX_VALUE-1-0);
 m1 = (double)(1)/(Integer.MAX_VALUE-1);
 m2 = 1.0/(Integer.MAX_VALUE-1);
 m2 = 0b0.0000000000000000000000000000001000000000000000000000000000001000000000000000000000000;

These slopes are small, and very close to each other but not quite identical. They do differ in fact, but the difference is calculated as 0b0.00000000000000000000000000000000000000000000000000000000000001 = 2^-62 when really it's slightly larger than that (it's 2^-62 + 2^-92 + ...). So, m1-m2 is either slightly too large or too small, which causes x to be either slightly too large or too small. So when the point of intersection is the endpoint of one of the segments, the algorithm determines the lines don't intersect, when in fact they do.

You can verify the binary floating point arithmetic using http://www.exploringbinary.com/binary-calculator/ . Do the calculation once to find out where the most significant digit is, then add on 53 more bits.

I also notice the parallel check needs wraparound protection, both on the subtractions and on the multiplications; with a little work I could find two segments that report as parallel when they're not, or report as non parallel when they are. If you cast everything to double you still run the risk of the 53 bits of precision not being enough (think about multiplying two 30-bit quantities). I think you may need to use java.math.BigDecimal or java.math.BigInteger here.

• baum: While these over/underflow issues do seem to be legitimate concerns (and answer my question, so I'll probably grant bounty), the funny thing was that the contest organizers said that all coordinate points would be between -1000 and 1000. If no one finds anything else, perhaps they were in the wrong and I was right.
• Snowbody: Don't award it to me yet, maybe someone will come along and spot an error with coordinate points in that range.

Answer 5

You're using a dangerous method to obtain the point of intersection. Your method works fine when one of the slopes is calculated as Infinite but (削除) it's numerically unstable if the two slopes differ by many orders of magnitude (more than the 15-17 decimal digits that a Java double supports). (削除ここまで) This can't happen because all the coordinates are integers. Even the most extreme points can't cause roundoff error.

(削除) The problem is with the slopes. When you calculate m1 and m2 the subtractions and divisions may not be exact, which introduces some error into the slopes themselves. This is then compounded when the two slopes are then subtracted from each other m1 - m2 in the quotient where you are calculating x. It's quite possible that this subtraction will remove all the signal and be left with only noise, which will make the calculated x completely invalid. (削除ここまで) While I wish this was true, because all the coordinates are ints the subtractions are exact, and the divisions don't introduce that much error.

(削除) How about the first segment being (0,1), (10^{17},10^{17}) and the second segment being (1,0), (10^{17},10^{17}). The lines are not parallel; they intersect at x == 1e17, y == 1e17. Since your test for parallel is correct, it will not print "parallel" -- so far so good. Unfortunately, due to the limits of Java's double, both m1 and m2 will be calculated as 1.0000000000. Your program will crash when it tries to compute x due to the divide by zero. (削除ここまで)

That doesn't work, outside range for int plus they'd be detected as parallel. No bounty for me.

• Snowbody: This answer was completely wrong (though it would be correct if the inputs had been doubles). Please downvote this one and upvote my other answer.
• Mathieu Guindon: If you ask for it so kindly... ;)
• user34073: @Snowbody You can always delete this one...
• Snowbody: Eh, I wanted to point out the hazards of widely varying slopes.