Contact
CoCalc Logo Icon
StoreFeaturesDocsShareSupport News Sign UpSign In
| Download

Here we provide code which can explicitly verify the sum of squares counterexample appearing in the article https://arxiv.org/abs/1909.00081.

Views: 373
License: MIT
Image: ubuntu2004
Kernel: SageMath 9.2

An SOS counterexample to an inequality of symmetric functions

This document verifies the fact that a certain symmetric polynomial can be written as a sum of squares. This provides a counterexample to a certain conjecture relating nonnegativity and the dominance (or majorization) partial order on partitions. See https://arxiv.org/abs/1909.00081 for more information. Our symmetry-adapted SDP returned a numerical matrix (with floating point entries). Since we had reduced the problem size using representation theory of the symmetric group, and also knowledge of the real zeros of our polynomial, we were able to successfully round the floating point entries into (exact) rational numbers.

We load this matrix AA, whose entries are rational numbers, as well as a permutation matrix which reorders the columns and rows so that we can apply naive LDLTLDL^T decomposition below. In order to apply LDLTLDL^T decomposition, it is sometimes required to reorder the rows and columns in order to avoid zero pivots.

P = load("counterexample-P.sobj") A = load("counterexample-A.sobj") B = P.T * A * P C = B.change_ring(QQ)

Below we apply LDLTLDL^T decomposition to extract pairs (d,)(d,\ell) where dd is the pivot entry and \ell is a vector. These pairs can be used to write C=dTC = \sum d \cdot \ell \ell^T .

def naiveLDLt(A): N = A.nrows() terms = [] # (d,l) pair of diagonal entry and vector for i in range(0,N): pivot = A[i,i] if pivot != 0: l = 1/pivot * vector(list(A.columns()[i])) terms.append((pivot,l)) A = A - pivot * l.outer_product(l) return terms terms = naiveLDLt(C)

Below we create a vector of monomials with total degree 88 in three variables. The permutation matrix PP loaded above reorders this basis appropriately to match the reordering of the matrix AA.

def create_monomials(n,d): varz = [var('x%s'%i) for i in range(1,n+1)] mons = [] for exps in IntegerVectors(d,n): prod = varz[0] - varz[0] + 1 for i,x in enumerate(varz): prod *= x^exps[i] mons.append(prod) return mons m = create_monomials(3,8) m_shuffled = P.T * vector(m) show(m_shuffled)
(x13x23x32,x13x22x33,x12x23x33,x14x22x32,x24x34,x14x24,x14x34,x14x23x3,x13x24x3,x13x2x34,x14x2x33,x1x24x33,x1x23x34,x12x24x32,x15x23,x15x33,x25x33,x12x22x34,x13x25,x13x35,x23x35,x12x25x3,x1x25x32,x15x22x3,x12x36,x12x2x35,x1x22x35,x15x2x32,x12x26,x26x32,x18,x1x2x36,x1x26x3,x16x2x3,x38,x28,x1x37,x27x3,x1x27,x16x32,x16x22,x17x2,x22x36,x2x37,x17x3)\renewcommand{\Bold}[1]{\mathbf{#1}}\left(x_{1}^{3} x_{2}^{3} x_{3}^{2},\,x_{1}^{3} x_{2}^{2} x_{3}^{3},\,x_{1}^{2} x_{2}^{3} x_{3}^{3},\,x_{1}^{4} x_{2}^{2} x_{3}^{2},\,x_{2}^{4} x_{3}^{4},\,x_{1}^{4} x_{2}^{4},\,x_{1}^{4} x_{3}^{4},\,x_{1}^{4} x_{2}^{3} x_{3},\,x_{1}^{3} x_{2}^{4} x_{3},\,x_{1}^{3} x_{2} x_{3}^{4},\,x_{1}^{4} x_{2} x_{3}^{3},\,x_{1} x_{2}^{4} x_{3}^{3},\,x_{1} x_{2}^{3} x_{3}^{4},\,x_{1}^{2} x_{2}^{4} x_{3}^{2},\,x_{1}^{5} x_{2}^{3},\,x_{1}^{5} x_{3}^{3},\,x_{2}^{5} x_{3}^{3},\,x_{1}^{2} x_{2}^{2} x_{3}^{4},\,x_{1}^{3} x_{2}^{5},\,x_{1}^{3} x_{3}^{5},\,x_{2}^{3} x_{3}^{5},\,x_{1}^{2} x_{2}^{5} x_{3},\,x_{1} x_{2}^{5} x_{3}^{2},\,x_{1}^{5} x_{2}^{2} x_{3},\,x_{1}^{2} x_{3}^{6},\,x_{1}^{2} x_{2} x_{3}^{5},\,x_{1} x_{2}^{2} x_{3}^{5},\,x_{1}^{5} x_{2} x_{3}^{2},\,x_{1}^{2} x_{2}^{6},\,x_{2}^{6} x_{3}^{2},\,x_{1}^{8},\,x_{1} x_{2} x_{3}^{6},\,x_{1} x_{2}^{6} x_{3},\,x_{1}^{6} x_{2} x_{3},\,x_{3}^{8},\,x_{2}^{8},\,x_{1} x_{3}^{7},\,x_{2}^{7} x_{3},\,x_{1} x_{2}^{7},\,x_{1}^{6} x_{3}^{2},\,x_{1}^{6} x_{2}^{2},\,x_{1}^{7} x_{2},\,x_{2}^{2} x_{3}^{6},\,x_{2} x_{3}^{7},\,x_{1}^{7} x_{3}\right)

Below we sum the squares. We use each vector \ell to create a polynomial m(x)\ell m(x), which we square and sum with positive coefficient dd. The resulting polynomial displayed below is therefore a sum of squares with positive coefficients.

okay = True for tup in terms: d,l = tup[0],tup[1] if d < 0: okay = False show(okay)
True\renewcommand{\Bold}[1]{\mathbf{#1}}\mathrm{True}
ans = [] for tup in terms: d,l = tup[0],tup[1] term = d * (l*m_shuffled)^2 term = term.expand() ans.append(term) res = sum(ans) show(res.simplify())
179450x116+11050x114x22+19450x112x24+1525x110x26+2315x18x28+1525x16x210+19450x14x212+11050x12x214+179450x216+11050x114x32164725x112x22x3281575x110x24x32+119450x18x26x32+119450x16x28x3281575x14x210x32164725x12x212x32+11050x214x32+19450x112x3481575x110x22x34114725x18x24x3411890x16x26x34114725x14x28x3481575x12x210x34+19450x212x34+1525x110x36+119450x18x22x3611890x16x24x3611890x14x26x36+119450x12x28x36+1525x210x36+2315x18x38+119450x16x22x38114725x14x24x38+119450x12x26x38+2315x28x38+1525x16x31081575x14x22x31081575x12x24x310+1525x26x310+19450x14x312164725x12x22x312+19450x24x312+11050x12x314+11050x22x314+179450x316\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{17}{9450} \, x_{1}^{16} + \frac{1}{1050} \, x_{1}^{14} x_{2}^{2} + \frac{1}{9450} \, x_{1}^{12} x_{2}^{4} + \frac{1}{525} \, x_{1}^{10} x_{2}^{6} + \frac{2}{315} \, x_{1}^{8} x_{2}^{8} + \frac{1}{525} \, x_{1}^{6} x_{2}^{10} + \frac{1}{9450} \, x_{1}^{4} x_{2}^{12} + \frac{1}{1050} \, x_{1}^{2} x_{2}^{14} + \frac{17}{9450} \, x_{2}^{16} + \frac{1}{1050} \, x_{1}^{14} x_{3}^{2} - \frac{16}{4725} \, x_{1}^{12} x_{2}^{2} x_{3}^{2} - \frac{8}{1575} \, x_{1}^{10} x_{2}^{4} x_{3}^{2} + \frac{11}{9450} \, x_{1}^{8} x_{2}^{6} x_{3}^{2} + \frac{11}{9450} \, x_{1}^{6} x_{2}^{8} x_{3}^{2} - \frac{8}{1575} \, x_{1}^{4} x_{2}^{10} x_{3}^{2} - \frac{16}{4725} \, x_{1}^{2} x_{2}^{12} x_{3}^{2} + \frac{1}{1050} \, x_{2}^{14} x_{3}^{2} + \frac{1}{9450} \, x_{1}^{12} x_{3}^{4} - \frac{8}{1575} \, x_{1}^{10} x_{2}^{2} x_{3}^{4} - \frac{11}{4725} \, x_{1}^{8} x_{2}^{4} x_{3}^{4} - \frac{1}{1890} \, x_{1}^{6} x_{2}^{6} x_{3}^{4} - \frac{11}{4725} \, x_{1}^{4} x_{2}^{8} x_{3}^{4} - \frac{8}{1575} \, x_{1}^{2} x_{2}^{10} x_{3}^{4} + \frac{1}{9450} \, x_{2}^{12} x_{3}^{4} + \frac{1}{525} \, x_{1}^{10} x_{3}^{6} + \frac{11}{9450} \, x_{1}^{8} x_{2}^{2} x_{3}^{6} - \frac{1}{1890} \, x_{1}^{6} x_{2}^{4} x_{3}^{6} - \frac{1}{1890} \, x_{1}^{4} x_{2}^{6} x_{3}^{6} + \frac{11}{9450} \, x_{1}^{2} x_{2}^{8} x_{3}^{6} + \frac{1}{525} \, x_{2}^{10} x_{3}^{6} + \frac{2}{315} \, x_{1}^{8} x_{3}^{8} + \frac{11}{9450} \, x_{1}^{6} x_{2}^{2} x_{3}^{8} - \frac{11}{4725} \, x_{1}^{4} x_{2}^{4} x_{3}^{8} + \frac{11}{9450} \, x_{1}^{2} x_{2}^{6} x_{3}^{8} + \frac{2}{315} \, x_{2}^{8} x_{3}^{8} + \frac{1}{525} \, x_{1}^{6} x_{3}^{10} - \frac{8}{1575} \, x_{1}^{4} x_{2}^{2} x_{3}^{10} - \frac{8}{1575} \, x_{1}^{2} x_{2}^{4} x_{3}^{10} + \frac{1}{525} \, x_{2}^{6} x_{3}^{10} + \frac{1}{9450} \, x_{1}^{4} x_{3}^{12} - \frac{16}{4725} \, x_{1}^{2} x_{2}^{2} x_{3}^{12} + \frac{1}{9450} \, x_{2}^{4} x_{3}^{12} + \frac{1}{1050} \, x_{1}^{2} x_{3}^{14} + \frac{1}{1050} \, x_{2}^{2} x_{3}^{14} + \frac{17}{9450} \, x_{3}^{16}

Below we create the symmetric polynomial H44H521H_{44} - H_{521}, evaluated at x12,x22,x32x_1^2, x_2^2, x_3^2 to check nonnegativity on the nonnegative orthant, as explained in https://arxiv.org/abs/1909.00081. The nonnegativity of this polynomial provides the counterexample to the conjecture.

def input_squares(p): # input is our polynomial p, like H_(2,1) - H_(1,1,1) # output will be inserting each variable squared, also making it start at x1, x2,... rather than x0,x1,... current_vars = p.args() # this is a tuple of variables x0, x1, x2,... etc. numvars = len(current_vars) new_vars = [var('x%i'%i) for i in range(1,numvars+1)] subz = {} for i,vr in enumerate(current_vars): subz[vr] = (new_vars[i])^2 ans = p.subs(subz) return ans def plugin_ones(p): # want to return the number/scalar that results from plugging in 1's to all variables current_vars = p.args() # this is a tuple of variables x0, x1, x2,... etc. numvars = len(current_vars) subz = {} for i,vr in enumerate(current_vars): subz[vr] = 1 ans = p.subs(subz) return ans def term_normalize(p): norm = plugin_ones(p) ans = p/norm return ans def create_difference(la,mu,n=3): h = SymmetricFunctions(QQ).h() # we will return the difference h_la - h_mu hla = h(la).expand(n) hmu = h(mu).expand(n) hla = input_squares(hla) hmu = input_squares(hmu) hla = term_normalize(hla) hmu = term_normalize(hmu) ans = hla - hmu return ans show(create_difference((4,4),(5,2,1)))
179450x116+11050x114x22+19450x112x24+1525x110x26+2315x18x28+1525x16x210+19450x14x212+11050x12x214+179450x216+11050x114x32164725x112x22x3281575x110x24x32+119450x18x26x32+119450x16x28x3281575x14x210x32164725x12x212x32+11050x214x32+19450x112x3481575x110x22x34114725x18x24x3411890x16x26x34114725x14x28x3481575x12x210x34+19450x212x34+1525x110x36+119450x18x22x3611890x16x24x3611890x14x26x36+119450x12x28x36+1525x210x36+2315x18x38+119450x16x22x38114725x14x24x38+119450x12x26x38+2315x28x38+1525x16x31081575x14x22x31081575x12x24x310+1525x26x310+19450x14x312164725x12x22x312+19450x24x312+11050x12x314+11050x22x314+179450x316\renewcommand{\Bold}[1]{\mathbf{#1}}\frac{17}{9450} \, x_{1}^{16} + \frac{1}{1050} \, x_{1}^{14} x_{2}^{2} + \frac{1}{9450} \, x_{1}^{12} x_{2}^{4} + \frac{1}{525} \, x_{1}^{10} x_{2}^{6} + \frac{2}{315} \, x_{1}^{8} x_{2}^{8} + \frac{1}{525} \, x_{1}^{6} x_{2}^{10} + \frac{1}{9450} \, x_{1}^{4} x_{2}^{12} + \frac{1}{1050} \, x_{1}^{2} x_{2}^{14} + \frac{17}{9450} \, x_{2}^{16} + \frac{1}{1050} \, x_{1}^{14} x_{3}^{2} - \frac{16}{4725} \, x_{1}^{12} x_{2}^{2} x_{3}^{2} - \frac{8}{1575} \, x_{1}^{10} x_{2}^{4} x_{3}^{2} + \frac{11}{9450} \, x_{1}^{8} x_{2}^{6} x_{3}^{2} + \frac{11}{9450} \, x_{1}^{6} x_{2}^{8} x_{3}^{2} - \frac{8}{1575} \, x_{1}^{4} x_{2}^{10} x_{3}^{2} - \frac{16}{4725} \, x_{1}^{2} x_{2}^{12} x_{3}^{2} + \frac{1}{1050} \, x_{2}^{14} x_{3}^{2} + \frac{1}{9450} \, x_{1}^{12} x_{3}^{4} - \frac{8}{1575} \, x_{1}^{10} x_{2}^{2} x_{3}^{4} - \frac{11}{4725} \, x_{1}^{8} x_{2}^{4} x_{3}^{4} - \frac{1}{1890} \, x_{1}^{6} x_{2}^{6} x_{3}^{4} - \frac{11}{4725} \, x_{1}^{4} x_{2}^{8} x_{3}^{4} - \frac{8}{1575} \, x_{1}^{2} x_{2}^{10} x_{3}^{4} + \frac{1}{9450} \, x_{2}^{12} x_{3}^{4} + \frac{1}{525} \, x_{1}^{10} x_{3}^{6} + \frac{11}{9450} \, x_{1}^{8} x_{2}^{2} x_{3}^{6} - \frac{1}{1890} \, x_{1}^{6} x_{2}^{4} x_{3}^{6} - \frac{1}{1890} \, x_{1}^{4} x_{2}^{6} x_{3}^{6} + \frac{11}{9450} \, x_{1}^{2} x_{2}^{8} x_{3}^{6} + \frac{1}{525} \, x_{2}^{10} x_{3}^{6} + \frac{2}{315} \, x_{1}^{8} x_{3}^{8} + \frac{11}{9450} \, x_{1}^{6} x_{2}^{2} x_{3}^{8} - \frac{11}{4725} \, x_{1}^{4} x_{2}^{4} x_{3}^{8} + \frac{11}{9450} \, x_{1}^{2} x_{2}^{6} x_{3}^{8} + \frac{2}{315} \, x_{2}^{8} x_{3}^{8} + \frac{1}{525} \, x_{1}^{6} x_{3}^{10} - \frac{8}{1575} \, x_{1}^{4} x_{2}^{2} x_{3}^{10} - \frac{8}{1575} \, x_{1}^{2} x_{2}^{4} x_{3}^{10} + \frac{1}{525} \, x_{2}^{6} x_{3}^{10} + \frac{1}{9450} \, x_{1}^{4} x_{3}^{12} - \frac{16}{4725} \, x_{1}^{2} x_{2}^{2} x_{3}^{12} + \frac{1}{9450} \, x_{2}^{4} x_{3}^{12} + \frac{1}{1050} \, x_{1}^{2} x_{3}^{14} + \frac{1}{1050} \, x_{2}^{2} x_{3}^{14} + \frac{17}{9450} \, x_{3}^{16}

Below we verify that our previous sum of squares reproduces the desired polynomial.

res - create_difference((4,4),(5,2,1))
0