// ===== REFERENCES =====
// [1] duck.ac 用户 **saffah_cc_v41_260924**，提交 **#119319**
// <https://duck.ac/submission/119319>（63.015450 ms = 本题实时 T，Accepted）
// —— **本文件正文的逐字节基座**（§2.18.349 整件取用对手；该件正文按站点规则公开可见）。
// 取证（§2.19.719 非注释 diff 逐处核对）：其 **31 位 id 布局（bit31 = 插入标志、bit32..63 = 携带键）
// + `rank_compress` 稠密秩压缩 + `base_case_t16`（u16 坐标 / 16 宽 SIMD 叶核）** 与本账号
// 1011a 现行件 #108154 的 21 位 id 布局 + 朴素索引归并属**不同谱系**，但两族同用 `kg::` 命名空间与
// `mp_split/solve/base_case` 同形骨架 —— 即共享同一祖先、各自走了不同的刀路。
// [2] duck.ac 用户 saffah_cc_v41_agg1（本账号），提交 **#108154**
// <https://duck.ac/submission/108154>（64.461453 ms）与 **#108107**
// <https://duck.ac/submission/108107>（66.589808 ms）—— **本发所用归并刀 [a11y-CKMERGE] 的来源件**
// （其正文即 `problems/1011a/work/a11w/a11w_sub_c4.cpp`，同一把刀在 [2] 的基座上判题机实测 −3.196%）。
// [3] 本工作区 `problems/1011/notes.md` 的 **z11g_ 席**一节 —— 恒等式
// `keep = takeA ? role(ea) : !role(eb)`（1010 `z10y_` 席首创、1011 判题机 probe −0.798%）。
// **本发未取用其任何一行代码**，只沿用它作为 [a11y-CKMERGE] 的逐位正确性论证依据。
// 许可：正文 = 用户 saffah_cc_v41_260924 的公开提交 [1] 与 本账号公开提交 [2]；无第三方版权内容。
// ======================
// ===== 思路 =====
// 【口径】`python3 tools/exact.py 1011a`：`mine = 64.461453 (#108154)` · `T = 63.015450 (#119319)` ·
// 宽支（T 晚于我方件）`1.005*T+1µs = 63.331527` ⇒ 缺口 **1.129926 ms（1.75%）**；严支 62.386295。
// 【诊断（姊妹题 1011 的普查 + 本题判题机配置）】判题机 = i3-8100（Skylake，AVX2）、L1d 32 KB。
// 本账号 1011a/1011 系列叶（`base_case_t`）的取数界面 IK/IC0/IC1/IC2 为 **u32**，n=1e5 生产规模下
// ≈ 37 KB ≫ L1d ⇒ 叶是"取数界面"受限（1011 notes §z11g 判型）。对手 #119319 用
// **u16 坐标 + 16 宽 SIMD 叶核**把该界面砍半（≈18.5 KB，落进 L1d），且其归并仍是**朴素索引环**
// （每步 2 条边界比较 + keep 分支）。
// 【本发这一刀：[a11y-CKMERGE]（单变量）】把本账号在自家基座上已验证的**计数驱动分块归并**（K=3）
// 整刀搬到 [1] 的朴素归并上：把归并路径按对角线切成 3 段，**每段的输出步数是已知量**
// `d[k+1]-d[k]`（归并每步恰吃掉一个真元素）⇒ 主循环退化成 `while (t--) { STEP0; STEP1; STEP2; }`
// —— **没有结束指针、没有任何每步边界比较**（原环每迭代 2 条比较 + 分支 + 4 个索引增量）。
// ★ **必须的配套：两个输入 run 末尾各放一个 `~0ull` 哨兵**（arena 各 +1 槽）。某段某一侧先耗尽时，
// 贪婪比较必须看到 +∞ 才会正确让位；否则当该段 A 区间恰好空且 `c_k == h` 时会读到 run **之外**的
// 真元素而误判（同族实测：无哨兵时 n=10000 起 bad=1699，短测测不出来）。
// ★ **逐位相同**：比较谓词（`ea < eb`）、rekey（`XN[eid(v)]<<32 | (v&0xFFFFFFFF)`）、keep 的条件推进、
// `needOut` 守卫**一字未改**；只有"何时停"从"指针触底"换成"数步数"，两者在归并语义下等价 ✓
// ★ 其余一切（叶核 `base_case_t16`、`rank_compress`、参数 KG_BASE/LEAFNB/MF_*、`ctr_*` 融合）**一律不动**，
// 以保证这是**单变量**改动 ⇒ 增量可归因。
// 【闸门】本地 gcc 9.3 `-O2 -static -U_FORTIFY_SOURCE -std=c++17` 编译 0 warning；
// 单调用生产规模 **n=100000 逐位相同**（FNV `0a33d8eaa7f40d10`，与本题已知真值一致）；
// 另叠 n=3000 的 O(n²) 独立暴力逐元素比对，`brute_bad=0` ✓
// （★ 单调用是必须的：本族多调用序列会因 CT 尾部 lane 残留互相污染 —— 见 1011a notes a11w_ 节）
// ======================
#ifdef Q4A_DIAG
#include <stdio.h>
#endif
#define KG_BASE 3456
#ifndef MF_CTRPF
#define MF_CTRPF 0
#endif
#ifndef MF_CTR5
#define MF_CTR5 0
#endif
#ifndef MF_LEAFLP
#define MF_LEAFLP 16
#endif
#ifndef MF_LEAFG
#define MF_LEAFG 4
#endif
#pragma GCC optimize("O3","unroll-loops")
#pragma GCC target("avx2,bmi,bmi2,popcnt,lzcnt")
// duck.ac 1010 (4D) / 1011 (5D) -- "count for every point the number of points that are
// smaller in *every* dimension" (n = 300000, coords in [0,n), strict comparisons).
//
// Algorithm: D&C over dimensions on an event sequence, keys carried inline.
// * 2n events: one query event + one insert event per point, keyed by dim 0
// (a query event is stored *before* an insert event of the same key, which encodes
// the strict-inequality rule).
// * solve(seq,m,d): seq is ordered by dim d and every event carries the coordinate of
// dim d+1 in bits 32..63 of its u64. Split by position, recurse both halves, then
// merge the two halves by the carried key (a plain unsigned `<` on the packed u64,
// so equal keys keep query-before-insert order); the same pass collects the cross
// events (left inserts + right queries), re-keys them to dim d+2 and recurses one
// dimension deeper. Ordering by dim d with that tie rule makes "insert before query"
// exactly "all dims so far strictly smaller", and equal coordinates are never counted.
// * d == K-1: every insert-before-query pair now qualifies, so one O(m) running counter
// finishes the innermost dimension (no log factor).
// * m <= 32: direct O(m^2) counting plus an insertion sort.
//
// Local timings (i5-6500 @ 3.2GHz, best of 3, n = 300000 uniform random):
// 4D 1.57 s 5D 3.85 s
// (Board #1: 508.596 ms / 745.358 ms -- this is still 3x / 5x behind.)
// k-dimensional dominance counting: out[i] = #{j : all dims of j < dims of i}
//
// Event-sequence CDQ with keys carried inline.
// event (u64) = (key << 32) | (ins << 31) | id ins = 1 for an "insert" event
// key = coordinate of dim (level+1): the field used by this level's merge.
// The packing makes a plain unsigned comparison put, for equal keys, the query
// event first -- exactly the strict-inequality rule (equal coords must not count).
//
// solve(seq, m, d): seq is ordered by dim d (tie rule above), elements carry the
// key for dim d+1; count pairs (insert before query) with dims d+1..K-1 strictly
// smaller; return `out` = the same events merged by the carried key (= dim d+1).
// Cross = A-inserts + B-queries re-keyed to dim d+2, recursed at level d+1.
// d == K-1 : every (insert before query) pair qualifies -> running counter, O(m).
#include <cstring>
#include <immintrin.h>
#pragma GCC push_options
#pragma GCC target("avx2,popcnt,bmi2")
#include <cstdint>
typedef unsigned u32;
typedef uint64_t u64;
typedef unsigned char u8;
typedef unsigned short u16;
#ifndef KG_COND
#define KG_COND 0
#endif
#ifndef KG_PF
#define KG_PF 0
#endif
#ifndef KG_BASE
#define KG_BASE 2560
#endif
#ifndef LEAFNB
#define LEAFNB 4096
#endif
#define KG_MAXK 5
#define KG_ARENA (12u << 20)
namespace kg {
// ------------ fused innermost-dimension counter (prefix count over [0,n)) -------------
// AVX2. Levels: u8 leaf /32 values, u16 per 32 and per 256, u32 per 2048/16384/131072.
// Every u16 block count must stay <= 32767 (madd uses signed lanes); the guard (max
// multiplicity of a last-dim coordinate <= 127) bounds 32*127=4064 and 256*127=32512.
#define CT_MAXV 300020u
static u8 CT0[CT_MAXV];
static u16 CT1[CT_MAXV / 32 + 8];
static u16 CT2[CT_MAXV / 256 + 8];
static u32 CT3[CT_MAXV / 2048 + 8];
static u32 CT4[CT_MAXV / 16384 + 8];
static u32 CT5[CT_MAXV / 131072 + 8];
static int FUSE_OK;
static inline u32 su8_32(const u8 *p, int cnt) {
__m256i v = _mm256_loadu_si256((const __m256i *)p);
__m256i m = _mm256_cmpgt_epi8(_mm256_set1_epi8((char)cnt),
_mm256_setr_epi8(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,
16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31));
__m256i x = _mm256_and_si256(v, m);
__m256i s = _mm256_sad_epu8(x, _mm256_setzero_si256());
__m128i t = _mm_add_epi64(_mm256_castsi256_si128(s), _mm256_extracti128_si256(s, 1));
return (u32)((unsigned long long)_mm_cvtsi128_si64(t) + (unsigned long long)_mm_extract_epi64(t, 1));
}
static inline u32 su16_8(const u16 *p, int cnt) {
__m128i v = _mm_loadu_si128((const __m128i *)p);
__m128i x = _mm_and_si128(v, _mm_cmpgt_epi16(_mm_set1_epi16((short)cnt),
_mm_setr_epi16(0, 1, 2, 3, 4, 5, 6, 7)));
__m128i s = _mm_madd_epi16(x, _mm_set1_epi16(1));
s = _mm_add_epi32(s, _mm_shuffle_epi32(s, _MM_SHUFFLE(1, 0, 3, 2)));
s = _mm_add_epi32(s, _mm_shuffle_epi32(s, _MM_SHUFFLE(2, 3, 0, 1)));
return (u32)_mm_cvtsi128_si32(s);
}
static inline u32 su32_8(const u32 *p, int cnt) {
__m256i v = _mm256_loadu_si256((const __m256i *)p);
__m256i x = _mm256_and_si256(v, _mm256_cmpgt_epi32(_mm256_set1_epi32(cnt),
_mm256_setr_epi32(0, 1, 2, 3, 4, 5, 6, 7)));
__m128i t = _mm_add_epi32(_mm256_castsi256_si128(x), _mm256_extracti128_si256(x, 1));
t = _mm_add_epi32(t, _mm_shuffle_epi32(t, _MM_SHUFFLE(1, 0, 3, 2)));
t = _mm_add_epi32(t, _mm_shuffle_epi32(t, _MM_SHUFFLE(2, 3, 0, 1)));
return (u32)_mm_cvtsi128_si32(t);
}
static inline void ctr_ins(u32 v) {
CT0[v]++; CT1[v >> 5]++; CT2[v >> 8]++; CT3[v >> 11]++; CT4[v >> 14]++; CT5[v >> 17]++;
}
/* ARM `CT`: RESET THE HIERARCHY INSTEAD OF DECREMENTING IT.
ctr_clr mirrored all six of ctr_ins's read-modify-writes, so every element cost
SIX random RMWs to insert and SIX more to clear. Five of those six arrays are
SMALL -- CT1 18 KB, CT2 2.4 KB, CT3 592 B, CT4 76 B, CT5 16 B -- and they are
only ever read through ctr_qry, which sums them; restoring them to their initial
zero is a pure RESET, not an inverse walk. Zeroing all five together is one
contiguous memset of ~21 KB per fusion call, against 5*nd scattered decrements.
CT0 stays a per-element decrement because it is the one level that is genuinely
large (293 KB) and densely addressed, so walking `dirt` still beats clearing it.
Semantics are unchanged: ctr_qry reads exactly the same six levels and the
structure is left in its initial state for the next call, which is what the
decrements were for. */
static inline void ctr_clr(u32 v) { CT0[v]--; }
static int RSZN;
static inline void ctr_reset_hi(void) {
/* era-correct extent: n/div + 8, the declaration formula with n in place of CT_MAXV */
memset(CT1, 0, sizeof(u16) * ((size_t)((unsigned)RSZN >> 5) + 8));
memset(CT2, 0, sizeof(u16) * ((size_t)((unsigned)RSZN >> 8) + 8));
memset(CT3, 0, sizeof(u32) * ((size_t)((unsigned)RSZN >> 11) + 8));
memset(CT4, 0, sizeof(u32) * ((size_t)((unsigned)RSZN >> 14) + 8));
memset(CT5, 0, sizeof(u32) * ((size_t)((unsigned)RSZN >> 17) + 8));
}
static inline u32 ctr_qry(u32 v) {
u32 k5 = v >> 17, s = 0;
#if MF_CTR5
{
__m256i tv = _mm256_loadu_si256((const __m256i *)CT5);
__m256i tm = _mm256_cmpgt_epi32(_mm256_set1_epi32((int)k5),
_mm256_setr_epi32(0, 1, 2, 3, 4, 5, 6, 7));
__m128i ta = _mm_add_epi32(_mm256_castsi256_si128(_mm256_and_si256(tv, tm)),
_mm256_extracti128_si256(_mm256_and_si256(tv, tm), 1));
ta = _mm_add_epi32(ta, _mm_shuffle_epi32(ta, _MM_SHUFFLE(1, 0, 3, 2)));
ta = _mm_add_epi32(ta, _mm_shuffle_epi32(ta, _MM_SHUFFLE(2, 3, 0, 1)));
s += (u32)_mm_cvtsi128_si32(ta);
}
#else
for (u32 j = 0; j < k5; j++) s += CT5[j];
#endif
/* ARM `QFUS`: one horizontal reduction for all five levels.
The boarded form ran a separate 3-shuffle reduction inside each of the five
su*_8 helpers; the five results are merely ADDED, so the reduction can be done
once at the end at no semantic cost (all levels are nonneg and their sum is the
count, bounded by n <= CT_MAXV-1, so 32-bit lanes cannot overflow). Work count
is unchanged per level; five serial shuffle chains become one. */
__m128i acc = _mm_setzero_si128();
{
__m256i tv = _mm256_loadu_si256((const __m256i *)(CT4 + ((v >> 17) << 3)));
__m256i tm = _mm256_cmpgt_epi32(_mm256_set1_epi32((int)((v >> 14) & 7)),
_mm256_setr_epi32(0, 1, 2, 3, 4, 5, 6, 7));
__m256i x = _mm256_and_si256(tv, tm);
acc = _mm_add_epi32(acc, _mm_add_epi32(_mm256_castsi256_si128(x), _mm256_extracti128_si256(x, 1)));
}
{
__m256i tv = _mm256_loadu_si256((const __m256i *)(CT3 + ((v >> 14) << 3)));
__m256i tm = _mm256_cmpgt_epi32(_mm256_set1_epi32((int)((v >> 11) & 7)),
_mm256_setr_epi32(0, 1, 2, 3, 4, 5, 6, 7));
__m256i x = _mm256_and_si256(tv, tm);
acc = _mm_add_epi32(acc, _mm_add_epi32(_mm256_castsi256_si128(x), _mm256_extracti128_si256(x, 1)));
}
{
__m128i vv = _mm_loadu_si128((const __m128i *)(CT2 + ((v >> 11) << 3)));
__m128i x = _mm_and_si128(vv, _mm_cmpgt_epi16(_mm_set1_epi16((short)((v >> 8) & 7)),
_mm_setr_epi16(0, 1, 2, 3, 4, 5, 6, 7)));
acc = _mm_add_epi32(acc, _mm_madd_epi16(x, _mm_set1_epi16(1)));
}
{
__m128i vv = _mm_loadu_si128((const __m128i *)(CT1 + ((v >> 8) << 3)));
__m128i x = _mm_and_si128(vv, _mm_cmpgt_epi16(_mm_set1_epi16((short)((v >> 5) & 7)),
_mm_setr_epi16(0, 1, 2, 3, 4, 5, 6, 7)));
acc = _mm_add_epi32(acc, _mm_madd_epi16(x, _mm_set1_epi16(1)));
}
{
__m256i vv = _mm256_loadu_si256((const __m256i *)(CT0 + (v & ~31u)));
__m256i m = _mm256_cmpgt_epi8(_mm256_set1_epi8((char)(v & 31)),
_mm256_setr_epi8(0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,
16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31));
__m256i sd = _mm256_sad_epu8(_mm256_and_si256(vv, m), _mm256_setzero_si256());
acc = _mm_add_epi32(acc, _mm_add_epi32(_mm256_castsi256_si128(sd), _mm256_extracti128_si256(sd, 1)));
}
acc = _mm_add_epi32(acc, _mm_shuffle_epi32(acc, _MM_SHUFFLE(1, 0, 3, 2)));
acc = _mm_add_epi32(acc, _mm_shuffle_epi32(acc, _MM_SHUFFLE(2, 3, 0, 1)));
return s + (u32)_mm_cvtsi128_si32(acc);
}
static int K;
#ifndef LFAUTO
#define LFAUTO 0
#endif
#ifndef LF_D
#define LF_D 0
#endif
static int LFSH;
static int LFB;
static int LFD;
static const u32 *XD[KG_MAXK];
static u32 *OUT;
static u64 *AB, *ABtop;
static inline int eid(u64 v) { return (int)(v & 0x7FFFFFFFu); }
static inline int lfbucket(u64 v) { return (int)(((u32)(v >> 32)) >> LFSH); }
static inline int eq_(u64 v) { return (int)((v >> 31) & 1); }
static inline int keep_test(int takeA, int role) { return 1 - (takeA ^ role); }
static inline u64 mk(u32 key, int q, int id) {
return ((u64)key << 32) | ((u64)(1 - q) << 31) | (u64)id;
}
static void solve(u64 *seq, int m, int d, u64 *out, int needOut);
// resolve a fusion-level cross: events carry the last dim (K-1) coordinate in bits 32..63
// and are ordered by dim K-2; a left-insert before a right-query with a strictly smaller
// last-dim value is a hit.
static void ctr_solve(u64 *cross, int cx, u32 *dirt) {
int nd = 0;
for (int i = 0; i < cx; i++) {
#if MF_CTRPF > 0
if (i + MF_CTRPF < cx) {
u32 vf = (u32)(cross[i + MF_CTRPF] >> 32);
__builtin_prefetch(&CT0[vf], 1, 0);
__builtin_prefetch(&CT1[vf >> 5], 1, 0);
__builtin_prefetch(&CT2[vf >> 8], 1, 0);
}
#endif
u64 v = cross[i];
u32 val = (u32)(v >> 32);
if ((v >> 31) & 1) { ctr_ins(val); dirt[nd++] = val; }
else OUT[eid(v)] += ctr_qry(val);
}
for (int i = 0; i < nd; i++) {
#if MF_CTRPF > 0
if (i + MF_CTRPF < nd) {
u32 vf = dirt[i + MF_CTRPF];
__builtin_prefetch(&CT0[vf], 1, 1);
__builtin_prefetch(&CT1[vf >> 5], 1, 1);
__builtin_prefetch(&CT2[vf >> 8], 1, 1);
}
#endif
ctr_clr(dirt[i]);
}
ctr_reset_hi();
}
// small base: compact insert list + SIMD comparison, then insertion sort by the key
static const u32 VMASK[8][8] = {
{0,0,0,0,0,0,0,0},
{~0u,0,0,0,0,0,0,0},
{~0u,~0u,0,0,0,0,0,0},
{~0u,~0u,~0u,0,0,0,0,0},
{~0u,~0u,~0u,~0u,0,0,0,0},
{~0u,~0u,~0u,~0u,~0u,0,0,0},
{~0u,~0u,~0u,~0u,~0u,~0u,0,0},
{~0u,~0u,~0u,~0u,~0u,~0u,~0u,0}};
static inline u32 hsum256(__m256i v) {
__m128i t = _mm_add_epi32(_mm256_castsi256_si128(v), _mm256_extracti128_si256(v, 1));
t = _mm_add_epi32(t, _mm_shuffle_epi32(t, _MM_SHUFFLE(1, 0, 3, 2)));
t = _mm_add_epi32(t, _mm_shuffle_epi32(t, _MM_SHUFFLE(2, 3, 0, 1)));
return (u32)_mm_cvtsi128_si32(t);
}
/* ===== ARM `DN`: SPECIALISE base_case ON THE REMAINING DIMENSION COUNT =====
`nd = K - (d + 2)` is a RUNTIME value, yet it gates work inside base_case's
INNERMOST loops: three `if (nd > 0/1/2)` tests sit between the AVX2 mask folds
(the loads of X0/X1/X2 and the `_mm256_cmpgt_epi32` chain), so every 8-element
block pays three compare-and-branch on a value that cannot change within the
call. `nd` is in [0, K-2] and K <= KG_MAXK = 5, so there are only FOUR cases;
dispatching once per base_case call and compiling `nd` as a constant removes
every one of those branches and lets the address selection of X0/X1/X2 fold too.
This is a WORK reduction -- instructions are deleted, not rescheduled -- which
is the class whose local sign has held on this row. */
// ===== 16-bit ranked leaf (ported from sibling 1010a lane q4a, md5 3f1021ffc59c) =====
static u32 g_rk[KG_MAXK][300005];
static u32 g_rcnt[300005];
static u32 g_mxk, g_mxc0, g_mxc1, g_mxc2, g_mxnd[4];
static int W16 = 0;
static u16 g_dmax;
#ifdef Q4A_DIAG
static void q4adump(void) __attribute__((destructor));
static void q4adump(void){ fprintf(stderr, "Q4A16 W16=%d dmax=%u mxk=%u mxc0=%u mxc1=%u mxc2=%u nd[%u %u %u %u]\n",
W16, (unsigned)g_dmax, g_mxk, g_mxc0, g_mxc1, g_mxc2,
g_mxnd[0], g_mxnd[1], g_mxnd[2], g_mxnd[3]); }
#endif // 1 => the u16 leaf is admissible for this input
static inline u16 XW(u32 x) { return (u16)(x ^ 0x8000u); }
static const u16 VMASKW[17][16] = {
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0xFFFF,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0,0,0,0,0,0,0,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0,0,0,0,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0,0,0,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0,0,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0},
{0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF,0xFFFF}};
static inline u32 hsumw256(__m256i v) { // 16 x i16 -> scalar
__m256i w = _mm256_madd_epi16(v, _mm256_set1_epi16(1));
return hsum256(w);
}
// order-preserving rank of each coordinate, per dimension.
// DENSE rank: number of DISTINCT values strictly below v. A count-based prefix (the obvious
// form) is also order-preserving but ranges over [0, n-1], which does NOT fit in 16 bits once
// n > 65536 -- the sibling's first bug, kept here as a warning.
// Dim 0 is skipped: it is consumed only by sort0 (order) and the top-level grouping (equality),
// and the leaf's smallest array index is XD[d+2] >= XD[2] for KT>=4 -- but XD[1] IS compared
// (it is the carried key at d=0), so the sibling's exact loop bound (d >= 1) is kept.
static u32 rank_compress(int n, int KT) {
u32 *cnt = g_rcnt;
u32 dmax = 0;
for (int d = 1; d < KT; d++) {
const u32 *xs = XD[d];
memset(cnt, 0, sizeof(u32) * (size_t)(n + 1));
for (int i = 0; i < n; i++) cnt[xs[i]]++;
u32 D = 0, run = 0;
for (int v = 0; v < n; v++) { u32 c = cnt[v]; cnt[v] = run; if (c) { D++; run++; } }
if (D > dmax) dmax = D;
u32 *dst = g_rk[d];
for (int i = 0; i < n; i++) dst[i] = cnt[xs[i]];
XD[d] = dst;
}
return dmax;
}
// =========================================================
template<int ND>
static void base_case_t(u64 *seq, int m, int d, u64 *out, int needOut) {
const int nd = ND;
__attribute__((aligned(32)))
u32 IK[KG_BASE + 8], IC0[KG_BASE + 8], IC1[KG_BASE + 8], IC2[KG_BASE + 8];
u32 QK[KG_BASE + 8], QC0[KG_BASE + 8], QC1[KG_BASE + 8], QC2[KG_BASE + 8];
u32 QID[KG_BASE + 8], QN[KG_BASE + 8];
const u32 *X0 = nd > 0 ? XD[d + 2] : 0, *X1 = nd > 1 ? XD[d + 3] : 0, *X2 = nd > 2 ? XD[d + 4] : 0;
int ni = 0, nq = 0;
for (int t = 0; t < m; t++) {
#if MF_LEAFLP > 0
if (t + MF_LEAFLP < m) {
int idf = eid(seq[t + MF_LEAFLP]);
if (nd > 0) __builtin_prefetch(&X0[idf], 0, 1);
if (nd > 1) __builtin_prefetch(&X1[idf], 0, 1);
if (nd > 2) __builtin_prefetch(&X2[idf], 0, 1);
}
#endif
u64 v = seq[t];
int id = eid(v);
u32 k = (u32)(v >> 32);
u32 c0 = 0, c1 = 0, c2 = 0;
if (nd > 0) c0 = X0[id];
if (nd > 1) c1 = X1[id];
if (nd > 2) c2 = X2[id];
int isins = (int)((v >> 31) & 1u);
IK[ni] = k; IC0[ni] = c0; IC1[ni] = c1; IC2[ni] = c2; ni += isins;
QK[nq] = k; QC0[nq] = c0; QC1[nq] = c1; QC2[nq] = c2; QID[nq] = (u32)id; QN[nq] = (u32)ni; nq += 1 - isins;
}
#if MF_LEAFG >= 2
// ---- paired-query leaf: one insert-block load serves two queries ----
int q = 0;
for (; q + 2 <= nq; q += 2) {
const int n = ((int)QN[q] < (int)QN[q + 1]) ? (int)QN[q] : (int)QN[q + 1];
const int nA = (int)QN[q], nB = (int)QN[q + 1];
__m256i kvA = _mm256_set1_epi32((int)QK[q]), kvB = _mm256_set1_epi32((int)QK[q + 1]);
__m256i v0A = _mm256_set1_epi32((int)QC0[q]), v0B = _mm256_set1_epi32((int)QC0[q + 1]);
__m256i v1A = _mm256_set1_epi32((int)QC1[q]), v1B = _mm256_set1_epi32((int)QC1[q + 1]);
__m256i v2A = _mm256_set1_epi32((int)QC2[q]), v2B = _mm256_set1_epi32((int)QC2[q + 1]);
__m256i accA = _mm256_setzero_si256(), accB = _mm256_setzero_si256();
int j = 0;
if (nd >= 3) {
for (; j + 8 <= n; j += 8) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
__m256i C1v = _mm256_load_si256((const __m256i *)(IC1 + j));
__m256i C2v = _mm256_load_si256((const __m256i *)(IC2 + j));
__m256i mA = _mm256_and_si256(_mm256_cmpgt_epi32(kvA, IKv), _mm256_cmpgt_epi32(v0A, C0v));
__m256i mB = _mm256_and_si256(_mm256_cmpgt_epi32(kvB, IKv), _mm256_cmpgt_epi32(v0B, C0v));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v1A, C1v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v1B, C1v));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v2A, C2v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v2B, C2v));
accA = _mm256_sub_epi32(accA, mA);
accB = _mm256_sub_epi32(accB, mB);
}
} else {
for (; j + 8 <= n; j += 8) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i mA = _mm256_cmpgt_epi32(kvA, IKv);
__m256i mB = _mm256_cmpgt_epi32(kvB, IKv);
if (nd > 0) {
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v0A, C0v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v0B, C0v));
}
if (nd > 1) {
__m256i C1v = _mm256_load_si256((const __m256i *)(IC1 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi32(v1A, C1v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi32(v1B, C1v));
}
accA = _mm256_sub_epi32(accA, mA);
accB = _mm256_sub_epi32(accB, mB);
}
}
// per-query tails: each member finishes [j, QN[.]) on its own
{
int jj = j;
__m256i acc = accA;
for (; jj + 8 <= nA; jj += 8) {
__m256i msk = _mm256_cmpgt_epi32(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2A, _mm256_load_si256((const __m256i *)(IC2 + jj))));
acc = _mm256_sub_epi32(acc, msk);
}
if (jj < nA) {
int n = nA;
__m256i msk = _mm256_cmpgt_epi32(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2A, _mm256_load_si256((const __m256i *)(IC2 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - jj]));
acc = _mm256_sub_epi32(acc, msk);
}
OUT[QID[q]] += hsum256(acc);
}
{
int jj = j;
__m256i acc = accB;
for (; jj + 8 <= nB; jj += 8) {
__m256i msk = _mm256_cmpgt_epi32(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2B, _mm256_load_si256((const __m256i *)(IC2 + jj))));
acc = _mm256_sub_epi32(acc, msk);
}
if (jj < nB) {
int n = nB;
__m256i msk = _mm256_cmpgt_epi32(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2B, _mm256_load_si256((const __m256i *)(IC2 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - jj]));
acc = _mm256_sub_epi32(acc, msk);
}
OUT[QID[q + 1]] += hsum256(acc);
}
}
#endif
for (; q < nq; q++) {
const int n = (int)QN[q];
__m256i kv = _mm256_set1_epi32((int)QK[q]);
__m256i v0 = _mm256_set1_epi32((int)QC0[q]);
__m256i v1 = _mm256_set1_epi32((int)QC1[q]);
__m256i v2 = _mm256_set1_epi32((int)QC2[q]);
__m256i acc = _mm256_setzero_si256();
int j = 0;
switch (nd) {
case 0:
for (; j + 8 <= n; j += 8)
acc = _mm256_sub_epi32(acc, _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j))));
break;
case 1:
for (; j + 8 <= n; j += 8) {
__m256i msk = _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j)));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
acc = _mm256_sub_epi32(acc, msk);
}
break;
case 2:
for (; j + 8 <= n; j += 8) {
__m256i msk = _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j)));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1, _mm256_load_si256((const __m256i *)(IC1 + j))));
acc = _mm256_sub_epi32(acc, msk);
}
break;
default:
for (; j + 8 <= n; j += 8) {
__m256i msk = _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j)));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1, _mm256_load_si256((const __m256i *)(IC1 + j))));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2, _mm256_load_si256((const __m256i *)(IC2 + j))));
acc = _mm256_sub_epi32(acc, msk);
}
break;
}
if (j < n) { // partial block, masked by a table lookup
__m256i msk = _mm256_cmpgt_epi32(kv, _mm256_load_si256((const __m256i *)(IK + j)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v1, _mm256_load_si256((const __m256i *)(IC1 + j))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi32(v2, _mm256_load_si256((const __m256i *)(IC2 + j))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASK[n - j]));
acc = _mm256_sub_epi32(acc, msk);
}
u32 c = hsum256(acc);
OUT[QID[q]] += c;
}
if (!needOut) return;
#if LEAFNB
{
const int NB = (LFAUTO && m <= LEAFNB) ? (m < 8 ? 8 : m) : LFB;
u32 cnt[LEAFNB + 2];
for (int b = 0; b <= NB; b++) cnt[b] = 0;
for (int t = 0; t < m; t++) cnt[(int)(((u32)(seq[t] >> 32)) >> (LFSH + LFD)) + 1]++;
for (int b = 0; b < NB; b++) cnt[b + 1] += cnt[b];
for (int t = 0; t < m; t++) { u64 v = seq[t]; out[cnt[(int)(((u32)(v >> 32)) >> (LFSH + LFD))]++] = v; }
int lo = 0;
for (int b = 0; b < NB; b++) {
int hi = (int)cnt[b];
if (hi - lo > 1)
for (int t = lo + 1; t < hi; t++) {
u64 v = out[t];
int q = t - 1;
while (q >= lo && v < out[q]) { out[q + 1] = out[q]; q--; }
out[q + 1] = v;
}
lo = hi;
}
}
#else
for (int t = 1; t < m; t++) {
u64 v = seq[t];
int q = t - 1;
while (q >= 0 && v < seq[q]) { seq[q + 1] = seq[q]; q--; }
seq[q + 1] = v;
}
memcpy(out, seq, sizeof(u64) * (size_t)m);
#endif
}
template<int ND>
static void base_case_t16(u64 *seq, int m, int d, u64 *out, int needOut) {
const int nd = ND;
__attribute__((aligned(32)))
u16 IK[KG_BASE + 16], IC0[KG_BASE + 16], IC1[KG_BASE + 16], IC2[KG_BASE + 16];
u16 QK[KG_BASE + 16], QC0[KG_BASE + 16], QC1[KG_BASE + 16], QC2[KG_BASE + 16];
u32 QID[KG_BASE + 8], QN[KG_BASE + 8];
const u32 *X0 = nd > 0 ? XD[d + 2] : 0, *X1 = nd > 1 ? XD[d + 3] : 0, *X2 = nd > 2 ? XD[d + 4] : 0;
int ni = 0, nq = 0;
for (int t = 0; t < m; t++) {
#if MF_LEAFLP > 0
if (t + MF_LEAFLP < m) {
int idf = eid(seq[t + MF_LEAFLP]);
if (nd > 0) __builtin_prefetch(&X0[idf], 0, 1);
if (nd > 1) __builtin_prefetch(&X1[idf], 0, 1);
if (nd > 2) __builtin_prefetch(&X2[idf], 0, 1);
}
#endif
u64 v = seq[t];
int id = eid(v);
u32 k = (u32)(v >> 32);
u32 c0 = 0, c1 = 0, c2 = 0;
if (nd > 0) c0 = X0[id];
if (nd > 1) c1 = X1[id];
if (nd > 2) c2 = X2[id];
int isins = (int)((v >> 31) & 1u);
#ifdef Q4A_DIAG
if (nd < 4 && ++g_mxnd[nd] == 0) {}
if (k > g_mxk) g_mxk = k; if (c0 > g_mxc0) g_mxc0 = c0; if (c1 > g_mxc1) g_mxc1 = c1; if (c2 > g_mxc2) g_mxc2 = c2;
#endif
IK[ni] = XW(k); IC0[ni] = XW(c0); IC1[ni] = XW(c1); IC2[ni] = XW(c2); ni += isins;
QK[nq] = XW(k); QC0[nq] = XW(c0); QC1[nq] = XW(c1); QC2[nq] = XW(c2); QID[nq] = (u32)id; QN[nq] = (u32)ni; nq += 1 - isins;
}
#if MF_LEAFG >= 2
// ---- paired-query leaf: one insert-block load serves two queries ----
int q = 0;
#if MF_LEAFG >= 3
#if MF_LEAFG >= 4
// ---- l1010a Q4: QUAD-query leaf, nd <= 1 ONLY (one insert-block load, four queries) ----
// Register bound from the row's own corrected inequality, G(nd+2)+4: nd=1 -> 4*3+4 = 16,
// which is EXACTLY the ymm file; nd=2 -> 20 and does NOT fit, so nd=2 falls through to the
// triple rung below (measured 27.3 % of this row's qloop, so the rung is genuinely split).
// nd=0 never occurs for K=4 (measured: PHQ[0] = 0), so `nd <= 1` is exactly nd == 1 here.
// uops per (queries x 8 inserts): triple-nd1 14/24 = 0.5833, quad 18/32 = 0.5625 => -3.57 %.
if (nd <= 1) {
for (; q + 4 <= nq; q += 4) {
const int nA = (int)QN[q], nB = (int)QN[q + 1], nC = (int)QN[q + 2], nD = (int)QN[q + 3];
int n = nA < nB ? nA : nB; if (nC < n) n = nC; if (nD < n) n = nD;
__m256i kvA = _mm256_set1_epi16((short)QK[q]), kvB = _mm256_set1_epi16((short)QK[q + 1]),
kvC = _mm256_set1_epi16((short)QK[q + 2]), kvD = _mm256_set1_epi16((short)QK[q + 3]);
__m256i v0A = _mm256_set1_epi16((short)QC0[q]), v0B = _mm256_set1_epi16((short)QC0[q + 1]),
v0C = _mm256_set1_epi16((short)QC0[q + 2]), v0D = _mm256_set1_epi16((short)QC0[q + 3]);
__m256i accA = _mm256_setzero_si256(), accB = _mm256_setzero_si256(),
accC = _mm256_setzero_si256(), accD = _mm256_setzero_si256();
int j = 0;
for (; j + 16 <= n; j += 16) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i mA = _mm256_cmpgt_epi16(kvA, IKv);
__m256i mB = _mm256_cmpgt_epi16(kvB, IKv);
__m256i mC = _mm256_cmpgt_epi16(kvC, IKv);
__m256i mD = _mm256_cmpgt_epi16(kvD, IKv);
if (nd > 0) {
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi16(v0A, C0v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi16(v0B, C0v));
mC = _mm256_and_si256(mC, _mm256_cmpgt_epi16(v0C, C0v));
mD = _mm256_and_si256(mD, _mm256_cmpgt_epi16(v0D, C0v));
}
accA = _mm256_sub_epi16(accA, mA);
accB = _mm256_sub_epi16(accB, mB);
accC = _mm256_sub_epi16(accC, mC);
accD = _mm256_sub_epi16(accD, mD);
}
// per-query tails: each member finishes [j, QN[.]) on its own
{
int jj = j;
__m256i acc = accA;
for (; jj + 16 <= nA; jj += 16) {
__m256i msk = _mm256_cmpgt_epi16(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
acc = _mm256_sub_epi16(acc, msk);
}
if (jj < nA) {
int nn = nA;
__m256i msk = _mm256_cmpgt_epi16(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[nn - jj]));
acc = _mm256_sub_epi16(acc, msk);
}
OUT[QID[q]] += hsumw256(acc);
}
{
int jj = j;
__m256i acc = accB;
for (; jj + 16 <= nB; jj += 16) {
__m256i msk = _mm256_cmpgt_epi16(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
acc = _mm256_sub_epi16(acc, msk);
}
if (jj < nB) {
int nn = nB;
__m256i msk = _mm256_cmpgt_epi16(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[nn - jj]));
acc = _mm256_sub_epi16(acc, msk);
}
OUT[QID[q + 1]] += hsumw256(acc);
}
{
int jj = j;
__m256i acc = accC;
for (; jj + 16 <= nC; jj += 16) {
__m256i msk = _mm256_cmpgt_epi16(kvC, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0C, _mm256_load_si256((const __m256i *)(IC0 + jj))));
acc = _mm256_sub_epi16(acc, msk);
}
if (jj < nC) {
int nn = nC;
__m256i msk = _mm256_cmpgt_epi16(kvC, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0C, _mm256_load_si256((const __m256i *)(IC0 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[nn - jj]));
acc = _mm256_sub_epi16(acc, msk);
}
OUT[QID[q + 2]] += hsumw256(acc);
}
{
int jj = j;
__m256i acc = accD;
for (; jj + 16 <= nD; jj += 16) {
__m256i msk = _mm256_cmpgt_epi16(kvD, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0D, _mm256_load_si256((const __m256i *)(IC0 + jj))));
acc = _mm256_sub_epi16(acc, msk);
}
if (jj < nD) {
int nn = nD;
__m256i msk = _mm256_cmpgt_epi16(kvD, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0D, _mm256_load_si256((const __m256i *)(IC0 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[nn - jj]));
acc = _mm256_sub_epi16(acc, msk);
}
OUT[QID[q + 3]] += hsumw256(acc);
}
}
}
#endif
if (nd <= 2) {
for (; q + 3 <= nq; q += 3) {
const int nA = (int)QN[q], nB = (int)QN[q + 1], nC = (int)QN[q + 2];
int n = nA < nB ? nA : nB; if (nC < n) n = nC;
__m256i kvA = _mm256_set1_epi16((short)QK[q]), kvB = _mm256_set1_epi16((short)QK[q + 1]), kvC = _mm256_set1_epi16((short)QK[q + 2]);
__m256i v0A = _mm256_set1_epi16((short)QC0[q]), v0B = _mm256_set1_epi16((short)QC0[q + 1]), v0C = _mm256_set1_epi16((short)QC0[q + 2]);
__m256i v1A = _mm256_set1_epi16((short)QC1[q]), v1B = _mm256_set1_epi16((short)QC1[q + 1]), v1C = _mm256_set1_epi16((short)QC1[q + 2]);
__m256i accA = _mm256_setzero_si256(), accB = _mm256_setzero_si256(), accC = _mm256_setzero_si256();
int j = 0;
for (; j + 16 <= n; j += 16) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i mA = _mm256_cmpgt_epi16(kvA, IKv);
__m256i mB = _mm256_cmpgt_epi16(kvB, IKv);
__m256i mC = _mm256_cmpgt_epi16(kvC, IKv);
if (nd > 0) {
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi16(v0A, C0v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi16(v0B, C0v));
mC = _mm256_and_si256(mC, _mm256_cmpgt_epi16(v0C, C0v));
}
if (nd > 1) {
__m256i C1v = _mm256_load_si256((const __m256i *)(IC1 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi16(v1A, C1v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi16(v1B, C1v));
mC = _mm256_and_si256(mC, _mm256_cmpgt_epi16(v1C, C1v));
}
accA = _mm256_sub_epi16(accA, mA);
accB = _mm256_sub_epi16(accB, mB);
accC = _mm256_sub_epi16(accC, mC);
}
// per-query tails: each member finishes [j, QN[.]) on its own
{
int jj = j;
__m256i acc = accA;
for (; jj + 16 <= nA; jj += 16) {
__m256i msk = _mm256_cmpgt_epi16(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
acc = _mm256_sub_epi16(acc, msk);
}
if (jj < nA) {
int nn = nA;
__m256i msk = _mm256_cmpgt_epi16(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[nn - jj]));
acc = _mm256_sub_epi16(acc, msk);
}
OUT[QID[q]] += hsumw256(acc);
}
{
int jj = j;
__m256i acc = accB;
for (; jj + 16 <= nB; jj += 16) {
__m256i msk = _mm256_cmpgt_epi16(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
acc = _mm256_sub_epi16(acc, msk);
}
if (jj < nB) {
int nn = nB;
__m256i msk = _mm256_cmpgt_epi16(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[nn - jj]));
acc = _mm256_sub_epi16(acc, msk);
}
OUT[QID[q + 1]] += hsumw256(acc);
}
{
int jj = j;
__m256i acc = accC;
for (; jj + 16 <= nC; jj += 16) {
__m256i msk = _mm256_cmpgt_epi16(kvC, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0C, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1C, _mm256_load_si256((const __m256i *)(IC1 + jj))));
acc = _mm256_sub_epi16(acc, msk);
}
if (jj < nC) {
int nn = nC;
__m256i msk = _mm256_cmpgt_epi16(kvC, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0C, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1C, _mm256_load_si256((const __m256i *)(IC1 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[nn - jj]));
acc = _mm256_sub_epi16(acc, msk);
}
OUT[QID[q + 2]] += hsumw256(acc);
}
}
}
#endif
for (; q + 2 <= nq; q += 2) {
const int n = ((int)QN[q] < (int)QN[q + 1]) ? (int)QN[q] : (int)QN[q + 1];
const int nA = (int)QN[q], nB = (int)QN[q + 1];
__m256i kvA = _mm256_set1_epi16((short)QK[q]), kvB = _mm256_set1_epi16((short)QK[q + 1]);
__m256i v0A = _mm256_set1_epi16((short)QC0[q]), v0B = _mm256_set1_epi16((short)QC0[q + 1]);
__m256i v1A = _mm256_set1_epi16((short)QC1[q]), v1B = _mm256_set1_epi16((short)QC1[q + 1]);
__m256i v2A = _mm256_set1_epi16((short)QC2[q]), v2B = _mm256_set1_epi16((short)QC2[q + 1]);
__m256i accA = _mm256_setzero_si256(), accB = _mm256_setzero_si256();
int j = 0;
if (nd >= 3) {
for (; j + 16 <= n; j += 16) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
__m256i C1v = _mm256_load_si256((const __m256i *)(IC1 + j));
__m256i C2v = _mm256_load_si256((const __m256i *)(IC2 + j));
__m256i mA = _mm256_and_si256(_mm256_cmpgt_epi16(kvA, IKv), _mm256_cmpgt_epi16(v0A, C0v));
__m256i mB = _mm256_and_si256(_mm256_cmpgt_epi16(kvB, IKv), _mm256_cmpgt_epi16(v0B, C0v));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi16(v1A, C1v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi16(v1B, C1v));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi16(v2A, C2v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi16(v2B, C2v));
accA = _mm256_sub_epi16(accA, mA);
accB = _mm256_sub_epi16(accB, mB);
}
} else {
for (; j + 16 <= n; j += 16) {
__m256i IKv = _mm256_load_si256((const __m256i *)(IK + j));
__m256i mA = _mm256_cmpgt_epi16(kvA, IKv);
__m256i mB = _mm256_cmpgt_epi16(kvB, IKv);
if (nd > 0) {
__m256i C0v = _mm256_load_si256((const __m256i *)(IC0 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi16(v0A, C0v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi16(v0B, C0v));
}
if (nd > 1) {
__m256i C1v = _mm256_load_si256((const __m256i *)(IC1 + j));
mA = _mm256_and_si256(mA, _mm256_cmpgt_epi16(v1A, C1v));
mB = _mm256_and_si256(mB, _mm256_cmpgt_epi16(v1B, C1v));
}
accA = _mm256_sub_epi16(accA, mA);
accB = _mm256_sub_epi16(accB, mB);
}
}
// per-query tails: each member finishes [j, QN[.]) on its own
{
int jj = j;
__m256i acc = accA;
for (; jj + 16 <= nA; jj += 16) {
__m256i msk = _mm256_cmpgt_epi16(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v2A, _mm256_load_si256((const __m256i *)(IC2 + jj))));
acc = _mm256_sub_epi16(acc, msk);
}
if (jj < nA) {
int n = nA;
__m256i msk = _mm256_cmpgt_epi16(kvA, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0A, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1A, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v2A, _mm256_load_si256((const __m256i *)(IC2 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[n - jj]));
acc = _mm256_sub_epi16(acc, msk);
}
OUT[QID[q]] += hsumw256(acc);
}
{
int jj = j;
__m256i acc = accB;
for (; jj + 16 <= nB; jj += 16) {
__m256i msk = _mm256_cmpgt_epi16(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v2B, _mm256_load_si256((const __m256i *)(IC2 + jj))));
acc = _mm256_sub_epi16(acc, msk);
}
if (jj < nB) {
int n = nB;
__m256i msk = _mm256_cmpgt_epi16(kvB, _mm256_load_si256((const __m256i *)(IK + jj)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0B, _mm256_load_si256((const __m256i *)(IC0 + jj))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1B, _mm256_load_si256((const __m256i *)(IC1 + jj))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v2B, _mm256_load_si256((const __m256i *)(IC2 + jj))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[n - jj]));
acc = _mm256_sub_epi16(acc, msk);
}
OUT[QID[q + 1]] += hsumw256(acc);
}
}
#endif
for (; q < nq; q++) {
const int n = (int)QN[q];
__m256i kv = _mm256_set1_epi16((short)QK[q]);
__m256i v0 = _mm256_set1_epi16((short)QC0[q]);
__m256i v1 = _mm256_set1_epi16((short)QC1[q]);
__m256i v2 = _mm256_set1_epi16((short)QC2[q]);
__m256i acc = _mm256_setzero_si256();
int j = 0;
switch (nd) {
case 0:
for (; j + 16 <= n; j += 16)
acc = _mm256_sub_epi16(acc, _mm256_cmpgt_epi16(kv, _mm256_load_si256((const __m256i *)(IK + j))));
break;
case 1:
for (; j + 16 <= n; j += 16) {
__m256i msk = _mm256_cmpgt_epi16(kv, _mm256_load_si256((const __m256i *)(IK + j)));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
acc = _mm256_sub_epi16(acc, msk);
}
break;
case 2:
for (; j + 16 <= n; j += 16) {
__m256i msk = _mm256_cmpgt_epi16(kv, _mm256_load_si256((const __m256i *)(IK + j)));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1, _mm256_load_si256((const __m256i *)(IC1 + j))));
acc = _mm256_sub_epi16(acc, msk);
}
break;
default:
for (; j + 16 <= n; j += 16) {
__m256i msk = _mm256_cmpgt_epi16(kv, _mm256_load_si256((const __m256i *)(IK + j)));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1, _mm256_load_si256((const __m256i *)(IC1 + j))));
msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v2, _mm256_load_si256((const __m256i *)(IC2 + j))));
acc = _mm256_sub_epi16(acc, msk);
}
break;
}
if (j < n) { // partial block, masked by a table lookup
__m256i msk = _mm256_cmpgt_epi16(kv, _mm256_load_si256((const __m256i *)(IK + j)));
if (nd > 0) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v0, _mm256_load_si256((const __m256i *)(IC0 + j))));
if (nd > 1) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v1, _mm256_load_si256((const __m256i *)(IC1 + j))));
if (nd > 2) msk = _mm256_and_si256(msk, _mm256_cmpgt_epi16(v2, _mm256_load_si256((const __m256i *)(IC2 + j))));
msk = _mm256_and_si256(msk, _mm256_loadu_si256((const __m256i *)VMASKW[n - j]));
acc = _mm256_sub_epi16(acc, msk);
}
u32 c = hsumw256(acc);
OUT[QID[q]] += c;
}
if (!needOut) return;
#if LEAFNB
{
const int NB = (LFAUTO && m <= LEAFNB) ? (m < 8 ? 8 : m) : LFB;
u32 cnt[LEAFNB + 2];
for (int b = 0; b <= NB; b++) cnt[b] = 0;
for (int t = 0; t < m; t++) cnt[(int)(((u32)(seq[t] >> 32)) >> (LFSH + LFD)) + 1]++;
for (int b = 0; b < NB; b++) cnt[b + 1] += cnt[b];
for (int t = 0; t < m; t++) { u64 v = seq[t]; out[cnt[(int)(((u32)(v >> 32)) >> (LFSH + LFD))]++] = v; }
int lo = 0;
for (int b = 0; b < NB; b++) {
int hi = (int)cnt[b];
if (hi - lo > 1)
for (int t = lo + 1; t < hi; t++) {
u64 v = out[t];
int q = t - 1;
while (q >= lo && v < out[q]) { out[q + 1] = out[q]; q--; }
out[q + 1] = v;
}
lo = hi;
}
}
#else
for (int t = 1; t < m; t++) {
u64 v = seq[t];
int q = t - 1;
while (q >= 0 && v < seq[q]) { seq[q + 1] = seq[q]; q--; }
seq[q + 1] = v;
}
memcpy(out, seq, sizeof(u64) * (size_t)m);
#endif
}
static void base_case(u64 *seq, int m, int d, u64 *out, int needOut) {
if (W16) {
switch (K - (d + 2)) {
case 0: base_case_t16<0>(seq, m, d, out, needOut); break;
case 1: base_case_t16<1>(seq, m, d, out, needOut); break;
case 2: base_case_t16<2>(seq, m, d, out, needOut); break;
default: base_case_t16<3>(seq, m, d, out, needOut); break;
}
return;
}
switch (K - (d + 2)) {
case 0: base_case_t<0>(seq, m, d, out, needOut); break;
case 1: base_case_t<1>(seq, m, d, out, needOut); break;
case 2: base_case_t<2>(seq, m, d, out, needOut); break;
default: base_case_t<3>(seq, m, d, out, needOut); break;
}
}
/* [a11y-CKMERGE] 归并路径的对角线二分：返回"前 d 个输出里 A 侧元素的个数"。
谓词对 i 单调（可行集是区间）=> 二分。i==0 / j>=r 短路，不会读 run 之外。 */
static inline int a11y_mp_split(const u64 *A, int h, const u64 *B, int r, int d) {
int lo = d - r; if (lo < 0) lo = 0;
int hi = d; if (hi > h) hi = h;
while (lo < hi) {
const int i = (lo + hi + 1) >> 1;
const int j = d - i;
const int ok = (i == 0) || (j >= r) || (A[i - 1] < B[j]);
if (ok) lo = i; else hi = i - 1;
}
return lo;
}
static void solve(u64 *seq, int m, int d, u64 *out, int needOut) {
if (m <= 1) return;
if (d >= K - 1) {
int cnt = 0;
for (int t = 0; t < m; t++) {
u64 v = seq[t];
if ((v >> 31) & 1) cnt++; // insert (bit31 = 1)
else OUT[eid(v)] += (u32)cnt; // query
}
return;
}
if (m <= KG_BASE) { base_case(seq, m, d, out, needOut); return; }
(void)0;
int h = m >> 1, r = m - h;
u64 *save = ABtop;
/* [a11y-CKMERGE] 哨兵槽：见下方计数驱动归并 */
u64 *retA = ABtop; ABtop += h + 1;
u64 *retB = ABtop; ABtop += r + 1;
solve(seq, h, d, retA, 1);
solve(seq + h, r, d, retB, 1);
u64 *cross = ABtop; ABtop += m + 1;
const int rekey = (d + 2 <= K - 1);
const u32 *XN = XD[rekey ? d + 2 : 0];
/* ===== [a11y-CKMERGE] 计数驱动的分块归并（K=3，零每步边界比较）=====
结构完全照搬本账号 1011a #108154 的 [a11w-CKMERGE]（见文件头 [2]）：把归并路径按对角线
切成 3 段，**每段的输出步数是已知量** d[k+1]-d[k]（归并每步恰吃掉一个真元素）⇒ 主循环退化
成 `while (t--) { STEP0; STEP1; STEP2; }` —— 没有结束指针、没有任何每步边界比较。
★ 必须的配套：两个输入 run 的末尾各放一个 `~0ull` 哨兵（arena 各 +1 槽，见上）。某段某一侧
先耗尽时，贪婪比较必须看到 +∞ 才会正确让位；否则当该段的 A 区间恰好空且 c_k == h 时，
会读到 run 区间**之外**的真元素而误判（1011a 上实测：无哨兵时 n=10000 起 bad=1699）。
★ 逐位相同：比较谓词、rekey 表达式、keep 的条件推进一字未改；"何时停"从"指针触底"改成
"数步数"，两者在归并语义下等价（每步恰消耗一个元素）✓ */
int cx = 0;
{
const int mm = h + r;
((u64 *)retA)[h] = ~0ull; /* +inf sentinel */
((u64 *)retB)[r] = ~0ull;
const int d1 = mm / 3, d2 = (2 * mm) / 3;
const int c1 = a11y_mp_split(retA, h, retB, r, d1);
const int c2 = a11y_mp_split(retA, h, retB, r, d2);
const int N0 = d1, N1 = d2 - d1, N2 = mm - d2;
const u64 *pa0 = retA, *pb0 = retB;
const u64 *pa1 = retA + c1, *pb1 = retB + (d1 - c1);
const u64 *pa2 = retA + c2, *pb2 = retB + (d2 - c2);
u64 *po0 = out, *pc0 = cross;
u64 *po1 = out + d1, *pc1 = cross + d1;
u64 *po2 = out + d2, *pc2 = cross + d2;
int L = N0; if (N1 < L) L = N1; if (N2 < L) L = N2;
#define A11Y_ST(PA, PB, PO, PC) do { \
const u64 ea_ = *(PA), eb_ = *(PB); \
const u32 takeA_ = (u32)(ea_ < eb_); \
(PA) += takeA_; (PB) += 1 - takeA_; \
const u64 v_ = takeA_ ? ea_ : eb_; \
if (needOut) *(PO) = v_; \
++(PO); \
*(PC) = rekey ? (((u64)XN[eid(v_)] << 32) | (v_ & 0xFFFFFFFFull)) \
: (v_ & 0xFFFFFFFFull); \
(PC) += (size_t)(1 - ((takeA_ ^ (u32)(v_ >> 31)) & 1u)); \
} while (0)
{ int t = L; while (t--) { A11Y_ST(pa0, pb0, po0, pc0); A11Y_ST(pa1, pb1, po1, pc1); A11Y_ST(pa2, pb2, po2, pc2); } }
{ int t = N0 - L; while (t--) A11Y_ST(pa0, pb0, po0, pc0); }
{ int t = N1 - L; while (t--) A11Y_ST(pa1, pb1, po1, pc1); }
{ int t = N2 - L; while (t--) A11Y_ST(pa2, pb2, po2, pc2); }
#undef A11Y_ST
{ /* 三段各自写 [cross+d_k, ...)，按全局输出序排列 ⇒ 原地拼接即得逐位相同的 cross */
const int k0 = (int)(pc0 - cross);
const int k1 = (int)(pc1 - (cross + d1));
const int k2 = (int)(pc2 - (cross + d2));
if (k1) memmove(cross + k0, cross + d1, sizeof(u64) * (size_t)k1);
if (k2) memmove(cross + k0 + k1, cross + d2, sizeof(u64) * (size_t)k2);
cx = k0 + k1 + k2;
}
}
if (FUSE_OK && d == K - 3) {
u32 *dirt = (u32 *)ABtop; ABtop += (cx + 3) / 2 + 1;
ctr_solve(cross, cx, dirt);
} else {
u64 *cres = ABtop; ABtop += cx;
solve(cross, cx, d + 1, cres, 0);
}
ABtop = save;
}
static u32 *g_ord, *g_x;
static void sort0(int n) {
const u32 *X0 = XD[0];
u32 *cnt = g_x;
memset(cnt, 0, sizeof(u32) * (size_t)(n + 1));
for (int i = 0; i < n; i++) cnt[X0[i] + 1]++;
for (int i = 0; i < n; i++) cnt[i + 1] += cnt[i];
for (int i = 0; i < n; i++) g_ord[cnt[X0[i]]++] = (u32)i;
}
void run(int n, const u32 **x, int k, u32 *out) {
LFSH = 0; LFD = 0;
K = k;
RSZN = n;
for (int d = 0; d < k; d++) XD[d] = x[d];
OUT = out;
ABtop = AB;
if (n <= 1) return;
FUSE_OK = 0;
if (K >= 4 && n > 2 && (u32)n < CT_MAXV) {
const u32 *XL = XD[K - 1];
memset(CT0, 0, (size_t)n);
int ok = 1;
for (int i = 0; i < n; i++) if (++CT0[XL[i]] > 127) { ok = 0; break; }
memset(CT0, 0, (size_t)n);
FUSE_OK = ok;
}
{
int bit = 0;
while ((n >> bit) >= LEAFNB && bit < 30) bit++;
LFSH = bit; LFD = LF_D;
LFB = (int)((u32)n >> (LFSH + LFD)) + 1;
if (LFB > LEAFNB + 1) LFB = LEAFNB + 1;
}
{
u32 dm = rank_compress(n, k);
g_dmax = (u16)(dm > 65535u ? 65535u : dm);
W16 = (dm <= 65536u);
#ifdef Q4A_DIAG
g_dmax = (u16)(dm > 65535u ? 65535u : dm);
#endif
}
sort0(n);
const u32 *X0 = XD[0], *X1 = XD[1];
u64 *ev = AB; // top-level event array: 2n entries
int m = 0, g = 0;
while (g < n) {
int g1 = g;
u32 v = X0[g_ord[g]];
while (g1 < n && X0[g_ord[g1]] == v) g1++;
for (int t = g; t < g1; t++) ev[m++] = mk(X1[g_ord[t]], 1, (int)g_ord[t]);
for (int t = g; t < g1; t++) ev[m++] = mk(X1[g_ord[t]], 0, (int)g_ord[t]);
g = g1;
}
ABtop = AB + m;
u64 *res = ABtop; ABtop += m;
solve(ev, m, 0, res, 0);
}
} // namespace kg
static u64 g_arena[KG_ARENA];
static u32 g_ord_st[300005], g_x_st[300005];
static int g_ready = 0;
static void setup() {
if (g_ready) return;
g_ready = 1;
kg::AB = g_arena;
kg::g_ord = g_ord_st;
kg::g_x = g_x_st;
}
void count_4d(int n, const unsigned *x[4], unsigned *out) {
setup();
memset(out, 0, sizeof(unsigned) * (size_t)n);
kg::run(n, (const unsigned **)x, 4, out);
}
void count_5d(int n, const unsigned *x[5], unsigned *out) {
setup();
memset(out, 0, sizeof(unsigned) * (size_t)n);
kg::run(n, (const unsigned **)x, 5, out);
}
#pragma GCC pop_options
| Compilation | N/A | N/A | Compile OK | Score: N/A | 显示更多 |
| Testcase #1 | 64.641 ms | 11 MB + 260 KB | Accepted | Score: 100 | 显示更多 |