# Writing a matrix as a sum of matrices

Given a matrix

M = matrix(ZZ, [[1, 1, 1, 1], [1, -1, 1, -1], [1, 1, -1, -1], [1, -1, -1, 1]])


write

M = A1 + A2 + A3 + A4 + A5 + A6


where each Ai is either a 0,1 or 0,-1 matrix. What are all such possible combinations of Ai's?

Ai being 0,1 or 0,-1 matrix means that 1 and -1 can't be in the same Ai.

edit retag close merge delete

Sort by » oldest newest most voted def sum_list(length, total):
"""
INPUT:

- length -- how many terms
- total -- desired total

Return list of lists, each of which has length entries chosen
from {0, 1, -1} and adds up to total.
"""
C = IntegerListsLex(length=length, n=length+total, min_part=0, max_part=2)
return [[a-1 for a in L] for L in C]


(It's only two lines long; everything between the """ lines is a comment.)

Sample use to get all lists of length 4 that add up to 1:

sage: sum_list(4, 1)
[[1, 1, 0, -1],
[1, 1, -1, 0],
[1, 0, 1, -1],
[1, 0, 0, 0],
[1, 0, -1, 1],
[1, -1, 1, 0],
[1, -1, 0, 1],
[0, 1, 1, -1],
[0, 1, 0, 0],
[0, 1, -1, 1],
[0, 0, 1, 0],
[0, 0, 0, 1],
[0, -1, 1, 1],
[-1, 1, 1, 0],
[-1, 1, 0, 1],
[-1, 0, 1, 1]]


Now try to build up matrices using sum_list(LENGTH, 1) and sum_list(LENGTH, -1) based on the entries of your given matrix.

The following is a draft, and it could probably be sped up, but it seems to work. It treats a solution (A,B) as different from (B,A): the order matters.

def sum_list(length, total):
"""
INPUT:

- length -- how many terms
- total -- desired total

Return list of lists, each of which has length entries chosen
from {0, 1, -1} and adds up to total.
"""
C = IntegerListsLex(length=length, n=length+total, min_part=0, max_part=2)
return [[a-1 for a in L] for L in C]

def mat_list(mat, length):
"""
INPUT:

- mat -- matrix, assumed to have entries in {1, -1}
- length -- how many terms

Return list of lists, each of which has length entries and
adds up to mat, and each entry is either a (0, 1) or a (0,
-1) matrix.
"""
sum_plus = sum_list(length, 1)
sum_minus = sum_list(length, -1)
old_attempts = []
for x in mat.list():
if x == 1:
sums = sum_plus
elif x == -1:
sums = sum_minus
else:
raise ValueError('each entry of matrix should be 1 or -1')
if not old_attempts: # initial step
new_attempts = [[[s] for s in S] for S in sums]
else:
new_attempts = []
for mats in old_attempts:
NEW = []
for S in sums:
new_M = [M + [s] for (M, s) in zip(mats, S)]
# keep only entries that correspond to (0,1) or (0,-1) matrices
if all(max(entries) - min(entries) < 2 for entries in new_M):
NEW.append([M + [s] for (M, s) in zip(mats, S)])
new_attempts.extend(NEW)
old_attempts = new_attempts

matrices = [[matrix(mat.base_ring(), mat.nrows(), mat.ncols(), entries) for entries in L] for L in old_attempts]
return matrices


With this code:

sage: mat = matrix(ZZ, [[1, 1, 1, 1], [1, -1, 1, -1], [1, 1, -1, -1], [1, -1, -1, 1]])
sage: %time L = mat_list(mat, 3)
CPU times: user 8.04 s, sys: 141 ms, total: 8.18 s
Wall time: 8.18 s
sage: len(L)
179328


EDIT: here is a modification, in case you don't want to distinguish the same lists with a different ordering:

def mat_list_NEW(mat, length):
"""
INPUT:

- mat -- matrix, assumed to have entries in {1, -1}
- length -- how many terms

Return list of lists, each of which has length entries and
adds up to mat, and each entry is either a (0, 1) or a (0,
-1) matrix.
"""
sum_plus = sum_list(length, 1)
sum_minus = sum_list(length, -1)
old_attempts = []
for x in mat.list():
if x == 1:
sums = sum_plus
elif x == -1:
sums = sum_minus
else:
raise ValueError('each entry of matrix should be 1 or -1')
if not old_attempts: # initial step
new_attempts = set(tuple([(s,) for s in S]) for S in sums)
else:
new_attempts = set()
for mats in old_attempts:
NEW = set()
for S in sums:
new_M = tuple(sorted(M + (s,) for (M, s) in zip(mats, S)))
# keep only entries that correspond to (0,1) or (0,-1) matrices
if all(max(entries) - min(entries) < 2 for entries in new_M):
new_attempts.update(NEW)
old_attempts = new_attempts

matrices = [[matrix(mat.base_ring(), mat.nrows(), mat.ncols(), entries) for entries in L] for L in old_attempts]
return matrices


The differences: use the set structure so that duplicates are automatically removed, and sort the lists before adding them to the set, so that it's easy to recognize duplicates. Elements of a set have to be hashable so use tuples instead of lists. With this new version:

sage: mat = matrix(ZZ, [[1, 1, 1, 1], [1, -1, 1, -1], [1, 1, -1, -1], [1, -1, -1, 1]])
sage: %time L = mat_list_NEW(mat, 3)
CPU times: user 1.47 s, sys: 41 ms, total: 1.51 s
Wall time: 1.53 s
sage:  len(L)
29889


But the lengths of the lists are still going to get large as the number of summands increases...

more

Note that these lists get long pretty quickly: for a 3x3 matrix M3, len(mat_list(M3, 4)) is 2666412.

Sir John, your code is really very fast as compared to mine. Thanks for you time and effort. But we require unique summands in a combination. Also, all-zero and all non-zero summands are also not required for our purpose. See, if these can be done. In my earlier code(which I have to discard , as it is very slow), later I specified the number of non-negative and non-positive summands to reduce the search space( I have mentioned in the code below my earler code). See, if similar or any other thing can be done in your code to futher reduce the output time. (In fact I had also broken into symmetric and non-symmetric summands, but thought breaking into non-negative and non-positive was better as we know for example that breaking into 6 we would have 2,4 or 3,3 or 4,2 summands ).Thanks.

To prevent repetitions, you can change the main if statement to if (len(set(new_M)) == len(new_M) and all(max(entries) - min(entries) < 2 for entries in new_M)). To avoid matrices with entries all the same, you can change the definition of matrices to matrices = [[matrix(mat.base_ring(), mat.nrows(), mat.ncols(), entries) for entries in L if max(entries) > min(entries)] for L in old_attempts]. With the 4x4 case, it produces just under 20 thousand lists, instead of just under 30 thousand, and it takes about 2/3 the amount of time.

Sir, thanks again for your time and effort.

(1)Sir, as per your second code, if for example, we take the 2X2 matrix [[1,1],[1,-1]] with length 3, then we get 15 combinations, all with distinct matrices. But after changing the "if "statement we get only 9. So valid combintions are also removed.

(2) As per the new definition of "matrices", we are still getting the combinations with only the all-same-entry matrices being removed , whereas we need complete removal of such combinations.

I tried to figure out the things on my own, but as I am fairly new to sagemath, I could not. So, kindly look into it. Thanks.

I'm sorry, I have put in a lot of time on this, and I can't do any more. You can consider hiring someone to write your code for you.

Well... This problem (decomposing a matrix as the sum of $p$ matrices with your constraints) can be turned in the search of the ideal of some field defined by a polynomial system :

• each element of M will be translated as a sum of $p+1$ terms (the corresponding elements of the $p$ matrices minus the corresponding element of M) ;

• each $(0, 1)$ $r\times c$ matrix (generic term $m$) can be created by nullifying $rc$ $m(m-1)$ polynomials ;

• each $(0, -1)$ $r\times c$ matrix (generic term $m$) can be created by multiplying a $(0, 1)$ matrix by -1 ;

• therefore, all the matrices can be given the same form, by multiplying each of them by $s_k,\ k\in1\dots p$ with the additional constraints of nullifying $s_k^2-1$.

This can be implemented as follows :

def solve_system_for(M, Nmat=6):
"""
Decompose M as a sum of Nmat matrices of the same shape, each of them
being composed either of 1 and 0 OR -1 and 0
"""
Nr = M.nrows()
Nc = M.ncols()
Nel = Nr*Nc
# Names for signs of matrices
sn = ["s%d"%u for u in (1..Nmat)]
# Names for matrix elements
mn = [[["m%d%d%d"%(t,u,v)for v in (1..Nc)] for u in (1..Nr)]
for t in (1..Nmat)]
# Polynomial ring : must be a field containing -1, 0, 1,
# with arithmetic compatible with ZZ's ==> QQ
R1 = PolynomialRing(QQ, sn+flatten(mn))
R1.inject_variables()
# Vector of inteterminates related to matrice's signs
s = R1.gens()[0:Nmat]
# vector of indeterminates for matrices elements
a = R1.gens()[Nmat:]
# List of matrices of indeterminates
A = [s[u-1]*matrix(R1.gens()[Nmat+u*Nel:Nmat+(u+1)*Nel],
nrows=Nr, ncols=Nc)
for u in range(Nmat)]
# Matrix of equations to solve
Eq = sum(A)-M.change_ring(QQ)
# Constraints on s elements : s ∈ {-1, 1} <==> s in roots of s^2-1 in QQ
Cs = [u^2-1 for u in s]
# Constraints on a elements : a ∈ {O, 1} <==> a in roots of a*(a-1) in QQ
Ca = [u*(u-1) for u in a]
# Solutions ideal described by :
J1 = R1.ideal(flatten(Eq.list()) + Cs + Ca)
# The hard, long part :
from time import time as stime
t0 = stime()
D = J1.dimension()
t1 = stime()
if D==0: # Expected case on the basis of finite set of possible values
t2 = stime()
S = J1.variety()
t3 = stime()
return {"Time for dimension":t1-t0,
"Time for solution":t3-t2,
"Dimension":D,
"Solution(s)":S}
if D<0: # No solution.
return {"Time for dimension":t1-t0,
"Time for solution":t3-t2,
"Dimension":D,
"Solution(s)":[]}
if D>0:
# Infinite solutions. This shouldn't happen,
# but I learned the virtues of defensive programming...
return {"Time for dimension":t1-t0,
"Time for solution":t3-t2,
"Dimension":D,
"Parent ring":R1,
"Ideal to reprocess":J1}


But I have some serious doubts regarding the usefulness of this implementation :

• Decomposing your $M$ as the sum of two matrices involves 34 indeterminates and 50 polynomials ; determining the variety of its ideal takes 4.17 seconds (0,3 seconds to check its dimension, 4 s for finding the two elements of its variety).

• Decomposing it in the sum of 3 matrices involves 51 indeterminates and 67 polynomials, and hasn't finished in about 11 hours.

Combinatorial explosion is fast...

Unless my reasoning or my programming is|are faulty (which is of course possible), I do not realistically contemplate the use of this solution to decompose M in a sum of 6 matrices...

Hoping this non-help helps finding a better solution...

more

I asked this thing because a friend of mine made code for this but it is very time taking. For 2nd order matrices output comes out instantly. For 3rd order , output starts coming but it is taking a lot of time. For 4th order , no output is shown even after hours. Since for 2nd order, number of elements in V.List is 2^4, for 3rd order it is 2^9 and for 4th order it is 2^16 so it is time taking. Can any change be made to make this code run faster?(Please note that after writing the code it appeared differently in preview, I tried changing but could not. So indentation is required where B and C are mentioned.)

#Write your matrix
H = matrix(2,2, [1, 1,1, -1])

show("H=",H)

n = int(input("Enter the order of the matrix H : "))

V = GF(2)^(n^2)

R1=V.list()[1:]
#We are removing the all-zero list because not needed for our purpose

R2=R1[:-1]
#We are removing the all-ones list because not needed for our purpose

A = []
for v in  R2:
A.append(matrix(ZZ,n,n,v))
# Forming n by n matrices from  list in  R2. Contains all 0,1 matrices.

B=[]
for a in A:
b=-a                                           # Forming n by n matrices from A by taking their negative. Contains all 0,-1 matrices
B.append(b)

C = []
for a in A:
C.append(a)
for b in B:
C.append(b)         #Union of A and B

test_sum = int(input("Into how many matrices you want to break H :"))

test_comb = Combinations(C,test_sum)
# Here we are forming all "test_sum"-set combinations  from C".
# That is, if test_sum =3, then test_comb, gives all 3 set combinations  from C"

for t in test_comb:
if (sum(t)==H):
#Checking which combinations in  test_sum give the  sum = H and printing those

print("A combination is: ")
show(t)
print("\n")


## To reduce search space I had specified the number of non-negative and non-positive summands as shown below:

 #Write your matrix
H = matrix(2,2, [1, 1,1, -1])

show("H=",H)

n = int(input("Enter the order of the matrix H : "))

V = GF(2)^(n^2)
R1=V.list()[1:]
#We are removing the all-zero list because not needed for our purpose

R2=R1[:-1]
#We are removing the all-ones list because not needed for our purpose

A = []
for v in  R2:
A.append(matrix(ZZ,n,n,v))
# Forming n by n matrices from  list in  R2. Contains all 0,1 matrices.

B=[]
for a in A:
b=-a
B.append(b)
# Forming n by n matrices from A by taking their negative. Contains all 0,-1 matrices

test_sum = int(input("Into how many matrices you want to break H :"))

num_nonneg = int(input("How many non-negative matrices u want :"))
num_nonpos = (test_sum)-(num_nonneg)

nonneg_comb = Combinations(A,num_nonneg)
nonpos_comb = Combinations(B,num_nonpos)

for x in nonneg_comb:
for y in nonpos_comb:
if (sum(x+y)==H):

print("A combination is: ")
show(x+y)
print("\n")
print("No Combination exists")

more

I took the liberty to edit your "answer" in order to get somethng barely legible. I hope that I haven't introduced errors doing so (In Python, whitespane is syntactic...) : please check my edits.

If you have 4x4 matrices and 5 summands, then you will be checking 10^22 possible combinations in this code. No surprise that it's slow.