Classroom Examples of Robustness Problems in Geometric Computations
Classroom Examples of Robustness Problems in Geometric Computations Lutz Kettner, Kurt Mehlhorn MPI für Informatik, Saarbrücken Sylvain Pion INRIA, Sophia Antipolis Stefan Schirra Otto-von-Guericke-Universität, Magdeburg 3/4/2021 Chee Yap New York University
Overview • Motivation • 2 D Convex Hull Algorithm • Constructing Bad Input Examples • Analysis Classroom Examples ESA'04 2
Motivation • • The algorithms of computational geometry are designed for a machine model with exact real arithmetic. Substituting floating point arithmetic for the assumed real arithmetic may cause implementations to fail. Although this is well known, it is not common knowledge. There is no good material for the classroom to show the inadequacy of floating point arithmetic. Classroom Examples ESA'04 3
Example • Arrangements of circular arcs Seen in Antibes, France, during SGP’ 04 Classroom Examples ESA'04 4
Example • Arrangements of circular arcs, with double arithmetic Classroom Examples ESA'04 5
Example • Arrangements of circular arcs, with float arithmetic Classroom Examples ESA'04 6
Motivation • Let’s be more serious … Classroom Examples ESA'04 7
Intersection of Four Simple Solids • Rhino 3 D: § ((( s 1 ∩ s 2) ∩ c 1) → successful § ((( c 1 ∩ c 2) ∩ s 1) ∩ s 2) → ``Boolean operation failed'' Classroom Examples ESA'04 8
Intersection of Four Simple Solids output is a combinatorial object plus coordinates (not a point set) • Rhino 3 D: § ((( s 1 ∩ s 2) ∩ c 1) → successful § ((( c 1 ∩ c 2) ∩ s 1) ∩ s 2) → ``Boolean operation failed'' Classroom Examples ESA'04 9
Motivation • • Let’s be more serious … …, but that is hard to explain in all detail in class. • Let’s do a really simple 2 D convex hull algorithm. Classroom Examples ESA'04 10
2 D-Orientation of Three Points Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) Implemented with IEEE 784 double precision float_orient (p, q, r) Classroom Examples ESA'04 11
Convex Hulls in the Plane v 3 v 1 r v 2 • • • maintain current hull as a circular list L=(v 0, v 1, . . . , vk-1) of its extreme points in counter-clockwise order start with three non-collinear points in S. consider the remaining points r one by one. Classroom Examples ESA'04 12
Convex Hulls in the Plane v 4 v 1 v 3 v 2 • • • maintain current hull as a circular list L=(v 0, v 1, . . . , vk-1) of its extreme points in counter-clockwise order start with three non-collinear points in S. consider the remaining points r one by one. Classroom Examples ESA'04 13
Convex Hulls in the Plane v 4 v 1 v 3 v 2 • • • maintain current hull as a circular list L=(v 0, v 1, . . . , vk-1) of its extreme points in counter-clockwise order start with three non-collinear points in S. consider the remaining points r one by one. Classroom Examples ESA'04 14
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • • • maintain current hull as a circular list L=(v 0, v 1, . . . , vk-1) of its extreme points in counter-clockwise order start with three non-collinear points in S. consider the remaining points r one by one. Classroom Examples ESA'04 15
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) Classroom Examples ESA'04 16
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) Classroom Examples ESA'04 17
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) Classroom Examples ESA'04 18
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) [Property B] If r is outside CH, then the set of edges that are weakly visible (= orientation is non-negative) from r forms a non-empty consecutive subchain; so does the set of edges that are not weakly visible from r. Classroom Examples ESA'04 19
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) [Property B] If r is outside CH, then the set of edges that are weakly visible (= orientation is non-negative) from r forms a non-empty consecutive subchain; so does the set of edges that are not weakly visible from r. Classroom Examples ESA'04 20
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) [Property B] If r is outside CH, then the set of edges that are weakly visible (= orientation is non-negative) from r forms a non-empty consecutive subchain; so does the set of edges that are not weakly visible from r. Classroom Examples ESA'04 21
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) [Property B] If r is outside CH, then the set of edges that are weakly visible (= orientation is non-negative) from r forms a non-empty consecutive subchain; so does the set of edges that are not weakly visible from r. Classroom Examples ESA'04 22
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) [Property B] If r is outside CH, then the set of edges that are weakly visible (= orientation is non-negative) from r forms a non-empty consecutive subchain; so does the set of edges that are not weakly visible from r. Classroom Examples ESA'04 23
Convex Hulls in the Plane v 4 r v 1 v 3 v 2 • • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) [Property B] If r is outside CH, then the set of edges that are weakly visible (= orientation is non-negative) from r forms a non-empty consecutive subchain; so does the set of edges that are not weakly visible from r. Classroom Examples ESA'04 24
Convex Hulls in the Plane v 4 v 3 v 1 v 2 • • [Property A] A point r is outside CH iff there is an i such that the edge (vi, vi+1) is visible for r. (orientation(vi, vi+1, r) > 0) [Property B] If r is outside CH, then the set of edges that are weakly visible (= orientation is non-negative) from r forms a non-empty consecutive subchain; so does the set of edges that are not weakly visible from r. Classroom Examples ESA'04 25
Devil’s Advocate • Systematic construction of instances leading to violations of properties (A) and (B) when executed with double’s • and in all possible ways § § • • a point outside sees no edge a point inside sees an edge a point outside sees all edges a point outside sees a non-contiguous set of edges examples involve nearly collinear points, of course examples are realistic as many real-life instances contain collinear points (which become nearly collinear by conversion to double’s) Classroom Examples ESA'04 26
Global Consequences • A point outside sees no edge of the current hull. p 1=(7. 3000000194, 7. 3000000167) p 2=(24. 000000068, 24. 000000071) p 3=(24. 00000005, 24. 000000053) p 4=(0. 50000001621, 0. 50000001243) p 5=(8, 4) p 6=(4, 9) p 7=(15, 27) p 8=(26, 25) p 9=(19, 11) float_orient(p 1, p 2, p 3) > 0 float_orient(p 1, p 2, p 4) > 0 float_orient(p 2, p 3, p 4) > 0 float_orient(p 3, p 1, p 4) > 0 Classroom Examples ESA'04 27
Global Consequences • A point outside sees all edges of the current hull. p 1=(200, 49. 200000003) p 2=(100, 49. 600000001) p 3=(-233. 33333334, 50. 93333333) p 4=(166. 66666669, 49. 33333336) float_orient(p 1, p 2, p 3) > 0 float_orient(p 1, p 2, p 4) < 0 float_orient(p 2, p 3, p 4) < 0 float_orient(p 3, p 1, p 4) < 0 Classroom Examples ESA'04 28
Global Consequences • A point outside sees all edges of the current hull. p 1=(200, 49. 200000003) p 2=(100, 49. 600000001) p 3=(-233. 33333334, 50. 93333333) p 4=(166. 66666669, 49. 33333336) float_orient(p 1, p 2, p 3) > 0 float_orient(p 1, p 2, p 4) < 0 float_orient(p 2, p 3, p 4) < 0 float_orient(p 3, p 1, p 4) < 0 Depending on the implementation: • Algorithm does not terminate! • Algorithm crashes! Classroom Examples ESA'04 29
Global Consequences • A point inside sees an edge of the current hull. p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 6. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) (p 1, p 2, p 3, p 4) form a convex quadrilateral p 4 is truly inside this quadrilateral, but float_orient(p 4, p 1, p 5) > 0 p 4 p 1 p 5 p 2 Classroom Examples ESA'04 p 3 30
Global Consequences • A point inside sees an edge of the current hull. p 4 p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 6. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) (p 1, p 2, p 3, p 4) form a convex quadrilateral p 4 is truly inside this quadrilateral, but float_orient(p 4, p 1, p 5) > 0 p 5 p 1 p 4 p 2 p 1 p 5 p 2 Classroom Examples ESA'04 p 3 31
Global Consequences • A point outside sees a non-contiguous set of edges p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 6. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) p 6=(6, 6) p 4 p 1 p 5 p 6 p 2 Classroom Examples ESA'04 p 3 32
Global Consequences • p 4 A point outside sees a non-contiguous set of edges p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 6. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) p 6=(6, 6) p 5 p 1 p 4 p 2 p 1 p 5 p 6 p 2 Classroom Examples ESA'04 p 3 33
Global Consequences • p 4 A point outside sees a non-contiguous set of edges p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 6. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) p 6=(6, 6) p 5 p 1 p 4 p 2 p 1 p 5 float_orient(p 4, p 5, p 6) < 0 float_orient(p 5, p 1, p 6) > 0 float_orient(p 1, p 2, p 6) < 0 p 6 p 2 Classroom Examples ESA'04 p 3 34
Global Consequences • p 4 A point outside sees a non-contiguous set of edges p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 10. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) p 6=(6, 6) p 5 p 1 p 4 p 2 p 1 p 5 p 6 p 2 Classroom Examples ESA'04 p 3 35
Global Consequences • A point outside sees a non-contiguous set of edges p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 10. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) p 6=(6, 6) p 4 p 1 p 5 p 6 p 2 Classroom Examples ESA'04 p 3 36
Global Consequences • p 4 A point outside sees a non-contiguous set of edges p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 10. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) p 6=(6, 6) p 5 p 1 p 6 p 4 p 1 p 5 p 6 p 2 Classroom Examples ESA'04 p 3 37
Global Consequences Classroom Examples ESA'04 38
Global Consequences • p 4 A point outside sees a non-contiguous set of edges p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 6. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) p 6=(6, 6) p 5 p 1 p 4 p 2 p 1 p 5 p 6 p 2 Classroom Examples ESA'04 p 3 39
Global Consequences • A point outside sees a non-contiguous set of edges p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 6. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) p 6=(6, 6) p 4 p 1 p 5 p 6 p 2 Classroom Examples ESA'04 p 3 40
Global Consequences • p 4 A point outside sees a non-contiguous set of edges p 6 p 1=(24. 00000005, 24. 000000053) p 2=(24. 0, 6. 0) p 3=(54. 85, 6. 0) p 4= (54. 8500000357, 61. 000000121) p 5=(24. 000000068, 24. 000000071) p 6=(6, 6) p 5 p 1 p 4 p 6 p 1 p 5 p 6 p 2 Classroom Examples ESA'04 p 3 41
Global Consequences Classroom Examples ESA'04 42
2 D-Orientation of Three Points Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) Classroom Examples ESA'04 43
2 D-Orientation of Three Points Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) 256 x 256 pixel image red=pos. , yellow=0, blue=neg. orientation evaluated with double Classroom Examples ESA'04 44
2 D-Orientation of Three Points Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) 256 x 256 pixel image red=pos. , yellow=0, blue=neg. orientation evaluated with double Classroom Examples ESA'04 45
2 D-Orientation of Three Points Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) Keep y-coordinate fix Classroom Examples ESA'04 46
2 D-Orientation of Three Points Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) Keep y-coordinate fix: block sizes of 25 and 24 Classroom Examples ESA'04 47
2 D-Orientation of Three Points Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) 256 x 256 pixel image red=pos. , yellow=0, blue=neg. orientation evaluated with double Classroom Examples ESA'04 48
2 D-Orientation of Three Points Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) 256 x 256 pixel image red=pos. , yellow=0, blue=neg. orientation evaluated with double Classroom Examples ESA'04 49
2 D-Orientation of Three Points Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) 256 x 256 pixel image red=pos. , yellow=0, blue=neg. orientation evaluated with ext double Classroom Examples ESA'04 50
(Semi-) Systematic Search • a point outside sees no edge § p 1 = ( 0. 5, 0. 5) § p 2 = ( 7. 3000000194, 7. 3000000167) § p 3 = ( 24. 000000068, 24. 000000071) § p 4 = ( 24. 00000005, 24. 000000053) • • (p 2, p 3, p 4) form a counter-clockwise triangle Classification of (x(p 1) + x∙u, y(p 1) + y∙u) with respect to the edges (p 2, p 3) and (p 4, p 2). red = sees no edge, ocher = collinear with one, yellow = collinear with both, blue = sees one but not the other, green = sees both 51 Classroom Examples ESA'04
Choice of the Pivot? Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) p Classroom Examples ESA'04 52
Choice of the Pivot? Orientation( p, q, r) = sign((qx-px)(ry-py) – (qy-py)(rx-px)) p q r Classroom Examples ESA'04 53
Conclusion • Classroom examples Data sets and C++ programs available online: http: //www. mpi-sb. mpg. de/~kettner/proj/Non. Robust/ • Long version (available soon): § More analysis § Graham’s scan algorithm § 3 D Delaunay triangulation Classroom Examples ESA'04 54
- Slides: 54