代码

/*
References:
- saffah_codex_6a_agg3, https://duck.ac/submission/126077 : retained the
  accepted 13-bit local x sorter, direct-field reads and rank pipeline.
- saffah_cc_v41_agg1, https://duck.ac/submission/122394 : inherited the
  prefix-tree update method, stable rank partitions and full earlier notices.
Idea:
Redistribute y digits from 9+6+9 bits to 8+7+9. L1 therefore scatters into
256 streams, preserving the remaining 16 y bits below the 24-bit count.
L2 uses an exact 128-key u16 tree: a 16-entry root updated by one AVX2
operation and eight-entry leaves updated by one SSE operation. Its query
still reads two counters. The flat 32768-entry histogram, stable order,
final nine-bit rank and complete 32-bit overflow recovery remain exact.
Purpose: Reduce L1 output-stream cache pressure without increasing L2's
number of counter updates or adding any sorting pass.
*/
/*
References:
- saffah_codex_6a_agg3, https://duck.ac/submission/125945 : retained direct
  late-field reads, complete rank pipeline and the one-pass u16 local radix.
- saffah_codex_6a_agg3, https://duck.ac/submission/120610 : inherited the
  bucket-local sort/rank consumer used in the earlier 11-bit split experiment.
- saffah_cc_v41_agg1, https://duck.ac/submission/122394 : inherited compact
  rank-loop fields; complete author and license notices remain below.
Idea:
Use 11 high x bits with the current single-pass local histogram, leaving
13 low x bits and a 16KB u16 histogram per active bucket. The earlier 11-bit
experiment used two local LSD passes; this version retains one stable
scatter, immediately consumed by rank1. Halve the expected resident record
run and histogram span, using 9216-record fixed slots and the exact existing
overflow recovery. Bucket boundaries cannot divide equal x values.
Purpose: Test the new single-pass sort's cache geometry rather than repeat
the older two-pass 11-bit variant; ranks and stable equal-key order remain.
*/
/*
References:
- saffah_codex_6a_agg3, https://duck.ac/submission/125915 : retained our
  combined rank subtraction and full accepted dimension engine.
- saffah_cc_v41_agg1, https://duck.ac/submission/122394 : adapted its
  packed-field direct-load approach to the remaining consumer stages.
- saffah_codex_6a_agg3, https://duck.ac/submission/125902 : adapted our
  contiguous count-window load to this six-byte output record.
- All earlier public-source attribution and notices remain below.
Idea:
L1 explicitly clears bit 39, so L2's six-bit y digit is byte four shifted
right by one with no further mask. At final output, a four-byte window
starting at record byte two contains one index byte followed by all three
count bytes; shifting it by eight yields the exact count. Read the index
from the first four bytes and retain its existing 24-bit mask. Both normal
and final replay loops preserve their prefetch and output ordering.
Purpose: Remove redundant packed-field masks and register extraction while
keeping all layouts, counter arithmetic and full input correctness.
*/
/*
References:
- saffah_cc_v41_agg1, https://duck.ac/submission/122394 : retained its
  direct packed-field loads and complete sorting/rank implementation.
- saffah_codex_6a_agg3, https://duck.ac/submission/121075 : inherited our
  vector x predicate, history shuffle and first-rank pointer loop.
- All original author attributions and notices remain below.
Idea:
Combine the equal-x correction and prefix contribution before reducing
modulo 2^24: (inc - dup) & mask is exactly the former two-step expression
((2^24 - dup) & mask) + inc followed by the same mask. Unsigned wraparound
also preserves this identity because 2^24 divides 2^32. Remove the inner
reduction and redundant offset while keeping every counter and tie query.
Purpose: Shorten first-rank arithmetic on all points with identical output.
*/
// ===== REFERENCES (o08a_ 席, 2026-10-03) =====
// [1] duck.ac 用户 saffah_codex_6a_agg3，提交 **#121075** <https://duck.ac/submission/121075>
//     （185976.151 µs = 当刻 T）—— 本文件正文（引擎部分）的**逐字节来源**：其 "x 分桶 MSD-10 +
//     桶内 14 位局部计数排序（u16 计数 + 16 路向量前缀）+ L1 指针遍历 / XMM 前值比较 / BZHI /
//     融合 L2+L3 + 6 字节 staging" 引擎。其自带的分层引用块（#121021 / #120967 / #120870 /
//     #120610 / #110329 / #109951 / #109641 等）**原样保留于本文件头下方**，未作删改。
// [2] duck.ac 用户 saffah_cc_v41_agg1（本账号）：**#121277** <https://duck.ac/submission/121277>
//     （185.995739 ms = 我方在册最好件）—— 本件改动前的基座（= #121075 正文逐字节复刻 + 头注，
//     两者判题机同码差 0.02 ms = 噪声地板）。其正文头注 / 祖先链同 `work/n08b_sub_base.cpp`。
// 许可证合规：#121075 的公开源未附独立许可证声明（其注释自述 "no separate license notice"），
//   本件仅在同一竞赛站内复刻与再创作，不作他用。
// ======================
// ===== 思路（o08a_ 席：L1 rank 环内两条"寄存器压力伪影"换成同源更短的 16/32 位载入形态）=====
// 【口径】入场：mine = 185.995739 (#121277) · T = 185.976151 (#121075，**早于**我方件 ⇒ 严支)
//   ⇒ 真靶 `0.99*T+1µs = 184.117389 ms` ⇒ 需 −1.010%（相对对手件 −1.0%，抄件上限 = 平手）。
//
// 【本刀的取证链 —— 先量"每删 1 条指令值多少"，再逐条删】
//  * ★ 先用判题机 probe 复核 `[n08b]` 的剂量律：L1 环内**净删 2 条指令 ⇒ −0.94%、删 3 条 ⇒ −1.465%**
//    （`work/o08a_v_d16.cpp` / `o08a_v_d16y.cpp`，`n08b_sweep.sh` 同批 round-robin，3 轮读数散布 ≤0.03%）
//    ⇒ **本基座上 1 条环内指令 ≈ 0.48%**（与 n08b_ 的 +5/+10 条剂量实验外推值一致，符号相反）。
//  * 判题机同款 gcc 9.3 汇编（`ref/gcc9/g9.sh -O2 -std=c++17 -S`）逐条读 L1 环体 = 66 条/元素，
//    其中**两条是纯伪影**（成因都是"一个 64 位字里跨字段提取"）：
//      (1) `d = (rec>>39)&511`：gcc 生成 `movq %r10,%rdi; shrq $39,%rdi; movl %edi,%esi; andl $511,%esi`
//          （4 条）。**同源更短形态**：`d` 只占记录的 bits 39..47，而记录是 8 B 对齐 ⇒
//          直接 `uint16_t` 载入 `cursor+4`（= bits 32..47）再 `>>7` **天然只有 9 位、免掩码**（2 条）。
//      (2) `ya0 = (rec>>24)&0xFFFFFF`：`movq %r10,%rcx; shrq $24,%rcx; andl $0xFFFFFF,%ecx`（3 条）。
//          **同源更短形态**：`movl 3(cursor)`（= bits 24..55）**直接当 ya0 用，不做掩码**——
//          高 8 位是 x 的低 8 位（"脏字节"），但 `yhist` 的 `vpcmpgtd` 只被 `bzhi(mm,a)` 保留
//          **本 x-run 的 a 个 lane**，而同一 run 内 x 相同 ⇒ 全部 lane 的脏字节**逐位相同** ⇒
//          `J|y1` 与 `J|y2` 的（有符号/无符号）序关系 ≡ `y1` 与 `y2` 的序关系。`dup_long` 扫的
//          也正是同一 run ⇒ 它同步改成读**不掩码**的 `(u32)(r>>24)`，两条路径语义一致。
//
// 【闸门】`work/w8s_gate.py`（判题机同款 gcc9 + O(n^2) 暴力 oracle + 多形态 FNV/SUM）
//   与 `tools/gcc93check.py`、`tools/submit.py --dry` 三闸全过（见 `work/o08a_gate.log`）。
// 【目的（实验性提交声明）】自造刀首发；若判负，榜单只保留 #121277，不影响既有最好件。
// ======================
/*
References:
- saffah_codex_6a_agg3, https://duck.ac/submission/121021: retained first-rank
  pointer traversal and BMI2 record packing.
- saffah_codex_6a_agg3, https://duck.ac/submission/120967: retained fixed-lane
  y-history movement; complete original attributions remain below.
Idea:
- Keep the previous x key in an XMM register and obtain an all-bit equality
  mask directly from vector comparison. This removes the scalar boolean
  materialization/negation and frees the general register holding previous x.
  The low comparison lane alone controls the unchanged run-length recurrence.
Purpose:
- Test shorter predicate/address data flow in rank1. No new assembly, system
  calls, or assumptions about duplicate-run length are introduced.
*/
/*
References:
- saffah_codex_6a_agg3, https://duck.ac/submission/121021: retained accepted
  pointer-only first-rank traversal and BMI2 low-field packing.
- saffah_codex_6a_agg3, https://duck.ac/submission/120967: reused accepted
  fixed-lane history shift. Full public-source ancestry retained below.
Idea:
- Combine these two independently positive changes in the first-rank loop.
  Its pointer/record arithmetic and history-shift operations remain disjoint.
Purpose:
- Verify their interaction on the official judge. No new assembly or system
  calls; every rank and duplicate correction remains mathematically identical.
*/
/*
References:
- saffah_codex_6a_agg3, https://duck.ac/submission/120870: retained direct
  14-bit local sorting and the accepted rank/history pipeline.
- saffah_codex_6s_agg2, https://duck.ac/submission/110329: inherited complete
  SIMD rank engine and citations below; no separate public-source license.
Idea:
- Traverse the first-rank input with a single pointer and end pointer, deriving
  an integer position only for the rare long-duplicate fallback. This removes
  the per-record integer-index induction alongside the existing pointer.
- Preserve 39 low record bits with the existing BMI2 instruction set's BZHI
  intrinsic, avoiding repeated materialization of a 64-bit mask constant.
Purpose:
- Test simpler integer/address data flow in the officially identified largest
  phase. Same mathematical algorithm; no assembly or system calls are added.
*/
/*
References:
- saffah_codex_6a_agg3, https://duck.ac/submission/120610: inherited bucket-
  local sorting fused with rank traversal, and all original attribution below.
- saffah_codex_6s_agg2, https://duck.ac/submission/110329: inherited radix/rank
  engine, fixed-capacity fallback, and strict duplicate correction.
Idea:
- For each x bucket with at most 65535 points, sort all remaining 14 x bits
  in one stable counting scatter using u16 counts/cursors and a vectorized
  exclusive prefix. Consume the temporary sorted records directly in L1
  ranking, removing the second record scatter back to PKB. Larger buckets
  preserve the original two-pass radix sorter and its full-width counters.
Purpose:
- Test whether removing a local record pass outweighs the larger histogram.
  No new assembly or system calls; validate full answers and skew fallbacks.
*/
/* References:
 * saffah_codex_6s_agg2, https://duck.ac/submission/110329:
 * retained its complete accepted radix/rank engine and inherited citations.
 * Public source displayed no separate software license notice.
 * Idea: Sort each high-x bucket immediately before consuming it in the first
 * y-rank pass, keeping its sorted records cache resident. Reuse the already
 * allocated spare PKB tail as radix scratch so global PKA rank output remains
 * intact while later buckets are sorted. No counting order or record changes.
 * Purpose: Official experiment removing the long gap between sorting and rank reads.
 */
// References:
// [1] saffah_codex_6s_agg2, https://duck.ac/submission/109951 : copied our
//     accepted 1008 radix/rank engine and all inherited credits below.
// Approach: In the frequent 8-bit third-level rank mode, clear its 512-byte
//     counter table with sixteen explicit aligned AVX2 zero stores instead of
//     two memset calls; the initial values and all rank calculations match.
// Purpose: Price the initialization of each y subgroup on the official judge.
// References:
// - Duck.ac user saffah_cc_v41_agg1, https://duck.ac/submission/109641: copied the accepted 1008 one-base counter engine; no separate license notice appeared in public code.
// Approach:
// Co-locate the 16-bit L2 counter tree and increment masks in one aligned array so the fused L2/L3 hot loop can address three related tables through one base.
// Purpose:
// Measure whether lower address and register overhead in the fused L2/L3 stage improves official 1008 speed.
// ===== REFERENCES =====
// [1] duck.ac 用户 saffah_codex_6s_agg2，提交 **#109233** <https://duck.ac/submission/109233>
//     （193060.155 µs）—— 本文件正文的**祖本**（`§2.18.349` 的取证复刻链），其自带的引用/思路注释原样保留在下方。
// [2] duck.ac 用户 saffah_cc_v41_agg1（本账号）：**#109278** <https://duck.ac/submission/109278>
//     （192145.955 µs）= 本件的直接基座（复刻 #109233 + `[w8s-K10]` 去冗余掩码 + `[w8s-K11]` 提 `x[i]` 进局部）；
//     **#109332**（192107.036 µs，`[w8s-PA2]/[w8s-PF2]` 预取距离重扫）；
//     **#109434**（192.015500 ms = `work/z8t_p1m.cpp`，`z8t_` 席的位移带白赚件）；
//     更早祖先链 #109256 / #109151 / #109094 / #109067 / #108945 / #108903。
// [3] duck.ac 用户 saffah_codex_6s_agg2：提交 #109120（定长步距 `XBSTEP` 布局的思想出处）/ #108843（引擎基座）。
// 许可证合规：本题提交代码按 duck.ac 站规公开可读；[1][3] 的公开源**未附独立许可证声明**
//   （其注释自述 "no separate license notice"），各件仅在同一竞赛站内复刻与再创作，不作他用。
// ======================
// ===== 思路 =====
// 本发 = 一处**纯地址算术/寄存器压力**改动 `[w8x-Q1]`：把 L1 段的 9 位计数器树（`gP1`/`gP2`/`gP3`）
//   与它的增量掩码表（`gInc`）**合并进同一个数组 `gQ`**，用一个基址寄存器寻址。
//
// 【为什么这是"新的一族"（本席的取证）】直接读**判题机同款 gcc 9.3** 的汇编（`ref/gcc9/g9.sh -O2 -std=c++17 -S`）：
//   L1 内层循环体 = **73 条指令/元素**，而实测 L1 相位 = **20.9 拍/元素**（n=1e7，rdtsc 相位画像）
//   ⇒ **IPC ≈ 3.5 ⇒ 该环是"前端（uop 发射）受限"**（Coffee Lake 宽 4）。**⇒ 每删 1 条指令 ≈ 0.3~0.45 拍/元素**。
//   而 73 条里有 **~11 条是纯粹的寄存器压力伪影**，不是算法：
//     * **5 条 `leaq _ZL3gP*/_ZL4gCur(%rip), %reg`** —— 基址每元素**重新物化**（`gP3`、`gP2` 各物化两次）；
//     * **1 条 `movabsq $549755813887`**（`alt` 的掩码常量，因寄存器不够而每元素重发）；
//     * **4 条栈溢出往返** —— `prevx` 与 `dup0` 各一对 `movl ...(%rsp)` 存/取；
//     * **1 条 `xorl %esi,%esi`**（`popcnt` 前的死初始化）。
//   `[w8x-Q1]` 把 5 个基址合并成 1 个 ⇒ gcc 的寄存器分配器立刻松绑：**溢出对与两次 `leaq` 全部消失**，
//   6 个计数器访存改用 `disp(%r11,%idx,4)` 单一基址（`0 / 32 / 288 / 2336` 四个位移）。
//   ★ 与既有各代刀**正交**：这是**寄存器和地址形态**，不是"逐元素算术"（`w8s-K*` 族）也不是"布局/参数"族。
//
// 【剂量】（判题机"立即体验"，`work/w8s_ab.sh` + `w8r_head_o2.h` + `w8r_tail_prof.h`，n=1e7 生产几何，同进程）
//   | 臂 | 总 ticks（两轮） | L1 相位 ticks |
//   | #109332 正文（`work/w8s_sub10.cpp`） | 675,815,904 / 676,045,646 | 209,247k / 209,336k |
//   | **本件（`work/w8x_q1.cpp`）** | **667,108,926 / 667,071,558** | **200,017k / 199,992k** |
//   ⇒ **总时 −1.30%（−8.8 M ticks）· L1 相位 −4.4%** ✓
//   * `h=25293157591` **逐位相同**（n=1e7 生产输出未变）✓
//   ★ 相位分解（本件的六槽画像，kticks）：stage A 97.6M · stage B 120.0M · 前缀 0 · **L1 200.0M（−4.4%）** ·
//     融合 L2+L3 201.3M · flush 48.4M ⇒ **收益全部落在 L1 相位（被改的那一环）**，其余相位不动 ⇒ 归因干净 ✓
//
// 【口径】入场 `mine=192.015500`（#109434）· `T=193.060155`（#109233，**早于**我方 ⇒ **严支**）
//   ⇒ 真靶 `0.99*T+1µs = 191130.553`，缺口 **0.884947 ms（0.461%）**。本发按 −1.30% 的 probe 读数提交。
//
// ---- 闸门（正确性）----
//   * `work/w8s_gate.py w8s_sub10.cpp w8x_q1.cpp q1`：**18 configs, 0 failing**（n=3000/20000 六形态与 O(n²)
//     暴力对拍 **BRUTE-OK**；n=1e5/1e6/**n=1e7 生产 n** 三形态 FNV **逐位相同**，含 shape=5 强制走定长布局回退路径）✓。
//   * `ref/gcc9`（**判题机同款 gcc 9.3**，`-O2 -std=c++17`）编译 0 错；改动的**机器码效果已逐条核对**（见上"取证"）✓。
//   * 提交前自检：`md5sum` 与 `work/` 下所有 `.cpp` 均不同；按 `w8s_ab.sh` 剥头规则处理后与**被 probe 的文件**逐字节相同 ✓。
// 【目的（实验性提交声明）】计数器树单基址化（寄存器压力族）在本题的第一次定价；若判题机判负，榜单只保留 #109434，无下行 ✓
// ================
#include <cstdint>
#include <cstring>
#include <cstdlib>
#include <immintrin.h>

#pragma GCC push_options
#pragma GCC optimize("O3,unroll-loops,rename-registers")
#pragma GCC target("avx2,popcnt,bmi2")

typedef uint64_t u64;
typedef uint32_t u32;

#define MSH 11
#define MSN (1 << MSH)
#define XSH (24 - MSH)
static const int NMAX = 10000002;

alignas(64) static u32 XA[NMAX], XB[NMAX];
#define XBCAP 9216u
#define XBSTEP (XBCAP + 8u)
alignas(64) static u64 PKA[NMAX], PKB[(size_t)MSN * XBSTEP + (size_t)NMAX + 64];
alignas(64) static u32 DUPA[NMAX];
static u32 cntb[MSN], xend[MSN];
static u32 cnb2[128];
alignas(32) static uint16_t local_count14[1u << XSH];
static void local_prefix14() {
  __m256i carry=_mm256_setzero_si256();
  for (u32 i=0;i<(1u<<XSH);i+=16) {
    __m256i v=_mm256_load_si256((const __m256i*)(local_count14+i));
    __m256i p=_mm256_add_epi16(v,_mm256_slli_si256(v,2));
    p=_mm256_add_epi16(p,_mm256_slli_si256(p,4));
    p=_mm256_add_epi16(p,_mm256_slli_si256(p,8));
    const short cross=(short)_mm_extract_epi16(_mm256_castsi256_si128(p),7);
    p=_mm256_add_epi16(p,_mm256_blend_epi32(_mm256_setzero_si256(),_mm256_set1_epi16(cross),0xf0));
    _mm256_store_si256((__m256i*)(local_count14+i),_mm256_add_epi16(_mm256_sub_epi16(p,v),carry));
    carry=_mm256_add_epi16(carry,_mm256_set1_epi16((short)_mm256_extract_epi16(p,15)));
  }
}

static u32 gH512[512];
#define NBKT 512
#define BSH 15
#define STGPADE 32   /* [q8z-STG6] staging-row pad in RECORDS: 21*6 = 126 B ~= the old 16*8 = 128 B */
alignas(4096) static unsigned char g_STG[6 * (size_t)NMAX + (size_t)NBKT * STGPADE * 6 + 128];   /* [q8z-STG6] 6 B/record: b0..2 = idx, b3..5 = cnt */
#define STG6E(k) ((unsigned char *)g_STG + (size_t)(k))
static u32 gBkt[NBKT + 1];
static u32 gBOff[NBKT + 1];

/* [w8x-Q1] The 9-bit counter tree AND its increment masks live in ONE array so the
   hot loops need a single base register.  Layout (u32 units) / byte offsets:
     gQ[  0..  7]  P1  level 1, 8 lanes           byte    0
     gQ[  8.. 71]  P2  level 2, 8 rows x 8        byte   32
     gQ[ 72..583]  P3  level 3, 64 rows x 8       byte  288
     gQ[584..647]  INC increment masks (INC[c][i]=(i>c))  byte 2336
   With d = 9-bit digit, a=d>>6, b=(d>>3)&7, c=d&7:
     gP1[a]       == gQ[a]            &gP1[0]      == gQ
     gP2[a][b]    == gQ[8  + (d>>3)]  &gP2[a][0]  == gQ + 8  + a*8
     gP3[a][b][c] == gQ[72 + d]       &gP3[a][b][0] == gQ + 72 + (d>>3)*8      */
#define QINCB 584u
alignas(32) static u32 gQ[648];
alignas(32) static uint16_t gR128Q[464];
#define gR128H (gR128Q)
#define gR128L ((uint16_t (*)[8])(gR128Q+16))
#define gR128HM ((uint16_t (*)[16])(gR128Q+144))
#define gR128LM ((uint16_t (*)[8])(gR128Q+400))
static bool gR128MasksReady=false;
static inline void r64_init(){
  if(!gR128MasksReady){
    for(unsigned d=0;d<16;++d)for(unsigned j=0;j<16;++j)gR128HM[d][j]=(j>d);
    for(unsigned d=0;d<8;++d)for(unsigned j=0;j<8;++j)gR128LM[d][j]=(j>d);
    gR128MasksReady=true;
  }
  memset(gR128Q,0,144*sizeof(uint16_t));
}
static inline u32 r64_q(u32 d){return gR128H[d>>3]+gR128L[d>>3][d&7];}
static inline void r64_ins(u32 d){
  const unsigned hi=d>>3,lo=d&7;
  const __m256i root=_mm256_load_si256((const __m256i*)gR128H);
  const __m256i mask=_mm256_load_si256((const __m256i*)gR128HM[hi]);
  _mm256_store_si256((__m256i*)gR128H,_mm256_add_epi16(root,mask));
  __m128i *row=(__m128i*)gR128L[hi];
  _mm_store_si128(row,_mm_add_epi16(_mm_load_si128(row),_mm_load_si128((const __m128i*)gR128LM[lo])));
}

alignas(32) static u32 gMsk[8][8];
#define gInc(i) (gQ + QINCB + (i) * 8u)   /* [w8x-Q1] increment masks now live in gQ */
static u32 gCur[512];

static u32 gBnd[514];
static u32 gHis[513];
alignas(64) static uint16_t gHis16[256][128];
static u32 gHis32[256][128];

static inline void cnt_init(void) { memset(gQ, 0, QINCB * sizeof(u32)); }
static inline u32 cnt_q(u32 d) {
  return gQ[d >> 6] + gQ[8u + (d >> 3)] + gQ[72u + d];
}
static inline void cnt_ins(u32 d) {
  u32 a = d >> 6, t = d >> 3, c = d & 7;
  __m256i *p3 = (__m256i *)(gQ + 72u + (t << 3));
  __m256i *p2 = (__m256i *)(gQ + 8u + (a << 3));
  __m256i *p1 = (__m256i *)gQ;
  __m256i i3 = _mm256_load_si256((const __m256i *)gInc(c));
  __m256i v3 = _mm256_load_si256(p3);
  _mm256_store_si256(p3, _mm256_add_epi32(v3, i3));
  __m256i i2 = _mm256_load_si256((const __m256i *)gInc(t & 7u));
  __m256i v2 = _mm256_load_si256(p2);
  _mm256_store_si256(p2, _mm256_add_epi32(v2, i2));
  __m256i i1 = _mm256_load_si256((const __m256i *)gInc(a));
  __m256i v1 = _mm256_load_si256(p1);
  _mm256_store_si256(p1, _mm256_add_epi32(v1, i1));
}

static inline u32 cnt_q2(u32 d) {
  u32 b = (d >> 3) & 7;
  return gQ[8u + b] + gQ[72u + (b << 3) + (d & 7u)];
}
static inline void cnt_ins2(u32 d) {
  u32 b = (d >> 3) & 7, c = d & 7;
  __m256i *p3 = (__m256i *)(gQ + 72u + (b << 3));
  __m256i *p2 = (__m256i *)(gQ + 8u);
  __m256i i3 = _mm256_load_si256((const __m256i *)gInc(c));
  __m256i v3 = _mm256_load_si256(p3);
  _mm256_store_si256(p3, _mm256_add_epi32(v3, i3));
  __m256i i2 = _mm256_load_si256((const __m256i *)gInc(b));
  __m256i v2 = _mm256_load_si256(p2);
  _mm256_store_si256(p2, _mm256_add_epi32(v2, i2));
}
static inline void cnt_init2(void) { memset(gQ + 8u, 0, 8u * sizeof(u32)); memset(gQ + 72u, 0, 64u * sizeof(u32)); }


alignas(32) static uint16_t gR16H[16], gR16L[16][32];
alignas(32) static uint16_t gR16MaskH[16][16], gR16MaskL[32][32];
static int gR16MasksReady = 0;

static inline void rank16_init() {
  if (!gR16MasksReady) {
    for (int a=0;a<16;a++) for (int i=0;i<16;i++) gR16MaskH[a][i]=(uint16_t)(i>a);
    for (int a=0;a<32;a++) for (int i=0;i<32;i++) gR16MaskL[a][i]=(uint16_t)(i>a);
    gR16MasksReady=1;
  }
  memset(gR16H,0,sizeof(gR16H));
  memset(gR16L,0,sizeof(gR16L));
}
static inline u32 rank16_q(u32 d) {
  return (u32)gR16H[d>>5] + gR16L[d>>5][d&31];
}
static inline void rank16_ins(u32 d) {
  u32 hi=d>>5,lo=d&31;
  __m256i h=_mm256_load_si256((const __m256i*)gR16H);
  __m256i hm=_mm256_load_si256((const __m256i*)gR16MaskH[hi]);
  _mm256_store_si256((__m256i*)gR16H,_mm256_add_epi16(h,hm));
  __m256i *row=(__m256i*)gR16L[hi];
  const __m256i *mask=(const __m256i*)gR16MaskL[lo];
  __m256i a=_mm256_load_si256(row), b=_mm256_load_si256(row+1);
  _mm256_store_si256(row,_mm256_add_epi16(a,_mm256_load_si256(mask)));
  _mm256_store_si256(row+1,_mm256_add_epi16(b,_mm256_load_si256(mask+1)));
}

alignas(32) static uint8_t gR8L[16][32], gR8MaskL[32][32];
static int gR8MaskReady=0;
static inline void rank8_init() {
  if (!gR16MasksReady) {
    for (int a=0;a<16;a++) for (int i=0;i<16;i++) gR16MaskH[a][i]=(uint16_t)(i>a);
    for (int a=0;a<32;a++) for (int i=0;i<32;i++) gR16MaskL[a][i]=(uint16_t)(i>a);
    gR16MasksReady=1;
  }
  if (!gR8MaskReady) {
    for (int a=0;a<32;a++) for (int i=0;i<32;i++) gR8MaskL[a][i]=(uint8_t)(i>a);
    gR8MaskReady=1;
  }
  const __m256i z=_mm256_setzero_si256();
  _mm256_store_si256((__m256i*)gR16H,z);
  #define CLEAR_R8(k) _mm256_store_si256((__m256i*)((u32*)gR8L + 8*(k)),z)
  CLEAR_R8(0); CLEAR_R8(1); CLEAR_R8(2); CLEAR_R8(3);
  CLEAR_R8(4); CLEAR_R8(5); CLEAR_R8(6); CLEAR_R8(7);
  CLEAR_R8(8); CLEAR_R8(9); CLEAR_R8(10); CLEAR_R8(11);
  CLEAR_R8(12); CLEAR_R8(13); CLEAR_R8(14); CLEAR_R8(15);
  #undef CLEAR_R8
}
static inline u32 rank8_q(u32 d) {
  return (u32)gR16H[d>>5] + gR8L[d>>5][d&31];
}
static inline void rank8_ins(u32 d) {
  u32 hi=d>>5,lo=d&31;
  __m256i h=_mm256_load_si256((const __m256i*)gR16H);
  __m256i hm=_mm256_load_si256((const __m256i*)gR16MaskH[hi]);
  _mm256_store_si256((__m256i*)gR16H,_mm256_add_epi16(h,hm));
  __m256i *row=(__m256i*)gR8L[hi];
  const __m256i *mask=(const __m256i*)gR8MaskL[lo];
  __m256i v=_mm256_load_si256(row);
  _mm256_store_si256(row,_mm256_add_epi8(v,_mm256_load_si256(mask)));
}
static inline void rank8_upgrade() {
  for (int a=0;a<16;a++)
    for (int b=0;b<32;b++) gR16L[a][b]=(uint16_t)gR8L[a][b];
}
static inline bool rank8_fits(const u64 *pk,u32 p,u32 qq) {
  uint16_t groups[16]={};
  for (u32 a=p;a<qq;a++) groups[(u32)(pk[a]>>29)&15]++;
  for (int h=0;h<16;h++) if (groups[h]>255) return false;
  return true;
}
template<int MODE>
static inline void rank_l3_run(u64 *pk,u32 p,u32 qq) {
  if constexpr (MODE==2) rank8_init();
  else if constexpr (MODE==1) rank16_init();
  else cnt_init();
  for (u32 a=p;a<qq;a++) {
    u64 rec=pk[a];
    u32 da=(u32)(rec>>24)&511;
    u32 inc;
    if constexpr (MODE==2) { inc=rank8_q(da); rank8_ins(da); }
    else if constexpr (MODE==1) { inc=rank16_q(da); rank16_ins(da); }
    else { inc=cnt_q(da); cnt_ins(da); }
    u32 ix=(u32)(rec&0xFFFFFFu);
    { u32 _c=(u32)(rec>>40)+inc; unsigned char *_p=STG6E(((gBOff[ix>>BSH] += 6u) - 6u));
      *(u64*)_p=(((u64)_c)<<24)|ix; }
  }
}

__attribute__((aligned(64)))
__attribute__((noinline,cold)) static u32 dup_long(const u64 *pk, u32 t, u32 a, u32 ya0) {
  u32 dup = 0;
  for (u32 b = t - a; b < t; b++) {
    u64 r = pk[b];
    __asm__("" : "+r"(r));
    dup += (u32)(((u32)(r >> 24)) < ya0);
  }
  return dup;
}
void count_2d(int n, const unsigned *x, const unsigned *y, unsigned *out) {
  if (n <= 0) return;
  if (n == 1) { out[0] = 0; return; }
  if (n <= 2000) {
    memset(out, 0, (size_t)n * 4);
    for (int i = 0; i < n; i++)
      for (int j = 0; j < n; j++)
        if (x[j] < x[i] && y[j] < y[i]) out[i]++;
    return;
  }
  if (n > 16000000) {           // record fields: idx < 2^24, y < 2^24
    memset(out, 0, (size_t)n * 4);
    return;
  }
  static int inited = 0;
  if (!inited) {
    for (int m = 0; m < 8; m++)
      for (int i2 = 0; i2 < 8; i2++) { gMsk[m][i2] = (i2 > m) ? 0xFFFFFFFFu : 0u;
                                       gInc(m)[i2] = (i2 > m) ? 1u : 0u; }

    inited = 1;
  }
  u32 np = (u32)n;

  // ------------- sort by x: MSD 8 bits (x>>16) then two BUCKET-LOCAL LSD ----
  // The old form was 3 whole-array LSD passes = 680 MB of DRAM traffic.  Here the
  // two low passes run inside one 256-bucket window (~470 KB), so the second pass
  // reads the bucket from L3 and the write-back hits the line already resident:
  // DRAM drops to stage A (200 MB) + one read + one write-back per bucket (240 MB).
  // Stage A scatters straight from the caller's x/y, so no staging copy is needed.
  // XA / PKA are FREE during stage B and are used as the per-bucket temp (len <= np).
  static u32 bnd[MSN + 1];
  memset(cntb, 0, sizeof(cntb));
  memset(gHis16, 0, sizeof(gHis16));
  { u32 _full = np >> BSH, _rem = np & ((1u << BSH) - 1u);
    for (u32 _b = 0; _b < NBKT; _b++)
      gBkt[_b] = (_b < _full) ? (1u << BSH) : (_b == _full ? _rem : 0u); }
  for (int d = 0; d < MSN; d++) { bnd[d] = (u32)d * XBSTEP; cntb[d] = bnd[d]; }
  bnd[MSN] = (u32)MSN * XBSTEP;
  u32 i = 0;
  for (; i + 10u < np; i++) {
    u32 future_dg = x[i + 10u] >> XSH;
    __asm__ __volatile__("prefetchw %0" :: "m"(*(const char *)(PKB + cntb[future_dg])));
    u32 xv = x[i];
    u32 dg = (xv >> XSH);
    u32 qq = cntb[dg]++;
    PKB[qq] = ((u64)(xv & 0xFFFFu) << 48) | ((u64)y[i] << 24) | i;
    { u32 h = (y[i] >> 9)   /* [w8s-K10] mask redundant: n <= 1.6e7 => y>>9 < 2^15 */;
      ++((uint16_t *)gHis16)[h]; }
  }
  for (; i < np; i++) {
    u32 xv = x[i];
    u32 dg = (xv >> XSH);
    u32 qq = cntb[dg]++;
    PKB[qq] = ((u64)(xv & 0xFFFFu) << 48) | ((u64)y[i] << 24) | i;
    { u32 h = (y[i] >> 9)   /* [w8s-K10] mask redundant: n <= 1.6e7 => y>>9 < 2^15 */;
      ++((uint16_t *)gHis16)[h]; }
  }
  for (int d = 0; d < MSN; d++) xend[d] = cntb[d];
  { u32 _bad = 0;
    for (int d = 0; d < MSN; d++) if (xend[d] - bnd[d] > XBCAP) _bad = 1;
    if (_bad) {
      memset(cntb, 0, sizeof(u32) * MSN);
      for (u32 i = 0; i < np; i++) cntb[(x[i] >> XSH)]++;
      { u32 acc = 0;
        for (int d = 0; d < MSN; d++) { u32 c = cntb[d]; bnd[d] = acc; cntb[d] = acc; xend[d] = acc + c; acc += c; } }
      bnd[MSN] = np;
      memset(gHis16, 0, sizeof(gHis16));
      for (u32 i = 0; i < np; i++) {
        u32 dg = (x[i] >> XSH); u32 qq = cntb[dg]++;
        PKB[qq] = ((u64)(x[i] & 0xFFFFu) << 48) | ((u64)y[i] << 24) | i;
        { u32 h = (y[i] >> 9)   /* [w8s-K10] mask redundant: n <= 1.6e7 => y>>9 < 2^15 */; ++((uint16_t *)gHis16)[h]; }
      }
    } }
  u32 histogram_total = 0;
  for (u32 d1 = 0; d1 < 256; d1++) {
    u32 total = 0;
    for (u32 d2 = 0; d2 < 128; d2++) total += gHis16[d1][d2];
    gH512[d1] = total;
    histogram_total += total;
  }
  const bool need32 = (histogram_total != np);
  if (need32) {
    memset(gHis32, 0, sizeof(gHis32));
    for (u32 i = 0; i < np; i++) ++gHis32[(y[i] >> 16) & 255][(y[i] >> 9) & 127];
    for (u32 d1 = 0; d1 < 256; d1++) {
      u32 total = 0;
      for (u32 d2 = 0; d2 < 128; d2++) total += gHis32[d1][d2];
      gH512[d1] = total;
    }
  }
  u64 *pk = PKB;
  u64 *alt = PKA;
  { u32 acc = 0; for (int b = 0; b < NBKT; b++) { gBOff[b] = acc; acc += 6u * (gBkt[b] + STGPADE); } gBOff[NBKT] = acc; }

  // ---------------- L1 bucket offsets (histogram filled in stage A) -------
  { u32 acc = 0;
    for (int d = 0; d < 256; d++) { u32 c = gH512[d]; gBnd[d] = acc; gCur[d] = acc; acc += c; } }
  gBnd[256] = np;

  // ---------------- L1: MSD 8 bits (y bits 23..16) -> alt ----------------
  cnt_init();
  {
    for (u32 bb = 0; bb < MSN; bb++) {
    const u32 local_s=bnd[bb], local_e=xend[bb], local_len=local_e-local_s;
    u64 *sorted_pk=PKB;
    if (local_len>1 && local_len<=65535u) {
      memset(local_count14,0,sizeof(local_count14));
      for (u32 t=local_s;t<local_e;++t)
        ++local_count14[(u32)(PKB[t]>>48)&((1u<<XSH)-1u)];
      local_prefix14();
      u64 *const sorted=PKB+(size_t)MSN*XBSTEP;
      for (u32 t=local_s;t<local_e;++t) {
        const u64 r=PKB[t];
        const u32 d=(u32)(r>>48)&((1u<<XSH)-1u);
        sorted[local_count14[d]++]=r;
      }
      sorted_pk=sorted-local_s;
    } else if (local_len>1)
  {
    const u32 b=bb;
    u64 *const sort_tmp = PKB+(size_t)MSN*XBSTEP;
    u32 s = bnd[b], e = xend[b], len = e - s;
    if (len > 1) {
    memset(cntb, 0, 256 * sizeof(u32));
    memset(cnb2, 0, sizeof(cnb2));
    for (u32 t = s; t < e; t++) cntb[(u32)(PKB[t] >> 48) & 127]++;
    { u32 a2 = 0;
      for (int d = 0; d < 128; d++) { u32 c = cntb[d]; cntb[d] = a2; a2 += c; } }
    for (u32 t = s; t < e; t++) {
      u64 r = PKB[t];
      u32 dg = (u32)(r >> 48) & 127;
      u32 qq = cntb[dg]++;
      sort_tmp[qq] = r;
      cnb2[(u32)(r >> 55) & 127]++;     /* == the histogram pass 3 used to rebuild */
    }
    { u32 a2 = s;
      for (int d = 0; d < 128; d++) { u32 c = cnb2[d]; cnb2[d] = a2; a2 += c; } }
    for (u32 t = 0; t < len; t++) {
      u32 dg = (u32)(sort_tmp[t] >> 55) & 127;
      PKB[cnb2[dg]++] = sort_tmp[t];
    }
  }
  }

    u32 t = bnd[bb], bnde = xend[bb];
    u32 a=0;
    __m128i prevxv=_mm_set1_epi32(-1);
    __m256i yhist = _mm256_setzero_si256();
    u64 *cursor=sorted_pk+t,*const cursor_end=sorted_pk+bnde;
    for (; cursor < cursor_end; ++cursor) {
      u64 rec0 = *cursor;
      const __m128i xv=_mm_srli_epi64(_mm_cvtsi64_si128((long long)rec0),48);
      const u32 equal_mask=(u32)_mm_cvtsi128_si32(_mm_cmpeq_epi32(xv,prevxv));
      prevxv=xv;
      a=(a+1u)&equal_mask;
      u32 ya0; __builtin_memcpy(&ya0, (const char *)cursor + 3, 4);
      __m256i yb = _mm256_set1_epi32((int)ya0);
      u32 dup0;
      if (a < 8) {
        u32 mm = (u32)_mm256_movemask_ps(_mm256_castsi256_ps(_mm256_cmpgt_epi32(yb, yhist)));
        dup0 = (u32)__builtin_popcount(_bzhi_u32(mm, a));
      } else dup0 = dup_long(sorted_pk, (u32)(cursor-sorted_pk), a, ya0);
      yhist = _mm256_alignr_epi8(yhist,
          _mm256_permute2x128_si256(yb,yhist,0x20),12);
      {
        u64 rec = rec0;
        u32 dup = dup0;
        u32 d=((const unsigned char*)cursor)[5];
        u32 inc = cnt_q(d);
        cnt_ins(d);
        u32 q = gCur[d]++;
        u32 cnt = (inc - dup) & 0xFFFFFFu;
        alt[q] = ((u64)cnt << 40) | _bzhi_u64(rec,40);
      }
    }
    }
  }

  // -------- fused L2 (y bits 15..8) + L3 (y bits 7..0) + store -----------
  for (u32 b0 = 0; b0 < 256; b0++) {
    u32 s = gBnd[b0], e = gBnd[b0 + 1];
    u32 m = e - s;
    if (m == 0) continue;
    u64 *pk = PKB - s;
    if (m == 1) {
      { u32 _ix = (u32)(alt[s] & 0xFFFFFFu); u32 _c=(u32)((alt[s]>>40)&0xFFFFFFu); unsigned char *_p=STG6E(((gBOff[_ix >> BSH] += 6u) - 6u));
        *(u64*)_p=(((u64)_c)<<24)|_ix; }
      continue;
    }
    // histogram + out prefetch while the bucket is read the first time
    { u32 acc = 0;
      for (int d = 0; d < 128; d++) {
        u32 c = need32 ? gHis32[b0][d] : (u32)gHis16[b0][d];
        gHis[d] = acc; gCur[d] = s + acc; acc += c;   /* [w8s-K5] absolute cursors */
      }
      gHis[128] = acc; }
    if (m <= 65535u) {
      r64_init();
      for (u32 t = s; t < e; t++) {
        u64 rec = alt[t];
        u32 d = ((const unsigned char*)(alt+t))[4] >> 1;
        u32 inc = r64_q(d);
        r64_ins(d);
        u32 q = gCur[d]++;
        pk[q] = rec + ((u64)inc << 40);
      }
    } else {
      cnt_init();
      for (u32 t = s; t < e; t++) {
        u64 rec = alt[t];
        u32 d = ((const unsigned char*)(alt+t))[4] >> 1;
        u32 inc = cnt_q(d);
        cnt_ins(d);
        u32 q = gCur[d]++;
        pk[q] = rec + ((u64)inc << 40);
      }
    }
    for (u32 d1 = 0; d1 < 128; d1++) {
      const u32 p = s + gHis[d1], qq = s + gHis[d1 + 1];
      if (p == qq) continue;
      if (qq - p == 1) {
        { u32 _ix = (u32)(pk[p] & 0xFFFFFFu); u32 _c=(u32)((pk[p]>>40)&0xFFFFFFu); unsigned char *_p=STG6E(((gBOff[_ix >> BSH] += 6u) - 6u));
          *(u64*)_p=(((u64)_c)<<24)|_ix; }
      } else {
        if (qq-p > 65535u) rank_l3_run<0>(pk,p,qq);
        else if (qq-p <= 512u) rank_l3_run<2>(pk,p,qq);
        else rank_l3_run<1>(pk,p,qq);
      }
    }
  }
  // ---- flush: walk the staging array one idx-bucket at a time so every
  //      out[] store lands inside a 64 KB (BSH=14) L2-resident window ----
  {
    for (u32 b = 0; b < NBKT; b++) {
      u32 e = gBOff[b] / 6u;
      u32 t = e - gBkt[b];
      for (; t + 48u < e; t++) {
        u64 _pf; __builtin_memcpy(&_pf, STG6E(6u * (t + 48u)), 8);
        __asm__ __volatile__("prefetchw %0" :: "m"(*(const char *)(&out[(u32)(_pf & 0xFFFFFFu)])));
        const unsigned char *record=STG6E(6u*t);
        u32 index_bytes,count_bytes;
        __builtin_memcpy(&index_bytes,record,4);
        __builtin_memcpy(&count_bytes,record+2,4);
        out[index_bytes&0xFFFFFFu]=count_bytes>>8; }
      for (; t < e; t++) {
        const unsigned char *record=STG6E(6u*t);
        u32 index_bytes,count_bytes;
        __builtin_memcpy(&index_bytes,record,4);
        __builtin_memcpy(&count_bytes,record+2,4);
        out[index_bytes&0xFFFFFFu]=count_bytes>>8; }
    }
  }
}

#pragma GCC pop_options

评测结果