ASKSAGE: Sage Q&A Forum - Individual question feedhttp://ask.sagemath.org/questions/Q&A Forum for SageenCopyright Sage, 2010. Some rights reserved under creative commons license.Wed, 02 May 2018 08:56:23 -0500Speed of enumerating integer points in polytopehttp://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/ I'm trying to enumerate all integer points in some polytopes, but sage takes way too long.
This is a pretty minimal example of one of these polytopes.
<pre><code>q = MixedIntegerLinearProgram(maximization=False)
x = q.new_variable(integer=True, nonnegative=True)
q.add_constraint( x[0] == 1)
q.add_constraint(88 * x[2] == -306 * x[0] + 1 * x[1])
q.add_constraint( 1 * x[3] == 42636 * x[0] + 42 * x[1])
q.add_constraint(11 * x[4] == 30804 * x[0] + -42 * x[1])
P = q.polyhedron()
P
# A 1-dimensional polyhedron in QQ^5 defined as the convex hull of 2 vertices
P.Vrepresentation()
# (A vertex at (1, 34/7, 5134/7, 73440, 0), A vertex at (1, 0, 306, 55488, 1632))
</code></pre>
Now, I would like to enumerate the integral points in this (1-dimensional) polytope. Looking at the second coordinate of the two boundary points, it is clear there can be at most five such points (in fact, there are exactly five).
However, `P.integral_points()` seems to take forever.
Is there anything I can do to speed this up? I would like to do similar problems that are higher-dimensional, and have a lot more constraints (bounding hyperplanes).
Earlier, I used a MILP library (COIN) to solve the system, giving one integer point, then further constrain the system to exclude this solution, solve again, etc. This is fast enough for this particular example, but becomes unfeasible when the dimension of the polytope and/or the number of integer points increases.Fri, 10 Oct 2014 11:53:00 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/Comment by FrédéricC for <p>I'm trying to enumerate all integer points in some polytopes, but sage takes way too long.</p>
<p>This is a pretty minimal example of one of these polytopes.</p>
<pre><code>q = MixedIntegerLinearProgram(maximization=False)
x = q.new_variable(integer=True, nonnegative=True)
q.add_constraint( x[0] == 1)
q.add_constraint(88 * x[2] == -306 * x[0] + 1 * x[1])
q.add_constraint( 1 * x[3] == 42636 * x[0] + 42 * x[1])
q.add_constraint(11 * x[4] == 30804 * x[0] + -42 * x[1])
P = q.polyhedron()
P
# A 1-dimensional polyhedron in QQ^5 defined as the convex hull of 2 vertices
P.Vrepresentation()
# (A vertex at (1, 34/7, 5134/7, 73440, 0), A vertex at (1, 0, 306, 55488, 1632))
</code></pre>
<p>Now, I would like to enumerate the integral points in this (1-dimensional) polytope. Looking at the second coordinate of the two boundary points, it is clear there can be at most five such points (in fact, there are exactly five).</p>
<p>However, <code>P.integral_points()</code> seems to take forever.</p>
<p>Is there anything I can do to speed this up? I would like to do similar problems that are higher-dimensional, and have a lot more constraints (bounding hyperplanes).</p>
<p>Earlier, I used a MILP library (COIN) to solve the system, giving one integer point, then further constrain the system to exclude this solution, solve again, etc. This is fast enough for this particular example, but becomes unfeasible when the dimension of the polytope and/or the number of integer points increases.</p>
http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24491#post-id-24491Maybe you can find one point, shift it to zero and then work on the subspace spanned by the polytope ? That would make the ambient dimension smaller.Wed, 15 Oct 2014 13:48:42 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24491#post-id-24491Answer by Dima for <p>I'm trying to enumerate all integer points in some polytopes, but sage takes way too long.</p>
<p>This is a pretty minimal example of one of these polytopes.</p>
<pre><code>q = MixedIntegerLinearProgram(maximization=False)
x = q.new_variable(integer=True, nonnegative=True)
q.add_constraint( x[0] == 1)
q.add_constraint(88 * x[2] == -306 * x[0] + 1 * x[1])
q.add_constraint( 1 * x[3] == 42636 * x[0] + 42 * x[1])
q.add_constraint(11 * x[4] == 30804 * x[0] + -42 * x[1])
P = q.polyhedron()
P
# A 1-dimensional polyhedron in QQ^5 defined as the convex hull of 2 vertices
P.Vrepresentation()
# (A vertex at (1, 34/7, 5134/7, 73440, 0), A vertex at (1, 0, 306, 55488, 1632))
</code></pre>
<p>Now, I would like to enumerate the integral points in this (1-dimensional) polytope. Looking at the second coordinate of the two boundary points, it is clear there can be at most five such points (in fact, there are exactly five).</p>
<p>However, <code>P.integral_points()</code> seems to take forever.</p>
<p>Is there anything I can do to speed this up? I would like to do similar problems that are higher-dimensional, and have a lot more constraints (bounding hyperplanes).</p>
<p>Earlier, I used a MILP library (COIN) to solve the system, giving one integer point, then further constrain the system to exclude this solution, solve again, etc. This is fast enough for this particular example, but becomes unfeasible when the dimension of the polytope and/or the number of integer points increases.</p>
http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?answer=24463#post-id-24463try using Latte Integrale (an experimental package in Sage), called latte_int; it implements Barvinok's algorithm for integer point counting. There is however no interface between Latte and Sage developed, so you need to create files to be real into Latte... Sat, 11 Oct 2014 16:31:10 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?answer=24463#post-id-24463Comment by Arie for <p>try using Latte Integrale (an experimental package in Sage), called latte_int; it implements Barvinok's algorithm for integer point counting. There is however no interface between Latte and Sage developed, so you need to create files to be real into Latte... </p>
http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24479#post-id-24479I'm trying to do an intermediately big case, where the dimension is 8, and the polytope has 1804 vertices.
Latte Integrale has already been computing for over 20 minutes (just counting the points, for now), and there is no indication how long it will take. It just says: "431000 unimodular cones done", this number increasing by a 1000 every few seconds, but I don't know how many unimodular cones it must do.Tue, 14 Oct 2014 06:40:05 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24479#post-id-24479Comment by Dima for <p>try using Latte Integrale (an experimental package in Sage), called latte_int; it implements Barvinok's algorithm for integer point counting. There is however no interface between Latte and Sage developed, so you need to create files to be real into Latte... </p>
http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24478#post-id-24478it is the most compact way known, and this is what counts.IIRC, the rational function consists of summands, each of them a rational function corresponding to a simplex, simplices forming a kind of partition of the polytope.So one can probably only take summands corresponding to x^K with positive signs, if you don't mind getting the same point several times.Tue, 14 Oct 2014 02:16:57 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24478#post-id-24478Comment by Dima for <p>try using Latte Integrale (an experimental package in Sage), called latte_int; it implements Barvinok's algorithm for integer point counting. There is however no interface between Latte and Sage developed, so you need to create files to be real into Latte... </p>
http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24476#post-id-24476it does count them, in a way; you have to call count with the option --multivariate-generating-function; such a function "simplifies" to a (Laurent) polynomial, with exponent vectors of the monomials corresponding to the integer points.Mon, 13 Oct 2014 10:07:10 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24476#post-id-24476Comment by Dima for <p>try using Latte Integrale (an experimental package in Sage), called latte_int; it implements Barvinok's algorithm for integer point counting. There is however no interface between Latte and Sage developed, so you need to create files to be real into Latte... </p>
http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24484#post-id-24484IMHO internally it would still be computing the same generating function I mentioned, no matter you "only" count, or not.Tue, 14 Oct 2014 14:49:24 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24484#post-id-24484Comment by Arie for <p>try using Latte Integrale (an experimental package in Sage), called latte_int; it implements Barvinok's algorithm for integer point counting. There is however no interface between Latte and Sage developed, so you need to create files to be real into Latte... </p>
http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24474#post-id-24474Unfortunately, it seems Latte Integrale does not enumerate the lattice points, it only counts them. I need the coordinates of all the integer points, for further calculations.Mon, 13 Oct 2014 06:27:45 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24474#post-id-24474Comment by Arie for <p>try using Latte Integrale (an experimental package in Sage), called latte_int; it implements Barvinok's algorithm for integer point counting. There is however no interface between Latte and Sage developed, so you need to create files to be real into Latte... </p>
http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24477#post-id-24477I see. It does not sound like a very convenient way to represent the thousands of points in somewhat higher-dimensional spaces of my other cases. I'll give it a try, though. Thanks for the suggestion!Mon, 13 Oct 2014 11:29:33 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?comment=24477#post-id-24477Answer by jipilab for <p>I'm trying to enumerate all integer points in some polytopes, but sage takes way too long.</p>
<p>This is a pretty minimal example of one of these polytopes.</p>
<pre><code>q = MixedIntegerLinearProgram(maximization=False)
x = q.new_variable(integer=True, nonnegative=True)
q.add_constraint( x[0] == 1)
q.add_constraint(88 * x[2] == -306 * x[0] + 1 * x[1])
q.add_constraint( 1 * x[3] == 42636 * x[0] + 42 * x[1])
q.add_constraint(11 * x[4] == 30804 * x[0] + -42 * x[1])
P = q.polyhedron()
P
# A 1-dimensional polyhedron in QQ^5 defined as the convex hull of 2 vertices
P.Vrepresentation()
# (A vertex at (1, 34/7, 5134/7, 73440, 0), A vertex at (1, 0, 306, 55488, 1632))
</code></pre>
<p>Now, I would like to enumerate the integral points in this (1-dimensional) polytope. Looking at the second coordinate of the two boundary points, it is clear there can be at most five such points (in fact, there are exactly five).</p>
<p>However, <code>P.integral_points()</code> seems to take forever.</p>
<p>Is there anything I can do to speed this up? I would like to do similar problems that are higher-dimensional, and have a lot more constraints (bounding hyperplanes).</p>
<p>Earlier, I used a MILP library (COIN) to solve the system, giving one integer point, then further constrain the system to exclude this solution, solve again, etc. This is fast enough for this particular example, but becomes unfeasible when the dimension of the polytope and/or the number of integer points increases.</p>
http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?answer=42243#post-id-42243In Sage 8.1, I can reproduce your code as follows (notice that I added "ppl" as a solver since the default behavior changed between your post and this answer).
sage: q = MixedIntegerLinearProgram(maximization=False, solver="ppl")
....: x = q.new_variable(integer=True, nonnegative=True)
....:
....: q.add_constraint( x[0] == 1)
....: q.add_constraint(88 * x[2] == -306 * x[0] + 1 * x[1])
....: q.add_constraint( 1 * x[3] == 42636 * x[0] + 42 * x[1])
....: q.add_constraint(11 * x[4] == 30804 * x[0] + -42 * x[1])
In order to speed things up, I use `normaliz`.
(You can install it by typing the commands in a terminal:
sage -i normaliz
sage -i pynormaliz
)
and then use the backend `normaliz` for the polyhedron:
....: P = q.polyhedron(backend="normaliz")
....: P.integral_points()
....:
((1, 0, 306, 55488, 1632),
(1, 1, 394, 59184, 1296),
(1, 2, 482, 62880, 960),
(1, 3, 570, 66576, 624),
(1, 4, 658, 70272, 288))
which should appear instantaneously.Wed, 02 May 2018 08:56:23 -0500http://ask.sagemath.org/question/24458/speed-of-enumerating-integer-points-in-polytope/?answer=42243#post-id-42243