Trac #15521: Deprecations: default LP variables will become real inst…

…ead of nonnegative

Following the discussion on Sage-devel [1], default LP variables will
become "real" instead of nonnegative as it is the case at the moment.

The aim is to have 4 types available through `new_variable()` : binary,
integer, nonnegative, and real. Right now, nonnegative does not exist,
and "real" represents nonnegative variables.

What this ticket does:

- A warning has to be displayed when `new_variable()` is called without
arguments, as the default behaviour will change. Users are advised to
add `nonnegative=True` instead.

- A warning has to be displayed when `new_variable(real=True)` is
called, for it represents nonnegative variables at the moment and will
represent real variables later. Users are advised to switch to
`nonnegative=True` instead.

Nathann

[1] https://groups.google.com/d/msg/sage-devel/3vrPzUqFpM4/hKFp0RjV8poJ

URL: http://trac.sagemath.org/15521
Reported by: ncohen
Ticket author(s): Nathann Cohen
Reviewer(s): Benjamin Jones

src/doc/en/thematic_tutorials/linear_programming.rst

@@ -87,6 +87,8 @@ values**, see next section)::
 
     sage: p = MixedIntegerLinearProgram()
     sage: x, y, z = p['x'], p['y'], p['z']
+    doctest:839: DeprecationWarning: The default behaviour of new_variable() will soon change ! It will return 'real' variables instead of nonnegative ones. Please be explicit and call new_variable(nonnegative=True) instead.
+    See http://trac.sagemath.org/15521 for details.
 
 Next, we set the objective function
 
@@ -166,8 +168,8 @@ write
 
 ::
 
-    sage: p.add_constraint(x["I am a valid key"]
-    ...                    + x[("a",pi)] <= 3)
+    sage: p.add_constraint(x["I am a valid key"] +
+    ....:                  x[("a",pi)] <= 3)
 
 
 And because any immutable object can be used as a key, doubly indexed variables
@@ -180,7 +182,6 @@ And because any immutable object can be used as a key, doubly indexed variables
 
     sage: p.add_constraint(x[3,2] + x[5] == 6)
 
-
 Typed variables and bounds
 """"""""""""""""""""""""""
 

src/sage/coding/delsarte_bounds.py

@@ -86,7 +86,7 @@ def _delsarte_LP_building(n, d, d_star, q, isinteger,  solver, maxc = 0):
     from sage.numerical.mip import MixedIntegerLinearProgram
 
     p = MixedIntegerLinearProgram(maximization=True, solver=solver)
-    A = p.new_variable(integer=isinteger) # A>=0 is assumed
+    A = p.new_variable(integer=isinteger, nonnegative=not isinteger) # A>=0 is assumed
     p.set_objective(sum([A[r] for r in xrange(n+1)]))
     p.add_constraint(A[0]==1)
     for i in xrange(1,d):

src/sage/graphs/generic_graph.py

@@ -4825,13 +4825,13 @@ def steiner_tree(self,vertices, weighted = False, solver = None, verbose = 0):
         R = lambda (x,y) : (x,y) if x<y else (y,x)
 
         # edges used in the Steiner Tree
-        edges = p.new_variable()
+        edges = p.new_variable(binary=True)
 
         # relaxed edges to test for acyclicity
-        r_edges = p.new_variable()
+        r_edges = p.new_variable(nonnegative=True)
 
         # Whether a vertex is in the Steiner Tree
-        vertex = p.new_variable()
+        vertex = p.new_variable(binary = True)
         for v in g:
             for e in g.edges_incident(v, labels=False):
                 p.add_constraint(vertex[v] - edges[R(e)], min = 0)
@@ -4861,7 +4861,6 @@ def steiner_tree(self,vertices, weighted = False, solver = None, verbose = 0):
 
         p.set_objective(p.sum([w(e)*edges[R(e)] for e in g.edges(labels = False)]))
 
-        p.set_binary(edges)
         p.solve(log = verbose)
 
         edges = p.get_values(edges)
@@ -4971,15 +4970,15 @@ def edge_disjoint_spanning_trees(self,k, root=None, solver = None, verbose = 0):
         # The colors we can use
         colors = range(0,k)
 
-        # edges[j][e] is equal to one if and only if edge e belongs to color j
-        edges = p.new_variable()
+        # edges[j,e] is equal to one if and only if edge e belongs to color j
+        edges = p.new_variable(binary=True)
 
         if root is None:
             root = self.vertex_iterator().next()
 
         # r_edges is a relaxed variable grater than edges. It is used to
         # check the presence of cycles
-        r_edges = p.new_variable()
+        r_edges = p.new_variable(nonnegative=True)
 
         epsilon = 1/(3*(Integer(self.order())))
 
@@ -5049,9 +5048,6 @@ def edge_disjoint_spanning_trees(self,k, root=None, solver = None, verbose = 0):
         for j in colors:
             for v in self:
                 p.add_constraint(p.sum(r_edges[j,(u,v)] for u in self.neighbors(v)), max=1-epsilon)
-
-        p.set_binary(edges)
-
         try:
             p.solve(log = verbose)
 
@@ -5364,16 +5360,16 @@ def vertex_cut(self, s, t, value_only=True, vertices=False, solver=None, verbose
             value_only = False
 
         p = MixedIntegerLinearProgram(maximization=False, solver=solver)
-        b = p.new_variable()
-        v = p.new_variable()
+        b = p.new_variable(binary=True)
+        v = p.new_variable(binary=True)
 
         # Some vertices belong to part 1, some others to part 0
-        p.add_constraint(v[s], min=0, max=0)
-        p.add_constraint(v[t], min=1, max=1)
+        p.add_constraint(v[s] == 0)
+        p.add_constraint(v[t] == 1)
 
         # b indicates whether the vertices belong to the cut
-        p.add_constraint(b[s], min=0, max=0)
-        p.add_constraint(b[t], min=0, max=0)
+        p.add_constraint(b[s] ==0)
+        p.add_constraint(b[t] ==0)
 
         if g.is_directed():
 
@@ -5389,11 +5385,8 @@ def vertex_cut(self, s, t, value_only=True, vertices=False, solver=None, verbose
             # adjacent vertices belong to the same part except if one of them
             # belongs to the cut
             for (x,y) in g.edges(labels=None):
-                p.add_constraint(v[x] + b[y] - v[y],min=0)
-                p.add_constraint(v[y] + b[x] - v[x],min=0)
-
-        p.set_binary(b)
-        p.set_binary(v)
+                p.add_constraint(v[x] + b[y] >= v[y])
+                p.add_constraint(v[y] + b[x] >= v[x])
 
         if value_only:
             return Integer(round(p.solve(objective_only=True, log=verbose)))
@@ -5505,7 +5498,7 @@ def multiway_cut(self, vertices, value_only = False, use_edge_labels = False, so
         p = MixedIntegerLinearProgram(maximization = False, solver= solver)
 
         # height[c,v] represents the height of vertex v for commodity c
-        height = p.new_variable()
+        height = p.new_variable(nonnegative=True)
 
         # cut[e] represents whether e is in the cut
         cut = p.new_variable(binary = True)
@@ -5650,7 +5643,7 @@ def max_cut(self, value_only=True, use_edge_labels=False, vertices=False, solver
         p = MixedIntegerLinearProgram(maximization=True, solver=solver)
 
         in_set = p.new_variable(binary = True)
-        in_cut = p.new_variable()
+        in_cut = p.new_variable(binary = True)
 
         # A vertex has to be in some set
         for v in g:
@@ -5960,7 +5953,7 @@ def longest_path(self, s=None, t=None, use_edge_labels=False, algorithm="MILP",
 
         # relaxed version of the previous variable, to prevent the
         # creation of cycles
-        r_edge_used = p.new_variable()
+        r_edge_used = p.new_variable(nonnegative=True)
 
         # vertex_used[v] == 1 if vertex v is used
         vertex_used = p.new_variable(binary=True)
@@ -6444,8 +6437,8 @@ def traveling_salesman_problem(self, use_edge_labels = False, solver = None, con
 
         p = MixedIntegerLinearProgram(maximization = False, solver = solver)
 
-        f = p.new_variable()
-        r = p.new_variable()
+        f = p.new_variable(binary=True)
+        r = p.new_variable(nonnegative=True)
 
         eps = 1/(2*Integer(g.order()))
         x = g.vertex_iterator().next()
@@ -6514,8 +6507,6 @@ def traveling_salesman_problem(self, use_edge_labels = False, solver = None, con
         else:
             p.set_objective(None)
 
-        p.set_binary(f)
-
         try:
             obj = p.solve(log = verbose)
             f = p.get_values(f)
@@ -7011,7 +7002,7 @@ def flow(self, x, y, value_only=True, integer=False, use_edge_labels=True, verte
         from sage.numerical.mip import MixedIntegerLinearProgram
         g=self
         p=MixedIntegerLinearProgram(maximization=True, solver = solver)
-        flow=p.new_variable()
+        flow=p.new_variable(nonnegative=True)
 
         if use_edge_labels:
             from sage.rings.real_mpfr import RR
@@ -7334,7 +7325,7 @@ def multicommodity_flow(self, terminals, integer=True, use_edge_labels=False,ver
             set_terminals.add(t)
 
         # flow[i,(u,v)] is the flow of commodity i going from u to v
-        flow=p.new_variable()
+        flow=p.new_variable(nonnegative=True)
 
         # Whether to use edge labels
         if use_edge_labels:
@@ -7721,7 +7712,7 @@ def dominating_set(self, independent=False, value_only=False, solver=None, verbo
         from sage.numerical.mip import MixedIntegerLinearProgram
         g=self
         p=MixedIntegerLinearProgram(maximization=False, solver=solver)
-        b=p.new_variable()
+        b=p.new_variable(binary=True)
 
         # For any vertex v, one of its neighbors or v itself is in
         # the minimum dominating set
@@ -7737,8 +7728,6 @@ def dominating_set(self, independent=False, value_only=False, solver=None, verbo
         # Minimizes the number of vertices used
         p.set_objective(p.sum([b[v] for v in g.vertices()]))
 
-        p.set_integer(b)
-
         if value_only:
             return Integer(round(p.solve(objective_only=True, log=verbose)))
         else:
@@ -8111,7 +8100,7 @@ def vertex_connectivity(self, value_only=True, sets=False, solver=None, verbose=
         p = MixedIntegerLinearProgram(maximization=False, solver=solver)
 
         # Sets 0 and 2 are "real" sets while set 1 represents the cut
-        in_set = p.new_variable()
+        in_set = p.new_variable(binary=True)
 
         # A vertex has to be in some set
         for v in g:

src/sage/graphs/graph.py

@@ -2912,7 +2912,7 @@ def degree_constrained_subgraph(self, bounds=None, solver=None, verbose=0):
         from sage.numerical.mip import MixedIntegerLinearProgram, MIPSolverException
 
         p = MixedIntegerLinearProgram(maximization=False, solver=solver)
-        b = p.new_variable()
+        b = p.new_variable(binary=True)
 
         reorder = lambda x,y: (x,y) if x<y else (y,x)
 
@@ -2935,7 +2935,6 @@ def degree_constrained_subgraph(self, bounds=None, solver=None, verbose=0):
             p.add_constraint(p.sum([ b[reorder(x,y)]*weight(l) for x,y,l in self.edges_incident(v)]), min=minimum, max=maximum)
 
         p.set_objective(p.sum([ b[reorder(x,y)]*weight(l) for x,y,l in self.edge_iterator()]))
-        p.set_binary(b)
 
         try:
             p.solve(log=verbose)
@@ -3120,9 +3119,9 @@ def minimum_outdegree_orientation(self, use_edge_labels=False, solver=None, verb
         # and indicates whether the edge uv
         # with u<v goes from u to v ( equal to 0 )
         # or from v to u ( equal to 1)
-        orientation = p.new_variable(binary = True)
+        orientation = p.new_variable(binary=True)
 
-        degree = p.new_variable()
+        degree = p.new_variable(nonnegative=True)
 
         # Whether an edge adjacent to a vertex u counts
         # positively or negatively
@@ -3755,7 +3754,7 @@ def has_homomorphism_to(self, H, core = False, solver = None, verbose = 0):
         if core:
 
             # the value of m is one if the corresponding vertex of h is used.
-            m = p.new_variable()
+            m = p.new_variable(nonnegative=True)
             for uh in H:
                 for ug in self:
                     p.add_constraint(b[ug,uh] <= m[uh])
@@ -3833,7 +3832,7 @@ def fractional_chromatic_index(self, verbose_constraints = 0, verbose = 0):
         p = MixedIntegerLinearProgram(constraint_generation = True)
 
         # One variable per edge
-        r = p.new_variable()
+        r = p.new_variable(nonnegative=True)
         R = lambda x,y : r[x,y] if x<y else r[y,x]
 
         # We want to maximize the sum of weights on the edges
@@ -3946,8 +3945,8 @@ def maximum_average_degree(self, value_only=True, solver = None, verbose = 0):
 
         p = MixedIntegerLinearProgram(maximization=True, solver = solver)
 
-        d = p.new_variable()
-        one = p.new_variable()
+        d = p.new_variable(nonnegative=True)
+        one = p.new_variable(nonnegative=True)
 
         # Reorders u and v so that uv and vu are not considered
         # to be different edges
@@ -4059,7 +4058,7 @@ def independent_set_of_representatives(self, family, solver=None, verbose=0):
 
         # Boolean variable indicating whether the vertex
         # is the representative of some set
-        vertex_taken=p.new_variable()
+        vertex_taken=p.new_variable(binary=True)
 
         # Boolean variable in two dimension whose first
         # element is a vertex and whose second element
@@ -4091,8 +4090,6 @@ def independent_set_of_representatives(self, family, solver=None, verbose=0):
 
         p.set_objective(None)
 
-        p.set_binary(vertex_taken)
-
         try:
             p.solve(log=verbose)
         except Exception:
@@ -4232,7 +4229,7 @@ def minor(self, H, solver=None, verbose=0):
 
         # a tree  has no cycle
         epsilon = 1/(5*Integer(self.order()))
-        r_edges = p.new_variable()
+        r_edges = p.new_variable(nonnegative=True)
 
         for h in H:
             for u,v in self.edges(labels=None):
@@ -4243,7 +4240,7 @@ def minor(self, H, solver=None, verbose=0):
 
         # Once the representative sets are described, we must ensure
         # there are arcs corresponding to those of H between them
-        h_edges = p.new_variable()
+        h_edges = p.new_variable(nonnegative=True)
 
         for h1, h2 in H.edges(labels=None):
 
@@ -5448,7 +5445,7 @@ def vertex_cover(self, algorithm = "Cliquer", value_only = False,
 
             from sage.numerical.mip import MixedIntegerLinearProgram
             p = MixedIntegerLinearProgram(maximization=False, solver=solver)
-            b = p.new_variable()
+            b = p.new_variable(binary=True)
 
             # minimizes the number of vertices in the set
             p.set_objective(p.sum([b[v] for v in g.vertices()]))
@@ -5457,8 +5454,6 @@ def vertex_cover(self, algorithm = "Cliquer", value_only = False,
             for (u,v) in g.edges(labels=None):
                 p.add_constraint(b[u] + b[v], min=1)
 
-            p.set_binary(b)
-
             if value_only:
                 size_cover_g = p.solve(objective_only=True, log=verbosity)
             else:

src/sage/graphs/graph_coloring.py

@@ -646,11 +646,11 @@ def grundy_coloring(g, k, value_only = True, solver = None, verbose = 0):
     classes = range(k)
 
     # b[v,i] is set to 1 if and only if v is colored with i
-    b = p.new_variable()
+    b = p.new_variable(binary = True)
 
     # is_used[i] is set to 1 if and only if color [i] is used by some
     # vertex
-    is_used = p.new_variable()
+    is_used = p.new_variable(binary = True)
 
     # Each vertex is in exactly one class
     for v in g:
@@ -679,10 +679,6 @@ def grundy_coloring(g, k, value_only = True, solver = None, verbose = 0):
     for i in classes:
         p.add_constraint( p.sum( b[v,i] for v in g ) - is_used[i], min = 0)
 
-    # Both variables are binary
-    p.set_binary(b)
-    p.set_binary(is_used)
-
     # Trying to use as many colors as possible
     p.set_objective( p.sum( is_used[i] for i in classes ) )
 
@@ -1264,7 +1260,7 @@ class :class:`MixedIntegerLinearProgram
     c = p.new_variable(binary = True)
 
     # relaxed value
-    r = p.new_variable()
+    r = p.new_variable(nonnegative=True)
 
     E = lambda x,y : (x,y) if x<y else (y,x)
 
@@ -1467,7 +1463,7 @@ class :class:`MixedIntegerLinearProgram
     c = p.new_variable(binary = True)
 
     # relaxed value
-    r = p.new_variable()
+    r = p.new_variable(nonnegative=True)
 
     E = lambda x,y : (x,y) if x<y else (y,x)
 

src/sage/graphs/graph_decompositions/vertex_separation.pyx

@@ -816,10 +816,10 @@ def vertex_separation_MILP(G, integrality = False, solver = None, verbosity = 0)
     p = MixedIntegerLinearProgram( maximization = False, solver = solver )
 
     # Declaration of variables.
-    x = p.new_variable( binary = integrality)
-    u = p.new_variable( binary = integrality)
-    y = p.new_variable( binary = True)
-    z = p.new_variable( integer = True, dim = 1 )
+    x = p.new_variable(binary = integrality, nonnegative=bool(not integrality)) # at least one has to be set (#15221)
+    u = p.new_variable(binary = integrality, nonnegative=bool(not integrality))
+    y = p.new_variable(binary = True)
+    z = p.new_variable(integer = True)
 
     N = G.num_verts()
     V = G.vertices()

src/sage/numerical/backends/glpk_backend.pyx

@@ -386,6 +386,8 @@ cdef class GLPKBackend(GenericBackend):
 
             sage: p = MixedIntegerLinearProgram(solver='GLPK')
             sage: x,y = p[0], p[1]
+            doctest:839: DeprecationWarning: The default behaviour of new_variable() will soon change ! It will return 'real' variables instead of nonnegative ones. Please be explicit and call new_variable(nonnegative=True) instead.
+            See http://trac.sagemath.org/15521 for details.
             sage: p.add_constraint(2*x + 3*y, max = 6)
             sage: p.add_constraint(3*x + 2*y, max = 6)
             sage: p.set_objective(x + y + 7)