Combinatorial Generation: Subsets, Lexicographic Permutations (Narayana) & Gray Codes

1. In-Place State Machines vs. Recursive Materialization

When enumerating $2^n$ subsets or $n!$ permutations, materializing all states into a list requires $\Omega(n \cdot 2^n)$ or $\Omega(n \cdot n!)$ memory. For $n = 13$, $13! = 6{,}227{,}020{,}800$, which exceeds 75 GB of RAM. Instead, production combinatorial generators maintain a single mutable state vector $A[0..n-1]$ and implement an in-place next() transition that mutates $A$ in $O(1)$ amortized time.

Why is std::next_permutation $O(1)$ amortized per permutation? In half of all $n!$ states, $A[n-2] < A[n-1]$, so the algorithm compares 1 pair and swaps 2 elements. In $1/k!$ of states, the suffix of length $k$ is reversed, doing $O(k)$ work. Summing over all $k$ gives $\sum_{k=2}^{n} \frac{c \cdot k}{k!} = O(1)$ amortized comparisons and swaps per call.

2. Reflected Binary Gray Code & Gosper's Hack

Standard binary counting from $0$ to $2^n - 1$ flips up to $n$ bits on a single increment (e.g., $0111_2 \to 1000_2$), which means maintaining a running subset property (such as a knapsack weight or XOR sum) takes $O(1)$ amortized work only if we inspect the changed bits, or $O(1)$ worst-case work if we use a Binary Reflected Gray Code (BRGC).

• Closed-Form $i$-th Gray Code: $G(i) = i \oplus (i \gg 1)$.
• Bit Flipped at Step $i \to i+1$: The index of the single bit that changes between $G(i)$ and $G(i+1)$ is the number of trailing zeros of $i+1$, i.e., __builtin_ctz(i + 1).
• Inverse Gray Code ($G(i) \to i$): Repeatedly XOR with shifted versions: for (uint64_t mask = g >> 1; mask; mask >>= 1) g ^= mask; (or in $\log_2(64) = 6$ parallel prefix XOR shifts).
• Gosper's Hack ($inom{n}{k}$ Fixed-Popcount Subsets): Given a bitmask $x > 0$ with exactly $k$ set bits, the lexicographically next integer with $k$ set bits is computed in 4 ALU ops without branching.

3. Complexity Matrix

4. Reference C++20 Implementation

#include <algorithm>
#include <cstdint>
#include <vector>

// 1. Narayana's Lexicographic Next Permutation (handles duplicates seamlessly!)
template <typename T>
bool next_permutation_narayana(std::vector<T>& a) {
    const int n = static_cast<int>(a.size());
    int i = n - 2;
    // Strict '<' ensures duplicate elements are never redundantly swapped
    while (i >= 0 && !(a[i] < a[i + 1])) --i;
    if (i < 0) {
        std::reverse(a.begin(), a.end());
        return false; // Wrapped around to lexicographically smallest state
    }
    int j = n - 1;
    while (!(a[i] < a[j])) --j;
    std::swap(a[i], a[j]);
    std::reverse(a.begin() + i + 1, a.end());
    return true;
}

// 2. Gosper's Hack: Iterate over all n-bit masks with exactly k set bits in O(1) per mask
template <typename Callback>
void enumerate_k_subsets_gosper(int n, int k, Callback&& visit) {
    if (k == 0) { visit(0ULL); return; }
    if (k > n) return;
    std::uint64_t mask = (1ULL << k) - 1;
    const std::uint64_t limit = (n == 64) ? ~0ULL : (1ULL << n);
    while (mask < limit) {
        visit(mask);
        std::uint64_t c = mask & -mask;              // Lowest set bit (BLSI)
        std::uint64_t r = mask + c;                  // Carry propagation across lowest 1-run
        if (r == 0) break;                           // Overflow guard when n = 64
        mask = (((r ^ mask) >> 2) / c) | r;          // Or: (((r ^ mask) >> (__builtin_ctzll(c) + 2)) | r)
    }
}

// 3. Gray Code Incremental Subset Sum Evaluator in O(2^n) total additions (not O(n 2^n))
long long count_subsets_with_target_gray(const std::vector<long long>& weights, long long target) {
    const int n = static_cast<int>(weights.size());
    long long current_sum = 0;
    long long matches = (current_sum == target) ? 1 : 0;
    std::uint64_t prev_gray = 0;
    const std::uint64_t total = 1ULL << n;

    for (std::uint64_t i = 1; i < total; ++i) {
        std::uint64_t curr_gray = i ^ (i >> 1);
        int bit = __builtin_ctzll(i); // Exactly bit index that flipped between prev_gray and curr_gray
        if ((curr_gray >> bit) & 1ULL) {
            current_sum += weights[bit];
        } else {
            current_sum -= weights[bit];
        }
        if (current_sum == target) ++matches;
        prev_gray = curr_gray;
    }
    return matches;
}