Introduction

The intersection of two sets represented by sorted collections (like lists or arrays) can be done in linear time. If we label with k the size of one collection, and with n the size of another collection, then the complexity of intersection is O(n + k).

Below is a naive C++ implementation; the C++ standard library comes with std::set_intersection.

template <typename INSERTER>
void custom_set_intersection(const vec& A, const vec& B, INSERTER output) {

    size_t i = 0;
    size_t j = 0;

    while (i < A.size() && j < B.size()) {
        if (A[i] < B[j]) {
            i += 1;
        } else if (B[j] < A[i]) {
            j += 1;
        } else {
            // A[i] == B[j]
            *output++ = A[i];

            i += 1;
            j += 1;
        }
    }
}

If there is a huge difference between sizes of collections, then we might use two different approaches described below. I use terms "smaller" (k items) and "larger" collection/set (n items).

Binary search

In this approach we iterate over the smaller collection. For each item the binary search on the larger collection is performed. The complexity of this solution is O(k ⋅ logn).

A clear disadvantage is that at least the larger collection must support random access. It must be an array, a vector or something similar.

To make a practical implementation faster, we utilize the fact that both collections are sorted. In each iteration we are not looking for an equal element but for the first element which is greater or equal to the searched one. Thanks to that the searched range is narrowed after each iteration: the part of the larger collection containing elements smaller than the elements already tested is never touched again. It doesn't change the asymptotic complexity, though.

In C++ the function std::lower_bound implements this look up; it's a modification of the binary search algorithm.

Below is a C++ implementation.

template <typename INSERTER>
void binsearch_set_intersection(const vec& A, const vec& B, INSERTER output) {

    auto it = B.begin();

    for (const auto& a: A) {
        it = std::lower_bound(it, B.end(), a);
        if (it == B.end()) {
            // since there are no values greater or equal the current,
            // there also won't be any for subsequent value (all greater than current)
            return;
        }

        if (*it == a) {
            output = a;
            ++it;
        }
    }
}

SIMD algorithm

A SIMD approach is vectorized linear algorithm in which a single value from the smaller set is compared with vectors read from the larger set. Thanks to that the number of comparisons is smaller and we can faster skip parts of the larger collection.

In the example below we assume that the first element of the smaller collection is 9, and the head of the larger collection contains values [-5, 1, 2, 5, 7, 8, 9, 12].

• Fill the vector A with the first item from the smaller collection:

  +-----+-----+-----+-----+
  |  9  |  9  |  9  |  9  |
  +-----+-----+-----+-----+
  

• Load into another vector B a part of the larger collection:

  +-----+-----+-----+-----+
  | -5  |  1  |  2  |  5  |
  +-----+-----+-----+-----+
  

• Compare B < A:

  +-----+-----+-----+-----+
  | -1  | -1  | -1  | -1  |
  +-----+-----+-----+-----+
  

• If all items in B are smaller (all bits of the comparison vector are set), then load the next chunk into B:

  +-----+-----+-----+-----+
  |  7  |  8  |  9  | 12  |
  +-----+-----+-----+-----+
  

• And keep comparing with A; now B < A vector is equal:

  +-----+-----+-----+-----+
  | -1  | -1  |  0  |  0  |
  +-----+-----+-----+-----+
  

• The result vector contains some zero items. It means that B may contain the value present in A. If it's true (in the example it is), save this value in the output collection. In either case fill A with the next value from the smaller collection and again compare with B.

Below is a sample implementation that uses SSE instructions.

template <typename INSERTER>
void sse_set_intersection(const vec& A, const vec& B, INSERTER output) {

    __m128i a_rep;
    __m128i b;

    size_t ai = 0;
    size_t bi = 0;

    a_rep = _mm_set1_epi32(A[ai]);
    b     = _mm_loadu_si128(reinterpret_cast<const __m128i*>(&B[bi]));

    while (ai < A.size() && bi < B.size()) {
        const __m128i lt = _mm_cmplt_epi32(b, a_rep);
        const uint16_t mask = _mm_movemask_epi8(lt);
        if (mask == 0xffff) {
            // all elements in b are smaller, fetch the next chunk from B
            bi += 4;
            b = _mm_loadu_si128(reinterpret_cast<const __m128i*>(&B[bi]));
        } else {
            // there might be element equal to A[ai]
            // a simple linear search, as there're only 4 elements to search in
            for (size_t i = 0; i < 4; i++) {
                if (B[bi + i] == A[ai]) {
                    output = A[ai];
                    break;
                }
            }

            // fetch the next value from A
            ai += 1;
            a_rep = _mm_set1_epi32(A[ai]);
        }
    }
}

Evaluation

Experiment:

• The larger collection is a C++ std::vector<int32_t> that has got 1024 * 1024 random elements.
• The size of smaller table varies. It is always a subset of the larger collection.
• Intersection is done 1000 times and the minimum time is noted.

Full source code is available.

Conculsions

• When differences in size of sets is huge the presented approaches might be significantly faster.
• The algorithm that use binary search is faster, but the size ratio is strictly limited.
• Both SSE and AVX2 implementations perform better than linear search in circumstances when binary search become slower.