# Ticket #7301: trac_7301.2.patch

File trac_7301.2.patch, 10.0 KB (added by wdj, 13 years ago)

apply this patch only to 4.3.

• ## sage/combinat/integer_vector.py

# HG changeset patch
# User DavidJoyner <wdjoyner@gmail.com>
# Date 1262003958 18000
# Node ID 04bf06a1338ed82c40a1f7295cbec9a9a39ef91c
# Parent  21efb0b3fc474972b5c7f617d99173536a3d79d0
implementation of Gale-Ryser theorem constructions

diff -r 21efb0b3fc47 -r 04bf06a1338e sage/combinat/integer_vector.py
 a import functools def is_gale_ryser(r,s): r""" INPUT:  r, s - integer sequences $r=(r_1,...,r_m)$ and $s=(s_1,...,s_n)$, where $r_1\geq r_2 ... geq r_m\geq 0$, $s_1\geq s_2 ... geq s_n\geq 0$. Let $r^*$ denote the conjugate of $r$. OUTPUT: The boolean value of $s\preceq r^*$, where $\preceq$ is the ordering on sequences above. We say that a pair of such sequences $(r,s)$ is {\bf Gale-Ryser} if  $s\leq r^*$ (in which case the hypothesis to the Gale-Ryser Theorem is satisfied). EXAMPLES: sage: from sage.combinat.integer_vector import is_gale_ryser sage: is_gale_ryser([4,2,2],[3,3,1,1]) True sage: is_gale_ryser([4,2,1,1],[3,3,1,1]) True sage: is_gale_ryser([3,2,1,1],[3,3,1,1]) False REMARK: In the literature, what we are calling a Gale-Ryser sequence sometimes goes by the (rather generic-sounding) term realizable sequence''. """ from sage.combinat.partition import Partition r2 = [x for x in r if x>0] # remove extraneous 0's s2 = [x for x in s if x>0] rstar = Partition(r2).conjugate() m = len(rstar) if m > len(s2): return False for i in range(1,m): #print i, sum(s2[:i]), sum(rstar[:i]) if sum(rstar[:i]) < sum(s2[:i]): return False if sum(rstar) <> sum(s2): return False return True def slider01(A, t, k): r""" Assumes $A$ is a $(0,1)$-matrix. For each of the $t$ rows with highest row sums, this function returns a matrix $B$ which is the same as $A$ except that it has slid exactly one $1$ in $A$ over to the $k$-th column. INPUT:: A - an $m\times n$ (0,1) matrix t, k - integers satisfying $0 < t < m$, $0 < k < n$ OUTPUT:: An $m\times n$ (0,1) matrix, which is the same as $A$ except that it has exactly one $1$ in $A$ slid over to the $k$-th column. EXAMPLES:: sage: from sage.combinat.integer_vector import slider01 sage: A = matrix([[1,1,1,0],[1,1,1,0],[1,0,0,0],[1,0,0,0]]) sage: A [1 1 1 0] [1 1 1 0] [1 0 0 0] [1 0 0 0] sage: slider01(A, 1, 3) [1 1 1 0] [1 1 0 1] [1 0 0 0] [1 0 0 0] sage: slider01(A, 3, 3) [1 1 1 0] [1 1 1 0] [1 0 0 0] [0 0 0 1] """ from sage.matrix.constructor import matrix maxrows = [] m = len(A.rows()) n = len(A.columns()) B = matrix(m,n) rowsums = [sum(A.row(i)) for i in range(m)] maxrows = [rowsums[i] for i in range(t)] # t rows with max rowsums maxRows = [] # want these to correspond to elts of maxrows B = [list(A.row(i)) for i in range(m)] for i in range(m): if sum(A.row(i)) in maxrows: maxRows.append(i)  # >= t rows numbers with max rowsums for j in range(len(maxrows)): i = maxRows[-1] rw = A.row(i) # to make mutable # now we want to move the rightmost left 1 to the k-th column rwlist = list(rw) if not(0 in rwlist): N = n-1 else: N = rwlist.index(0) rwlist[int(N)-1] = 0 rwlist[k] = 1 B[i] = rwlist return matrix(B) def gale_ryser_theorem(p1, p2, algorithm="ryser"): r""" Returns the binary matrix given by the Gale-Ryser theorem. The Gale Ryser theorem asserts that if p_1,p_2 are two partitions of n of respective lengths k_1,k_2, then there is a binary k_1\times k_2 matrix M such that p_1 is the vector of row sums and p_2 is the vector of column sums of M, if and only if p_2 dominates p_1. INPUT:: - p1 -- the first partition of n (trailing 0's allowed) - p2 -- the second partition of n (trailing 0's allowed) OUTPUT:: - A binary matrix OPTION:: algorithm - "ryser" (default) implements the construction due to Ryser [1]. - "gale" (or anything other than "gale") implements the construction due to Gale [2]. ALGORITHM: (Gale:) The Gale-Ryser theorem asserts that the following Linear Program admits a solution. .. MATH:: \forall i&\sum_{j=1}^{k_2} b_{i,j}=p_{1,j}\\ \forall i&\sum_{j=1}^{k_1} b_{j,i}=p_{2,j}\\ &b_{i,j}\mbox{ is a binary variable} (Ryser:) The construction of an $m\times n$ matrix $A=A_{r,s}$, due to Ryser [1], is described as follows. The construction works if and only if have $s\preceq r^*$. * Construct the $m\times n$ matrix $B$ from $r$ by defining the $i$-th row of $B$ to be the vector whose first $r_i$ entries are $1$, and the remainder are $0$'s, $1\leq i\leq m$. This maximal matrix $B$ with row sum $r$ and ones left justified has column sum $r^{*}$. * Shift the last $1$ in certain rows of $B$ to column $n$ in order to achieve the sum $s_n$. Call this $B$ again. * The $1$'s in column n are to appear in those rows in which $A$ has the largest row sums, giving preference to the bottom-most positions in case of ties. * Note: When this step automatically "fixes" other columns, one must skip ahead to the first column index with a wrong sum in the step below. * Proceed inductively to construct columns $n-1$, ..., $2$, $1$. * Set $A = B$. Return $A$. EXAMPLES:: Computing the matrix for p_1=p_2=2+2+1 :: sage: from sage.combinat.integer_vector import gale_ryser_theorem sage: p1 = [2,2,1] sage: p2 = [2,2,1] sage: print gale_ryser_theorem(p1, p2, algorithm="gale") # Optional - requires GLPK or CBC [0 1 1] [1 1 0] [1 0 0] Or for a non-square matrix with p_1=3+3+2+1 and p_2=3+2+2+1+1 :: sage: p1 = [3,3,2,1] sage: p2 = [3,2,2,1,1] sage: print gale_ryser_theorem(p1, p2, algorithm="gale") # Optional - requires GLPK or CBC [1 0 1 1 0] [1 0 1 0 1] [1 1 0 0 0] [0 1 0 0 0] sage: p1 = [3,3,1,1] sage: p2 = [3,3,1,1] sage: gale_ryser_theorem(p1, p2) [1 1 1 0] [1 1 0 1] [1 0 0 0] [0 1 0 0] sage: p1 = [4,2,2] sage: p2 = [3,3,1,1] sage: gale_ryser_theorem(p1, p2) [1 1 1 1] [1 1 0 0] [1 1 0 0] sage: p1 = [4,2,2,0] sage: p2 = [3,3,1,1,0,0] sage: gale_ryser_theorem(p1, p2) [1 1 1 1 0 0] [1 1 0 0 0 0] [1 1 0 0 0 0] [0 0 0 0 0 0] REFERENCES:: [1] H. J. Ryser, Combinatorial Mathematics, Carus Monographs, MAA, 1963. [2] D. Gale, A theorem on flows in networks, Pacific J. Math. 7(1957)1073-1082. """ from sage.combinat.partition import Partition from sage.matrix.constructor import matrix if sum(p1) != sum(p2): raise ValueError("The two lists must sum to the same value.") if not Partition(p2).conjugate().dominates(Partition(p1)): raise ValueError("The domination criterion is not fulfilled.") if algorithm == "gale": from sage.numerical.mip import MixedIntegerLinearProgram p1 = Partition(p1) p2 = Partition(p2) k1=p1.length() k2=p2.length() p = MixedIntegerLinearProgram() b = p.new_variable(dim=2) for (i,c) in enumerate(p1): p.add_constraint(sum([b[i][j] for j in xrange(k2)]),min=c,max=c) for (i,c) in enumerate(p2): p.add_constraint(sum([b[j][i] for j in xrange(k1)]),min=c,max=c) p.set_objective(None) p.set_binary(b) p.solve() b = p.get_values(b) M = [[0]*k2 for i in xrange(k1)] for i in xrange(k1): for j in xrange(k2): M[i][j] = int(b[i][j]) return matrix(M) else: r = p1; s = p2 if not(is_gale_ryser(r,s)): return False rstar = Partition(r).conjugate() m = len(r) n = len(s) rowsA0 = [[0]*n]*m for j in range(m): if js[n-k]: break A0 = slider01(A0,s[n-k],n-k) return A0 def _default_function(l, default, i): """ EXAMPLES:: """ return infinity class IntegerVectors_nk(CombinatorialClass): def __init__(self, n, k): """