Bug in volume computation of polyhedron

[videlec]

sage: polytopes.dodecahedron(base_ring=RDF).volume()
Traceback (most recent call last):
...
ZeroDivisionError: input matrix must be nonsingular

The same error is triggered by

sage: D=polytopes.dodecahedron(base_ring=RDF)
sage: D.triangulate()

CC: @jplab @mo271 @tom111

Component: geometry

Keywords: bug

Author: Dima Pasechnik

Branch/Commit: b6758cd

Reviewer: Matthias Koeppe,Frédéric Chapoton

Issue created by migration from https://trac.sagemath.org/ticket/18214

[videlec]

Changed keywords from none to bug

[fchapoton]

comment:4

Here is a proposal for a solution. In case of inexact base ring, I introduce a non-zero epsilon constant, to take care of approximate equality of points. Maybe this could become a parameter.

---

New commits:

d409a4c

trac 18214 triangulation over inexact base rings

[fchapoton]

Branch: u/chapoton/18214

[fchapoton]

Commit: d409a4c

[fchapoton]

Author: Frédéric Chapoton

[mkoeppe]

comment:5

I think you also need to change method PointConfiguration.contained_simplex.

And that epsilon should indeed be customizable. And when PointConfiguration is called from Polyhedron_RDF, it should use the epsilons from Polyhedron_RDF._is_zero etc.

[mkoeppe]

Reviewer: Matthias Koeppe

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

11dc033

Merge branch 'u/chapoton/18214' in 7.5.b3

[sagetrac-git]

Changed commit from d409a4c to 11dc033

[sagetrac-git]

Branch pushed to git repo; I updated commit sha1. New commits:

f086a14	Merge branch 'u/chapoton/18214' in 8.2.b6
e042c42	trac 18214 adding abs tol to volume doctest

[sagetrac-git]

Changed commit from 11dc033 to e042c42

[fchapoton]

New commits:

f086a14	Merge branch 'u/chapoton/18214' in 8.2.b6
e042c42	trac 18214 adding abs tol to volume doctest

[vbraun]

comment:13

Is this still reproducable? I get

sage: polytopes.dodecahedron(base_ring=RDF).volume()
6.4520359589542045

[videlec]

comment:14

Replying to @vbraun:

Is this still reproducable? I get

sage: polytopes.dodecahedron(base_ring=RDF).volume()
6.4520359589542045

On my computer it is.

[fchapoton]

comment:15

patchbot petitbonum is seeing this as a random failure:

./sage -t --long --warn-long 54.5 src/sage/geometry/polyhedron/library.py
...
sage: ico.volume()
...
      File "/home/chapoton/sage/local/lib/python2.7/site-packages/sage/geometry/triangulation/point_configuration.py", line 1994, in facets_of_simplex
        normals = span.inverse().columns()
      File "sage/matrix/matrix2.pyx", line 8854, in sage.matrix.matrix2.Matrix.inverse (build/cythonized/sage/matrix/matrix2.c:71887)
        return ~self
      File "sage/matrix/matrix_double_dense.pyx", line 466, in sage.matrix.matrix_double_dense.Matrix_double_dense.__invert__ (build/cythonized/sage/matrix/matrix_double_dense.c:6026)
        raise ZeroDivisionError("input matrix must be nonsingular")
    ZeroDivisionError: input matrix must be nonsingular

[dimpase]

comment:16

Replying to @videlec:

Replying to @vbraun:

Is this still reproducable? I get

sage: polytopes.dodecahedron(base_ring=RDF).volume()
6.4520359589542045

On my computer it is.

Also reproducible on my Google CE host, with 8 Xeon CPUs:
$ uname -a
Linux linux 4.9.0-8-amd64 #1 SMP Debian 4.9.130-2 (2018-10-27) x86_64 GNU/Linux

[dimpase]

comment:17

with random test order I get more errors like this:

...
File "src/sage/geometry/polyhedron/library.py", line 1158, in sage.geometry.polyhedron.library.Polytopes.truncated_dodecahedron
Failed example:
    td = polytopes.truncated_dodecahedron(exact=False)
Expected:
    doctest:warning
    ...
    UserWarning: This polyhedron data is numerically complicated; cdd could not convert between the inexact V and H representation without loss of data. The resulting object might show inconsistencies.
Got:
    <BLANKLINE>
...
2 items had failures:
   1 of   9 in sage.geometry.polyhedron.library.Polytopes.truncated_dodecahedron
   1 of   7 in sage.geometry.polyhedron.library.Polytopes.truncated_icosidodecahedron
    [196 tests, 2 failures, 9.83 s]
----------------------------------------------------------------------
sage -t src/sage/geometry/polyhedron/library.py  # 2 doctests failed
----------------------------------------------------------------------
Total time for all tests: 10.0 seconds
    cpu time: 7.9 seconds
    cumulative wall time: 9.8 seconds
dimpase@linux:~/sage$ ./sage -tp --randorder 127 src/sage/geometry/polyhedron/library.py

[vbraun]

comment:18

Degenerate polyhedra over floating-point numbers aren't really supported; Usually you are lucky and numerical errors triangulate your faces. But if you are unlucky then numerical fuzz leads to inconsistent incidence relations where a point is on both (or neither) side of a plane. Just adding fuzzy zero is afaik not enough.

[dimpase]

comment:19

This seems to be state of the art: https://dl.acm.org/citation.cfm?id=3194656

Surely, naive triangulation approach is doomed to fail...
Should we just mark the corresponding tests are "known bug" and move on?

[dimpase]

comment:20

And here is the code: https://github.com/GeomScale/volume_approximation (it's C++).

[dimpase]

Changed branch from u/chapoton/18214 to u/dimpase/GQ

[dimpase]

Changed author from Frédéric Chapoton to Frédéric Chapoton, Dima Pasechnik

[dimpase]

Changed commit from e042c42 to b6758cd

[dimpase]

comment:21

see #27056 for a follow-up

---

New commits:

b6758cd

volume computation of fuzzy RDF polytope might fail

[fchapoton]

comment:22

ok

[fchapoton]

Changed author from Frédéric Chapoton, Dima Pasechnik to Dima Pasechnik

[fchapoton]

Changed reviewer from Matthias Koeppe to Matthias Koeppe,Frédéric Chapoton

[vbraun]

Changed branch from u/dimpase/GQ to b6758cd