代码

#include <algorithm>
#include <cstdint>
#include <cstdlib>
#include <cstring>
#include <limits>

#ifndef SORT_CHAINS
#define SORT_CHAINS 16
#endif

namespace fast_u32_detail {
constexpr unsigned TOP = 256, MID = 16, RADIX = 1024;
constexpr unsigned LIMIT = 1u << 20;
static_assert(std::numeric_limits<unsigned>::digits == 32,
              "This implementation requires 32-bit unsigned.");

// Prefix sums are exact. All holes belong to the bucket currently processed.
template<unsigned Bits, unsigned Shift, unsigned Chains>
void partition(unsigned* a, const unsigned* start) {
    constexpr unsigned K = 1u << Bits;
    alignas(64) unsigned next[K];
    std::memcpy(next, start, sizeof(next));
    for (unsigned k = 0; k < K; ++k) {
        const unsigned end = start[k + 1];
        if (end - start[k] == start[K]) return; // One nonempty bucket.
        if constexpr (Chains > 1) {
            if (end - next[k] >= Chains) {
                unsigned value[Chains], hole[Chains];
                for (unsigned r = 0; r < Chains; ++r) {
                    hole[r] = next[k]++;
                    value[r] = a[hole[r]];
                }
                while (end - next[k] >= Chains) {
                    #pragma GCC unroll 16
                    for (unsigned r = 0; r < Chains; ++r) {
                        unsigned d = (value[r] >> Shift) & (K - 1);
                        if (d != k) {
                            unsigned j = next[d]++;
                            unsigned x = a[j];
                            a[j] = value[r];
                            value[r] = x;
                        } else {
                            a[hole[r]] = value[r];
                            hole[r] = next[k]++;
                            value[r] = a[hole[r]];
                        }
                    }
                }
                for (unsigned r = 0; r < Chains; ++r) {
                    unsigned v = value[r], d = (v >> Shift) & (K - 1);
                    while (d != k) {
                        unsigned j = next[d]++;
                        unsigned x = a[j];
                        a[j] = v; v = x;
                        d = (v >> Shift) & (K - 1);
                    }
                    a[hole[r]] = v;
                }
            }
        }
        while (next[k] < end) {
            unsigned v = a[next[k]], d = (v >> Shift) & (K - 1);
            while (d != k) {
                unsigned j = next[d]++;
                unsigned x = a[j];
                a[j] = v; v = x;
                d = (v >> Shift) & (K - 1);
            }
            a[next[k]++] = v;
        }
    }
}

void counting(const unsigned* src, unsigned* dst, unsigned n,
              unsigned prefix, unsigned bits, unsigned* hist) {
    const unsigned size = 1u << bits, mask = size - 1;
    std::memset(hist, 0, size * sizeof(unsigned));
    for (unsigned i = 0; i < n; ++i) ++hist[src[i] & mask];
    for (unsigned k = 0; k < size; ++k) {
        unsigned c = hist[k];
        if (c) {
            std::fill_n(dst, c, prefix | k);
            dst += c;
        }
    }
}

// The upper 12 bits are fixed. Only the lower 20 bits need sorting.
// Stable low-10-bit scatter drops that digit; it is recovered from bucket order.
void finish(const unsigned* src, unsigned* dst, unsigned n,
            unsigned prefix, uint16_t* tmp, unsigned* dense) {
    if (!n) return;
    if (n <= 32) {
        if (src != dst) std::memcpy(dst, src, size_t(n) * sizeof(unsigned));
        std::sort(dst, dst + n);
        return;
    }
    if (n > LIMIT) {
        counting(src, dst, n, prefix, 20, dense);
        return;
    }
    alignas(64) unsigned low[RADIX] = {}, high[RADIX] = {};
    alignas(64) unsigned pl[RADIX], ph[RADIX];
    for (unsigned i = 0; i < n; ++i) {
        unsigned v = src[i];
        ++low[v & 1023];
        ++high[(v >> 10) & 1023];
    }
    unsigned x = 0, y = 0;
    for (unsigned k = 0; k < RADIX; ++k) {
        pl[k] = x; x += low[k];
        ph[k] = y; y += high[k];
    }
    for (unsigned i = 0; i < n; ++i) {
        unsigned v = src[i];
        tmp[pl[v & 1023]++] = uint16_t((v >> 10) & 1023);
    }
    unsigned at = 0;
    for (unsigned lo = 0; lo < RADIX; ++lo) {
        unsigned end = at + low[lo], common = prefix | lo;
        for (; at < end; ++at) {
            unsigned hi = tmp[at];
            dst[ph[hi]++] = common | (hi << 10);
        }
    }
}

struct Workspace {
    unsigned* scatter = nullptr;
    uint16_t* temp = nullptr;
    unsigned* dense = nullptr;
    ~Workspace() { std::free(scatter); std::free(temp); std::free(dense); }
};
}

void sort(unsigned* a, int n) {
    using namespace fast_u32_detail;
    if (n <= 1) return;
    if (n <= 32) { std::sort(a, a + n); return; }

    // Count the high 12 bits once; derive both partition levels from this scan.
    alignas(64) unsigned hist[2][4096] = {};
    unsigned different = 0, first = a[0];
    int i = 0;
    for (; i <= n - 4; i += 4) {
        unsigned x0 = a[i], x1 = a[i+1], x2 = a[i+2], x3 = a[i+3];
        ++hist[0][x0 >> 20]; ++hist[1][x1 >> 20];
        ++hist[0][x2 >> 20]; ++hist[1][x3 >> 20];
        different |= (x0 ^ first) | (x1 ^ first) | (x2 ^ first) | (x3 ^ first);
    }
    for (; i < n; ++i) {
        ++hist[0][a[i] >> 20]; different |= a[i] ^ first;
    }
    if (!different) return;

    Workspace w;
    if (different <= 65535u) {
        w.dense = static_cast<unsigned*>(std::malloc(65536u * sizeof(unsigned)));
        if (!w.dense) { std::sort(a, a + n); return; }
        counting(a, a, unsigned(n), first & 0xffff0000u, 16, w.dense);
        return;
    }

    alignas(64) unsigned start[4097], top[257];
    start[0] = 0;
    unsigned max_temp = 0; bool need_dense = false;
    for (unsigned k = 0; k < 4096; ++k) {
        unsigned c = hist[0][k] + hist[1][k];
        start[k+1] = start[k] + c;
        if (c > LIMIT) need_dense = true;
        else if (c > 32) max_temp = std::max(max_temp, c);
    }
    unsigned max_scatter = 0;
    for (unsigned k = 0; k <= TOP; ++k) top[k] = start[k * MID];
    for (unsigned k = 0; k < TOP; ++k) {
        unsigned c = top[k+1] - top[k];
        if (c > 32 && c <= LIMIT) max_scatter = std::max(max_scatter, c);
    }
    if (max_scatter) w.scatter = static_cast<unsigned*>(std::malloc(size_t(max_scatter) * sizeof(unsigned)));
    if (max_temp) w.temp = static_cast<uint16_t*>(std::malloc(size_t(max_temp) * sizeof(uint16_t)));
    if (need_dense) w.dense = static_cast<unsigned*>(std::malloc(size_t(LIMIT) * sizeof(unsigned)));
    if ((max_scatter && !w.scatter) || (max_temp && !w.temp) || (need_dense && !w.dense)) {
        std::sort(a, a + n); return;
    }

    partition<8, 24, SORT_CHAINS>(a, top);
    for (unsigned k = 0; k < TOP; ++k) {
        const unsigned base = top[k], c = top[k+1] - base;
        if (c <= 32) { std::sort(a + base, a + base + c); continue; }
        unsigned sub[17], p[16];
        for (unsigned j = 0; j <= MID; ++j) sub[j] = start[k * MID + j] - base;
        const unsigned* src;
        if (c <= LIMIT) {
            std::memcpy(p, sub, sizeof(p));
            for (unsigned j = 0; j < c; ++j) {
                unsigned v = a[base+j];
                w.scatter[p[(v >> 20) & 15]++] = v;
            }
            src = w.scatter;
        } else {
            partition<4, 20, 1>(a + base, sub);
            src = a + base;
        }
        for (unsigned j = 0; j < MID; ++j) {
            unsigned m = sub[j+1] - sub[j];
            finish(src + sub[j], a + base + sub[j], m,
                   ((k << 4) | j) << 20, w.temp, w.dense);
        }
    }
}

评测结果