Fast algorithm for checking if binary arrays can be rotated to not have an elementwise sum over 1

Question

Let's say I have a set of constant-length arrays containing only zeroes and ones. My goal is to find out whether, after any rotation of any of the arrays, the element-wise sums of the arrays will not exceed 1.

For instance, let's say I have the following three arrays: [1, 0, 0, 0], [1, 0, 1, 0], and [1, 0, 0, 0]. I can rotate the second array by one element and the third array by two elements to obtain the arrays [1, 0, 0, 0], [0, 1, 0, 1], [0, 0, 1, 0], the element-wise sum of which is [1, 1, 1, 1]. However, had I not applied the rotations, I would have gotten a sum of [3, 0, 1, 0], which does not fit my requirements as one of the elements (the 3) is greater than 1.

Now, my question is, what is a fast way to determine if this is possible for an arbitrary number of arrays? For instance, there is no way to rotate [1, 0, 0, 0], [1, 0, 1, 0], [1, 0, 1, 0] so that the elements of the sum do not exceed 1.

Current heuristics

Obviously, if the total sum of the arrays, which are, say, length n, exceeds n, then this is trivially impossible.
The best idea for an approach I can think of so far is taking two arrays, finding a way to merge them together, and inverting the result. Then, we take this result and the next array, and repeat this process. However, this method does not guarantee to find a solution if one exists.

My question is, short of trying every possible rotation, what would be a good algorithm for this problem?

Mandorla

How big is N for your task? – Balch

Quite small—the length of the arrays is less than 100. – Mandorla

Answer 1

You could reduce this problem to an exact cover problem and use one of the known algorithms for exact cover (Knuth's Algorithm X, integer linear programming, SAT solving as sascha reminds me, possibly others). The reduction involves creating all rotations of each input array and extending them with an indicator to ensure that exactly one rotation is chosen. For the instance [1, 0, 0, 0], [1, 0, 1, 0], [1, 0, 0, 0], for example, the exact cover instance is

[1, 0, 0, 0; 1, 0, 0]  # from [1, 0, 0, 0]
[0, 1, 0, 0; 1, 0, 0]
[0, 0, 1, 0; 1, 0, 0]
[0, 0, 0, 1; 1, 0, 0]
[1, 0, 1, 0; 0, 1, 0]  # from [1, 0, 1, 0]
[0, 1, 0, 1; 0, 1, 0]
[1, 0, 0, 0; 0, 0, 1]  # from [1, 0, 0, 0]
[0, 1, 0, 0; 0, 0, 1]
[0, 0, 1, 0; 0, 0, 1]
[0, 0, 0, 1; 0, 0, 1]
[1, 0, 0, 0; 0, 0, 0]  # extra columns to solve the impedance mismatch
[0, 1, 0, 0; 0, 0, 0]  # between zeros being allowed and exact cover
[0, 0, 1, 0; 0, 0, 0]
[0, 0, 0, 1; 0, 0, 0]

I have a feeling that your problem is NP-hard, which would mean that there's a reduction in the reverse direction too and hence no hope for an algorithm whose provable worst-case running time is subexponential.

EDIT: yes, this problem is NP-hard. There's an easy reduction from 3-partition that I'll demonstrate by example.

Suppose that the 3-partition instance is [20, 23, 25, 45, 27, 40]. Then we make a binary array

[1, ..(20 ones in total).., 1, 0, ..., 0]
[1, ..(23 ones in total).., 1, 0, ..., 0]
[1, ..(25 ones in total).., 1, 0, ..., 0]
[1, ..(45 ones in total).., 1, 0, ..., 0]
[1, ..(27 ones in total).., 1, 0, ..., 0]
[1, ..(40 ones in total).., 1, 0, ..., 0].

We're looking for a partition where each of the two parts sums to 90, so the final array is a "stencil"

[1, 0, ..(90 zeros in total).., 0, 1, 0, ..(90 zeros in total).., 0]

that enforces the 3-partition constraints.

Xerxes

I worry about your reduction because it looks like the binary arrays could have length exponential in the size of the original input, since numbers in the original input turn into array lengths. – Elegant

@Elegant Fortunately 3-partition is strongly NP-hard, so unary is OK – Xerxes

Answer 2

I'm still undecided on the question if that problem is in P or NP-hard. There surely is a lot of mathematical-structure to exploit.. See David's answer.

But until someone else comes with a nice solution, here something which will work in theory and might also do in practice.

The basic idea is: formulate it as SAT-problem and use very efficient solvers for this kind of combinatorial-problems.

The SAT-solver kind we use here are CDCL-based solvers which are complete and sound (they will find a feasible solution or proof there is none!).

Analysis (naive approach)

While SAT-solving is NP-hard in general, often instances with millions of variables and constraints can be solved. But this does not guarantee that this will work here. It's hard to say without testing and sadly, no test-data was provided.

The formulation is as simple as it gets:

• N * M binary-variables:
  • N marks the data-row; M the roation/shift value
• A: preprocess all possible pairwise conflicts
• B: add constraints that at least one configuration of each row is used
• C: add constraints forbidding conflicts

For N=100 rows of M=100 cols, there will be 4950 ordered pairs, each multiplied by M*M (pairwise rotation-combinations) each. So there are <= 4950 * 100 * 100 = 49.500.000 conflict-checks (which is even feasible in slow languages). This also is an upper-bound on the number of conflict-constraints.

Of course this might change if you got very sparse-data which allows N to grow while M is fixed (in feasible-instance world).

There are probably a lot of potential reductions possible.

The take-away message here:

• Preprocessing is a lot of work (asymptotics!) and the approach is based on the comment the length of the arrays is less than 100
• SAT-solving is very very fast in regards to propagation and if P or NP-hard, the kind of constraints we provided are very powerful in terms of propagation-efficiency
• Testing this empirically (on your data) is recommended !

Remark:

There is no <= per row constraint and in some cases maybe two configurations are selected. The solution-reconstruction code does not check for this case (but there is no theoretical problem -> just pick one => will be compatible).

Code

from pyCadical import PyCadical  # own wrapper; not much tested; @github
import itertools
import numpy as np

""" DATA """
data = np.array([[1, 0, 0, 0],
                 [1, 0, 1, 0],
                 [1, 0, 0, 0]])

""" Preprocessing """
N = len(data)
M = len(data[0])

conflicts = []
for i, j in itertools.combinations(range(N), 2):
    for rotA in range(M):
        for rotB in range(M):
            if np.amax(np.roll(data[i], rotA) + np.roll(data[j], rotB)) > 1:
                conflicts.append((i, j, rotA, rotB))
conflicts = np.array(conflicts)

""" SAT """
cad = PyCadical()
vars = np.arange(N*M, dtype=int).reshape(N,M) + 1

# at least one rotation chosen per element
for i in range(N):
    cad.add_clause(vars[i, :])  # row0rot0 OR row0rot1 OR ...

# forbid conflicts
for i, j, rotA, rotB in conflicts:
    cad.add_clause([-vars[i, rotA], -vars[j, rotB]])  # (not rowIrotA) or (not rowJrotB)

""" Solve """
cad.solve()

""" Read out solution """
sol = cad.get_sol_np().reshape(N, M)
chosen = np.where(sol > 0)

solution = []  # could be implemented vectorized
for i in range(N):
    solution.append(np.roll(data[i], chosen[1][i]))

print(np.array(solution))

Output

[[0 1 0 0]
 [1 0 1 0]
 [0 0 0 1]]

Balch

Answer 3

I will view each collection of bits as (sufficiently large) Integer.

Say I have such a collection of integers on n bits. Here is some Squeak Smalltalk code that show how to reduce a little bit the combinatorial:

SequenceableCollection>>canPreventsOverlapingBitByRotatingOver: n
    "Answer whether we can rotate my elements on n bits, such as to obtain non overlaping bits"
    | largestFirst nonNul nonSingletons smallest |

    "Exclude trivial case when there are more than n bits to dispatch in n holes"
    (self detectSum: #bitCount) > n ifTrue: [^false].

    "Exclude non interesting case of zero bits"
    nonNul := self reject: [:each | each = 0].

    "Among all possible rotations, keep the smallest"
    smallest := nonNul collect: [:each | each smallestAmongBitRotation: n].

    "Note that they all have least significant bit set to 1"
    [smallest allSatisfy: [:each | (each bitAnd: 1) = 1]] assert.

    "Bit singletons can occupy any hole, skip them"
    nonSingletons := smallest reject: [:each | each = 1].

    "Sort those with largest bitCount first, so as to accelerate detection of overlaping"
    largestFirst := nonSingletons sorted: #bitCount descending.

    "Now try rotations: all the shift must differ, otherwise the shifted LSB would overlap"
    ^largestFirst checkOverlapingBitRotated: n

Where we have the following utilities defined as:

SequenceableCollection>>checkOverlapingBitRotated: n
    "Answer true if the bits of my elements can be rotated on n bits so as to not overlap"
    ^self checkOverlapingBitRotatedBy: (1 << n - 1) among: n startingAt: 2 accum: self first

SequenceableCollection>>checkOverlapingBitRotatedBy: shiftMask among: n startingAt: index accum: accum
    index > self size ifTrue: [^true].
    (shiftMask bitClear: accum) bitsDo: [:bit |
        | shifted |
        shifted := (self at: index) bitRotate: bit lowBit - 1 among: n.
        ((accum bitAnd: shifted) = 0
            and: [self
                    checkOverlapingBitRotatedBy: shiftMask
                    among: n
                    startingAt: index + 1
                    accum: (accum bitOr: shifted)])
            ifTrue: [^true]].
    ^ false

This requires an additional explanation: every bit in the shiftMask indicates (the rank of) a possible shift. Since the accumulator already occupy a number of bits, and since the LSB of each element is 1, we cannot shift the remaining element to the bits already occupied by the accumulator. Thus we have to clear the bits occupied by the accumulator from the mask. This considerably reduce the combinations, and that's why it is beneficial to sort largest bitCount first.

Secondly, the guard (accum bitAnd: shifted) = 0 is cutting the recursion as soon as we can rather than generating useless combinations and testing infeasibility a posteriori.

We then have those small bit utilities:

Integer>>bitRotate: shift among: n
    "Rotate the n lowest bits of self, by shift places"
    "Rotate left if shift is positive."
    "Bits of rank higher than n are erased."
    | highMask lowMask r |
    (r := shift \\ n) = 0 ifTrue: [^self].
    lowMask := 1 << (n - r) - 1.
    highMask := 1 << n - 1 - lowMask.
    ^((self bitAnd: lowMask) << r)
        bitOr: ((self bitAnd: highMask) >> (n - r))

Integer>>smallestAmongBitRotation: n
    "Answer the smallest rotation of self on n bits"
    ^self
        bitRotate: ((1 to: n) detectMin: [:k | self bitRotate: k among: n])
        among: n

Integer>>bitsDo: aBlock
    "Evaluate aBlock will each individual non nul bit of self"
    | bits lowBit |
    bits := self.
    [bits = 0] whileFalse: [
        lowBit := (bits bitAnd: 0 - bits).
        aBlock value: lowBit.
        bits := bits - lowBit].

It works instantly on small collections like this:

| collec bitCount |
collec := #( 2r11 2r1001  2r1101 2r11011 2r1110111 2r11001101
       2r11010010111010 2r1011101110101011100011).
bitCount := collec detectSum: #bitCount.
(bitCount to: bitCount*2) detect:
    [:n | collec canPreventsOverlapingBitByRotatingOver: n].

would answer 52, that means that we need at least 52 bits to obtain a non overlaping combination, though the bitCount is only 44.

All is performed with simple bit operations and should scale well (once translated in a static language).

That's not the case of my 32 bits interpreter which creates boxed large integers, pressures the garbage collector, and takes a bit of time with 10 sets for a total of about 100 bits:

| collec bitCount |
collec := (1 to: 10) collect: [:i | (1 << 18 - 1) atRandom].
bitCount := collec detectSum: #bitCount.
bitCount ->
    [ collec canPreventsOverlapingBitByRotatingOver: bitCount + 10] timeToRun.

First try took 75 seconds for bitCount=88

More fair (sparse) bit distributions lead to a faster average time (and still horrible worst time):

| collec bitCount |
collec := (1 to: 15) collect: [:i |
    ((1 to: 4) collect: [:j | (1 to: 1<<100-1) atRandom])
        reduce: #bitAnd:].
bitCount := collec detectSum: #bitCount.
bitCount ->
    [ collec canPreventsOverlapingBitByRotatingOver: (bitCount + 10 max: 100)] timeToRun.

104->1083ms
88->24ms    
88->170ms
91->3294ms
103->31083ms

Sexy