代码

// 4D strict dominance via bitset blocking (for n where n^2/64 * (n/B) fits in memory).
#include <cstring>
#include <algorithm>
#include <cstdlib>

typedef unsigned u32;
typedef unsigned long long u64;

#define B 64

static int order[4][100010];
static u64 pref[4][(100000/64 + 4) * (100000/64 + 4)];  // will not use for large n
static int NW;

void count_4d(int n, const unsigned *x[4], unsigned *out) {
    NW = (n + 63) >> 6;
    int M = (n + B - 1) / B;  // number of blocks

    // ---- per dimension: histogram prefix (first_pos) + counting sort order ----
    static int cnt[100002];
    static int first_pos[4][100010];
    for (int d = 0; d < 4; d++) {
        const unsigned *xd = x[d];
        memset(cnt, 0, sizeof(int) * (n + 1));
        for (int i = 0; i < n; i++) cnt[xd[i]]++;
        int s = 0;
        for (int v = 0; v <= (int)n; v++) { int c = cnt[v]; cnt[v] = s; s += c; }
        for (int i = 0; i < n; i++) first_pos[d][i] = cnt[xd[i]];
        // counting sort: order[d]
        memset(cnt, 0, sizeof(int) * (n + 1));
        for (int i = 0; i < n; i++) cnt[xd[i]]++;
        for (int v = 1; v <= (int)n; v++) cnt[v] += cnt[v-1];
        for (int i = n - 1; i >= 0; i--) order[d][--cnt[xd[i]]] = i;
    }

    // ---- build prefix bitsets pref[d][b] ----
    // pref[d] layout: [b][NW]
    for (int d = 0; d < 4; d++) {
        u64 *p = pref[d];
        memset(p, 0, (size_t)(M + 1) * NW * sizeof(u64));
        for (int b = 0; b < M; b++) {
            u64 *cur = p + (size_t)(b + 1) * NW;
            u64 *prev = p + (size_t)b * NW;
            memcpy(cur, prev, (size_t)NW * sizeof(u64));
            int e = (b + 1) * B; if (e > n) e = n;
            for (int t = b * B; t < e; t++) {
                int pt = order[d][t];
                cur[pt >> 6] |= 1ULL << (pt & 63);
            }
        }
    }

    // ---- query each point ----
    static int seen[100010];
    static int stamp = 0;
    for (int j = 0; j < n; j++) {
        int b0 = first_pos[0][j] / B;
        int b1 = first_pos[1][j] / B;
        int b2 = first_pos[2][j] / B;
        int b3 = first_pos[3][j] / B;
        const u64 *p0 = pref[0] + (size_t)b0 * NW;
        const u64 *p1 = pref[1] + (size_t)b1 * NW;
        const u64 *p2 = pref[2] + (size_t)b2 * NW;
        const u64 *p3 = pref[3] + (size_t)b3 * NW;
        u32 c = 0;
        for (int w = 0; w < NW; w++) {
            c += (u32)__builtin_popcountll(p0[w] & p1[w] & p2[w] & p3[w]);
        }
        // partial blocks (points with coord < x[d][j] in the same block)
        stamp++;
        int fp0 = first_pos[0][j];
        for (int t = b0 * B; t < fp0; t++) {
            int q = order[0][t];
            if (x[1][q] < x[1][j] && x[2][q] < x[2][j] && x[3][q] < x[3][j]) {
                if (seen[q] != stamp) { seen[q] = stamp; c++; }
            }
        }
        int fp1 = first_pos[1][j];
        for (int t = b1 * B; t < fp1; t++) {
            int q = order[1][t];
            if (x[0][q] < x[0][j] && x[2][q] < x[2][j] && x[3][q] < x[3][j]) {
                if (seen[q] != stamp) { seen[q] = stamp; c++; }
            }
        }
        int fp2 = first_pos[2][j];
        for (int t = b2 * B; t < fp2; t++) {
            int q = order[2][t];
            if (x[0][q] < x[0][j] && x[1][q] < x[1][j] && x[3][q] < x[3][j]) {
                if (seen[q] != stamp) { seen[q] = stamp; c++; }
            }
        }
        int fp3 = first_pos[3][j];
        for (int t = b3 * B; t < fp3; t++) {
            int q = order[3][t];
            if (x[0][q] < x[0][j] && x[1][q] < x[1][j] && x[2][q] < x[2][j]) {
                if (seen[q] != stamp) { seen[q] = stamp; c++; }
            }
        }
        out[j] = c;
    }
}

评测结果