Find the k non-repeating elements in a list with "little" additional space

Question

The original problem statement is this one:

Given an array of 32bit unsigned integers in which every number appears exactly twice except three of them (which appear exactly once), find those three numbers in O(n) time using O(1) extra space. The input array is read-only. What if there are k exceptions instead of 3?

It's easy to solve this in Ο(1) time and Ο(1) space if you accept a very high constant factor because of the input restriction (the array can have at most 2³³ entries):

for i in lst:
    if sum(1 for j in lst if i == j) == 1:
        print i

So, for the sake of this question, let's drop the restriction in bit length and concentrate on the more general problem where the numbers can have up to m bits.

Generalizing an algorithm for k = 2, what I had in mind is the following:

1. XOR those numbers with a least significant bit of 1 and those with a 0 separately. If for both of the partitions, the resulting value is not zero, we know that we have partitioned the non-repeating numbers into two groups, each of which has at least one member
2. For each of those groups, try to partition it further by examining the second-least significant bit and so on

There is a special case to be considered, though. If after partitioning a group, the XOR values of one of the groups are both zero, we don't know whether one of the resulting sub-groups is empty or not. In this case my algorithm just leaves this bit out and continues with the next one, which is incorrect, for example it fails for the input [0,1,2,3,4,5,6].

Now the idea I had was to compute not only the XOR of the element, but also the XOR of the values after applying a certain function (I had chosen f(x) = 3x + 1 here). See Evgeny's answer below for a counter-example for this additional check.

Now although the below algorithm is not correct for k >= 7, I still include the implementation here to give you an idea:

def xor(seq):
  return reduce(lambda x, y: x ^ y, seq, 0)

def compute_xors(ary, mask, bits):
  a = xor(i for i in ary if i & mask == bits)
  b = xor(i * 3 + 1 for i in ary if i & mask == bits)
  return a if max(a, b) > 0 else None

def solve(ary, high = 0, mask = 0, bits = 0, old_xor = 0):
  for h in xrange(high, 32):
    hibit = 1 << h
    m = mask | hibit
    # partition the array into two groups
    x = compute_xors(ary, m, bits | hibit)
    y = compute_xors(ary, m, bits)
    if x is None or y is None:
      # at this point, we can't be sure if both groups are non-empty,
      # so we check the next bit
      continue
    mask |= hibit
    # we recurse if we are absolutely sure that we can find at least one
    # new value in both branches. This means that the number of recursions
    # is linear in k, rather then exponential.
    solve(ary, h + 1, mask, bits | hibit, x)
    solve(ary, h + 1, mask, bits, y)
    break
  else:
    # we couldn't find a partitioning bit, so we output (but 
    # this might be incorrect, see above!)
    print old_xor

# expects input of the form "10 1 1 2 3 4 2 5 6 7 10"
ary = map(int, raw_input().split())
solve(ary, old_xor=xor(ary))

From my analysis, this code has a worst-case time complexity of O(k * m² * n) where n is the number of input elements (XORing is O(m) and at most k partitioning operations can be successful) and space complexity O(m²) (because m is the maximum recursion depth and the temporary numbers can be of length m).

The question is of course if there is a correct, efficient approach with good asymptotic runtime (let's assume that k << n and m << n here for the sake of completeness), which also needs little additional space (for example, approaches that sort the input will not be accepted, because we'd need at least O(n) additional space for that, as we can't modify the input!).

EDIT: Now that the algorithm above is proven to be incorrect, it would of course be nice to see how it could be made correct, possibly by making it a bit less effient. Space complexity should be in o(n*m) (that is, sublinear in the total number of input bits). It would be okay to take k as an additional input if that makes the task easier.

Answer 1

I went offline and proved the original algorithm subject to the conjecture that the XOR tricks worked. As it happens, the XOR tricks don't work, but the following argument may still interest some people. (I re-did it in Haskell because I find proofs much easier when I have recursive functions instead of loops and I can use data structures. But for the Pythonistas in the audience I tried to use list comprehensions wherever possible.)

Compilable code at http://pastebin.com/BHCKGVaV.

Beautiful theory slain by an ugly fact

Problem: we're given a sequence of n nonzero 32-bit words in which every element is either singleton or doubleton:

• If a word appears exactly once, it is singleton.

• If a word appears exactly twice, it is doubleton.

• No word appears three or more times.

The problem is to find the singletons. If there are three singletons, we should use linear time and constant space. More generally, if there are k singletons, we should use O(k*n) time and O(k) space. The algorithm rests on an unproven conjecture about exclusive or.

We begin with these basics:

module Singleton where
import Data.Bits
import Data.List
import Data.Word
import Test.QuickCheck hiding ((.&.))

Key abstraction: Partial specification of a word

To tackle the problem I'm going to introduce an abstraction: to describe the least significant $w$ bits of a 32-bit word, I introduce a Spec:

data Spec = Spec { w :: Int, bits :: Word32 }
   deriving Show
width = w -- width of a Spec

A Spec matches a word if the least significant w bits are equal to bits. If w is zero, by definition all words match:

matches :: Spec -> Word32 -> Bool
matches spec word = width spec == 0 ||
                    ((word `shiftL` n) `shiftR` n) == bits spec
  where n = 32 - width spec

universalSpec = Spec { w = 0, bits = 0 }

Here are some claims about Specs:

• All words match the universalSpec, which has width 0

• If matches spec word and width spec == 32, then word == bits spec

Key idea: "extend" a partial specification

Here is the key idea of the algorithm: we can extend a Spec by adding another bit to the specification. Extending a Spec produces a list of two Specs

extend :: Spec -> [Spec]
extend spec = [ Spec { w = w', bits = bits spec .|. (bit `shiftL` width spec) }
              | bit <- [0, 1] ]
  where w' = width spec + 1

And here's the crucial claim: if spec matches word and if width spec is less than 32, then exactly one of the two specs from extend spec match word. The proof is by case analysis on the relevant bit of word. This claim is so important that I'm going to call it Lemma One Here's a test:

lemmaOne :: Spec -> Word32 -> Property
lemmaOne spec word =
  width spec < 32 && (spec `matches` word) ==> 
      isSingletonList [s | s <- extend spec, s `matches` word]

isSingletonList :: [a] -> Bool
isSingletonList [a] = True
isSingletonList _   = False

We're going to define a function which given a Spec and a sequence of 32-bit words, returns a list of the singleton words that match the spec. The function will take time proportional to the length of the input times the size of the answer times 32, and extra space proportional to the size of the answer times 32. Before we tackle the main functio, we define some constant-space XOR functions.

XOR ideas that are broken

Function xorWith f ws applies function f to every word in ws and returns the exclusive or of the result.

xorWith :: (Word32 -> Word32) -> [Word32] -> Word32
xorWith f ws = reduce xor 0 [f w | w <- ws]
  where reduce = foldl'

Thanks to stream fusion (see ICFP 2007), function xorWith takes constant space.

A list of nonzero words has a singleton if and only if either the exclusive or is nonzero, or if the exclusive or of 3 * w + 1 is nonzero. (The "if" direction is trivial. The "only if" direction is a conjecture that Evgeny Kluev has disproven; for a counterexample, see array testb below. I can make Evgeny's example work by adding a third function g, but obviously this situation calls for a proof, and I don't have one.)

hasSingleton :: [Word32] -> Bool
hasSingleton ws = xorWith id ws /= 0 || xorWith f ws /= 0 || xorWith g ws /= 0
  where f w = 3 * w + 1
        g w = 31 * w + 17

Efficient search for singletons

Our main function returns a list of all the singletons matching a spec.

singletonsMatching :: Spec -> [Word32] -> [Word32]
singletonsMatching spec words =
 if hasSingleton [w | w <- words, spec `matches` w] then
   if width spec == 32 then
     [bits spec]       
   else
     concat [singletonsMatching spec' words | spec' <- extend spec]
 else
   []

We'll prove its correctness by induction on the width of the spec.

• The base case is that spec has width 32. In this case, the list comprehension will give the list of words that are exactly equal to bits spec. Function hasSingleton will return True if and only if this list has exactly one element, which will be true exactly when bits spec is singleton in words.

• Now let's prove that if singletonsMatching is correct for with m+1, it is also correct for width m, where *m < 32$. (This is the opposite direction as usual for induction, but it doesn't matter.)

  Here is the part that is broken: for narrower widths, hasSingleton may return False even when given an array of singletons. This is tragic.

  Calling extend spec on a spec of width m returns two specs that have width $m+1$. By hypothesis, singletonsMatching is correct on these specs. To prove: that the result contains exactly those singletons that match spec. By Lemma One, any word that matches spec matches exactly one of the extended specs. By hypothesis, the recursive calls return exactly the singletons matching the extend specs. When we combine the results of these calls with concat, we get exactly the matching singletons, with no duplicates and no omissions.

Actually solving the problem is anticlimactic: the singletons are all the singletons that match the empty spec:

singletons :: [Word32] -> [Word32]
singletons words = singletonsMatching universalSpec words

Testing code

testa, testb :: [Word32]
testa = [10, 1, 1, 2, 3, 4, 2, 5, 6, 7, 10]
testb = [ 0x0000
        , 0x0010
        , 0x0100
        , 0x0110
        , 0x1000
        , 0x1010
        , 0x1100
        , 0x1110
        ]

Beyond this point, if you want to follow what's going on, you need to know QuickCheck.

Here's a random generator for specs:

instance Arbitrary Spec where
  arbitrary = do width <- choose (0, 32)
                 b <- arbitrary
                 return (randomSpec width b)
  shrink spec = [randomSpec w' (bits spec) | w' <- shrink (width spec)] ++
                [randomSpec (width spec) b | b  <- shrink (bits spec)]
randomSpec width bits = Spec { w = width, bits = mask bits }     
  where mask b = if width == 32 then b
                 else (b `shiftL` n) `shiftR` n
        n = 32 - width

Using this generator, we can test Lemma One using quickCheck lemmaOne.

We can test to see that any word claimed to be a singleton is in fact singleton:

singletonsAreSingleton nzwords = 
  not (hasTriple words) ==> all (`isSingleton` words) (singletons words)
  where isSingleton w words = isSingletonList [w' | w' <- words, w' == w]
        words = [w | NonZero w <- nzwords]

hasTriple :: [Word32] -> Bool
hasTriple words = hasTrip (sort words)
hasTrip (w1:w2:w3:ws) = (w1 == w2 && w2 == w3) || hasTrip (w2:w3:ws)
hasTrip _ = False

Here's another property that tests the fast singletons against a slower algorithm that uses sorting.

singletonsOK :: [NonZero Word32] -> Property
singletonsOK nzwords = not (hasTriple words) ==>
  sort (singletons words) == sort (slowSingletons words)
 where words = [w | NonZero w <- nzwords ]
       slowSingletons words = stripDoubletons (sort words)
       stripDoubletons (w1:w2:ws) | w1 == w2 = stripDoubletons ws
                                  | otherwise = w1 : stripDoubletons (w2:ws)
       stripDoubletons as = as

Answer 2

Disproof of algorithm in OP for k >= 7

This algorithm uses the possibility to recursively split a set of k unique values into two groups using value of a single bit when at least one of these groups is XORed to a nonzero value. For example, the following numbers

01000
00001
10001

may be split into

01000

and

00001
10001

using the value of the least significant bit.

If properly implemented, this works for k <= 6. But this approach fails for k = 8 and k = 7. Let's assume m = 4 and use 8 even numbers from 0 to 14:

0000
0010
0100
0110
1000
1010
1100
1110

Each bit, except the least significant one, has exactly 4 nonzero values. If we try to partition this set, because of this symmetry, we'll always get a subset with 2 or 4 or 0 nonzero values. XOR of these subsets is always 0. Which does not allow algorithm to make any split, so else part just prints the XOR of all these unique values (a single zero).

3x + 1 trick does not help: it only shuffles these 8 values and toggles the least significant bit.

Exactly the same arguments are applicable for k = 7 if we remove the first (all-zero) value from the above subset.

Since any group of unique values may be split into a group of 7 or 8 values and some other group, this algorithm also fails for k > 8.

---

Probabilistic algorithm

It is possible not to invent a completely new algorithm, but instead modify the algorithm in OP, making it work for any input values.

Each time the algorithm accesses an element of the input array, it should apply some transformation function to this element: y=transform(x). This transformed value y may be used exactly as x was used in the original algorithm - for partitioning the sets and XORing the values.

Initially transform(x)=x (unmodified original algorithm). If after this step we have less than k results (some of the results are several unique values XORed), we change transform to some hash function and repeat computations. This should be repeated (each time with different hash function) until we get exactly k values.

If these k values are obtained on the first step of the algorithm (without hashing), these values are our result. Otherwise we should scan the array once more, computing hash of each value and reporting those values, that match one of k hashes.

Each subsequent step of computations with different hash function may be performed either on the original set of k values or (better) separately on each of the subsets, found on previous step.

To obtain different hash function for each step of the algorithm, you can use Universal hashing. One necessary property for the hash function is reversibility - original value should be (in theory) reconstructible from the hash value. This is needed to avoid hashing of several "unique" values to the same hash value. Since using any reversible m-bit hash function has not much chances to solve the problem of "counterexample", hash values should be longer than m bits. One simple example of such hash function is concatenation of the original value and some one-way hash function of this value.

If k is not very large, it is not likely that we get a set of data similar to that counter-example. (I have no proof that there are no other "bad" data patterns, with different structure, but let's hope they are also not very probable). In this case average time complexity is not much larger than O(k * m² * n).

---

Other improvements for the original algorithm

• While computing the XOR of all the (yet unpartitioned) values it is reasonable to check for a unique zero value in the array. If there is one, just decrement k.
• On each recursion step we cannot always know the exact size of each partition. But we know if it is odd or even: each split on a non-zero bit gives odd-sized subset, the other subset's parity is "toggled" parity of the original subset.
• On the latest recursion steps, when the only non-split subset is of size 1, we may skip the search for splitting bit and report the result immediately (this is an optimization for very small k).
• If we get an odd-sized subset after some split (and if we don't know for sure its size is 1), scan the array and try to find a unique value, equal to XOR of this subset.
• There is no need to iterate through every bit to split even-sized set. Just use any non-zero bit of its XORed values. XORing one of the resulting subsets may produce zero, but this split is still valid because we have odd number of "ones" for this splitting bit but even set size. This also means that any split, that produces even-sized subset which is non-zero when XORed, is a valid split, even if the remaining subset XORs to zero.
• You shouldn't continue splitting bit search on each recursion (like solve(ary, h + 1...). Instead you should restart search from the beginning. It is possible to split the set on bit 31, and have the only splitting possibility for one of the resulting subsets on bit 0.
• You shouldn't scan the whole array twice (so second y = compute_xors(ary, m, bits) is not needed). You already have XOR of the whole set and XOR of a subset where the splitting bit is non-zero. Which means you can compute y immediately: y = x ^ old_xor.

---

Proof of algorithm in OP for k = 3

This is a proof not for the actual program in OP, but for its idea. The actual program currently rejects any split when one of the resulting subsets is zero. See suggested improvements for the cases when we may accept some of such splits. So the following proof may be applied to that program only after if x is None or y is None is changed to some condition that takes into account parity of the subset sizes or after a preprocessing step is added to exclude unique zero element from the array.

We have 3 different numbers. They should be different in at least 2 bit positions (if they are different in only one bit, the third number must be equal to one of the others). Loop in the solve function finds leftmost of these bit positions and partitions these 3 numbers into two subsets (of a single number and of 2 distinct numbers). The 2-number subset has equal bits in this bit position, but the numbers still should be different, so there should be one more splitting bit position (obviously, to the right of the first one). Second recursion step easily splits this 2-number subset into two single numbers. Trick with i * 3 + 1 is redundant here: it only doubles complexity of the algorithm.

Here is an illustration for the first split in a set of 3 numbers:

 2  1
*b**yzvw
*b**xzvw
*a**xzvw

We have a loop that iterates through every bit position and computes XOR of the whole words, but separately, one XOR value (A) for true bits in given position, other XOR value (B) for false bit. If number A has zero bit in this position, A contains XOR of some even-sized subset of values, if non-zero - odd-sized subset. The same is true for B. We are interested only in the even-sized subset. It may contain either 0 or 2 values.

While there is no difference in the bit values (bits z, v, w), we have A=B=0, which means we cannot split our numbers on these bits. But we have 3 non-equal numbers, which means at some position (1) we should have different bits (x and y). One of them (x) can be found in two of our numbers (even-sized subset!), other (y) - in one number. Let's look at the XOR of values in this even-sized subset. From A and B select value (C), containing bit 0 at position 1. But C is just a XOR of two non-equal values. They are equal at bit position 1, so they must differ in at least one more bit position (position 2, bits a and b). So C != 0 and it corresponds to the even-sized subset. This split is valid because we can split this even-sized subset further either by very simple algorithm or by next recursion of this algorithm.

If there are no unique zero elements in the array, this proof may be simplified. We always split unique numbers into 2 subsets - one with 2 elements (and it cannot XOR to zero because elements are different), other with one element (non-zero by definition). So the original program with little pre-processing should work properly.

Complexity is O(m² * n). If you apply the improvements I suggested earlier, the expected number of times this algorithm scans the array is m / 3 + 2. Because the first splitting bit position is expected to be m / 3, a single scan is needed to deal with 2-element subset, every 1-element subset does not need any array scans, and one more scan is needed initially (outside of solve method).

---

Proof of algorithm in OP for k = 4 .. 6

Here we assume that all the suggested improvements to the original algorithm are applied.

k=4 and k=5: Since there is at least one position with different bits, this set of numbers can be split in such a way that one of the subsets has size 1 or 2. If subset's size is 1, it is non-zero (we have no zero unique values). If subset's size is 2, we have XOR of two different numbers, which is non-zero. So in both cases the split is valid.

k=6: If XOR of the whole set is non-zero, we can split this set by any position where this XOR has non-zero bit. Otherwise we have even number of non-zero bit in each position. Since there is at least one position with different bits, this position splits the set into subsets of sizes 2 and 4. Subset of size 2 has always non-zero XOR because it contains 2 different numbers. Again, in both cases we have the valid split.

---

Deterministic algorithm

Disproof for k >= 7 shows the pattern where original algorithm does not work: we have a subset of size greater than 2 and at each bit position we have even number of non-zero bits. But we can always find a pair of positions where non-zero bits overlap in single number. In other words, it is always possible to find a pair of positions in the subset of size 3 or 4 with non-zero XOR of all bits in the subset in both positions. This suggest us to use an additional split-position: iterate through bit positions with two separate pointers, group all numbers in the array into two subsets where one subset has both non-zero bits in these positions, and other - all the remaining numbers. This increases the worst case complexity my m, but allows more values for k. Once there is no more possible to obtain a subset of size less than 5, add the third "splitting pointer", and so on. Each time k doubles, we may need an additional "splitting pointer", which increases the worst case complexity my m once more.

This might be considered as a sketch of a proof for the following algorithm:

1. Use original (improved) algorithm to find zero or more unique values and zero or more non-splittable subsets. Stop when there are no more non-splittable subsets.
2. For any of these non-splittable subsets, try to split it while increasing the number of "splitting pointers". When split is found, continue with step 1.

Worst case complexity is O(k * m² * n * m^{max(0, floor(log(floor(k/4))))}), which may be approximated by O(k * n * m^log(k)) = O(k * n * k^log(m)).

Expected run time of this algorithm for small k is a little bit worse than for probabilistic algorithm, but still not much larger than O(k * m² * n).