代码

// 1001 v14: [8][8][16] MSD + THREE byte-plane scratch + DIRECT scatter into a[].
// (leaner version of v10: 3 separate byte planes instead of packed 3-byte)
// Pass 1: hist top-8 -> B[257].
// Pass 2: stream a[]; k=a[i]>>24; p=curs[k]++; s0[p]=lo, s1[p]=mid, s2[p]=hi16 (3 byte streams/bucket).
// Pass 3: per bucket, hist byte2 (read s2 region) -> 256 sub-boundaries; reconstruct u32=(k<<24)|(g<<16)|low16
//         and store directly into a[] at [B[k]+sub_start[g]+cursor[g]++] (no swap).
// Pass 4: per byte2-sub-bucket, 2x8-bit LSD sort of low-16 in a[] (L1-resident).
#include <string.h>
typedef unsigned int u32;
typedef unsigned char u8;
typedef unsigned short u16;

#ifndef N
#define N 100000000
#endif

#define SUBMAX 8192

static u8 s0[N] __attribute__((aligned(64)));
static u8 s1[N] __attribute__((aligned(64)));
static u8 s2[N] __attribute__((aligned(64)));
static u32 cnt[256];
static u32 B[257];
static u32 curs[256];      // per-bucket element cursor
static u32 c2[256];        // byte2 histogram within one bucket
static u32 sstart[256];    // sub-bucket start (element offset within bucket)
static u32 scurs[256];     // sub-bucket cursor
static u32 h[256];         // low-byte histogram (pass 4)
static u32 h2[256];        // high-byte histogram (pass 4)
static u32 tmpbuf[SUBMAX] __attribute__((aligned(64)));

// Stable 2x8-bit LSD sort of low-16 bits of m u32 values in base[] (top-16 fixed).
static void sort_low16(u32 *base, u32 m) {
    int i, j;
    if (m <= 1) return;
    for (j = 0; j < 256; j++) h[j] = 0;
    int m8 = (int)m & ~7;
    for (i = 0; i < m8; i += 8) {
        h[base[i]&255]++; h[base[i+1]&255]++; h[base[i+2]&255]++; h[base[i+3]&255]++;
        h[base[i+4]&255]++; h[base[i+5]&255]++; h[base[i+6]&255]++; h[base[i+7]&255]++;
    }
    for (; i < (int)m; i++) h[base[i] & 255]++;
    for (j = 1; j < 256; j++) h[j] += h[j-1];
    int rm = (int)m & ~7;
    for (i = (int)m - 1; i >= rm; i--) { u32 x = base[i]; tmpbuf[--h[x & 255]] = x; }
    for (i = rm - 1; i >= 0; i -= 8) {
        u32 x0=base[i],x1=base[i-1],x2=base[i-2],x3=base[i-3],x4=base[i-4],x5=base[i-5],x6=base[i-6],x7=base[i-7];
        tmpbuf[--h[x0&255]]=x0; tmpbuf[--h[x1&255]]=x1; tmpbuf[--h[x2&255]]=x2; tmpbuf[--h[x3&255]]=x3;
        tmpbuf[--h[x4&255]]=x4; tmpbuf[--h[x5&255]]=x5; tmpbuf[--h[x6&255]]=x6; tmpbuf[--h[x7&255]]=x7;
    }
    for (j = 0; j < 256; j++) h2[j] = 0;
    for (i = 0; i < m8; i += 8) {
        h2[(tmpbuf[i]>>8)&255]++; h2[(tmpbuf[i+1]>>8)&255]++; h2[(tmpbuf[i+2]>>8)&255]++; h2[(tmpbuf[i+3]>>8)&255]++;
        h2[(tmpbuf[i+4]>>8)&255]++; h2[(tmpbuf[i+5]>>8)&255]++; h2[(tmpbuf[i+6]>>8)&255]++; h2[(tmpbuf[i+7]>>8)&255]++;
    }
    for (; i < (int)m; i++) h2[(tmpbuf[i] >> 8) & 255]++;
    for (j = 1; j < 256; j++) h2[j] += h2[j-1];
    for (i = (int)m - 1; i >= rm; i--) { u32 x = tmpbuf[i]; base[--h2[(x >> 8) & 255]] = x; }
    for (i = rm - 1; i >= 0; i -= 8) {
        u32 x0=tmpbuf[i],x1=tmpbuf[i-1],x2=tmpbuf[i-2],x3=tmpbuf[i-3],x4=tmpbuf[i-4],x5=tmpbuf[i-5],x6=tmpbuf[i-6],x7=tmpbuf[i-7];
        base[--h2[(x0>>8)&255]]=x0; base[--h2[(x1>>8)&255]]=x1; base[--h2[(x2>>8)&255]]=x2; base[--h2[(x3>>8)&255]]=x3;
        base[--h2[(x4>>8)&255]]=x4; base[--h2[(x5>>8)&255]]=x5; base[--h2[(x6>>8)&255]]=x6; base[--h2[(x7>>8)&255]]=x7;
    }
}

void sort(unsigned *aa, int n) {
    u32 *a = (u32*)aa;
    int i, j;
    (void)n;

    // ---- Pass 1: histogram top-8 ----
    for (j = 0; j < 256; j++) cnt[j] = 0;
    int n8 = N & ~7;
    for (i = 0; i < n8; i += 8) {
        cnt[a[i]>>24]++; cnt[a[i+1]>>24]++; cnt[a[i+2]>>24]++; cnt[a[i+3]>>24]++;
        cnt[a[i+4]>>24]++; cnt[a[i+5]>>24]++; cnt[a[i+6]>>24]++; cnt[a[i+7]>>24]++;
    }
    for (; i < N; i++) cnt[a[i] >> 24]++;
    { u32 acc = 0; B[0] = 0; for (j = 0; j < 256; j++) { acc += cnt[j]; B[j+1] = acc; } }

    // ---- Pass 2: compress-scatter to 3 byte planes ----
    for (j = 0; j < 256; j++) curs[j] = B[j];
    for (i = 0; i < n8; i += 8) {
        u32 x0=a[i+0],k0=x0>>24,p0=curs[k0]++; s0[p0]=(u8)x0; s1[p0]=(u8)(x0>>8); s2[p0]=(u8)(x0>>16);
        u32 x1=a[i+1],k1=x1>>24,p1=curs[k1]++; s0[p1]=(u8)x1; s1[p1]=(u8)(x1>>8); s2[p1]=(u8)(x1>>16);
        u32 x2=a[i+2],k2=x2>>24,p2=curs[k2]++; s0[p2]=(u8)x2; s1[p2]=(u8)(x2>>8); s2[p2]=(u8)(x2>>16);
        u32 x3=a[i+3],k3=x3>>24,p3=curs[k3]++; s0[p3]=(u8)x3; s1[p3]=(u8)(x3>>8); s2[p3]=(u8)(x3>>16);
        u32 x4=a[i+4],k4=x4>>24,p4=curs[k4]++; s0[p4]=(u8)x4; s1[p4]=(u8)(x4>>8); s2[p4]=(u8)(x4>>16);
        u32 x5=a[i+5],k5=x5>>24,p5=curs[k5]++; s0[p5]=(u8)x5; s1[p5]=(u8)(x5>>8); s2[p5]=(u8)(x5>>16);
        u32 x6=a[i+6],k6=x6>>24,p6=curs[k6]++; s0[p6]=(u8)x6; s1[p6]=(u8)(x6>>8); s2[p6]=(u8)(x6>>16);
        u32 x7=a[i+7],k7=x7>>24,p7=curs[k7]++; s0[p7]=(u8)x7; s1[p7]=(u8)(x7>>8); s2[p7]=(u8)(x7>>16);
    }
    for (; i < N; i++) {
        u32 x = a[i]; u32 k = x >> 24; u32 p = curs[k]++;
        s0[p] = (u8)x; s1[p] = (u8)(x >> 8); s2[p] = (u8)(x >> 16);
    }

    // ---- Pass 3 + 4: per top-8 bucket ----
    for (int k = 0; k < 256; k++) {
        u32 m = cnt[k];
        if (m == 0) continue;
        u32 base = B[k];
        int mb8 = (int)m & ~7;

        // (a) histogram byte2
        for (j = 0; j < 256; j++) c2[j] = 0;
        for (i = 0; i < mb8; i += 8) {
            c2[s2[base+i]]++; c2[s2[base+i+1]]++; c2[s2[base+i+2]]++; c2[s2[base+i+3]]++;
            c2[s2[base+i+4]]++; c2[s2[base+i+5]]++; c2[s2[base+i+6]]++; c2[s2[base+i+7]]++;
        }
        for (; i < (int)m; i++) c2[s2[base + i]]++;
        { u32 acc = 0; for (j = 0; j < 256; j++) { sstart[j] = acc; scurs[j] = acc; acc += c2[j]; } }

        // (b) reconstruct + scatter directly into a[]
        u32 kk = (u32)k << 24;
        for (i = 0; i < mb8; i += 8) {
            u32 g0=s2[base+i]; u32 p0=scurs[g0]++; a[base+p0] = kk | ((u32)g0<<16) | ((u32)s1[base+i]<<8) | (u32)s0[base+i];
            u32 g1=s2[base+i+1]; u32 p1=scurs[g1]++; a[base+p1] = kk | ((u32)g1<<16) | ((u32)s1[base+i+1]<<8) | (u32)s0[base+i+1];
            u32 g2=s2[base+i+2]; u32 p2=scurs[g2]++; a[base+p2] = kk | ((u32)g2<<16) | ((u32)s1[base+i+2]<<8) | (u32)s0[base+i+2];
            u32 g3=s2[base+i+3]; u32 p3=scurs[g3]++; a[base+p3] = kk | ((u32)g3<<16) | ((u32)s1[base+i+3]<<8) | (u32)s0[base+i+3];
            u32 g4=s2[base+i+4]; u32 p4=scurs[g4]++; a[base+p4] = kk | ((u32)g4<<16) | ((u32)s1[base+i+4]<<8) | (u32)s0[base+i+4];
            u32 g5=s2[base+i+5]; u32 p5=scurs[g5]++; a[base+p5] = kk | ((u32)g5<<16) | ((u32)s1[base+i+5]<<8) | (u32)s0[base+i+5];
            u32 g6=s2[base+i+6]; u32 p6=scurs[g6]++; a[base+p6] = kk | ((u32)g6<<16) | ((u32)s1[base+i+6]<<8) | (u32)s0[base+i+6];
            u32 g7=s2[base+i+7]; u32 p7=scurs[g7]++; a[base+p7] = kk | ((u32)g7<<16) | ((u32)s1[base+i+7]<<8) | (u32)s0[base+i+7];
        }
        for (; i < (int)m; i++) {
            u32 g = s2[base+i];
            u32 p = scurs[g]++;
            a[base+p] = kk | ((u32)g<<16) | ((u32)s1[base+i]<<8) | (u32)s0[base+i];
        }

        // (c) pass 4: sort low-16 per sub-bucket
        for (j = 0; j < 256; j++) {
            u32 c = c2[j];
            if (c > 0) sort_low16(a + base + sstart[j], c);
        }
    }
}

评测结果