Upgrade cddlib to fix wrong intersection of polytopes

[videlec]

As reported in this sage-devel thread, the intersection of polytopes is sometimes very wrong

a = Polyhedron([[0, -1, 1], [1, -1, 1], [1, 1, -1]])
b = Polyhedron([[0.0, -0.5, 1.5], [1.0, -0.5, 1.5], [1.0, 1.5, -0.5]])
c = a.intersection(b)

Here c should have been empty but is equal to b.

We fix it by upgrading to cddlib 0.94m, installing symlinks so that consumers of cddlib that have not been updated to work with the new header file locations (#29413) continue to work.

Upstream: Fixed upstream, in a later stable release.

CC: @jplab @mkoeppe @dimpase @vbraun

Component: geometry

Author: Jonathan Kliem, Matthias Koeppe

Branch/Commit: 167c832

Reviewer: Matthias Koeppe, Dima Pasechnik

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

[kliem]

comment:2

Looks like a cdd bug.

sage: a = Polyhedron([[0, -1, 1], [1, -1, 1], [1, 1, -1]])
sage: b = Polyhedron([[0.0, -0.5, 1.5], [1.0, -0.5, 1.5], [1.0, 1.5, -0.5]])
sage: a = a.change_ring(RDF)
sage: new_eqns = a.equations() + b.equations()
sage: new_ieqs = a.inequalities() + b.inequalities()
sage: Polyhedron(eqns=new_eqns, ieqs=new_ieqs[1:2]+new_ieqs[4:5], verbose=True)
---- CDD input -----
H-representation
linearity 2 1 2
begin
 5 4 real
 0.0 0.0 1.0 1.0
 -1.0 0.0 1.0 1.0
 -1.0 2.0 -1.0 0.0
 -1.0 4.0 -2.0 0.0
 1 0 0 0
end

---- CDD output -----
Canonicalize the matrix.
Implicit linearity rows are:

Redundant rows are:1 3 5 

Nonredundant representation:
The new row positions are as follows (orig:new).
Each redundant row has the new number 0.
Each deleted duplicated row has a number nagative of the row that
represents its equivalence class.
 1:0 2:1 3:0 4:2 5:0
H-representation
linearity 1  1
begin
 2 4 real
 -1  0  1  1
 -1  4 -2  0
end

size = 5 x 4
Number Type = real

---- CDD input -----
Canonicalize the matrix.
Implicit linearity rows are:

Redundant rows are:1 3 5 

Nonredundant representation:
The new row positions are as follows (orig:new).
Each redundant row has the new number 0.
Each deleted duplicated row has a number nagative of the row that
represents its equivalence class.
 1:0 2:1 3:0 4:2 5:0
H-representation
linearity 1  1
begin
 2 4 real
 -1  0  1  1
 -1  4 -2  0
end

---- CDD output -----
The second representation:
V-representation
linearity 1  3
begin
 3 4 real
  0  1  0  0
  1  7.500000000E-01  1  0
  0 -5.000000000E-01 -1  1
end

Vertex incidence
ecd_file: Incidence of generators and inequalities
begin
  3    3
 1 -2 : 2 
 2 -2 : 3 
 3 -3 : 
end

Vertex adjacency
ead_file: Adjacency of generators
begin
  3    3
 1 1 : 2 
 2 1 : 1 
 3 0 : 
end

The first (input) representation
H-representation
linearity 1  1
begin
 2 4 real
 -1  0  1  1
 -1  4 -2  0
end

Facet incidence
icd_file: Incidence of inequalities and generators
begin
  3    3
 1 -3 : 
 2 -2 : 1 
 3 -2 : 2 
end

Facet adjacency
iad_file: Adjacency of inequalities
begin
  3    3
 1 -2 : 1 
 2 -2 : 2 
 3 -2 : 3 
end

size = 2 x 4
Number Type = real

---- CDD input -----
The second representation:
V-representation
linearity 1  3
begin
 3 4 real
  0  1  0  0
  1  7.500000000E-01  1  0
  0 -5.000000000E-01 -1  1
end
---- CDD output -----
The second representation:
H-representation
linearity 1  3
begin
 3 4 real
  1  0  0  0
 -1  4 -2  0
 -1  0  1  1
end

size = 3 x 4
Number Type = real

A 2-dimensional polyhedron in RDF^3 defined as the convex hull of 1 vertex, 1 ray, 1 line

For some reason cdd thinks that row 1 is redundant. If you drop an inequality, everything is fine.

[sagetrac-tmonteil]

comment:3

Could someone please report it upstream (i do not have a github account) ?

[sagetrac-tmonteil]

Upstream: Not yet reported upstream; Will do shortly.

[kliem]

comment:4

cddlib/cddlib#45

[kliem]

Changed upstream from Not yet reported upstream; Will do shortly. to Reported upstream. No feedback yet.

[dimpase]

comment:5

Note that with data in comment:2 one has

sage: Polyhedron(eqns=new_eqns)                                                                                                      
The empty polyhedron in RDF^3

so the emptiness should have been detected earlier.

[kliem]

comment:6

Are you proposing that we test whether the equations are infeasible ourselves as a workaround?

Sounds like a workaround that should work for us.

Is cdd specifically excluding this case?

[kliem]

comment:7

Something like

new_eqns_matrix = matrix(new_eqns)
feasible_directions = new_eqns_matrix.right_kernel_matrix()
if feasible_directions.column(0).is_zero():
    initialize_empty_polyhedron()

should work.

[dimpase]

comment:8

On the other hand, with RDF data and the internal representation used by cdd, it is not so hard to come up with (ugly) examples on which it can fail to produce a meaningful answer.

I'd say that computing c like in the ticket description should throw an error,
as a is exact and b is not.
Or at the very least print a warning like

sage: a = Polyhedron([[0, -1, 1], [1, -1, 1], [1, 1, -1]],base_ring=RDF) 
....: b = Polyhedron([[0.0, -0.5, 1.5], [1.0, -0.5, 1.5], [1.0, 1.5, -0.5]])                                                         
sage: a.intersection(b)                                                                                                              
/mnt/opt/Sage/sage-dev/local/lib/python3.8/site-packages/sage/geometry/polyhedron/backend_cdd.py:151: 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.
  warn("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.")
A 2-dimensional polyhedron in RDF^3 defined as the convex hull of 3 vertices

[dimpase]

comment:9

it's also weird that no warning is printed in this case:

sage: a = Polyhedron([[0, -1, 1], [1, -1, 1], [1, 1, -1]])
....: b = Polyhedron([[0.0, -0.5, 1.5], [1.0, -0.5, 1.5], [1.0, 1.5, -0.5]])      
....: b.intersection(a)                                  
A 2-dimensional polyhedron in RDF^3 defined as the convex hull of 3 vertices

as if seeing the exact a to intersect b with convinces cdd that everything should be fine.

[kliem]

comment:10

The behavior of coercing base rings is the correct and expected behaviour, I would say.

The failure of cdd has nothing to do with inexactness:

sage: a = Polyhedron([[0, -1, 1], [1, -1, 1], [1, 1, -1]], backend='cdd')                                                                                                           
sage: b = Polyhedron([[0, -1/2, 3/2], [1, -1/2, 3/2], [1, 3/2, -1/2]], backend='cdd')                                                                                               
sage: a.intersection(b)                                                                                                                                                             
A 2-dimensional polyhedron in QQ^3 defined as the convex hull of 3 vertices

[kliem]

comment:11

cdd for some reason determines that one of the equations is redundant. After this the calculation is correct. For this reason, no warning is printed. Because converting from Vrepresention to Hrepresentation works.

The output of cdd is numerical consistent, but incorrect.

[dimpase]

comment:12

Replying to @kliem:

The behavior of coercing base rings is the correct and expected behaviour, I would say.

The failure of cdd has nothing to do with inexactness:

OK, I see, you're right - I thought it's an inexactness issue, such as computing rank of an RDF matrix, sorry.

[kliem]

Changed upstream from Reported upstream. No feedback yet. to Reported upstream. Developers acknowledge bug.

[dimpase]

comment:15

The fix is in
https://github.com/cddlib/cddlib/releases/tag/0.94m

let's upgrade

[dimpase]

Changed upstream from Reported upstream. Developers acknowledge bug. to Fixed upstream, in a later stable release.

[kliem]

comment:17

How important do we consider this bug. Should we reject older system installs?

[kliem]

Changed upstream from Fixed upstream, in a later stable release. to Reported upstream. Developers acknowledge bug.

[kliem]

Changed upstream from Reported upstream. Developers acknowledge bug. to Fixed upstream, in a later stable release.

[kliem]

comment:18

Ok, it warned me because of a conflict, but I did not expect it to change the description back.

[kliem]

comment:19

There is a bug in latte_int with the new header location that needs fixing yet:

latte-int/latte#15

[kliem]

comment:20

@mkoeppe

There is already a patch for latte, but I'm incapable of creating a patch for latte-int not latte-integral.

[mkoeppe]

Changed work issues from rebase on dependencies, -I$SAGE_LOCAL/include/cddlib for gfan, topcom configure to none

[mkoeppe]

comment:36

Do you know how to test for the cddlib bug in spkg-configure.m4?

[mkoeppe]

Author: Jonathan Kliem, Matthias Koeppe

[kliem]

comment:37

Replying to @sagetrac-git:

Branch pushed to git repo; I updated commit sha1. This was a forced push. New commits:

308274d	update cddlib
ed3ed2a	build/pkgs/cddlib/spkg-install.in: Install symlinks for cddlib/*.h

I'm a bit confused on why you deleted the work issues. Did you forget to push something?

[kliem]

comment:38

Create a file

test.ine

with

H-representation
linearity 2 1 2
begin
 5 4 real
 0.0 0.0 1.0 1.0
 -1.0 0.0 1.0 1.0
 -1.0 2.0 -1.0 0.0
 -1.0 4.0 -2.0 0.0
 1 0 0 0
end

and run cddexec --redcheck <test.ine.

The line Redundant rows are should not contain a 1.

[mkoeppe]

comment:39

Replying to @kliem:

I'm a bit confused on why you deleted the work issues. Did you forget to push something?

I thought the result of our discussion was that by installing the symlinks, we do not have to fix anything in the other packages.

[mkoeppe]

comment:40

Or is there a package that no longer finds a prefix-less system cdd?

[kliem]

comment:41

I see.

Ok, then technically the latte upgrade isn't a decency anymore.

[sagetrac-git]

Changed commit from ed3ed2a to 167c832

[sagetrac-git]

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

167c832

build/pkgs/cddlib/spkg-configure.m4: Add test for cddexec --redcheck bug

[mkoeppe]

Changed dependencies from #31482, #25993 to none

[mkoeppe]

comment:46

OK... the new configure check works correctly on homebrew-macos-standard (result: yes - but then the package is rejected because the headers are in the new place) and on ubuntu-focal-standard (result: no).

[mkoeppe]

Reviewer: https://github.com/mkoeppe/sage/actions/runs/666698487

[mkoeppe]

comment:49

Let's please get this fix into 9.3

[mkoeppe]

Changed reviewer from https://github.com/mkoeppe/sage/actions/runs/666698487 to Matthias Koeppe, ...

[dimpase]

Changed reviewer from Matthias Koeppe, ... to Matthias Koeppe, Dima Pasechnik

[dimpase]

comment:51

ok

[kliem]

comment:52

Thank you.

[mkoeppe]

comment:53

Thanks!

[mkoeppe]

comment:54

Setting priority to blocker to bring this ticket to the attention of the release bot.

[vbraun]

Changed branch from public/30319 to 167c832