Opened 5 years ago
Closed 5 years ago
#21270 closed defect (fixed)
Polyhedron RDF face lattice bug / intersection of polyhedra
Reported by:  mkoeppe  Owned by:  moritz 

Priority:  major  Milestone:  sage7.4 
Component:  geometry  Keywords:  polyhedron, days84 
Cc:  chapoton, moritz, jipilab  Merged in:  
Authors:  Moritz Firsching  Reviewers:  JeanPhilippe Labbé 
Report Upstream:  N/A  Work issues:  
Branch:  736f8a3 (Commits, GitHub, GitLab)  Commit:  736f8a3f3dcba14b43f527e9ab3a1a82f40280fb 
Dependencies:  Stopgaps: 
Description (last modified by )
As noted first on #20570, current Sage is not able to plot some LPs with irrational data. This appears to be due to a bug in the Sage polyhedron code for RDF data (which InteractiveLP uses when the data are not rational).
sage: poly = polytopes.regular_polygon(7) sage: lp, x = poly.to_linear_program(solver='InteractiveLP', return_variable=True) sage: lp.set_objective(x[0] + x[1]) sage: b = lp.get_backend() sage: P = b.interactive_lp_problem() sage: p = P.plot() ### ERROR sage: p.show()
Change History (29)
comment:1 Changed 5 years ago by
 Component changed from PLEASE CHANGE to geometry
comment:2 Changed 5 years ago by
 Cc jipilab added
comment:3 Changed 5 years ago by
 Cc moritz added
comment:4 Changed 5 years ago by
 Keywords polyhedron days84 added
 Owner changed from (none) to moritz
comment:5 Changed 5 years ago by
 Summary changed from Polyhedron RDF plotting bug to Polyhedron RDF face lattice bug
comment:6 Changed 5 years ago by
Apperently something goes wrong when taking the intersection of a polyhedron in RDF and one in QQ:
Q plays the role of box in the function plot_feasible_set
of interactive_simplex_method.py
.
sage: L=[[1349811595/3385908, 488029882/7345133], ....: [1349811595/3385908, 1349811595/5078862], ....: [732044823/7345133, 1349811595/5078862], ....: [732044823/7345133, 488029882/7345133]] ....: sage: sage: Q=Polyhedron(L) sage: P=Polyhedron(polytopes.regular_polygon(7).vertices_list(), base_ring=RDF) sage: P.intersection(Q) A 1dimensional polyhedron in RDF^2 defined as the convex hull of 7 vertices sage: Q.intersection(P) The empty polyhedron in RDF^2
comment:7 Changed 5 years ago by
 Summary changed from Polyhedron RDF face lattice bug to Polyhedron RDF face lattice bug / intersection of polyhedra
Here is a smaller example of the failing intersection between QQ and RDF polyhedra:
sage: Q=Polyhedron(ieqs=[[499999,1000000],[1499999,1000000]]) sage: P=Polyhedron(ieqs=[[0,1.0],[1.0,1.0]], base_ring=RDF) sage: Q.intersection(P) The empty polyhedron in RDF^1 sage: P.intersection(Q) A 1dimensional polyhedron in RDF^1 defined as the convex hull of 2 vertices sage: Q=Polyhedron(ieqs=[[4999999,10000000],[14999999,10000000]]) sage: P.intersection(Q)
The last command then raises
ValueError: *Error: Numerical inconsistency is found. Use the GMP exact arithmetic.
, which is better than a wrong result, I guess. In any case one should perhaps do the conversions of the base ring more explicitly.
The relevant function in base.py reads:
new_ieqs = self.inequalities() + other.inequalities() new_eqns = self.equations() + other.equations() parent = self.parent() try: return parent.element_class(parent, None, [new_ieqs, new_eqns]) except TypeError as msg: if self.base_ring() is ZZ: parent = parent.base_extend(QQ) return parent.element_class(parent, None, [new_ieqs, new_eqns]) else: raise TypeError(msg)
So in the example above we have in the new_ieqs variable:
(An inequality (1.0) x + 0.0 >= 0, An inequality (1.0) x + 1.0 >= 0, An inequality (10000000) x + 14999999 >= 0, An inequality (10000000) x  4999999 >= 0)
Which leads to the unexpected behaviour
comment:8 Changed 5 years ago by
 Branch set to u/moritz/21270
comment:9 Changed 5 years ago by
 Commit set to 83d5e8f9fc68d4dbaf09ed87db216a5cd3c24922
This commit fixes the instances of the problem, that I added.
The original thing is still not working, but this is now another bug (more later).
New commits:
83d5e8f  normalizing for rational input, when converting polyhedra to RDF

comment:10 Changed 5 years ago by
 Commit changed from 83d5e8f9fc68d4dbaf09ed87db216a5cd3c24922 to 3e1ebc0c5746c5ee68e296c344ad1464f3d7e4e8
Branch pushed to git repo; I updated commit sha1. New commits:
3e1ebc0  use the correct base ring

comment:11 Changed 5 years ago by
Now it works..
comment:12 Changed 5 years ago by
 Description modified (diff)
comment:13 Changed 5 years ago by
I submitted a fix to this bug.
Explanation:
What went wrong was was happens if a polyhedron in QQ
is intersection with a polyhedron in RDF
; basically it takes the union of the inequalities. The inequalities are typically stored as linear inequality with coefficients in ZZ
. Those coefficients can get very large and cause numerical trouble when converted to RDF
. Thus it is advisable to convert the inequalities to RDF
only after some renormalization. We choose to divide all terms of the inequality by the absolute values of the coefficient with the largest absolute value. (In case all the coefficients are zero, we leave it as is.)
This helped with the examples above. In the original question there was another problem: in the method _solve
of the class InteractiveLPProblem(
in the interactive_simplex_method.py
file, taking the base ring to be something other than the base ring of the feasible region (which might have been coerced to something else) caused trouble. I changed this by simply setting
R = F.base_ring()
Then I was able to produce the plot.
comment:14 Changed 5 years ago by
 Description modified (diff)
comment:15 Changed 5 years ago by
By the way a better plotting result is obtained by
sage: p = P.plot(xmin=3, xmax=8, ymin=2, ymax=4) sage: p.show(dpi=500)
But giving an explicit bounding box circumvents the first problem..
comment:16 Changed 5 years ago by
 Commit changed from 3e1ebc0c5746c5ee68e296c344ad1464f3d7e4e8 to 85bb29e255626dce99c63987158f5982000a0525
Branch pushed to git repo; I updated commit sha1. New commits:
85bb29e  added doctests

comment:17 Changed 5 years ago by
 Status changed from new to needs_review
comment:18 Changed 5 years ago by
comment:19 Changed 5 years ago by
There seems to be a merge conflict. Could you rebase it?
comment:20 Changed 5 years ago by
 Status changed from needs_review to needs_work
comment:21 Changed 5 years ago by
 Commit changed from 85bb29e255626dce99c63987158f5982000a0525 to 455e83e035fc7acb18d5ad565ae10a77d68c7d2d
comment:22 Changed 5 years ago by
 Status changed from needs_work to needs_review
This should now merge cleanly. (There was just an issue with some trailing whitespace that I had removed in a line that was changed in rc0.)
comment:23 Changed 5 years ago by
Hi Moritz,
There are a couple of missing spaces between binary operator.
Let's see what the bot thinks.
comment:24 Changed 5 years ago by
Hi Moritz,
The bot is happy and it looks good to me. If you correct the 3 spacing around operators in the parent.py
, you can set it to positive review on my behalf.
comment:25 Changed 5 years ago by
 Reviewers set to JeanPhilippe Labbé
comment:26 Changed 5 years ago by
 Commit changed from 455e83e035fc7acb18d5ad565ae10a77d68c7d2d to 736f8a3f3dcba14b43f527e9ab3a1a82f40280fb
comment:27 Changed 5 years ago by
Thank for your review, JP!
comment:28 Changed 5 years ago by
 Status changed from needs_review to positive_review
comment:29 Changed 5 years ago by
 Branch changed from u/moritz/21270 to 736f8a3f3dcba14b43f527e9ab3a1a82f40280fb
 Resolution set to fixed
 Status changed from positive_review to closed
misleading comment deleted