# Ticket #11279: trac_11279.patch

File trac_11279.patch, 21.8 KB (added by ncohen, 9 years ago)
• ## doc/en/reference/graphs.rst

# HG changeset patch
# User Nathann Cohen <nathann.cohen@gmail.com>
# Date 1304282108 -7200
# Node ID 21fa710c545fdd76a157e1677a77d1f6e78e617d
# Parent  e79daf00e92ab5227088289de977a30ae6beef7e
trac 11279 -- Convexity in graphs

diff --git a/doc/en/reference/graphs.rst b/doc/en/reference/graphs.rst
 a sage/graphs/pq_trees sage/graphs/matchpoly sage/graphs/graph_decompositions/vertex_separation sage/graphs/convexity_properties sage/graphs/graph_latex sage/graphs/graph_list
• ## module_list.py

diff --git a/module_list.py b/module_list.py
 a Extension('sage.graphs.graph_decompositions.vertex_separation', sources = ['sage/graphs/graph_decompositions/vertex_separation.pyx']), Extension('sage.graphs.convexity_properties', sources = ['sage/graphs/convexity_properties.pyx']), Extension('sage.graphs.generic_graph_pyx', sources = ['sage/graphs/generic_graph_pyx.pyx'], libraries = ['gmp']),
• ## new file sage/graphs/convexity_properties.pxd

diff --git a/sage/graphs/convexity_properties.pxd b/sage/graphs/convexity_properties.pxd
new file mode 100644
 - include "../misc/bitset_pxd.pxi" cdef class ConvexityProperties: cdef int _n cdef list _list_integers_to_vertices cdef dict _dict_vertices_to_integers cdef bitset_t * _cache_hull_pairs cdef list _vertices_to_integers(self, vertices) cdef list _integers_to_vertices(self, integers) cdef _bitset_convex_hull(self, bitset_t hull) cpdef hull(self, list vertices) cdef _greedy_increase(self, bitset_t bs) cpdef hull_number(self, value_only = *, verbose = *)
• ## new file sage/graphs/convexity_properties.pyx

diff --git a/sage/graphs/convexity_properties.pyx b/sage/graphs/convexity_properties.pyx
new file mode 100644
 - r""" Convexity properties of graphs This class gathers the algorithms related to convexity in a graph. AUTHOR: Nathann Cohen """ include "../misc/bitset.pxi" from sage.numerical.backends.generic_backend cimport GenericBackend from sage.numerical.backends.generic_backend import get_solver cdef class ConvexityProperties: r""" This class gathers the algorithms related to convexity in a graph. **Definitions** A set S \subseteq V(G) of vertices is said to be convex if for all u,v\in S the set S contains all the vertices located on a shortest path between u and v. Alternatively, a set S is said to be convex if the distances satisfy \forall u,v\in S, \forall w\in V\backslash S : d_{G}(u,w) + d_{G}(w,v) > d_{G}(u,v). The convex hull h(S) of a set S of vertices is defined as the smallest convex set containing S. It is a closure operator, as trivially S\subseteq h(S) and h(h(S)) = h(S). **What this class contains** As operations on convex sets generally involve the computation of distances between vertices, this class' purpose is to cache that information so that computing the convex hulls of several different sets of vertices does not imply recomputing several times the distances between the vertices. In order to compute the convex hull of a set S it is possible to write the following algorithm. *For any pair u,v of elements in the set S, and for any vertex w* *outside of it, add w to S if d_{G}(u,w) + d_{G}(w,v) = d_{G}(u,v).* *When no vertex can be added anymore, the set S is convex* The distances are not actually that relevant. The same algorithm can be implemented by remembering for each pair u, v of vertices the list of elements w satisfying the condition, and this is precisely what this class remembers, encoded as bitsets to make storage and union operations more efficient. .. NOTE:: * This class is useful if you compute the convex hulls of many sets in the same graph, or if you want to compute the hull number itself as it involves many calls to :meth:hull * Using this class on non-conected graphs is a waste of space and efficiency ! If your graph is disconnected, the best for you is to deal independently with each connected component, whatever you are doing. **Possible improvements** When computing a convex set, all the pairs of elements belonging to the set S are enumerated several times. * There should be a smart way to avoid enumerating pairs of vertices which have already been tested. The cost of each of them is not very high, so keeping track of those which have been tested already may be too expensive to gain any efficiency. * The ordering in which they are visited is currently purely lexicographic, while there is a Poset structure to exploit. In particular, when two vertices u, v are far apart and generate a set h(\{u,v\}) of vertices, all the pairs of vertices u', v'\in h(\{u,v\}) satisfy h(\{u',v'\}) \subseteq h(\{u,v\}), and so it is useless to test the pair u', v' when both u and v where present. * The information cached is for any pair u,v of vertices the list of elements z with d_{G}(u,w) + d_{G}(w,v) = d_{G}(u,v). This is not in general equal to h(\{u,v\}) ! Nothing says these recommandations will actually lead to any actual improvements. There are just some ideas remembered while writing this code. Trying to optimize may well lead to lost in efficiency on many instances. EXAMPLE:: sage: from sage.graphs.convexity_properties import ConvexityProperties sage: g = graphs.PetersenGraph() sage: CP = ConvexityProperties(g) sage: CP.hull([1,3]) [1, 2, 3] sage: CP.hull_number() 3 TESTS:: sage: ConvexityProperties(digraphs.Circuit(5)) Traceback (most recent call last): ... ValueError: This is currenly implemented for Graphs only.Only minor updates are needed if you want to makeit support DiGraphs too. """ def __init__(self, G): r""" Constructor EXAMPLE:: sage: from sage.graphs.convexity_properties import ConvexityProperties sage: g = graphs.PetersenGraph() sage: ConvexityProperties(g) """ from sage.graphs.digraph import DiGraph if isinstance(G, DiGraph): raise ValueError("This is currenly implemented for Graphs only."+ "Only minor updates are needed if you want to make"+ "it support DiGraphs too.") # Cached number of vertices cdef int n = G.order() self._n = n cdef int i = 0 cdef int j,k # Temporary variables cdef dict d_i cdef dict d_j cdef int d_ij self._dict_vertices_to_integers = {} self._list_integers_to_vertices = [] # Remembering integers instead of the labels, and building dictionaries # in both directions. for v in G: self._dict_vertices_to_integers[v] = i self._list_integers_to_vertices.append(v) i = i + 1 # Computation of distances between all pairs. Costly. cdef dict distances = G.distance_all_pairs() # _cache_hull_pairs[u*n + v] is a bitset whose 1 bits are the vertices located on a shortest path from vertex u to v # # Note that  u < v self._cache_hull_pairs = sage_malloc(((n*(n-1))>>1)*sizeof(bitset_t)) cdef bitset_t * p_bitset = self._cache_hull_pairs # Filling the cache # # The p_bitset variable iterates over the successive elements of the cache # # For any pair i,j of vertices (i """ cdef bitset_t * p_bitset = self._cache_hull_pairs cdef int i for 0 <= i < ((self._n*(self._n-1))>>1): bitset_free(p_bitset[0]) p_bitset = p_bitset + 1 sage_free(self._cache_hull_pairs) cdef list _vertices_to_integers(self, vertices): r""" Converts a list of vertices to a list of integers with the cached data. """ cdef list answer = [] for v in v: answer.append(self._dict_vertices_to_integers[v]) return answer cdef list _integers_to_vertices(self, integers): r""" Converts a list of integers to a list of vertices with the cached data. """ cdef list answer = [] for v in integers: answer.append(self._list_integers_to_vertices[v]) return answer cdef _bitset_convex_hull(self, bitset_t hull): r""" Computes the convex hull of a list of vertices given as a bitset. (this method returns nothing and modifies the input) """ cdef int count cdef int tmp_count cdef int i,j cdef bitset_t * p_bitset # Current size of the set count = bitset_len(hull) while True: # Iterating over all the elements in the cache p_bitset = self._cache_hull_pairs # For any vertex i for 0 <= i < self._n-1: # If i is not in the current set, we skip it ! if not bitset_in(hull, i): p_bitset = p_bitset + (self._n-1-i) continue # If it is, we iterate over all the elements j for i < j < self._n: # If both i and j are inside, we add all the (cached) # vertices on a shortest ij-path if bitset_in(hull, j): bitset_union(hull, hull, p_bitset[0]) # Next bitset ! p_bitset = p_bitset + 1 tmp_count = bitset_len(hull) # If we added nothing new during the previous loop, our set is # convex ! if tmp_count == count: return # Otherwise, update and back to the loop count = tmp_count cpdef hull(self, list vertices): r""" Returns the convex hull of a set of vertices. INPUT: * vertices -- A list of vertices. EXAMPLE:: sage: from sage.graphs.convexity_properties import ConvexityProperties sage: g = graphs.PetersenGraph() sage: CP = ConvexityProperties(g) sage: CP.hull([1,3]) [1, 2, 3] """ cdef bitset_t bs bitset_init(bs, self._n) bitset_set_first_n(bs, 0) for v in vertices: bitset_add(bs, self._dict_vertices_to_integers[v]) self._bitset_convex_hull(bs) #cdef list answer = bitset_list(bs) cdef list answer = self._integers_to_vertices(bitset_list(bs)) bitset_free(bs) return answer cdef _greedy_increase(self, bitset_t bs): r""" Given a bitset whose hull is not the whole set, greedily add vertices and stop before its hull is the whole set. NOTE: * Counting the bits at each turn is not the best way... """ cdef bitset_t tmp bitset_init(tmp, self._n) for 0<= i < self._n: if not bitset_in(bs, i): bitset_copy(tmp, bs) bitset_add(tmp, i) self._bitset_convex_hull(tmp) if bitset_len(tmp) < self._n: bitset_add(bs, i) cpdef hull_number(self, value_only = True, verbose = False): r""" Computes the hull number and a correspondin generating set. The hull number hn(G) of a graph G is the cardinality of a smallest set of vertices S such that h(S)=V(G). INPUT: * value_only (boolean) -- whether to return only the hull number (default) or a minimum set whose convex hull is the whole graph. * verbose (boolean) -- whether to display information on the LP. **COMPLEXITY:** This problem is NP-Hard [CHZ02]_, but seems to be of the "nice" kind. Update this comment if you fall on hard instances :-) **ALGORITHM:** This is solved by linear programming. As the function h(S) associating to each set S its convex hull is a closure operator, it is clear that any set S_G of vertices such that h(S_G)=V(G) must satisfy S_G \not \subseteq C for any *proper* convex set C \subsetneq V(G). The following formulation is hence correct .. MATH:: \text{Minimize :}& \sum_{v\in G}b_v\\ \text{Such that :}&\\ &\forall C\subsetneq V(G)\text{ a proper convex set }\\ &\sum_{v\in V(G)\backslash C} b_v \geq 1 Of course, the number of convex sets -- and so the number of constraints -- can be huge, and hard to enumerate, so at first an incomplete formulation is solved (it is missing some constraints). If the answer returned by the LP solver is a set S generating the whole graph, then it is optimal and so is returned. Otherwise, the constraint corresponding to the set h(S) can be added to the LP, which makes the answer S infeasible, and another solution computed. This being said, simply adding the constraint corresponding to h(S) is a bit slow, as these sets can be large (and the corresponding constrait a bit weak). To improve it a bit, before being added, the set h(S) is "greedily enriched" to a set S' with vertices for as long as h(S')\neq V(G). This way, we obtain a set S' with h(S)\subseteq h(S')\subsetneq V(G), and the constraint corresponding to h(S') -- which is stronger than the one corresponding to h(S) -- is added. This can actually be seen as a hitting set problem on the complement of convex sets. EXAMPLE: The Hull number of Petersen's graph:: sage: from sage.graphs.convexity_properties import ConvexityProperties sage: g = graphs.PetersenGraph() sage: CP = ConvexityProperties(g) sage: CP.hull_number() 3 sage: generating_set = CP.hull_number(value_only = False) sage: CP.hull(generating_set) [0, 1, 2, 3, 4, 5, 6, 7, 8, 9] REFERENCE: .. [CHZ02] F. Harary, E. Loukakis, C. Tsouros The geodetic number of a graph Mathematical and computer modelling vol. 17 n11 pp.89--95, 1993 """ cdef int i cdef list constraint # temporary variable to add constraints to the LP if self._n <= 2: if value_only: return self._n else: return self._list_integers_to_vertices cdef GenericBackend p = get_solver(constraint_generation = True) # Minimization p.set_sense(False) # We have exactly n binary variables, all of them with a coefficient of # 1 in the objective function p.add_variables(self._n, 0, None, True, False, False, 1, None) # We know that at least 2 vertices are required to cover the whole graph p.add_linear_constraint(zip(range(self._n), [1]*self._n), 2, None) # The set of vertices generated by the current LP solution cdef bitset_t current_hull bitset_init(current_hull, self._n) # Which is at first empty bitset_set_first_n(current_hull,1) while True: # Greedily increase it to obtain a better constraint self._greedy_increase(current_hull) if verbose: print "Adding a constraint corresponding to convex set ", print bitset_list(current_hull) # Building the corresponding constraint constraint = [] for 0 <= i < self._n: if not bitset_in(current_hull, i): constraint.append((i,1)) p.add_linear_constraint(constraint, 1, None) p.solve() # Computing the current solution's convex hull bitset_set_first_n(current_hull,0) for 0 <= i < self._n: if p.get_variable_value(i) > .5: bitset_add(current_hull, i) self._bitset_convex_hull(current_hull) # Are we done ? if bitset_len(current_hull) == self._n: break bitset_free(current_hull) if value_only: return p.get_objective_value() constraint = [] for 0 <= i < self._n: if p.get_variable_value(i) > .5: constraint.append(i) return self._integers_to_vertices(constraint)
• ## sage/graphs/graph.py

diff --git a/sage/graphs/graph.py b/sage/graphs/graph.py
 a from matchpoly import matching_polynomial return matching_polynomial(self, complement=complement, name=name) ### Convexity def convexity_properties(self): r""" Returns a ConvexityProperties objet corresponding to self. This object contains the methods related to convexity in graphs (convex hull, hull number) and caches useful information so that it becomes comparatively cheaper to compute the convex hull of many different sets of the same graph. .. SEEALSO:: In order to know what can be done through this object, please refer to module :mod:sage.graphs.convexity_properties. .. NOTE:: If you want to compute many convex hulls, keep this object in memory ! When it is created, it builds a table of useful information to compute convex hulls. As a result :: sage: g = graphs.PetersenGraph() sage: g.convexity_properties().hull([1, 3]) [1, 2, 3] sage: g.convexity_properties().hull([3, 7]) [2, 3, 7] Is a terrible waste of computations, while :: sage: g = graphs.PetersenGraph() sage: CP = g.convexity_properties() sage: CP.hull([1, 3]) [1, 2, 3] sage: CP.hull([3, 7]) [2, 3, 7] Makes perfect sense. """ from sage.graphs.convexity_properties import ConvexityProperties return ConvexityProperties(self) ### Centrality def centrality_betweenness(self, normalized=True):
• ## sage/numerical/backends/generic_backend.pxd

diff --git a/sage/numerical/backends/generic_backend.pxd b/sage/numerical/backends/generic_backend.pxd
 a cpdef bint is_variable_binary(self, int) cpdef bint is_variable_integer(self, int) cpdef bint is_variable_continuous(self, int) cpdef problem_name(self, char * name = *) cpdef row_bounds(self, int index) cpdef col_bounds(self, int index) cpdef col_name(self, int index) cpdef variable_upper_bound(self, int index, value = *) cpdef variable_lower_bound(self, int index, value = *) cpdef GenericBackend get_solver(constraint_generation = ?, solver = ?)
• ## sage/numerical/backends/glpk_backend.pyx

diff --git a/sage/numerical/backends/glpk_backend.pyx b/sage/numerical/backends/glpk_backend.pyx
 a glp_set_col_kind(self.lp, n_var - i, GLP_IV) if obj: self.bjective_coefficient(n_var - i - 1, obj) self.objective_coefficient(n_var - i - 1, obj) if names is not None: glp_set_col_name(self.lp, n_var, names[number - i - 1])