Given a list of numbers, find all matrices such that each column and row sum up to 264

Question

Let's say I have a list of 16 numbers. With these 16 numbers I can create different 4x4 matrices. I'd like to find all 4x4 matrices where each element in the list is used once, and where the sum of each row and each colum is equal to 264.

First I find all combination of elements within the list that sum up to 264

numbers = [11, 16, 18, 19, 61, 66, 68, 69, 81, 86, 88, 89, 91, 96, 98, 99]

candidates = []
result = [x for x in itertools.combinations(numbers, 4) if sum(x) == 264]

result becomes a list where each element is a list with 4 elements, where the sum of the 4 elements = 264. I think of these as my rows. Then I'd like to take all permutations of my rows, since addition is commutative.

for i in range(0, len(result)):
    candidates.append(list(itertools.permutations(result[i])))

Now given all my possible rows where the sum is 264. I'd like to choose all combinations of 4 rows, such that every column's sum is 264.

test = []
for i in range(0, len(candidates)):
    test = test + candidates[i]
result2 = [x for x in itertools.combinations(test, 4) if list(map(add, x[0], list(map(add, x[1], list( map(add, x[2], x[3])))))) == [264, 264, 264, 264]]

Is there a faster/better way? The last part, finding all combinations of 4 rows, takes a lot of time and computer power.

Leporine

Isn't result which you correctly say is all the rows, also all the columns? Pretty sure there's something you can do with that, but I can't quite see it right now. – Mixture

@Mixture Yes in a sense. But if a+b+c+d=264, then b+a+c+d=264 as well and thus I need all permutations. So actually candidates holds all rows and columns with all its permutations. So is there a clever way to choose elementa from candidates... cant see it either – Leporine

For anyone else wondering, there's 30 combinations that add up to 264, and 30^4 is 810,000 permutations to check. – Postdiluvian

Answer 1

This is a kind of constraint satisfaction problem; there are sixteen variables each with the same domain, eight constraints about their sums, and one constraint that they should all have different values from the domain.

There are potentially a large number of solutions, so any algorithm which generates a bigger set of candidates and then checks which candidates really are solutions is probably inefficient by a large factor, since the true solutions are likely to be a very low proportion of your candidates. A backtracking search is generally better, since it allows partial candidates to be rejected when they violate any constraint, potentially eliminating many complete candidates without having to generate them all in the first place.

Rather than write your own backtracking search algorithm, you can use an existing constraint-solver such as the python-constraint library. Here's an example:

numbers = [11, 16, 18, 19, 61, 66, 68, 69, 81, 86, 88, 89, 91, 96, 98, 99]
target = 264

from constraint import *

problem = Problem()
problem.addVariables(range(16), numbers)

for i in range(4):
    # column i
    v = [ i + 4*j for j in range(4) ]
    problem.addConstraint(ExactSumConstraint(target), v)
    # row i
    v = [ 4*i + j for j in range(4) ]
    problem.addConstraint(ExactSumConstraint(target), v)

problem.addConstraint(AllDifferentConstraint())

Example:

>>> problem.getSolution()
{0: 99, 1: 88, 2: 66, 3: 11, 4: 16, 5: 61, 6: 89, 7: 98, 8: 81, 9: 96, 10: 18, 11: 69, 12: 68, 13: 19, 14: 91, 15: 86}
>>> import itertools
>>> for s in itertools.islice(problem.getSolutionIter(), 10):
...     print(s)
... 
{0: 99, 1: 68, 2: 81, 3: 16, 4: 66, 5: 91, 6: 18, 7: 89, 8: 88, 9: 19, 10: 96, 11: 61, 12: 11, 13: 86, 14: 69, 15: 98}
{0: 99, 1: 68, 2: 81, 3: 16, 4: 66, 5: 91, 6: 18, 7: 89, 8: 11, 9: 86, 10: 69, 11: 98, 12: 88, 13: 19, 14: 96, 15: 61}
{0: 99, 1: 68, 2: 81, 3: 16, 4: 18, 5: 89, 6: 66, 7: 91, 8: 86, 9: 11, 10: 98, 11: 69, 12: 61, 13: 96, 14: 19, 15: 88}
{0: 99, 1: 68, 2: 81, 3: 16, 4: 18, 5: 89, 6: 66, 7: 91, 8: 61, 9: 96, 10: 19, 11: 88, 12: 86, 13: 11, 14: 98, 15: 69}
{0: 99, 1: 68, 2: 81, 3: 16, 4: 11, 5: 86, 6: 69, 7: 98, 8: 66, 9: 91, 10: 18, 11: 89, 12: 88, 13: 19, 14: 96, 15: 61}
{0: 99, 1: 68, 2: 81, 3: 16, 4: 11, 5: 86, 6: 69, 7: 98, 8: 88, 9: 19, 10: 96, 11: 61, 12: 66, 13: 91, 14: 18, 15: 89}
{0: 99, 1: 68, 2: 81, 3: 16, 4: 61, 5: 96, 6: 19, 7: 88, 8: 18, 9: 89, 10: 66, 11: 91, 12: 86, 13: 11, 14: 98, 15: 69}
{0: 99, 1: 68, 2: 81, 3: 16, 4: 61, 5: 96, 6: 19, 7: 88, 8: 86, 9: 11, 10: 98, 11: 69, 12: 18, 13: 89, 14: 66, 15: 91}
{0: 99, 1: 68, 2: 81, 3: 16, 4: 88, 5: 19, 6: 96, 7: 61, 8: 11, 9: 86, 10: 69, 11: 98, 12: 66, 13: 91, 14: 18, 15: 89}
{0: 99, 1: 68, 2: 81, 3: 16, 4: 88, 5: 19, 6: 96, 7: 61, 8: 66, 9: 91, 10: 18, 11: 89, 12: 11, 13: 86, 14: 69, 15: 98}

That's the first ten solutions. The problem.getSolutions() method returns a list containing all of them, but this takes quite a bit of time to run (about 2 minutes on my machine) because there are 6,912 of them to find.

One issue is that each solution has many symmetrical counterparts; you can permute the rows, and permute the columns, and take the transpose. It is possible to eliminate symmetries by adding more constraints, so that you just get one solution from each symmetry class. This makes the search more feasible:

# permute rows/cols so that lowest element is in top-left corner
m = min(numbers)
problem.addConstraint(InSetConstraint([m]), [0])

from operator import lt as less_than

for i in range(3):
    # permute columns so first row is in order
    problem.addConstraint(less_than, [i, i+1])
    # permute rows so first column is in order
    problem.addConstraint(less_than, [4*i, 4*i + 4])

# break transpose symmetry by requiring grid[0,1] < grid[1,0]
problem.addConstraint(less_than, [1, 4])

This breaks all symmetries, so now it returns 6,912 / (4! * 4! * 2) = 6 solutions in about 0.2 seconds.

Shaky

@StefanPochmann These 7 constraints disallow permutation of rows, permutation of columns and the mirror symmetry. Do there exist other symmetries in this problem? – Troop

@StefanPochmann The symmetry group has size 4! * 4! * 2, and the symmetry constraints reduce the number of solutions by that factor exactly, so it can't be outputting more than one solution from any symmetry class. Or do you mean a more general proof? Your first solution doesn't obey the "minimum number is in top-left corner" constraint. – Shaky

@StefanPochmann So, you mean that also the 0 element should be fixed (to 11) to avoid that it ends up somewhere else in the grid? – Troop

@Troop Breaking row and column symmetries separately doesn't mean that all symmetries are necessarily broken; just that the remaining symmetry group doesn't include any symmetries that only act on rows, or only act on columns. There can still be symmetries which act on both rows and columns, which the "min in top left" constraint is needed to eliminate. I think for problems where all the numbers are distinct, these constraints should break all symmetries, but I'm not 100% sure of that off the top of my head. Breaking symmetries in CSPs is surprisingly difficult. – Shaky

Ok I believe it now. Given any solution, permuting the smallest number to top-left also determines what else is in the first row, so your first-row constraints determines the columns permutation. Same with rows permutation. – Instrumentation

How do you know that there is 6912 number of solutions? – Leporine

@Leporine because I ran the solver for the required two minutes and it found 6,912. – Shaky

Alright! I believed that you had some non-practical argument. And I couldnt find such an argument, hence the question. Thanks! – Leporine

Answer 2

Here is an approach using z3py, Python's version of the Z3 SAT/SMT solver. Note that every permutation of rows and/or columns as well as mirroring gives an additional solution. Together, each primitive solution leads to 24*24*2 equivalent solutions.

Adding constraints to force an order, should enable to find all primitive solutions. If there are no mistakes, the following program finds all 6 of them. So, all together there should be 6*24*24*2 = 6912 solutions.

from z3 import Solver, BitVec, Or, Distinct, sat

numbers = [11, 16, 18, 19, 61, 66, 68, 69, 81, 86, 88, 89, 91, 96, 98, 99]

# X is a table to store the 16 variables for the solution
X = [BitVec(f'x{i}{j}', 16) for i in range(4) for j in range(4)]
s = Solver()
for x in X:
    s.add(Or([x == n for n in numbers]))  # all X[i] should be one of the given numbers

# constraints to avoid reordered solutions
s.add(X[0] == 11)
s.add(X[0] < X[1])
s.add(X[1] < X[2])
s.add(X[2] < X[3])
s.add(X[1] < X[4])
s.add(X[4] < X[8])
s.add(X[8] < X[12])

# all X[i] have to be distinct
s.add(Distinct(X))
for i in range(4):
    # all rows and all columns need to sum to 264
    s.add(sum([X[4*i+j] for j in range(4)]) == 264)
    s.add(sum([X[4*j+i] for j in range(4)]) == 264)

# start solving
res = s.check()

while res == sat:
    m = s.model()
    # show the solution
    for i in range(4):
        print([m[X[i*4+j]] for j in range(4)])
    print()

    # add the just found solution as a constraint so it doesn't get outputted again
    s.add(Or([X[i] != m[X[i]].as_long() for i in range(16)]))

    # solve again to find different solutions
    res = s.check()

Output:

[11, 68, 89, 96]
[69, 16, 91, 88]
[86, 99, 18, 61]
[98, 81, 66, 19]

[11, 68, 86, 99]
[69, 16, 98, 81]
[88, 91, 19, 66]
[96, 89, 61, 18]

[11, 66, 89, 98]
[69, 18, 91, 86]
[88, 99, 16, 61]
[96, 81, 68, 19]

[11, 66, 88, 99]
[68, 19, 91, 86]
[89, 98, 16, 61]
[96, 81, 69, 18]

[11, 66, 88, 99]
[69, 18, 96, 81]
[86, 91, 19, 68]
[98, 89, 61, 16]

[11, 66, 89, 98]
[68, 19, 96, 81]
[86, 91, 18, 69]
[99, 88, 61, 16]

Troop

It takes about 1.8 seconds – Troop

To confirm, I get six solutions too. – Shaky