Imprecision in Qhull

Up: Home page for Qhull
Up: Qhull manual: Table of Contents
To: Programs • Options • Output • Formats • Geomview • Print • Qhull • Precision • Trace • Functions
To: Qhull imprecision: Table of Contents (please wait while loading)

# Imprecision in Qhull

This section of the Qhull manual discusses the problems caused by coplanar points and why Qhull uses options 'C-0' or 'Qx' by default. If you ignore precision issues with option 'Q0', the output from Qhull can be arbitrarily bad. Qhull avoids precision problems if you merge facets (the default) or joggle the input ('QJ').

Use option 'Tv' to verify the output from Qhull. It verifies that adjacent facets are clearly convex. It verifies that all points are on or below all facets.

Qhull automatically tests for convexity if it detects precision errors while constructing the hull.

Copyright © 1995-2015 C.B. Barber

## » Qhull imprecision: Table of Contents

- Precision problems
- Merged facets or joggled input
- Delaunay triangulations
- Halfspace intersection/a>
- Merged facets
- How Qhull merges facets
- Limitations of merged facets
- Joggled input
- Exact arithmetic
- Approximating a convex hull

## » Precision problems

Since Qhull uses floating point arithmetic, roundoff error occurs with each calculation. This causes problems for geometric algorithms. Other floating point codes for convex hulls, Delaunay triangulations, and Voronoi diagrams also suffer from these problems. Qhull handles most of them.

There are several kinds of precision errors:

- Representation error occurs when there are not enough digits to represent a number, e.g., 1/3.
- Measurement error occurs when the input coordinates are from measurements.
- Roundoff error occurs when a calculation is rounded to a fixed number of digits, e.g., a floating point calculation.
- Approximation error occurs when the user wants an approximate result because the exact result contains too much detail.

Under imprecision, calculations may return erroneous results. For example, roundoff error can turn a small, positive number into a small, negative number. See Milenkovic ['93] for a discussion of strict robust geometry. Qhull does not meet Milenkovic's criterion for accuracy. Qhull's error bound is empirical instead of theoretical.

Qhull 1.0 checked for precision errors but did not handle them. The output could contain concave facets, facets with inverted orientation ("flipped" facets), more than two facets adjacent to a ridge, and two facets with exactly the same set of vertices.

Qhull 2.4 and later automatically handles errors due to machine round-off. Option 'C-0' or 'Qx' is set by default. In 5-d and higher, the output is clearly convex but an input point could be outside of the hull. This may be corrected by using option 'C-0', but then the output may contain wide facets.

Qhull 2.5 and later provides option 'QJ' to joggled input. Each input coordinate is modified by a small, random quantity. If a precision error occurs, a larger modification is tried. When no precision errors occur, Qhull is done.

Qhull 3.1 and later provides option 'Qt' for triangulated output. This removes the need for joggled input ('QJ'). Non-simplicial facets are triangulated. The facets may have zero area. Triangulated output is particularly useful for Delaunay triangulations.

By handling round-off errors, Qhull can provide a variety of output formats. For example, it can return the halfspace that defines each facet ('n'). The halfspaces include roundoff error. If the halfspaces were exact, their intersection would return the original extreme points. With imprecise halfspaces and exact arithmetic, nearly incident points may be returned for an original extreme point. By handling roundoff error, Qhull returns one intersection point for each of the original extreme points. Qhull may split or merge an extreme point, but this appears to be unlikely.

The following pipe implements the identity function for extreme points (with roundoff):

> qconvex FV n | qhalf Fp

Bernd Gartner published his Miniball algorithm ["Fast and robust smallest enclosing balls", Algorithms - ESA '99, LNCS 1643]. It uses floating point arithmetic and a carefully designed primitive operation. It is practical to 20-D or higher, and identifies at least two points on the convex hull of the input set. Like Qhull, it is an incremental algorithm that processes points furthest from the intermediate result and ignores points that are close to the intermediate result.

## » Merged facets or joggled input

This section discusses the choice between merged facets and joggled input. By default, Qhull uses merged facets to handle precision problems. With option 'QJ', the input is joggled. See examples of joggled input and triangulated output.

- Use merged facets (the default) when you want non-simplicial output (e.g., the faces of a cube).
- Use merged facets and triangulated output ('Qt') when you want simplicial output and coplanar facets (e.g., triangles for a Delaunay triangulation).
- Use joggled input ('QJ') when you need clearly-convex, simplicial output.

The choice between merged facets and joggled input depends on the application. Both run about the same speed. Joggled input may be faster if the initial joggle is sufficiently large to avoid precision errors.

Most applications should used merged facets with triangulated output.

Use merged facets (the default, 'C-0') or triangulated output ('Qt') if

- Your application supports non-simplicial facets, or it allows degenerate, simplicial facets (option 'Qt').
- You do not want the input modified.
- You want to set additional options for approximating the hull.
- You use single precision arithmetic (realT).

Use joggled input ('QJ') if

- Your application needs clearly convex, simplicial output
- Your application supports perturbed input points and narrow triangles.
- Seven significant digits is sufficient accuracy.

You may use both techniques or combine joggle with post-merging ('Cn').

Other researchers have used techniques similar to joggled input. Sullivan and Beichel [ref?] randomly perturb the input before computing the Delaunay triangulation. Corkum and Wyllie [news://comp.graphics, 1990] randomly rotate a polytope before testing point inclusion. Edelsbrunner and Mucke [Symp. Comp. Geo., 1988] and Yap [J. Comp. Sys. Sci., 1990] symbolically perturb the input to remove singularities.

Merged facets ('C-0') handles precision problems directly. If a precision problem occurs, Qhull merges one of the offending facets into one of its neighbors. Since all precision problems in Qhull are associated with one or more facets, Qhull will either fix the problem or attempt to merge the last remaining facets.

## » Delaunay triangulations

Programs that use Delaunay triangulations traditionally assume a triangulated input. By default, qdelaunay merges regions with cocircular or cospherical input sites. If you want a simplicial triangulation use triangulated output ('Qt') or joggled input ('QJ').

For Delaunay triangulations, triangulated output should produce good results. All points are within roundoff error of a paraboloid. If two points are nearly incident, one will be a coplanar point. So all points are clearly separated and convex. If qhull reports deleted vertices, the triangulation may contain serious precision faults. Deleted vertices are reported in the summary ('s', 'Fs'

You should use option 'Qbb' with Delaunay triangulations. It scales the last coordinate and may reduce roundoff error. It is automatically set for qdelaunay, qvoronoi, and option 'QJ'.

Edelsbrunner, H, Geometry and Topology for Mesh Generation, Cambridge University Press, 2001. Good mathematical treatise on Delaunay triangulation and mesh generation for 2-d and 3-d surfaces. The chapter on surface simplification is particularly interesting. It is similar to facet merging in Qhull.

Veron and Leon published an algorithm for shape preserving polyhedral simplification with bounded error [Computers and Graphics, 22.5:565-585, 1998]. It remove nodes using front propagation and multiple remeshing.

## » Halfspace intersection

The identity pipe for Qhull reveals some precision questions for halfspace intersections. The identity pipe creates the convex hull of a set of points and intersects the facets' hyperplanes. It should return the input points, but narrow distributions may drop points while offset distributions may add points. It may be better to normalize the input set about the origin. For example, compare the first results with the later two results: [T. Abraham]

> rbox 100 s t | tee r | qconvex FV n | qhalf Fp | cat - r | /bin/sort -n | tail
> rbox 100 L1e5 t | tee r | qconvex FV n | qhalf Fp | cat - r | /bin/sort -n | tail
> rbox 100 s O10 t | tee r | qconvex FV n | qhalf Fp | cat - r | /bin/sort -n | tail

## » Merged facets

Qhull detects precision problems when computing distances. A precision problem occurs if the distance computation is less than the maximum roundoff error. Qhull treats the result of a hyperplane computation as if it were exact.

Qhull handles precision problems by merging non-convex facets. The result of merging two facets is a thick facet defined by an inner plane and an outer plane. The inner and outer planes are offsets from the facet's hyperplane. The inner plane is clearly below the facet's vertices. At the end of Qhull, the outer planes are clearly above all input points. Any exact convex hull must lie between the inner and outer planes.

Qhull tests for convexity by computing a point for each facet. This point is called the facet's centrum. It is the arithmetic center of the facet's vertices projected to the facet's hyperplane. For simplicial facets with d vertices, the centrum is equivalent to the centroid or center of gravity.

Two neighboring facets are convex if each centrum is clearly below the other hyperplane. The 'Cn' or 'C-n' options sets the centrum's radius to n . A centrum is clearly below a hyperplane if the computed distance from the centrum to the hyperplane is greater than the centrum's radius plus two maximum roundoff errors. Two are required because the centrum can be the maximum roundoff error above its hyperplane and the distance computation can be high by the maximum roundoff error.

With the 'C-n' or 'A-n ' options, Qhull merges non-convex facets while constructing the hull. The remaining facets are clearly convex. With the 'Qx ' option, Qhull merges coplanar facets after constructing the hull. While constructing the hull, it merges coplanar horizon facets, flipped facets, concave facets and duplicated ridges. With 'Qx', coplanar points may be missed, but it appears to be unlikely.

If the user sets the 'An' or 'A-n' option, then all neighboring facets are clearly convex and the maximum (exact) cosine of an angle is n.

If 'C-0' or 'Qx' is used without other precision options (default), Qhull tests vertices instead of centrums for adjacent simplices. In 3-d, if simplex abc is adjacent to simplex bcd, Qhull tests that vertex a is clearly below simplex bcd , and vertex d is clearly below simplex abc. When building the hull, Qhull tests vertices if the horizon is simplicial and no merges occur.

## » How Qhull merges facets

If two facets are not clearly convex, then Qhull removes one or the other facet by merging the facet into a neighbor. It selects the merge which minimizes the distance from the neighboring hyperplane to the facet's vertices. Qhull also performs merges when a facet has fewer than d neighbors (called a degenerate facet), when a facet's vertices are included in a neighboring facet's vertices (called a redundant facet), when a facet's orientation is flipped, or when a ridge occurs between more than two facets.

Qhull performs merges in a series of passes sorted by merge angle. Each pass merges those facets which haven't already been merged in that pass. After a pass, Qhull checks for redundant vertices. For example, if a vertex has only two neighbors in 3-d, the vertex is redundant and Qhull merges it into an adjacent vertex.

Merging two simplicial facets creates a non-simplicial facet of d+1 vertices. Additional merges create larger facets. When merging facet A into facet B, Qhull retains facet B's hyperplane. It merges the vertices, neighbors, and ridges of both facets. It recomputes the centrum if a wide merge has not occurred (qh_WIDEcoplanar) and the number of extra vertices is smaller than a constant (qh_MAXnewcentrum).

## » Limitations of merged facets

- Uneven dimensions -- If one coordinate has a larger absolute value than other coordinates, it may dominate the effect of roundoff errors on distance computations. You may use option 'QbB' to scale points to the unit cube. For Delaunay triangulations and Voronoi diagrams, qdelaunay and qvoronoi always set option 'Qbb'. It scales the last coordinate to [0,m] where m is the maximum width of the other coordinates. Option 'Qbb' is needed for Delaunay triangulations of integer coordinates and nearly cocircular points.

  For example, compare

              rbox 1000 W0 t | qconvex Qb2:-1e-14B2:1e-14

  with

              rbox 1000 W0 t | qconvex

  The distributions are the same but the first is compressed to a 2e-14 slab.

- Post-merging of coplanar facets -- In 5-d and higher, option 'Qx' (default) delays merging of coplanar facets until post-merging. This may allow "dents" to occur in the intermediate convex hulls. A point may be poorly partitioned and force a poor approximation. See option 'Qx' for further discussion.

  This is difficult to produce in 5-d and higher. Option 'Q6' turns off merging of concave facets. This is similar to 'Qx'. It may lead to serious precision errors, for example,

              rbox 10000 W1e-13  | qhull Q6  Tv

- Maximum facet width -- Qhull reports the maximum outer plane and inner planes (if more than roundoff error apart). There is no upper bound for either figure. This is an area for further research. Qhull does a good job of post-merging in all dimensions. Qhull does a good job of pre-merging in 2-d, 3-d, and 4-d. With the 'Qx' option, it does a good job in higher dimensions. In 5-d and higher, Qhull does poorly at detecting redundant vertices.

  In the summary ('s'), look at the ratio between the maximum facet width and the maximum width of a single merge, e.g., "(3.4x)". Qhull usually reports a ratio of four or lower in 3-d and six or lower in 4-d. If it reports a ratio greater than 10, this may indicate an implementation error. Narrow distributions (see following) may produce wide facets.

  For example, if special processing for narrow distributions is turned off ('Q10'), qhull may produce a wide facet:

               rbox 1000 L100000 s G1e-16 t1002074964 | qhull Tv Q10

- Narrow distribution -- In 3-d, a narrow distribution may result in a poor approximation. For example, if you do not use qdelaunay nor option 'Qbb', the furthest-site Delaunay triangulation of nearly cocircular points may produce a poor approximation:

               rbox s 5000 W1e-13 D2 t1002151341 | qhull d Qt
               rbox 1000 s W1e-13 t1002231672 | qhull d Tv

  During construction of the hull, a point may be above two facets with opposite orientations that span the input set. Even though the point may be nearly coplanar with both facets, and can be distant from the precise convex hull of the input sites. Additional facets leave the point distant from a facet. To fix this problem, add option 'Qbb' (it scales the last coordinate). Option 'Qbb' is automatically set for qdelaunay and qvoronoi.

  Qhull generates a warning if the initial simplex is narrow. For narrow distributions, Qhull changes how it processes coplanar points -- it does not make a point coplanar until the hull is finished. Use option 'Q10' to try Qhull without special processing for narrow distributions. For example, special processing is needed for:

               rbox 1000 L100000 s G1e-16 t1002074964 | qhull Tv Q10

  You may turn off the warning message by reducing qh_WARNnarrow in `user.h` or by setting option 'Pp'.

  Similar problems occur for distributions with a large flat facet surrounded with many small facet at a sharp angle to the large facet. Qhull 3.1 fixes most of these problems, but a poor approximation can occur. A point may be left outside of the convex hull ('Tv'). Examples include the furthest-site Delaunay triangulation of nearly cocircular points plus the origin, and the convex hull of a cone of nearly cocircular points. The width of the band is 10^-13.

              rbox s 1000 W1e-13 P0 D2 t996799242 | qhull d Tv
              rbox 1000 s Z1 G1e-13 t1002152123 | qhull Tv
              rbox 1000 s Z1 G1e-13 t1002231668 | qhull Tv

- Quadratic running time -- If the output contains large, non-simplicial facets, the running time for Qhull may be quadratic in the size of the triangulated output. For example, `rbox 1000 s W1e-13 c G2 | qhull d` is 4 times faster for 500 points. The convex hull contains two large nearly spherical facets and many nearly coplanar facets. Each new point retriangulates the spherical facet and repartitions the remaining points into all of the nearly coplanar facets. In this case, quadratic running time is avoided if you use qdelaunay, add option 'Qbb', or add the origin ('P0') to the input.

- Nearly coincident points on an edge -- Multiple, nearly coincident points on an edge with non-simplicial facets may lead to wide facets or quadratic running time. For example, either the convex hull in 4-D of a narrow lens with coincident points or the Delaunay triangulation of nearly coincident points may lead to very wide facets (e.g., 2267021951.3x).

  For Delaunay triangulations an edge occurs between the corresponding upper and lower convex hull. Points on this edge correspond to vertices of the convex hull of the input sites. After multiple facet merges, four facets may share the same, duplicate ridge and must be merged. Some of these facets may be long and narrow, leading to a very wide merged facet. Starting with Qhull 2015.1, an error is reported.

  Duplicate ridges occur when the horizon facets for a new point is "pinched". In a duplicate ridge, a subridge (e.g., a line segment in 3-d) is shared by two horizon facets. At least two of its vertices are nearly coincident. Qhull 2015.1 adds the 'Cn,r,m' option to rbox. This option generates nearly coincident points. For example, every point of the following distributions has a nearly coincident point within a 1e-13 ball. Substantially smaller or larger balls do not lead to pinched horizons.

              rbox L100 1000 D4 s C1,1e-13 t | qhull
              rbox 75 C1,1E-13 t | qhull d

  For Delaunay triangulations, a bounding box may alleviate this error (e.g., `rbox 500 C1,1E-13 t c G1 | qhull d`). A later release of qhull will avoid pinched horizons by merging duplicate subridges. A subridge is merged by merging adjacent vertices.

- Facet with zero-area -- It is possible for a zero-area facet to be convex with its neighbors. This can occur if the hyperplanes of neighboring facets are above the facet's centrum, and the facet's hyperplane is above the neighboring centrums. Qhull computes the facet's hyperplane so that it passes through the facet's vertices. The vertices can be collinear.

- No more facets -- Qhull reports an error if there are d+1 facets left and two of the facets are not clearly convex. This typically occurs when the convexity constraints are too strong or the input points are degenerate. The former is more likely in 5-d and higher -- especially with option 'C-n'.

- Deleted cone -- Lots of merging can end up deleting all of the new facets for a point. This is a rare event that has only been seen while debugging the code.

- Triangulated output leads to precision problems -- With sufficient merging, the ridges of a non-simplicial facet may have serious topological and geometric problems. A ridge may be between more than two neighboring facets. If so, their triangulation ('Qt') will fail since two facets have the same vertex set. Furthermore, a triangulated facet may have flipped orientation compared to its neighbors.

The triangulation process detects degenerate facets with only two neighbors. These are marked degenerate. They have zero area.

- Coplanar points -- Option 'Qc' is determined by qh_check_maxout() after constructing the hull. Qhull needs to retain all possible coplanar points in the facets' coplanar sets. This depends on qh_RATIOnearInside in `user.h.` Furthermore, the cutoff for a coplanar point is arbitrarily set at the minimum vertex. If coplanar points are important to your application, remove the interior points by hand (set 'Qc Qi') or make qh_RATIOnearInside sufficiently large.

- Maximum roundoff error -- Qhull computes the maximum roundoff error from the maximum coordinates of the point set. Usually the maximum roundoff error is a reasonable choice for all distance computations. The maximum roundoff error could be computed separately for each point or for each distance computation. This is expensive and it conflicts with option 'C-n'.

- All flipped or upper Delaunay -- When a lot of merging occurs for Delaunay triangulations, a new point may lead to no good facets. For example, try a strong convexity constraint:

              rbox 1000 s t993602376 | qdelaunay C-1e-3

## » Joggled input

Joggled input is a simple work-around for precision problems in computational geometry ["joggle: to shake or jar slightly," Amer. Heritage Dictionary]. Other names are jostled input or random perturbation. Qhull joggles the input by modifying each coordinate by a small random quantity. If a precision problem occurs, Qhull joggles the input with a larger quantity and the algorithm is restarted. This process continues until no precision problems occur. Unless all inputs incur precision problems, Qhull will terminate. Qhull adjusts the inner and outer planes to account for the joggled input.

Neither joggle nor merged facets has an upper bound for the width of the output facets, but both methods work well in practice. Joggled input is easier to justify. Precision errors occur when the points are nearly singular. For example, four points may be coplanar or three points may be collinear. Consider a line and an incident point. A precision error occurs if the point is within some epsilon of the line. Now joggle the point away from the line by a small, uniformly distributed, random quantity. If the point is changed by more than epsilon, the precision error is avoided. The probability of this event depends on the maximum joggle. Once the maximum joggle is larger than epsilon, doubling the maximum joggle will halve the probability of a precision error.

With actual data, an analysis would need to account for each point changing independently and other computations. It is easier to determine the probabilities empirically ('TRn') . For example, consider computing the convex hull of the unit cube centered on the origin. The arithmetic has 16 significant decimal digits.

> Convex hull of unit cube

> joggle | error prob.
> --- | ---
> 1.0e-15 | 0.983
> 2.0e-15 | 0.830
> 4.0e-15 | 0.561
> 8.0e-15 | 0.325
> 1.6e-14 | 0.185
> 3.2e-14 | 0.099
> 6.4e-14 | 0.051
> 1.3e-13 | 0.025
> 2.6e-13 | 0.010
> 5.1e-13 | 0.004
> 1.0e-12 | 0.002
> 2.0e-12 | 0.001

A larger joggle is needed for multiple points. Since the number of potential singularities increases, the probability of one or more precision errors increases. Here is an example.

> Convex hull of 1000 points on unit cube

> joggle | error prob.
> --- | ---
> 1.0e-12 | 0.870
> 2.0e-12 | 0.700
> 4.0e-12 | 0.450
> 8.0e-12 | 0.250
> 1.6e-11 | 0.110
> 3.2e-11 | 0.065
> 6.4e-11 | 0.030
> 1.3e-10 | 0.010
> 2.6e-10 | 0.008
> 5.1e-09 | 0.003

Other distributions behave similarly. No distribution should behave significantly worse. In Euclidean space, the probability measure of all singularities is zero. With floating point numbers, the probability of a singularity is non-zero. With sufficient digits, the probability of a singularity is extremely small for random data. For a sufficiently large joggle, all data is nearly random data.

Qhull uses an initial joggle of 30,000 times the maximum roundoff error for a distance computation. This avoids most potential singularities. If a failure occurs, Qhull retries at the initial joggle (in case bad luck occurred). If it occurs again, Qhull increases the joggle by ten-fold and tries again. This process repeats until the joggle is a hundredth of the width of the input points. Qhull reports an error after 100 attempts. This should never happen with double-precision arithmetic. Once the probability of success is non-zero, the probability of success increases about ten-fold at each iteration. The probability of repeated failures becomes extremely small.

Merged facets produces a significantly better approximation. Empirically, the maximum separation between inner and outer facets is about 30 times the maximum roundoff error for a distance computation. This is about 2,000 times better than joggled input. Most applications though will not notice the difference.

## » Exact arithmetic

Exact arithmetic may be used instead of floating point. Singularities such as coplanar points can either be handled directly or the input can be symbolically perturbed. Using exact arithmetic is slower than using floating point arithmetic and the output may take more space. Chaining a sequence of operations increases the time and space required. Some operations are difficult to do.

Clarkson's hull program and Shewchuk's triangle program are practical implementations of exact arithmetic.

Clarkson limits the input precision to about fifteen digits. This reduces the number of nearly singular computations. When a determinant is nearly singular, he uses exact arithmetic to compute a precise result.

## » Approximating a convex hull

Qhull may be used for approximating a convex hull. This is particularly valuable in 5-d and higher where hulls can be immense. You can use 'Qx C-n' to merge facets as the hull is being constructed. Then use 'Cn' and/or 'An' to merge small facets during post-processing. You can print the n largest facets with option 'PAn'. You can print facets whose area is at least n with option 'PFn'. You can output the outer planes and an interior point with 'FV Fo' and then compute their intersection with 'qhalf'.

To approximate a convex hull in 6-d and higher, use post-merging with 'Wn' (e.g., qhull W1e-1 C1e-2 TF2000). Pre-merging with a convexity constraint (e.g., qhull Qx C-1e-2) often produces a poor approximation or terminates with a simplex. Option 'QbB' may help to spread out the data.

You will need to experiment to determine a satisfactory set of options. Use rbox to generate test sets quickly and Geomview to view the results. You will probably want to write your own driver for Qhull using the Qhull library. For example, you could select the largest facet in each quadrant.

Up: Home page for Qhull
Up: Qhull manual: Table of Contents
To: Programs • Options • Output • Formats • Geomview • Print • Qhull • Precision • Trace • Functions
To: Qhull imprecision: Table of Contents

The Geometry Center Home Page

Comments to: qhull@qhull.org
Created: Sept. 25, 1995 --- Last modified: see top