Polyhedron_normaliz._triangulate_normaliz should not use NmzResult directly

[mkoeppe]

_triangulate_normaliz calls ._cone_generators, which calls NmzResult directly instead of going through ._nmz_result.

This relies on the persistence of the rational and number field handlers in PyNormaliz, which is subject to change in future PyNormaliz versions.

[gh-kliem]

This has nothing to do with the future. It is already a bug.

sage: P = polytopes.cube(backend='normaliz')*AA(2).sqrt()                                                                                                                                                                                                                    
sage: P._triangulate_normaliz()                                                                                                                                                                                                                                              
---------------------------------------------------------------------------
NormalizInterfaceError                    Traceback (most recent call last)
<ipython-input-7-ce2270188ea9> in <module>
----> 1 P._triangulate_normaliz()

/srv/public/kliem/sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/backend_normaliz.py in _triangulate_normaliz(self)
   1583         data.pop('inhom_equations', None)
   1584         data.pop('inhom_inequalities', None)
-> 1585         cone = self._make_normaliz_cone(data)
   1586 
   1587         nmz_triangulation = self._nmz_result(cone, "Triangulation")

/srv/public/kliem/sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/backend_normaliz.py in _make_normaliz_cone(data, verbose)
   1110         if verbose:
   1111             print("# Calling PyNormaliz.NmzCone(**{})".format(data))
-> 1112         cone = PyNormaliz.NmzCone(**data)
   1113         assert cone, "NmzCone(**{}) did not return a cone".format(data)
   1114         return cone

NormalizInterfaceError: When parsing vertices: coefficient in matrix must be PyFloat, PyInt, PyLong, Sequence, must be able to be converted to a valid gmp input string

[gh-kliem]

I guess normaliz can handle this by now:

sage: P._nmz_result(P._normaliz_cone, "Triangulation")                                                                                                                                                                                                                       
[[[0, 1, 2, 4], 16*a],
 [[1, 2, 4, 3], 16*a],
 [[1, 3, 4, 5], 16*a],
 [[3, 5, 6, 7], 16*a],
 [[6, 2, 4, 3], 16*a],
 [[6, 3, 4, 5], 16*a]]

This is exactly what we want.

[gh-kliem]

New commits:

​71f169f	fix normaliz triangulation
​bb90807	documentation

[git]

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

​10b5312

typo

[gh-kliem]

Is homogenous correct in any context?

git grep homogenous reveals quite something. Could be on the list of common typos that are discovered by this tox testing framework (and should of course be a typo ticket before that, but is annoying as one needs to deprecate at least one function).

[Winfried]

I think the correct spelling is "homogeneous".

If there is any problem on the Normaliz/PyNormaliz? side, let me know.

[gh-kliem]

Normaliz works all right.

Is there a way to ignore the dehomogenization for Triangulation of non-compact polyhedra? I currently recreate the cone from precomputed data without specifying the dehomogenization to achieve this.

[Winfried]

I think you want to get a triangulation (of the homogenization) even if the input defines an unbounded polyhedron. At present forbidden. It would be no real problem to do do this computation, but it requires a change in the source code of Normaliz. I keep it in mind.

[mkoeppe]

Is this ready for review? The method _cone_generators of the ticket description is unchanged

[gh-kliem]

I guess the static method _cone_generators can just be removed. What do you think about it? It's not just anywhere else as much as I can see.

The correct way to do it anyway is by self._nmz_result(self._normaliz_cone, "Generators"). I wasn't sure about how to resolve this issue for a static method.

[mkoeppe]

If it's not used any more, yes, please remove it

[git]

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

​b45d775

remove unused private method that incorrectly uses the API

[Winfried]

I have changed the Normaliz code so that triangulations can be computed for arbitrary input (just pushed to GitHub/master; will be part of version 3.8.9). It is clear that this notion must be interpreted with care if the input defines an unbounded polyhedron. A corresponding remark has been added to Normaliz.pdf (p. 121 in the current version).

I am not sure whether the computation in which the triangulations come up, are genuinely inhomogeneous, or the inhomogeneity is due to the implementation of the Sage polyhedron class. As mentioned earlier, it would be really good if the Sage cone class could use Normaliz as a backend.

[Winfried]

IMPORTANT: Triangulations are defined with respect to an ordered set of generators. The only reliable interactive way to connect the triangulation with the generators is to ask for them AFTER the computation of the triangulation since this computation may have changed the order (or even the set). Normaliz knows this problem and therefore stores the generators after the computation of the triangulation for future reference. So with disciplined use there is no problem.

But often generators are already defined before the computation of the triangulation, and if users ask for them too early, chaos can ensue. For this reason I will introduce the matrix TriangulationGenerators? that stores the ordered set of generators used in building the triangulation. The name should help to avoid confusion.

[mkoeppe]

I've opened #30569 for the update to the upcoming version and left a note there regarding the triangulations

[gh-kliem]

Replying to Winfried:

IMPORTANT: Triangulations are defined with respect to an ordered set of generators. The only reliable interactive way to connect the triangulation with the generators is to ask for them AFTER the computation of the triangulation since this computation may have changed the order (or even the set). Normaliz knows this problem and therefore stores the generators after the computation of the triangulation for future reference. So with disciplined use there is no problem.

But often generators are already defined before the computation of the triangulation, and if users ask for them too early, chaos can ensue. For this reason I will introduce the matrix TriangulationGenerators? that stores the ordered set of generators used in building the triangulation. The name should help to avoid confusion.

From my point of view, this means that this ticket should go back to needs work.

[gh-kliem]

@Winfried: Could you maybe provide an example where the order changes or point me towards such a direction.

Or is it just in principal a wrong assumption for the API to make and this could go wrong at some point.

If I initialize with precomputed data and random order of vertices, calling Triangulation doesn't change the order (at least for the examples I tried).

[Winfried]

Look at

>>> import PyNormaliz
>>> from PyNormaliz import *
>>> vert = [[2, -1, -1, -2, 1], [2, -1, 0, -2, 1], [2, -1, 1, -2, 1], [2, 0, -1, -2, 1], [2, 0, 0, -2, 1], [2, 0, 1, -2, 1], [2, 1, -1, -2, 1], [2, 1, 0, -2, 1], [2, 1, 1, -2, 1]]
>>> C= Cone(vertices = vert)
>>> C.Generators()
[[2, -1, -1, -2, 1], [2, -1, 0, -2, 1], [2, -1, 1, -2, 1], [2, 0, -1, -2, 1], [2, 0, 0, -2, 1], [2, 0, 1, -2, 1], [2, 1, -1, -2, 1], [2, 1, 0, -2, 1], [2, 1, 1, -2, 1]]
>>> C.Triangulation()
[[[0, 1, 2], 1], [[0, 1, 3], 1], [[0, 2, 4], 1], [[0, 3, 4], 1], [[1, 2, 5], 1], [[1, 3, 6], 1], [[2, 4, 7], 1], [[3, 4, 8], 1]]
>>> C.Generators()
[[2, 0, 0, -2, 1], [2, -1, 0, -2, 1], [2, 0, -1, -2, 1], [2, 0, 1, -2, 1], [2, 1, 0, -2, 1], [2, -1, -1, -2, 1], [2, -1, 1, -2, 1], [2, 1, -1, -2, 1], [2, 1, 1, -2, 1]]
>>>

I would not bet that this impossible with precomputed generators. I stumbled on this problem when the last Travis run of PyNormaliz? failed because the triangulation had changed relative to the reference. In the Normaliz workflow the triangulation can be computed at two places, and this place has changed with the last commit in order to make triangulations available also for unbounded polyhedra.

If the triangulation is computed together with the Ehrhart series (==> old place), then we get a different result:

>>> D= Cone(vertices = vert)
>>> D.Compute("EhrhartSeries", "Triangulation")
True
>>> D.Generators()
[[2, -1, -1, -2, 1], [2, -1, 0, -2, 1], [2, -1, 1, -2, 1], [2, 0, -1, -2, 1], [2, 0, 0, -2, 1], [2, 0, 1, -2, 1], [2, 1, -1, -2, 1], [2, 1, 0, -2, 1], [2, 1, 1, -2, 1]]
>>> D.Triangulation()
[[[0, 1, 3], 1], [[1, 2, 3], 1], [[2, 3, 4], 1], [[2, 4, 5], 1], [[3, 4, 6], 1], [[4, 5, 6], 1], [[5, 6, 7], 1], [[5, 7, 8], 1]]
>>> D.Generators()
[[2, -1, -1, -2, 1], [2, -1, 0, -2, 1], [2, -1, 1, -2, 1], [2, 0, -1, -2, 1], [2, 0, 0, -2, 1], [2, 0, 1, -2, 1], [2, 1, -1, -2, 1], [2, 1, 0, -2, 1], [2, 1, 1, -2, 1]]
>>>

In this case the order has not changed. Both results are correct --- triangulations are not uniquely determined and there is no canonical order of the generators. The next days I will try to find out what is really going on. Nevertheless I think it is really better to introduce TriangulationGenereors? that are not defined before the triangulation has been computed.

[gh-kliem]

Thanks. You were faster than me. I just wanted to post that I stumped upon this while trying to get a doctest for algebraic polyhedra (that finding the correct indices also works for them).

The current branch produces:

sage: C1 = Polyhedron(rays=[[0,0,1],[1,0,AA(2).sqrt()],[0,1,1],[1,1,1]], backend='normaliz')                                                                                        
sage: C1.triangulate(engine='normaliz')                                                                                                                                             
(<0,1,2>, <0,2,3>)
sage: list(a.ambient_Vrepresentation() for a in C1.face_generator())                                                                                                                
[(A vertex at (0, 0, 0),
  A ray in the direction (0, 0, 1),
  A ray in the direction (0, 1, 1),
  A ray in the direction (0.7071067811865475?, 0, 1),
  A ray in the direction (1, 1, 1)),
 (),
 (A vertex at (0, 0, 0),
  A ray in the direction (0, 0, 1),
  A ray in the direction (0, 1, 1)),
 (A vertex at (0, 0, 0),
  A ray in the direction (0, 0, 1),
  A ray in the direction (0.7071067811865475?, 0, 1)),
 (A vertex at (0, 0, 0),
  A ray in the direction (0, 1, 1),
  A ray in the direction (1, 1, 1)),
 (A vertex at (0, 0, 0),
  A ray in the direction (0.7071067811865475?, 0, 1),
  A ray in the direction (1, 1, 1)),
 (A vertex at (0, 0, 0), A ray in the direction (0, 0, 1)),
 (A vertex at (0, 0, 0), A ray in the direction (0, 1, 1)),
 (A vertex at (0, 0, 0),),
 (A vertex at (0, 0, 0), A ray in the direction (0.7071067811865475?, 0, 1)),
 (A vertex at (0, 0, 0), A ray in the direction (1, 1, 1))]
sage: C1.rays()                                                                                                                                                                     
(A ray in the direction (0, 0, 1),
 A ray in the direction (0, 1, 1),
 A ray in the direction (0.7071067811865475?, 0, 1),
 A ray in the direction (1, 1, 1))

This is wrong (this is why I printed all the extra stuff).

Checking the indices I get: (<0,1,3>, <0,2,3>).

So accidentally I stumbled upon an example myself.

[git]

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

​63f1b79

account for maybe changed order of the generators

[gh-kliem]

It's not pretty, but it seems to work.

[mkoeppe]

My guess is that this workaround might fail if this not full dimensional.

But I think this is good enough as a hot fix for Sage 9.2 while we wait for the new Normaliz.

[Winfried]

To day I have changed the master branch of Normliz as follows:

-- Generators are no longer accessible from the outside

-- they are replaced by TriangulationGenerators?

-- TriangulationGeneretors? are set immediately after a triangulation has been computed

It is possible to compute several triangulations successively, for example an ordinary one first and then a unimodular one. I think it is clear that one must pull the TriangulationGenerators? again after the second computation.

One could go a step further by bundling the Triangulation and the TriangulationGenerators? so that they are returned together. If you think this is necessary, please let me know.

PyNormaliz? has not yet been updated. I think only changes in the tests are necessary.

[vbraun]

**********************************************************************
File "src/sage/geometry/polyhedron/backend_normaliz.py", line 1574, in sage.geometry.polyhedron.backend_normaliz.Polyhedron_normaliz._triangulate_normaliz
Failed example:
    C1 = Polyhedron(rays=[[0,0,1],[1,0,AA(2).sqrt()],[0,1,1],[1,1,1]], backend='normaliz')
Exception raised:
    Traceback (most recent call last):
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/doctest/forker.py", line 720, in _run
        self.compile_and_execute(example, compiler, test.globs)
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/doctest/forker.py", line 1145, in compile_and_execute
        exec(compiled, globs)
      File "<doctest sage.geometry.polyhedron.backend_normaliz.Polyhedron_normaliz._triangulate_normaliz[0]>", line 1, in <module>
        C1 = Polyhedron(rays=[[Integer(0),Integer(0),Integer(1)],[Integer(1),Integer(0),AA(Integer(2)).sqrt()],[Integer(0),Integer(1),Integer(1)],[Integer(1),Integer(1),Integer(1)]], backend='normaliz')
      File "sage/misc/lazy_import.pyx", line 353, in sage.misc.lazy_import.LazyImport.__call__ (build/cythonized/sage/misc/lazy_import.c:3736)
        return self.get_object()(*args, **kwds)
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/constructor.py", line 662, in Polyhedron
        return parent(Vrep, Hrep, convert=convert, verbose=verbose)
      File "sage/structure/parent.pyx", line 902, in sage.structure.parent.Parent.__call__ (build/cythonized/sage/structure/parent.c:9364)
        return mor._call_with_args(x, args, kwds)
      File "sage/structure/coerce_maps.pyx", line 180, in sage.structure.coerce_maps.DefaultConvertMap_unique._call_with_args (build/cythonized/sage/structure/coerce_maps.c:5154)
        raise
      File "sage/structure/coerce_maps.pyx", line 175, in sage.structure.coerce_maps.DefaultConvertMap_unique._call_with_args (build/cythonized/sage/structure/coerce_maps.c:5042)
        return C._element_constructor(x, *args, **kwds)
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/parent.py", line 598, in _element_constructor_
        return self.element_class(self, Vrep, Hrep, **kwds)
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/backend_normaliz.py", line 238, in __init__
        Polyhedron_base.__init__(self, parent, Vrep, Hrep, **kwds)
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/base.py", line 249, in __init__
        self._init_from_Vrepresentation(vertices, rays, lines, **kwds)
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/backend_normaliz.py", line 565, in _init_from_Vrepresentation
        self._init_from_normaliz_data(data, normaliz_field=normaliz_field, verbose=verbose)
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/backend_normaliz.py", line 387, in _init_from_normaliz_data
        cone = self._cone_from_normaliz_data(data, verbose)
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/geometry/polyhedron/backend_normaliz.py", line 419, in _cone_from_normaliz_data
        PythonModule("PyNormaliz", spkg="pynormaliz").require()
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/features/__init__.py", line 191, in require
        raise FeatureNotPresentError(self, presence.reason, presence.resolution)
    sage.features.FeatureNotPresentError: PyNormaliz is not available.
    Failed to import `PyNormaliz`.
    To install PyNormaliz you can try to run 'sage -i pynormaliz'.
**********************************************************************
File "src/sage/geometry/polyhedron/backend_normaliz.py", line 1575, in sage.geometry.polyhedron.backend_normaliz.Polyhedron_normaliz._triangulate_normaliz
Failed example:
    C1._triangulate_normaliz()
Exception raised:
    Traceback (most recent call last):
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/doctest/forker.py", line 720, in _run
        self.compile_and_execute(example, compiler, test.globs)
      File "/home/release/Sage/local/lib/python3.8/site-packages/sage/doctest/forker.py", line 1145, in compile_and_execute
        exec(compiled, globs)
      File "<doctest sage.geometry.polyhedron.backend_normaliz.Polyhedron_normaliz._triangulate_normaliz[1]>", line 1, in <module>
        C1._triangulate_normaliz()
    NameError: name 'C1' is not defined
**********************************************************************

[git]

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

​8d934b2

optional flag

[gh-kliem]

Thanks.