6.838 March 16 Meeting
|1. Robust Predicates|
See Shewchuk's Robust Predicates page.
See: Notes on IEEE 754 Floating Point Arithmetic (postscript)
252 numbers between tick marks. Smallest number > 1 is
1 + 2-52.
"Island" of high precision around 0.
fl(1+2-53) = 1
Predicates implemented using floating point might give wrong answers:
get unexpected or inconsistent state of a program.
Why? Cannot distinguish the two cases below if a = 1:
| orientation | = (a + e) - a
Example: a = 1, e = 2-53: | orientation | = 0 !
Not invariant w.r.t transformations (a = 0, 1).
GEOMETRIC = NUMERIC + COMBINATORIAL
Geometric algorithm : primitive tests on continuous data --> combinatorial output [Ex. convex hull]
There is a consistency condition between numerical and combinatorial data.
Can geometric algorithms be implemented using floating point arithmetic? At all?
Assumption: input is exact.
Numbers can be represented implicitly.
Evaluate all predicates exactly or up to the required precision (ex. sign = 1bit).
What kind of numbers [subgroups of R] do we need?
Given integer input, what is the smallest class of numbers sufficient for
See survey by Dube and Yap (postscript) for survey and references.
Numerical quantities are computed to sufficient precision so that the geometric predicates return correct values. (Underlying combinatorial structure is mathematically exact.)
Geometric predicates are expressions over input coordinates.
How to speed up exact arithmetic? Avoid computing all digits.
Must specify desired precision for the expressions: ex. 1 bit = sign.
First idea: evaluate expression in fp, then exactly if needed [Fortune and van Wyk].
Second idea: divide expression into subexpressions, evaluate incrementally.
Rational Numbers, finitely representable.
Priest & Shewchuk:
non-overlapping, floating point terms.
Idea: allow roundoff error, but account for it after the fact by using more terms.
Example: (3 mantissa bit floating point numbers):
101.0011 = 101 + 11. 2-4
Example: 101 + 0.0011 = 101 + roundoff error terms
Operations: +, -, *, /
Highest order term is a crude approximation.
Evaluate expression using MST first, then add additional terms as needed.
Algebraic Numbers. Real/Expr, LEDA
Example: Need to compare an expression to 0.
Need precision bound for root (ex. 1 bit).
1. Propagate precision bound top down (from root to leafs).
Example: a+b: if precision bound for (a+b) is 0.02, then precision bound for a and b is 0.01.
At each leaf: procedure to extract an approximate value to satisfy precision bound at the leaf.
2. Propagate value up.
Changing precision bound at root may force (incremental) recomputation.
To how many bits do we need to evaluate an expression to get at least 1 significant bit?
Example: Graham's problem:
a - b = 0 ? (a=b to 40 digits!!)
Algebraic numbers are roots of polynomials (characteristic polynomial.)
How close is an expression to 0? --> Root bounds.
For each node in expression tree maintain
Expression tree keeps track of required precision, which is propagated down the tree to know how precisely a given subexpression (or number) has to be computed.
Disadvantages: bit size grows, speed about 4 to 10 times slower than fp.
Robustness: A geometric algorithm (using f.p.) is robust if
Stability: An algorithm is stable if it is robust and the required perturbation is small.
Applications: line arrangements, 2D Delaunay triangulation (get
"almost" delaunay property), convex hull.
Difficulty: careful error analysis of predicates, careful design of algorithm.
Predicates : independent and dependent.
Numbers are represented by an interval
Interval arithmetic allows computations to be performed on intervals, such that the result of the computation provides a guaranteed enclosure for the correct mathematical result over the domain defined by the the interval operands.
Basic arithmetic operations , defined for the the intervals X and Y as all possible results
unless and .
Example: [a,b] + [c,d] = [a+c, b+d]
Implementation rounds each interval endpoint outwards.
Predicates can return YES, NO, MAYBE (it is not possible to tell within precision available.)
Properties: (inclusion isotomy)
about 3 times slower.
(thick lines) [qhull]
Fatten geometric objects, use tolerancing zones (conservative estimates.)
Predicates can return "don't know".
(polylines replace lines)
Example: Line arrangement construction.
Every intersection between two lines is "snapped" to a grid point.
(line coeffs rounded).
Goal: preserve combinatorial structure.
Fact: rounding a collection of disjoint simple polygons while preserving combinatorial structure is NP-complete.
Tolerances / epsilons useful in practice, but no proofs of correctness.
Consider an expression = 0 if |expression| < epsilon.
Tweak until your program works :)
Partity Algorithm: is a point inside or outside a polygon?
Voronoi Diagram: co-circular points.
Convex Hull: co-planar points. Is a point inside? on a face? on an edge? on a point?
Def.: A degeneracy occurs when arbitrarily small perturbation of the input can result in qualitatively different geometric structures. It is in some sense a discontinuity.
Inherent and algorithm-induced degeneracies.
Ex.: In divide and conquer algorithm for convex hull, colinear points on convex hull are an inherent degeneracy, collinear points inside are an induced degeneracy.
What is general position? (Depends on predicates: no predicates return 0).
Algorithm working for generic input --> algorithm working for general input.
Let a* be a non-degenerate instance of a problem.
Given a possibly degenerate input instance a, the instance (a + e a*) is non-degenerate for infinitesimal e.
Trivial perturbation: (ex. sorting)
(xj- xi) ? 0. Want >0 or <0
if xi = xj, return "xi < xj" iff i<j
Let e = very small number, with (xi+ei) != xj for any i,j.
( (xj+ ej) - (xi+ei) ) < 0 if xi = xj and i<j,
General comparisons: xij are variables (i=point id, j=coordinate).
Predicates = polynomials in x's and e's.
work with `high' probability
Never compute the values.
Compare x's, then successively higher degrees of e.
|Test query point if on edge|
|Convex hull||Merge adjacent faces if coplanar
Remove faces subset of larger faces
|Voronoi diagram||Merge Voronoi vertices generated by co-circular points
Remove 0 length edges