15-211 Fundamental Data Structures and Algorithms

1. 15-211Fundamental Data Structures and Algorithms Aleks Nanevski April 22, 2004 Computational Geometry

2. Announcements • HW7 the clock is ticking ... • Quiz 3: Today, after class • Final Exam on Tuesday May 4, 5:30 pm • review session April 29

3. The Basics

4. Computational Geometry Lots of applications. - computer graphics - CAD/CAM - motion planning Could be done in dimension D, but we will stick to D=2. D=3 is already much harder.

5. Airplane wing triangulation

6. The Basics What are the basic objects in plane geometry, how do we represent them as data structures? - point two floats - line two points - ray two points, constraints - line segment two points, constraints and triangles, rectangles, polygons, ... Never mind curves such as circles and ellipses.

7. Disclaimer There are two annoying issues one has to confront in CompGeo algorithms. - Degenerate Cases Points are collinear, lines parallel, … Dealing with degeneracy is a huge pain, not to mention boring. In this lecture we ignore the degenerate cases. - Floating point numbers are not reals. Floating point arithmetic introduces rounding errors and hence leads to accuracy problems.

8. High Precision Arithmetic If 32 or 64 bits are not enough, use 128, or 512, or 1024, or … Problems: We have to be able to analyze how many bits are required. Slows down computation significantly. See GNU gmp.

9. Interval Arithmetic Replace “real number” x by an interval xL x  xU Problems: Arithmetic becomes quite a bit more complicated. Slows down computation significantly.

10. Symbolic Computation Use (arbitrary precision) rationals and whatever algebraic operations that are necessary symbolically: compute with the formal expressions, not their numerical values. Sqrt[2] Problems: Manipulating these expressions is complicated. Slows down computation significantly.

11. Warm-Up Let's figure out how to check if a point P lies on a line L. P = (p1,p2) L = (A,B) two-points representation Need P = A + l(B - A) where l is a real parameter. P B P A

12. Linear Equations The vector equation P = A + l (B – A) = (1- l) A + l B is just short-hand for two linear equations: p1 = a1 + l (b1 – a1) p2 = a2 + l (b2 – a2) Will have a common solution if P lies on L.

13. Subroutines Solving systems of equations such as p1 = a1 + l (b1 – a1) p2 = a2 + l (b2 – a2) is best left to specialized software: send all linear equation subproblems to an equation solver - that will presumably handle all the tricky cases properly and will produce reasonable accuracy.

14. Rays and Line Segments What if L =[A,B) is the ray from A through B? No problem. Now we need P = A + l (B - A) where l ≥ 0. Likewise, if L = [A,B] is the line segment from A to B we need 0 l 1. Both l and (1 – l) must be non-negative.

15. Intersections How do we check if two lines intersect? A + l (B - A) = C + k (D – C) Two linear equations, two unknowns l, k. Add constraints on l, k for intersection of rays, line segments. B C A D

16. General Placement Here is a slightly harder problem: Which side of L is the point P on? P P left on B P right A Warning: orientation of the line matters.

17. left on B right A Left/Right Turns March from A to B, then perform a left/right turn.

18. Calculus: determinants Given vectors X=(x1, x2) and Y=(y1, y2), the determinant | x1 x2 | | y1 y2 | expresses the (signed) area of the parallelogram spanned by X and Y. = x1 y2 – x2 y1 (y1,y2) Can be generalized to higher dimension. (x1,x2)

19. Determinants and left/right turns Let P=(p1,p2), A=(a1,a2), B=(b1, b2) If the area of the triangle ABP is 0: then P is on the line (A, B) >0: then P is to the left of (A, B) <0: then P is to the right of (A,B) The triangle ABP is spanned by vectors B-A and P-A, so its area is half of the determinant |b1-a1, b2-a2| |p1-a1, p2-a2| P = (p1,p2) D = B = (b1,b2) Note: we only look for the sign of the ABP area, so no need to divide D by 2. A = (a1,a2)

20. Some Applications

21. Membership How do we check if a point X belongs to some region R in the plane (triangle, rectangle, polygon, ... ). R X R X X X

22. Membership Clearly the difficulty of a membership query depends a lot on the complexity of the region. For some regions it is not even clear that there is an inside and an outside.

23. Simple Polygons A polygon is simple if its lines do not intersect except at the vertices of the polygon (and then only two). By the Jordan curve theorem any simple polygon partitions the plane into inside and outside.

24. Point Membership The proof of the JCT provides a membership test. Let P be a simple polygon and X a point. Draw a line through X and count the intersections of that line with the edges of P.

25. Odd One can show that X lies inside of P iff the number of intersections (on either side) is odd. Note that there are bothersome degenerate cases: X lies on an edge, an edge lies on the line through X. Maybe wiggle the line a little.

26. Algorithm Theorem: One can test in linear time whether a point lies inside a simple polygon. Linear here refers to the number of edges of the polygon. And, of course, we assume that these edges are given as a list of points, say, in counterclockwise order.

27. Convexity Here is one type of simple region. A region R in the plane is convex if for all A, B in R: the line segment [A,B] is a subset of R. In particular convex polygons: boundary straight line segments.

28. Membership in Convex Polygon Traverse the boundary of the convex polygon in counterclockwise order, check that the point in question is always to the left of the corresponding line. This is clearly linear in the number of boundary points.

29. Intersection of line segments • Suppose we are given a collection of line segments and want to determine the intersections between them. • The brute force algorithms is (n2): consider all pairs of segments and check each pair for intersection. • Bad if there are only a few intersections.

30. Intersection of line segments Suppose there are s intersections. We would like an algorithm with running time O( (n+s) ??? ) where ??? is hopefully small. That way we only pay the quadratic overhead when there actually is a quadratic number of intersections.

31. Line Sweep Algorithm • Idea: use a sweep line that moves from left to right. • At each point, sweep line defines an ordering on segments (by decreasing y-coordinate). The ordering changes as the sweep progresses. • Check for intersection only if two segments are consecutive in this ordering. If two segments intersect, there will be a sweep line for which they are consecutive. Ordering for this sweep line is a > b > c b a c

32. Events • No need to consider all possible sweep lines • Just look at: • - segment edges • - intersections • These are the only • times when segment • ordering changes.

33. Events • What is an appropriate data structure to keep the positions of the sweep line? • Hint: remember the • Car dodger

34. Events • What is an appropriate data structure to keep the positions of the sweep line? • Answer: priority queue! • segment edges queued at initialization time • - intersections queued dynamically

35. Sweep line operations The sweep line has to keep track of all the segments that currently intersect it. So we need operations insert(s) delete(s) where s is a given segment. When a new segment is entered, we have to check if it interacts with any of the neighbouring segments. So we need additional operations above(s) below(s) which return the immediate neighbours above/below on the sweep line.

36. Sweep line operations What is an appropriate data structure for these operations? insert(s) delete(s) above(s) below(s)

37. Sweep line operations What is an appropriate data structure for these operations? insert(s) delete(s) above(s) below(s) Answer: some kind of dictionary, like balanced binary trees. Then each operation can be done in O(log n) time.

38. Sweep algorithm • Insert all segment endpoints into the queue Q • While Q not empty, extract next event E • If E = segment left point, insert segment into SL dictionary based on its y-coordinate. Test for intersection with segment immediately above and below. • If E = segment right point, delete from SL dictionary. Test for intersection the segments immediately preceding and succeeding this one. • If E = intersection, swap the intersecting segments in the SL dictionary. Test the new upper segment for intersection with predecessor. Test the new lower segment for intersection with successor. Theorem: Can compute all s intersections of n line segments in O((n+s) log n) time steps.

39. Testing Simplicity How do we check if a polygon is simple? Use the line segment intersection algorithm on the edges.

40. Polygon Intersection How do we check if two simple polygons P and Q intersect? Note: two cases! P P Q Q

41. Polygon Intersection The first case can be handled with LS intersection testing. How? For the second case, pick any vertex in Q and check whether it lies in (the interior of) P. If yes, then Q is inside P (assuming there were no excessive intersections in the first case). Total running time: O(n log n).

42. Convex Hulls

43. A Hull Operation Suppose P is a set of points in the plane. The convex hull of P is the least set of points Q such that: - P is a subset of Q, - Q is convex. Written CH(P). This pattern should look very familiar by now (reachability in graphs, equivalence relations, ...)

44. Convex Combinations Abstractly it is easy to describe CH(P): is the set of all points X =  i Pi where 0 i and  i = 1. X is a convex combination of the Pi. So the convex combinations of A and B are the line segment [A,B]. And the convex combinations of three points (in general position) form a triangle. And so on.

45. So? Geometrically this is nice, but computationally this characterization of the convex hull is not too useful. We want an algorithm that takes as input a simple polygon P and returns as output the simple polygon CH(P). CH(P) P

46. Extremal Points Point X in region R is extremal iff X is not a convex combination of other points in R. For a polygon, the extremal points make up the CH.

47. An Observation So let's deal with a set of points P rather than regions. Lemma: Point X in P is not extremal iff X lies in the interior of a triangle spanned by three other points in P.

48. An Algorithm Hence we can find the convex hull of a simple polygon: for all n points we eliminate non-extremal points by trying all possible triangles using the membership test for convex polygons. In the end we sort the remaining extremal points in counterclockwise order. Unfortunately, this is O(n4). Could it be faster? When do we get n4?

49. Next Time Clearly, O(n4) is way too slow to be of any practical use. Next time we will see a number of simple algorithms that push this down to O(n log n). This turns out to be optimal: sorting can be reduced to convex hull. Note, though, that in dimension 3 and higher things get much more messy.