代码

#define _GNU_SOURCE
#pragma GCC target("avx2")
#include <immintrin.h>
#include <sched.h>

#define N (1U << 27)
#define BUCKETS 256U
#define BUFFER_SIZE 16U
#define MAX_THREADS 4U
#define STACK_SIZE 65536U

typedef struct __attribute__((aligned(4096))) {
    unsigned high_count[BUCKETS];
    unsigned position[BUCKETS];
    unsigned local_count[3][BUCKETS];
    unsigned fill[BUCKETS];
    unsigned buffer[BUCKETS][BUFFER_SIZE] __attribute__((aligned(64)));
    unsigned id;
} Worker;

static unsigned work[N] __attribute__((aligned(4096)));
static unsigned starts[BUCKETS + 1];
static Worker workers[MAX_THREADS];
static unsigned char stacks[MAX_THREADS - 1][STACK_SIZE]
    __attribute__((aligned(64)));
static unsigned count_done[MAX_THREADS] __attribute__((aligned(64)));
static unsigned scatter_done[MAX_THREADS] __attribute__((aligned(64)));
static unsigned local_done[MAX_THREADS] __attribute__((aligned(64)));
static unsigned active_threads, start_flag, scatter_flag, local_flag;
static unsigned next_bucket;
static unsigned *input_array;

static inline void wait_for(const unsigned *flag) {
    while (!__atomic_load_n(flag, __ATOMIC_ACQUIRE)) _mm_pause();
}

static inline void wait_all(const unsigned *flags) {
    unsigned threads = __atomic_load_n(&active_threads, __ATOMIC_ACQUIRE);
    for (unsigned i = 0; i < threads; ++i) wait_for(flags + i);
}

static inline void prefix_from(unsigned *count, unsigned base) {
    unsigned sum = base;
    for (unsigned i = 0; i < BUCKETS; ++i) {
        unsigned value = count[i]; count[i] = sum; sum += value;
    }
}

static __attribute__((always_inline)) inline void flush_nt(
        unsigned *destination, const unsigned *source) {
    _mm256_stream_si256((__m256i *)(destination + 0),
                        _mm256_load_si256((const __m256i *)(source + 0)));
    _mm256_stream_si256((__m256i *)(destination + 8),
                        _mm256_load_si256((const __m256i *)(source + 8)));
}

static __attribute__((always_inline)) inline void radix_pass_range(
        Worker *worker, const unsigned *source, unsigned *destination,
        unsigned *position, unsigned shift, unsigned begin, unsigned end) {
    unsigned *fill = worker->fill;
    for (unsigned i = begin; i < end; ++i) {
        unsigned value = source[i];
        unsigned digit = (value >> shift) & 255U;
        unsigned output = position[digit];
        if (output & 7U) {
            destination[output] = value;
            position[digit] = output + 1;
        } else {
            unsigned slot = fill[digit];
            worker->buffer[digit][slot] = value;
            if (slot + 1 == BUFFER_SIZE) {
                flush_nt(destination + output, worker->buffer[digit]);
                position[digit] = output + BUFFER_SIZE;
                fill[digit] = 0;
            } else fill[digit] = slot + 1;
        }
    }
    for (unsigned digit = 0; digit < BUCKETS; ++digit) {
        unsigned amount = fill[digit], output = position[digit];
        for (unsigned j = 0; j < amount; ++j)
            destination[output + j] = worker->buffer[digit][j];
        position[digit] = output + amount;
        fill[digit] = 0;
    }
    _mm_sfence();
}

static int worker_main(void *argument) {
    Worker *worker = argument;
    wait_for(&start_flag);
    unsigned id = worker->id;
    unsigned threads = __atomic_load_n(&active_threads, __ATOMIC_ACQUIRE);
    unsigned begin = (unsigned)(((unsigned long long)N * id) / threads);
    unsigned end = (unsigned)(((unsigned long long)N * (id + 1)) / threads);

    for (unsigned i = begin; i < end; ++i)
        ++worker->high_count[input_array[i] >> 24];
    __atomic_store_n(count_done + id, 1, __ATOMIC_RELEASE);
    if (id == 0) {
        wait_all(count_done);
        unsigned sum = 0;
        for (unsigned digit = 0; digit < BUCKETS; ++digit) {
            starts[digit] = sum;
            for (unsigned thread = 0; thread < threads; ++thread) {
                workers[thread].position[digit] = sum;
                sum += workers[thread].high_count[digit];
            }
        }
        starts[BUCKETS] = N;
        __atomic_store_n(&scatter_flag, 1, __ATOMIC_RELEASE);
    } else wait_for(&scatter_flag);

    radix_pass_range(worker, input_array, work, worker->position, 24,
                     begin, end);
    __atomic_store_n(scatter_done + id, 1, __ATOMIC_RELEASE);
    if (id == 0) {
        wait_all(scatter_done);
        next_bucket = 0;
        __atomic_store_n(&local_flag, 1, __ATOMIC_RELEASE);
    } else wait_for(&local_flag);

    for (;;) {
        unsigned bucket = __atomic_fetch_add(&next_bucket, 1,
                                              __ATOMIC_RELAXED);
        if (bucket >= BUCKETS) break;
        begin = starts[bucket]; end = starts[bucket + 1];
        __builtin_memset(worker->local_count, 0,
                         sizeof(worker->local_count));
        for (unsigned i = begin; i < end; ++i) {
            unsigned value = work[i];
            ++worker->local_count[0][(unsigned char)value];
            ++worker->local_count[1][(unsigned char)(value >> 8)];
            ++worker->local_count[2][(unsigned char)(value >> 16)];
        }
        prefix_from(worker->local_count[0], begin);
        prefix_from(worker->local_count[1], begin);
        prefix_from(worker->local_count[2], begin);
        radix_pass_range(worker, work, input_array, worker->local_count[0],
                         0, begin, end);
        radix_pass_range(worker, input_array, work, worker->local_count[1],
                         8, begin, end);
        radix_pass_range(worker, work, input_array, worker->local_count[2],
                         16, begin, end);
    }
    __atomic_store_n(local_done + id, 1, __ATOMIC_RELEASE);
    if (id == 0) wait_all(local_done);
    _mm256_zeroupper();
    return 0;
}

void sort(unsigned *a, int n) {
    (void)n;
    input_array = a;
    workers[0].id = 0;
    unsigned threads = 1;
    const int flags = CLONE_VM | CLONE_FS | CLONE_FILES | CLONE_SIGHAND |
                      CLONE_THREAD | CLONE_SYSVSEM;
    for (unsigned i = 1; i < MAX_THREADS; ++i) {
        workers[i].id = i;
        void *top = stacks[i - 1] + STACK_SIZE;
        if (clone(worker_main, top, flags, &workers[i]) < 0) break;
        ++threads;
    }
    __atomic_store_n(&active_threads, threads, __ATOMIC_RELEASE);
    __atomic_store_n(&start_flag, 1, __ATOMIC_RELEASE);
    worker_main(&workers[0]);
}

评测结果