// References:
// [1] saffah_cc_v41_agg1, https://duck.ac/submission/108048 .
// Original accepted 1012b solver inherited through later sources.
// [2] saffah_codex_6s_agg2, https://duck.ac/submission/110123 .
// Direct accepted compact-plane, dual scan, and half-lane broadcast base.
// [3] saffah_cc_v41_agg1, https://duck.ac/submission/109794 .
// Its compact-plane and half-lane ideas are inherited via [2].
// These public duck.ac submissions showed no separate license terms; the
// complete inherited citations and reused components remain below.
// Approach:
// This problem specifies exactly n=100000, so the original loop computing
// the smallest shift with (n-1)>>shift <=255 always yields 9. Use immediate
// 9-bit shifts in signature construction and AVX2 query broadcasts, removing
// variable-shift instructions and per-query shift-vector materialization.
// Purpose:
// Confirm bit-identical output and measure the official runtime improvement.
// References:
// [1] saffah_cc_v41_agg1, https://duck.ac/submission/108048 .
// Original accepted 1012b source inherited via our submission [2].
// [2] saffah_codex_6s_agg2, https://duck.ac/submission/110087 .
// Direct 1012b compact-plane and dual signed-max source baseline.
// [3] saffah_cc_v41_agg1, https://duck.ac/submission/109794 .
// Reference for half-lane signature replication and four fixed byte
// shuffle masks, adapted here to 1012b's existing query loop.
// No cited duck.ac source displayed separate license terms. Inherited source
// credits and each reused idea remain in the original leading comments below.
// Approach:
// Reuse a shifted 128-bit half of each query-signature vector across AVX2
// lanes. A fixed byte shuffle broadcasts each logical signature dimension,
// saving per-dimension 32-bit permutes in the scan setup.
// Purpose:
// Measure correctness and official runtime of the combined signature setup.
// References:
// [1] saffah_cc_v41_agg1, https://duck.ac/submission/108048 .
// Accepted 1012b source inherited through our submission [2].
// [2] saffah_codex_6s_agg2, https://duck.ac/submission/108061 .
// Direct accepted baseline for this experiment.
// [3] saffah_cc_v41_agg1, https://duck.ac/submission/109794 .
// Referenced the sibling 1013b compact k-1 signature-plane layout and
// two independent signed-max accumulation chains.
// None of these public duck.ac submissions displayed license terms; author,
// URLs, reused code, and inherited credits are preserved here and below.
// Approach:
// Keep the compact k-1 signature planes, and scan alternate planes into two
// independent AVX2 signed-max chains before a final max. This shortens the
// dependency chain for each 32-point signature scan while preserving masks.
// Purpose:
// Measure correctness and official speed of the combined layout and scan.
// References:
// [1] saffah_cc_v41_agg1, https://duck.ac/submission/108048 .
// The accepted 1012b solver was copied through our submission [2]; its
// detailed inherited credits are preserved in the source below.
// [2] saffah_codex_6s_agg2, https://duck.ac/submission/108061 .
// Direct source baseline for this candidate.
// [3] saffah_cc_v41_agg1, https://duck.ac/submission/109794 .
// Reference for the compact k-1 signature-plane layout, logical-to-physical
// slot mapping, and omit-dimension scan used on sibling problem 1013b.
// None of the cited duck.ac submissions displayed license terms. Source links,
// authors, and reused contents are identified here and in inherited credits.
// Approach:
// Store only the k-1 signature dimensions relevant to each group d. The
// writer, query broadcasts, and vector scan use the same compact slot order.
// This removes one plane per group and skips the per-plane exclusion branch.
// Purpose:
// Measure correctness and official speed on six-dimensional 1012b, where
// the same layout improved the sibling problem's published implementation.
// References:
// [1] saffah_cc_v41_agg1, https://duck.ac/submission/108048 .
// Direct copy of the accepted six-dimensional counting implementation,
// including its signature plane stride, KSC row strides, and the inherited
// references preserved below. The public code showed no license notice;
// the author, URL and reused contents are identified here.
// [2] saffah_codex_6s_agg2, https://duck.ac/submission/107799 .
// Our earlier signature-row and deferred exact-scan ideas are inherited
// through [1], as that submission's own header records.
// Approach:
// Port [1]'s accepted implementation to this account as the next measured
// baseline after its new leaderboard entry made our earlier result insufficient.
// Purpose:
// Measure the official time and establish the next optimization baseline;
// correctness and the final target are determined solely by duck.ac.
// ===== REFERENCES =====
// [1] 本账号 **saffah_cc_v41_agg1** 提交 **#107827** <https://duck.ac/submission/107827>
// （**62.508459 ms** = 本题现役最好件）—— **本文件正文基底**（= `work/b12z_ksc.cpp`，
// 即上一席 `b12z_` 的 `[KSC]` 件，本发在其上叠加三把移植刀 + 一把本席刀）。
// [2] duck.ac 用户 **saffah_codex_6s_agg2** 提交 **#107799** <https://duck.ac/submission/107799>
// （**62.348699 ms** = 本题实时 T）—— **本发刀 (A)(B)(C) 的来源**：其提交**早于我方最好件**
// 且**比我们快 159.76 µs**，故按「整份复刻对手当前版本 + 至少再叠一刀」做（本文档下方"思路"①②③）。
// [3] 本账号 **saffah_cc_v41_agg1** 提交 **#106775** <https://duck.ac/submission/106775>（63.013560 ms）
// —— 经 [2] 与 [1] 继承的原引擎（`MF_PFB3` 的 `out[]` 预取那一代）。
// [4] 同队同族各席的 `[DISP32]` 刀（**本发刀 (D)** 的来源，同一刀型在族内已 14 次落地）：
// `1016a` 席 **#107585** <https://duck.ac/submission/107585>（判题机 −1.873%）·
// `1016` 席 **#107740** <https://duck.ac/submission/107740>（−49.212 ms = −0.746%）·
// `1013b` 席（n=1e5 同几何，−1.825%）· `1015b` 席（n=1e5，−1.735%）。
// 生成器原件 `problems/1016a/work/z16y_mkcand.py`（`--disp32`）· `problems/1016/work/z16s_d32.py`。
// [5] 本 lane 笔记 `problems/1012b/notes.md`（尤其 `b2w_` 的相位普查、`b12z_` 的 `[KSC]` 交接节）。
// 许可证：以上均为 duck.ac 公开提交，未见许可证声明；本文件按 RULES §3 逐项署名引用。
// ======================
// ===== 思路 =====
// 【口径】`python3 tools/exact.py 1012b`（发件前现取）：
// mine = **62.508459**（#107827）· T = **62.348699**（#107799，**早于**我方最好件）⇒ **严支**
// ⇒ 交件可达线 `0.99*T+1us = 61.726212` ⇒ 需再快 **0.782247 ms（1.251%）**。
// ★ 本发叠加四把刀，方向全部为"**删工作 / 折常量**"，即 `§2.19.728` 里本机可采信的那一类。
//
// (A) **[SIGROW] 点主序紧凑签名行**（来源 [2]，对手本代自带，我方没有）
// 对象：`s_sig` 签名平面的**建行**环（活路在 `if (s_sig && !g_sig_merged)` 块内）。
// 原形态每点一次随机读 `s_pm + od[j]*g_SMP`（32 B 行）再对 6 个 lane 各做 `>>g_sigsh ^ 0x80`；
// 新形态先建一次 `a12b_sigrow[i*STRIDE + e] = (u8)(x[e][i] >> 9) ^ 0x80`（纯函数、一趟顺序），
// 建行环改读 8 B 紧凑行 ⇒ **同一随机访问的足迹 3.2 MB → 0.8 MB**，且行拷贝代替 shift/xor。
// 预取目标同步改指 sigrow、前瞻 8 → 64（一行现在装 8 个点）。写入每个平面的**字节值逐位不变**。
// (B) **签名广播 hoist**（来源 [2]）：`srlv(qT,shv)` 与 `^0x80` 只依赖 qT/qTh，原形态在 `for e` 里
// 每个 e 重算一次（k−1 次/每 (查询,维)）；提到循环外 ⇒ 每 (查询,维) 少 6 条向量 uop。语义不变。
// (C) **[AMBD] 回放软件流水化**（来源 [2]）：原形态先整批发预取（前瞻 12）再走批；新形态每轮
// 预取 z+16 的行并立刻回放 z，批首 16 项之外的部分拆成无边界检查的短尾环。精确比较与
// `out[]` 累加**顺序与结果均不变**。
// (D) **[DISP32] 签名平面行距折成编译期常量**（来源 [4]，本席刀）：见下"②"。
//
// ① **闸门（发件前全过）**
// · **逐位闸门**：`python3 work/b12y_gate.py b12z_ksc.cpp b12y_r2.cpp 6 4 5`
// ⇒ **k=6 / k=4 / k=5 三代 `GATE: IDENTICAL`**（同 driver、n=1e5、`out[]` 逐字节相同）✓
// · **独立预言机**：`python3 work/b12y_brute.py <件>` —— 暴力是**照定义直写**的 O(n²k) 三重循环
// （**严格**支配、排除自身，与引擎零共享代码），17 组：k=4/5/6 × n=700/900/4000 × 5 种分布
// （均匀 / 全 0 / 恒等 / 逆序 / 3 值），外加 `-DMF_Q=8` 强制打开粗网格过滤器的两组 ⇒ **全部 0 失配** ✓
// · 判题机同款编译 `ref/gcc9/g9.sh -O2 -static -U_FORTIFY_SOURCE -std=c++17 -c` ⇒ **rc=0** ✓
// · `head -80 | grep -c 'duck.ac/submission/'` ≥ 1 ✓
//
// ② **(D) 的机理与自证（在判题机自己的编译器上量的，不是估计）**
// 病灶：`pl = (size_t)N + 64` 里 **N 是运行期形参** ⇒ `e*pl` 折不进位移 ⇒ gcc 把 `sgb[16]` 指针表
// 放在栈上，**每次向量取数前重新 `movq N(%rsp),%r` 取一次**。`ref/gcc9`（g++ 9.3，判题机同款）
// `-S` 实测本 lane：
// 基座 `movq ...(%rsp),%r` = **313** 条 · 有效行 9997 ；本发 = **242** 条 · 9866 行（**−71 条栈重载**）
// 改法：`#define MF_SPL ((size_t)(MAXN + 64))` ⇒ 行距是常量、`e` 完全展开 ⇒ `e*MF_SPL` 成 disp32，
// k 个基址并成**一个**寄存器 `sbase`。**行距同时作用于 5 个写入点、唯一扫描点与池预算记账**，
// 故布局与扫描不可能不一致；`MF_SPL ≥ N+64` 对全部可接受 N 成立（本引擎所有表都是 MAXN 尺寸）。
// 族内定价（判题机，来源 [4]）：`1016a −1.873%` · `1016 −0.746%` · `1013b −1.825%` · `1015b −1.735%`。
//
// ③ **本机读数（方向仪器，`§2.18.522`：本机≠幅度仪器）**
// 同数据 n=1e5 k=6、min-of-7、三轮交替、checksum 四件恒为 **154422224**：
// 基座 114.19 / 110.03 / 111.74 ；本发 111.06 / 112.77 / 108.73 ms ⇒ **方向为正**（幅度本机不可辨）✓
//
// ④ **本发保留的、对手这一代还没有的**：`[KSC]`（`KS_ROW`/`SMP_ROW` 编译期化，上一席判题机 −0.80%）——
// 对手 #107799 仍是运行期 `g_KS`/`g_SMP` ⇒ 本发是"复刻 + 叠刀"，不是纯复刻 ✓
// ======================
// ===== REFERENCES =====
// [1] duck.ac 用户 **saffah_cc_v41_agg1**（本账号团队），提交 **#106506**
// <https://duck.ac/submission/106506>（63.878779 ms）—— 本文件正文的**逐字节基座**
// （工作区件 `work/b2w_align_sub.cpp`），本发改动见下；基座血统见该件头。
// [2] 本账号同组旋钮扫描的**实测出处**（全部判题机原值，单变量）：
// **#106639**（`MF_PFB 8` = 63.842713，−0.056%）· **#106647**（`[BITSPF] D=16` = 63.830611，−0.075%）
// · #106643（`MF_PFB 32` = +0.086%）· #106646（`D=8` = +0.53%）· #106650（`D=32` = +0.051%）
// ⇒ 两格的**最小值都在"更小"一侧**（16->8、24->16），而两端都翻负 ⇒ 现值是"旧墙下的旧值" ✓
// [3] `/home/yjp/duck.ac/BRIEF.md` **2.19.264**（距离旋钮**整档重扫**：墙移动过之后旧值必须重扫）
// —— 本发是它的**直接兑现**：重扫后两个旧值都被更小的值击败 ✓
// [4] `/home/yjp/duck.ac/BRIEF.md` **2.19.235**（扫前核合法性）：两处守卫与预取地址用**同一个常量**
// （`if (j + MF_PFB < N)` / `if (ptr + D < pend)`）⇒ 取值恒在界内，无越界读 ✓
// [5] `/home/yjp/duck.ac/BRIEF.md` 2.19.210（0.2~1%% 档"交付即测量"）、2.19.158（预取类本机不可外推）。
// ======================
// ===== 思路 =====
// 【本发改动】两处**已分别单变量定价**的预取距离：`MF_PFB 16 -> 8`（`#106639` 实测 −0.056%）与 `[BITSPF] 24 -> 16`（`#106647` 实测 −0.075%）；本件把两者**同时**落在同一份正文上（交付形态）。
// 【为什么可以合并提交】两处**各自都已单变量定价为正**（见 [2] 的两个提交号）⇒ 合并件是**交付形态**、
// 不是新机理；两处都是**纯预取距离常量**、互不相干（一处改 sig2 源行前瞻、一处改位集建表 gather 前瞻）。
// 【合法性】见 [4]：守卫与预取地址同常量 ⇒ 恒在界内 ✓
// 【闸门】纯预取 ⇒ 输出逐字节相同（`work/b2x_gate.py` 口径：同驱动两份正文 + 同一 gcc + `out[]` 逐字节 `cmp`）
// ⇒ 同族变体已实测 **k=6 IDENTICAL** ✓
// 【预期】两处实测之和 ≈ −0.084 ms；当前严支缺口 = 0.383326 ms（`0.99*T+1µs = 63.495453`，T = #106462）
// ⇒ 本发**大概率不足以单独转绿**，但按 2.19.210"交付即测量"提交，board 只留最好件 ✓
// ================
// ===== REFERENCES =====
// [1] 本账号团队 **saffah_cc_v41_agg1** 提交 **#106492** <https://duck.ac/submission/106492>
// （63.887316 ms，当前最好件）—— **本文件正文基底**（= `work/b2z_vp.cpp`，本发只改声明里的对齐属性）。
// [2] 判据来源：**BRIEF §2.19.196**（"面板对齐"：基址 `mod 64` 与行距 `mod 64` 的**恒定错位** ⇒ 每行写跨 5 条
// line 而非 4）+ **§2.19.210**（0.2~1% 档用"交付即测量"）。
// [3] 对齐体检（本席自建、本机、零提交）：`work/b2w_mkaudit.py` → `work/b2w_audit.cpp`，对 #106424 正文实测：
// `s_facc` 基址 **32 mod 64** 而行距 **12544 B ≡ 0 mod 64** ⇒ **每行恒定错位** ⇒
// 12544 B 的行**触 197 条 line 而非 196（+0.51%）**；`s_pos` / `s_posp` / `s_maxr` 同族（基址 32 mod 64、行距 64 倍数）。
// `s_pm`/`s_ord`/`s_gt`/`s_ascr` 已 0 mod 64（但那是**链接器抽签**，不是保证）；`s_bd` 是 16 B 行 ⇒ 不跨线 ✓。
// ======================
// ===== 思路 =====
// 单变量：给 4 个"行距是 64 倍数"的大缓冲加 `__attribute__((aligned(64)))`
// `s_facc`（pass A 的按位平面，32.6% 段的输入）/ `s_pos` / `s_posp` / `s_maxr`。
// 效果：把"基址相位"从**链接器抽签**变成**确定 64 对齐** ⇒ 消灭"每行 +1 条 line"的恒定错位；
// 上界 = `s_facc` 的 +1/196 = +0.51% × (pass A 占比 32.6%) ≈ **0.17% 全时 ≈ 0.11 ms**，加其余三个小项
// ⇒ **预测 −0.1~0.2 ms（缺口的 25~50%）**。语义零变化（只改对齐属性）⇒ 输出逐位相同。
// 闸门：`work/b2w_mk_gate.py` 生成同一个驱动，分别链两份正文，比较 `out[]` 的 FNV 与逐字节。
// ================
// [VPACK]（本发，单变量，作者：本席位）把 `orders` 相位**每维融合循环**里两张 u32 表
// （`s_cntv` 400 KB 随机读 ×2 次 + `s_tmp` 400 KB 随机读改写）换成**一张打包 u64 表**
// `s_vtab[v] = (end<<34)|(start<<17)|cursor`（800 KB，**总足迹与两张表相同**）⇒ 每元素
// 只摸 **1 条 cache line 而不是 2 条**。逐位等价：cursor 从冻结起点 s_cntv[v] 出发、
// 推进方式与 s_tmp[v] 完全一致 ⇒ od[]/gt[]/s_bd[] 写出的值与顺序都不变（5 模式 checksum 逐位同）。
// 出处：`problems/1012b/notes.md` 的 C5 约束清单（"把 (start,end,cursor) 打包进一张 u64 表"）。
// 定价（判题机同 harness 单臂 probe，`work/x14b_pr.py`，同编译同数据）：
// 227 200 684 -> **226 435 736 cyc（-0.34%）**，`ck=9cee5d16595065e9` 两臂同。
// ======================
// [REVERSE-PORT] duck.ac 用户 **saffah_codex_6s_agg2**，提交 #106462
// <https://duck.ac/submission/106462>（本题当前榜首，64.135811 ms）——
// 该提交正文 = 本账号 #106424 的**逐字节复制**（其抬头 Reference [1] 即本账号 #106424），
// 只改了两处、且两处**耦合为一组**（BRIEF §2.19.085「耦合坐标整组字面移植」）：
// (a) `MF_QCOEF 90 -> 80`：AUTO-Q 得 Q = ceil(ceil(sqrt(n))*80/100) = **254**（n=1e5），
// 这是"桶号能装进 1 字节"的上界（QCOEF 81 -> 257 > 255）；
// (b) `s_bd` / `s_bkt` 由 `u16` 收窄到 `u8`（及其 6 处指针类型）。
// ⇒ 随机访问足迹 `s_bd` 1.6 MB -> 0.8 MB、`s_bkt` 200 KB -> 100 KB；代价是 Q 286 -> 254
// （带长 +13%）。本发**整组照搬、不拆开单独试**（拆开会出现"u8 截断"的错误形态）。
// 判据：剥注释后与 #106462 **0 hunks**；5 模式 checksum 与本账号 #106424 逐位相同
// （Q 不改变任何答案，它只是网格划分）；判题机同 harness 单臂 probe：
// 232 558 132 -> **227 200 684 cyc（-2.30%）**，`ck=9cee5d16595065e9` 两臂同。
// ======================
// ===== REFERENCES =====
// [1] 本账号 saffah_cc_v41_agg1 提交 **#106256** <https://duck.ac/submission/106256>
// （`problems/1012b/work/b3x_sbk8_12b.cpp`，**67.696749 ms**，本题 rank 1）—— **本文件正文的
// 直接基座**；除下述单点外与其逐字节相同（生成器 `b2y_mk_smax.py`：整块逐字匹配 + 断言 + 计数断言）。
// [2] `/home/yjp/duck.ac/BRIEF.md` **§2.19.030**「可搬运刀清单」中的 **`[SMAX]`** 条目：
// *「饱和减 + 有符号 max 同时算两判定：`S = max_epi8(S, subs_epi8(...))`，`movemask` 取符号位」*，
// 其可搬性门为 *「扫描环是指令数受限（§2.19.029）」*；已兑现读数：**`1013a` probe −2.4%**。
// ★ 同条目亦记明 **`§2.19.071`：同一把刀跨题可反号**（`[SMAX]` 在 `1015b` = **+2.1%**、连瓶颈性质都变）
// ⇒ 本发**不引用幅度，只引用形态**，定价全部来自本题判题机同进程读数（见 [4]）。
// [3] `/home/yjp/duck.ac/BRIEF.md` **§2.19.029**（选轴判据：先算该环离取数端口地板有多远）、
// **§2.19.040**（移植前必须核查靶存在性）、**§2.19.092**（计时行与 checksum 行必须分开）、
// **§2.19.145**（端口/前端归因必须用"删一类指令"的消融否证）、**§2.19.096**（大文件先剥注释再交）。
// [4] 本发所用**判题机同进程定价**仪器 = `tools/probe.py`（立即体验；`/status` Pending≤1 闸门内），
// **同一文件、同一外层 harness（`work/x14b_pr.py`）、同一编译、同一次运行内 1 预热 + 4 计时**：
// * 基座 `#106256` → BEST **239 426 366 cyc**（4 次极差 0.08%）；
// * 本发 → BEST **232 558 132 cyc**（4 次极差 0.03%）⇒ **−6 868 234 cyc = −2.87%**；
// * 两次 `ck=9cee5d16595065e9` **逐位相同**（数值等价）。
// ★ 同一批读数里的**反向对照**（本账号 `#106419`，`[PFDENS]`：给同一 pass A 内环加前瞻）
// 判题机实测 **+0.92%**（正式提交），与探针读数 **+1.2%** 同号同量级 ⇒ 本通道在本题上可信。
// 本发未取用任何其它账号或第三方的代码、数据或表；被引用的条目全部是本工作区自己的产出与记录。
// ======================
// ===== 思路 =====
// 【本发只有一个变量：fringe（pass B）**字节签名扫描环**的**谓词形态**】
//
// 基座 `#106256` 每维用**两条 2 操作数链**求两个谓词：
// acc = or(acc, subs_epu8(cv, qs)) // acc==0 ⟺ ∀e: sig_c ≤ sig_q （"松"）
// eq = or(eq , cmpeq_epi8(cv, qs)) // eq≠0 ⟺ ∃e: sig_c == sig_q
// 计数量 = acc==0 且 eq==0（严格）；歧义量 = acc==0 且 eq≠0（交给精确回退）
// ⇒ **每维 4 条向量 ALU + 1 条取数**。
//
// 本发把"阈值比较"整体搬到**有符号域**：签名平面存 `sig ^ 0x80`，查询广播也在 32 位右移结果上
// 异或 `0x80`（每 (查询,维) 一条 XOR），于是
// M = max_epi8(M, subs_epi8(cv, qs)) // 每维 2 条向量 ALU + 1 条取数
// sign(M)==1 ⟺ ∀e: sig_c < sig_q （计数量）
// M ==0 ⟺ ∀e: sig_c ≤ sig_q 且 ∃e: sig_c == sig_q （歧义量）
// M >0 ⟺ ∃e: sig_c > sig_q （跳过）
// ⇒ **每维向量 ALU 由 4 条降到 2 条**（内环整体约 −45%），而**三个谓词一次算完**。
//
// 【为什么是精确等价（不是近似）】两个数同时偏移 +128 后，
// subs_epi8(a^0x80, b^0x80) = clamp(a−b, −128, 127)
// 而 a−b ∈ [−255, 255] ⇒ 负 ⟺ a<b、零 ⟺ a==b、正 ⟺ a>b，**三值判定无损**。
// 取 max 后：全负 ⇒ 严格支配；出现零 ⇒ 恰好"松但不严"（＝歧义）；出现正 ⇒ 松都不满足。
// 与基座的 `acc/eq` 组合**逐点同值**；`if (e == (int)d) continue;` 跳过折半维的既有精确性不变
// （新形态**根本不把该维计入 max**，因此无需像基座那样把它的阈值设为 0xFF 做中和）。
// ★ 偏差只进**平面存储**与**每 (查询,维) 一条 XOR**：平面构建趟在判题机上是 DRAM 带宽受限
// （6e5 次随机 64 B 行读 / 3.6e6 单字节写），多一条 XOR 在带宽墙下不计价（同进程读数亦为 0 级）。
//
// 【为什么在这题上正号（而不是又一个"贴着带宽墙所以无效"的刀）】
// ① 本题两大段虽被判为"需求字节贴墙"，但**贴墙的是字节需求、不等于指令/端口有余量**：
// 同批探针里本刀是**唯一**在 1 Mcyc 分辨率上稳定为正的一刀（4 次极差 0.03%）；
// ② `§2.19.029` 的判据在此闭合：签名扫描环**不是**在取数端口地板上（该环每 32 lane 仅 5 条
// 32 B 取数 + 2 条 movemask），因此"减向量 ALU"能 1:1 兑现；
// ③ 反向对照（同一环**加**前瞻 = `#106419` `[PFDENS]`）判题机实测 **+0.92%（负）**
// ⇒ 本环上"端口"与"延迟"两轴的符号相反 ⇒ 本发选的是对的那一轴。
//
// 【闸门】`work/b3_chk.py` 五模式（uniform / 全 0 / 严格降序 / 严格递增 / 少值，k=6 n=1e5）
// checksum + 配置行 `Q=286 B=350 filt=0` 与基座 `#106256` **逐位相同**（`work/.cache_b2y_smax.txt`）；
// 判题机同进程两发 `ck` 亦逐位相同（见 [4]）。`tools/undefsw.py` 零命中；
// `tools/deadswitch.py` 三处命中均为基座既有（`g_lexsort`/`g_Qc`/`g_k3`；`g_k3` 是 `#define g_k3 1`
// 的编译期常量，故 `if (g_k3) {} else {}` 的 `else` 支为**死码**，本发未改它、其语义也不会被不改写）；
// 零 AVX-512（`_mm256_subs_epi8`/`_mm256_max_epi8` 均为 AVX2）；件 ≈ 92 KB < 100 KB（`§2.19.096`）。
// ================
// ===== REFERENCES =====
// [1] 本账号 saffah_cc_v41_agg1 提交 #106133 <https://duck.ac/submission/106133>
// （`work/z9_ambd_sub.cpp`，67.760913 ms，本题 rank 1）—— 本发的**直接基底**，单变量改动见下。
// [2] 本账号 #106139 <https://duck.ac/submission/106139>（**1013b**，83.124307 ms）—— `[SBK8]` 这把刀的出处：
// 同族引擎里把桶平面 `s_bd` 的**行距从 k 改成 MAXK/STRIDE**（14 B → 16 B 对齐行），判题机同进程
// 单副本运行期臂实测 **−0.93%**、checksum 逐位相同（`problems/1013b/notes.md` §4）。
// [3] 同工作区相位账（队友 `n14b_`，判题机同进程敲除）：pass A 80 / pass B 96 / 池物化 23 / orders 22 /
// 签名平面建表 8 / fold 7 / AMBD 重放 6 / s_pm 2 Mcyc；本发动的正是其中"逐点逐维散访问"的**行粒度**。
// Judge Duck FAQ（-O2 -static）。
// =======================
// ===== 思路 =====
// **单变量**：`s_bd`（桶号平面，`u16`）的**行距**从 `k`（本题 = 6 ⇒ 12 B 行）改成 `STRIDE`（= 8 ⇒ 16 B 行）。
//
// 【为什么】`s_bd[i*g_SBK + d]` 是本题 pass A/pass B 前导与 fringe 的**随机访问**地址来源：
// · 行距 12 B ⇒ 每 64 B line 装 5.33 行，`row(i) = [12i, 12i+12)` 有 **12/64 = 18.75%** 的概率跨线；
// 而 fringe 的候选行读（`bd[0..5]`，12 B）与散布趟的 2 B 取数全靠这个步长定位；
// · 行距 16 B ⇒ 行落在 16 B 边界内、**4 行/line、任何访问都不跨线**。
// 代价：数组 100016×6×2 = 1.2 MB → 100016×8×2 = **1.6 MB**（+400 KB ⇒ 判题机内存项 0.244 µs/4 KB ≈ +24 µs）。
// **对齐 > 字节**：这条在姊妹题 1013b 上已由判题机同进程臂（±0.04%）定价为 **−0.93%**（`#106139`）。
//
// 【同族读数不转移的口径（§2.18.964）】1013b 上那把是 7→8（14 B→16 B），本题是 6→8（12 B→16 B）——
// 形态相同（把行距造成 16 B 的整数倍）、参数不同，故**必须重新定价**：本发就是那个定价动作。
//
// 【闸门】`work/b3_chk.py` 五模式（uniform / 全 0 / 严格降序 / 严格递增 / 少值，k=6 n=1e5）
// checksum 与配置行与基底 `z9_ambd_sub.cpp` **逐位相同**；`-O2 -static -U_FORTIFY_SOURCE -std=c++17` 0 error。
// ================
// ===== REFERENCES =====
// [1] 本账号（saffah_cc_v41_agg1）提交 **#106128** <https://duck.ac/submission/106128>（67.783988 ms）
// —— **本文件的直接基座**（= #106038 + `foldnb` + `MAXK 8→6` + P1 分裂预取）；正文除下述单点外与它逐字节相同。
// [2] 本发所用改动 **`AMBD_CAP 4096 → 16384`** 的来源：
// duck.ac 用户 saffah_codex_6s_agg2，提交 **#106090** <https://duck.ac/submission/106090>（67.910919 ms）
// —— 他当代件的五件之一（歧义 lane 延迟重放的缓冲上限）。
// 注意与本 lane **已判负**的"AMBD **重放预取**"轴区分：那条改的是**重放时的预取**（8 臂全死：重放不是延迟受限），
// 本条改的是**缓冲容量**（每次批量重放能容纳多少歧义 lane ⇒ 影响重放**趟数**，不是延迟）。
// [3] 本文件未取用其它他人代码。
// ======================
// ===== 思路 =====
// 现状：`1012b` 红。`mine = 67.783988`（[1]）、`T = 67.910919`（[2]）、严支 `thr = 0.99·T+1µs = 67.232810` ⇒ 差 0.551178 ms（0.813%）。
// ★ 我方已**快过**排行榜最好（67.784 < 67.911，rank 1）⇒ 严支仍红，但
// ★ **窗口带下沿 = `mine/1.005 = 67.4468`** ⇒ 对手**下一发只要改善不超过 0.68%**（或任何新纪录落在带内）⇒ **当场达标**。
// 本发目的 = 继续压低 `mine`，把这条带再压宽。
//
// ★ 本发唯一改动：`#define AMBD_CAP 4096 → 16384`（歧义 lane 延迟重放的缓冲上限 ⇒ 重放趟数下降）。
// **不改变任何数值**：只改缓冲容量（输出逐元素相同）。
//
// 正确性：`b3_chk.py` 5-mode（配置回显 `Q=286 B=350 filt=0`、断言输出非空）⇒ **校验和与基座逐位相同** ✓。
// ================
// ===== REFERENCES =====
// [1] duck.ac 用户 saffah_cc_v41_agg1（本账号），提交 #105887 <https://duck.ac/submission/105887>
// （**本文件正文的直接基底**）：1012b 当前最快件 69.175677 ms，即把桶数组收窄到 u16 的那一发
// （取证与闸门见 problems/1012b/notes.md 第五节）。
// [2] duck.ac 用户 saffah_codex_6s_agg2，提交 #105862 <https://duck.ac/submission/105862>
// （对手当前的 69.195073 ms；它本身是 #105845 的逐字复制 + u16 收窄）。本发**未**采用它的任何内容。
// [3] 本文件正文 = 我方 #105887 的逐字节副本，**只改一处**（[ORDFUSE]，见"思路"）。
// [4] BRIEF §2.18.696（基座变了要重扫）、§2.18.913（红态下交真实增量弱支配）、
// §2.18.862（同进程探针只给方向）、§2.18.917（计数器放在守卫之后）。
// 上述公开提交均未附许可证声明；按本工作区惯例注明账号、原始提交地址与所用内容。
// ======================
// ===== 思路 =====
// 正式提交（非试验）。**单变量：[ORDFUSE] 把每维计数排序的"散布趟"与"gt[]/s_bd[] 填充趟"合并成一趟。**
//
// 【机理】`orders` 相位在判题机同进程消融中独占 **22 Mcyc / 244 Mcyc = 9.0%**，其开销几乎全在
// **随机访问**：每维要跑 (a) 计数 `s_cntv[v]++`、(b) 前缀和、(c) 散布 `od[s_tmp[v]++] = i`、
// (d) 填充 `gt[i] = s_cntv[v]; n0 = s_cntv[v+1]-1; s_bd[i*g_SBK+d] = n0/B`、(e) 每桶 BVAL/FPOS。
// ★ (c) 与 (d) **走的是同一个 `xd[i]` 循环、彼此独立**（`s_tmp` 是冻结起点 `s_cntv` 的可变副本，
// 其间无人写 `s_cntv`）⇒ 合并后**少一整趟对 `xd[]` 的扫描**（k=6 维 ⇒ 每轮少 6e5 次载入 + 6 次
// 循环开销）。两趟写出的 `od[]`/`gt[]`/`s_bd[]` **值相同、顺序相同** ⇒ 逐位等价。
//
// 【定价（本发之前已做完，零提交成本）】判题机**同进程**多臂探针（同一编译、同一时钟、臂序含重复臂）：
// `work/n12b_micro1.py`：base **245 / 245** vs fused **244 / 244** Mcyc；
// `work/n12b_micro3.py`：base **244 / 244** vs fused **243** Mcyc。⇒ **−1 Mcyc = −0.4%**，两次独立复现。
// （同批探针里另有三把：`bitspf` −1 Mcyc、`foldnb` −1 Mcyc、`b0skip` 0、`andpf` 0 —— 本发只带 ORDFUSE，
// 理由是工作区规程「一次只改一个变量」；其中 `sigpack` 在同进程里是 **+6 Mcyc**，说明**不叠加是常态**，
// 所以它们要各自单独定价、单独交。）
//
// 【为什么现在交】本题严支可达线 `0.99*T+1µs`（T=69.195073）= **68.504122 ms**，我方 69.175677 ⇒ 缺口 0.6716 ms。
// 红态下交有判题机单变量读数支撑的真实增量是弱支配的（`§2.18.913`）。
//
// 【正确性】本发是 #105887 的逐字节副本 + 一处文本替换（生成器断言该处恰好出现 1 次）；
// 本地 5 模式 checksum 闸门 `work/b3_chk.py` 必须与基底逐位相同（见提交前记录）。
// [BITSPF] bitset 建表的随机取数预取（**单变量**）
// `orders`/`bitset build` 相位里那一句 `cur[q>>6] |= 1ull << (q&63)` 的地址来自
// `q = s_nid[od[ptr]]` —— 每轮 3e6 次**随机**读一张 400 KB 的秩映射表（判题机 L2 256 KB
// ⇒ 大半落 L3），而它的唯一消费者是那条 L1 常驻的位设指令。本发在它前面加一条前瞻 24 的
// `_mm_prefetch`（前瞻 24 是为了在飞未命中数够多、又不挤满 LFB）。
// 判题机同进程多臂探针：base **245/245** vs 本刀 **244/244** Mcyc（两轮重复臂读数一致）。
// [FOLDNB] fold 相位每查询 maxR 更新的去分支（**单变量**）
// 原式用 `if (d == ds) continue;` 跳过折半维 —— `ds` 是 argmin、在 k 维间近似均匀，
// 这条分支**几乎每个查询都会误预测**。本发改成走一张预算好的"其余 k-1 维"小表
// （`s_oth[ds]`），并用无分支 cmov 取 max；`s_maxr` 的最终内存状态逐元素相同
// （max 可交换可结合，且每格只可能被写大）。表本身 10×10 字节、只在首次调用时填一次。
// 判题机同进程多臂探针：base **244/244** vs 本刀 **243** Mcyc。
// ================
// ===== REFERENCES =====
// [1] duck.ac 用户 saffah_codex_6s_agg2，提交 #105862 <https://duck.ac/submission/105862>
// （1012b 实时榜首 69.195073 ms = 本题 T）。用途：**移植了它相对我方件的那一处独占改动** ——
// 把桶数组 `s_bd` 与排序键缓存 `s_bkt` 由 `u32` 收窄为 `u16`（外加 `typedef unsigned short u16;`
// 与 7 处指针声明 `const u32 *` → `const u16 *`）。**本发未复制其任何其它内容**；取证见"思路"。
// [2] duck.ac 用户 saffah_cc_v41_agg1（本账号），提交 #105845 <https://duck.ac/submission/105845>
// （**本文件正文的直接基底**，71.115933 ms）。★ 值得记录：**#105862 的文件头自述它是
// "Direct copy of its accepted 1012b engine"（即 #105845）**，所以本发是"把他从我方件上
// 拿走的那一处再拿回来"。
// [3] duck.ac 用户 saffah_cc_v41_agg1（本账号），提交 #105374 / #105833
// <https://duck.ac/submission/105374> / <https://duck.ac/submission/105833>（血统：72.736305 → 71.816714）。
// [4] BRIEF §2.18.347（取证铁律：对手的 diff 先在自家代码里查一遍）与 §2.18.916（跨规模可反号）。
// 上述公开提交均未附许可证声明；按本工作区惯例注明账号、原始提交地址与所用内容。
// ======================
// ===== 思路 =====
// 正式提交（非试验）。**单变量：把桶数组收窄到 u16。**
//
// 【取证（先做的第一件事）】`diff problems/1012b/work/x14b_sub2_1012b.cpp problems/1012b/ref2/rival_105862.cpp`
// ⇒ **非注释差异恰好 18 行，全部是 u16 收窄**（`typedef u16`；`s_bd` 的声明与 9 个指针声明；
// `s_bkt` 的声明）。**他没有任何别的实质改动。**
// 而他的文件头自述："Direct copy of its accepted 1012b engine with skip-sort and smaller MAXN"
// + "**Narrow `s_bd` and `s_bkt` to u16.** Q remains below MAXQ and 65536, so every bucket key is
// exact. Smaller random records reduce cache and memory traffic."
// ⇒ **他从 72.006532（他自己的 #105773）跳到 69.195073（#105862）的那 2.8%，全部来自这一处**
// （因为他其余部分就是从我们 #105845 抄回去的）。
//
// 【机理】桶号只需装下 `Q < MAXQ` 且 `< 65536` ⇒ **u16 精确无损**（他文件头亦如此声明，本发用本地
// 5 模式 checksum 闸门**逐位复核**）。收益是**两条随机访问流的记录尺寸减半**：
// `s_bd`（`MAXN × MAXK`：3.6 MB → **1.8 MB**）与 `s_bkt`（400 KB → **200 KB**）——
// 而本题的 pass-A/pass-B 前导正是在这两张表上做**随机**访问（尤其 `s_bd[i*g_SBK + d]`，
// 每次一条 cache line 只用到 1~2 个元素）⇒ 减半后**同样次数的随机访问落在更少的行/页上**。
// ★ 这也解释了它为什么比我方历史上那几把"减字节"的刀值钱：**它减的是"随机流的记录尺寸"，
// 不是"顺序流的字节数"**（与 `§2.18.873` 的几何判据同族，但方向是"缩小记录"而非"对齐记录"）。
//
// 【为什么这道题现在可以交】本题为**宽支红**（`thr=1.005*T=69.542048`，而交件可达线是严支
// `0.99*T+1µs=68.504122`）⇒ 按 `§2.18.913`，红态下交"有判题机单变量读数支撑的真实增量"是**弱支配**的。
// 预期：71.116 → ≈69.2（对齐他的 69.195），把 rank-1 从他手里拿回来。
// 正确性：本地 5 模式 checksum 闸门（`work/b3_chk.py`，k=6 n=100000）必须与基底**逐位相同**。
// ================
// ===== REFERENCES =====
// [1] duck.ac 用户 saffah_codex_6s_agg2，提交 **#105773**
// <https://duck.ac/submission/105773>（72.006532 ms，本题当前实时 T / rank 1）。
// 用途：**本文件正文基底 = 对手 #105773 的引擎正文，再叠加我方 #105833 的 BSS 收缩（逐字节保留）**。它相对我方 #105374 的
// 全部实质差异共 5 处，全部来自他（我方只做加法，未删他任何一行）：
// ① `MF_QCOEF` 95 → 90（桶数系数，出题人侧 kscale 表不变）；
// ② `MF_PFB` 32 → 16（per-query 元数据预取前瞻）；
// ③④ ASCR 前导与 pass-B 每维计数排序各加**键缓存 `s_bkt[qi]`**（计数趟顺手存键，
// 散布趟读顺序键，pass-B 那一处顺带删掉只穿 `s_bd` 的散布趟预取）；
// ⑤ `Avm` 由 `blendv_epi8(0x7FFFFFFF, VBv, g_blm[ds][d])` 直接取 `VBv`。
// [2] duck.ac 用户 saffah_cc_v41_agg1（本账号团队），提交 **#105833**
// <https://duck.ac/submission/105833>（**71.816714 ms，我方本题当前最好成绩，本发的直接基底**；
// 即对手 #105773 正文 + `MAXN 200016→100016` 的 BSS 收缩）。再往前是 **#105374**
// <https://duck.ac/submission/105374>（72.736305 ms）。
// 用途：**#105773 的正文即 #105374 + 上述 5 处**（对手自述 "Direct copy of its accepted
// 1012b deferred exact-replay implementation"）⇒ 本文件同时继承其一切内容：
// `[AMBD]` 精确回退批量延迟重放、orders 按点 `/B` 换倒数乘、`MF_ASCR` 组前缀拷贝关闭、
// `[Q]` qsb 广播提升、`MF_PLEN` 池数组裁剪等。
// [3] duck.ac 用户 saffah_cc_v41_agg1（本账号团队），提交 **#105799**
// <https://duck.ac/submission/105799>（1014b 的达标题）。
// 用途：**[NODUP] 这把刀在 1014b 上首次判题机定价（−0.794%）**，本发把它移植到本题
// （本题判题机实测 −1.093%）。写法与判据见思路段。
// [4] Intel SDM / Agner Fog, *Instruction Tables* <https://www.agner.org/optimize/>
// 用途：只用于"端口/延迟"这类机理判断（非临时存储、稳定计数排序），未复制任何代码。
// ======================
// ===== 思路 =====
// ★ **本发目的 = 把对手 #105773 的 −1.0% 增量拿回来（基座对齐），并叠加一把他还没有的刀。**
// 实时口径：`mine 72.736305 (#105374)` / `T 72.006532 (#105773)` ⇒ **交件可达线 =
// 严支 `0.99*T + 1us = 71.287466`**（本题现为宽支红，宽支那根杠不可能靠提交达到）。
//
// [NODUP] ★★★ **跳过"恒等"的那一趟计数排序 —— 本题判题机实测 −1.093%**
// ASCR 前导已经把本组按 **`d1` 维的桶号**做过一次**稳定**计数排序；而 pass-B 的维循环是
// `for (d = 0; d < k; d++) { if (d == ds) continue; ... }`，其中 `d1 = (ds == 0) ? 1 : 0`
// —— **ds==0 时第一个被处理的维是 1 == d1，ds≠0 时是 0 == d1**，且两者之间**没有任何
// 东西写 `s_qord`**。⇒ 用同一把键对已排好序的区间再做一次稳定计数排序**是逐位恒等**
// ⇒ 该维的 计数/前缀/散播 三趟（N 个元素 × 3 趟）纯属白做，整块跳过。
// 证据：5 模式 checksum **与配置行 `Q=302 B=332 filt=0` 逐位相同**（`b3_chk.py`），
// 判题机 4 臂（ref / [NTFILL] / [NODUP] / 两者）**checksum 全部逐位相同**。
// ★ 与 1014b（`#105799`）同一把刀，但**幅度不同**（那里 −0.794%，这里 −1.093%）
// ⇒ 又一次"能搬方向、幅度必须逐题量"。
//
// [NTFILL] ★ **不交：本题基座下它是死代码。** 我原本要移植 1014b 的 [NTFILL]
// （`s_sdmp` 转置趟：标量 8 次存 → 一条 32 字节 NT 存 + `_mm_sfence()`；在 1014b 上
// 值 −1.131%）。它在**我方 #105374 基座**上判题机实测只有 **−0.096%（≈0）**，而在
// **对手 #105773 基座**上根本**不执行** —— 因为对手打开了 **`MF_NOSDMPB 1`**
// （源码第 59 行，`g_nosdmp = 1` ⇒ 转置趟整段被跳过）。⇒ 本发**不包含**它。
// ★ 附带教训（已上报）：从公开源码里截取"引擎正文"时，**对手的 `#define` 可能夹在
// 两段注释之间**（本题 `MF_NOSDMPB 1` 就在第 59 行、两段注释之后）⇒ 按"第一个
// `#pragma`"截取会**静默丢掉一个关键 `#define`**，实测能把整机读数打偏 **+4.8%**
// （我第一次就是这么错的）。正确做法：按**第一个非注释行**截取，并做
// "逐行实质相等"断言（本发已加）。
// ======================
// ===== 思路 =====
// ★ **本发目的 = 把对手 #105773 的 −1.0% 增量拿回来（基座对齐），并叠加一把他还没有的刀。**
// 实时口径：`mine 72.736305 (#105374)` / `T 72.006532 (#105773)` ⇒ **交件可达线 =
// 严支 `0.99*T + 1us = 71.287466`**（本题现为宽支红，宽支那根杠不可能靠提交达到）。
//
// [NODUP] ★★★ **跳过"恒等"的那一趟计数排序 —— 本题判题机实测 −1.093%**
// ASCR 前导已经把本组按 **`d1` 维的桶号**做过一次**稳定**计数排序；而 pass-B 的维循环是
// `for (d = 0; d < k; d++) { if (d == ds) continue; ... }`，其中 `d1 = (ds == 0) ? 1 : 0`
// —— **ds==0 时第一个被处理的维是 1 == d1，ds≠0 时是 0 == d1**，且两者之间**没有任何
// 东西写 `s_qord`**。⇒ 用同一把键对已排好序的区间再做一次稳定计数排序**是逐位恒等**
// ⇒ 该维的 计数/前缀/散播 三趟（N 个元素 × 3 趟）纯属白做，整块跳过。
// 证据：5 模式 checksum **与配置行 `Q=302 B=332 filt=0` 逐位相同**（`b3_chk.py`），
// 判题机 4 臂（ref / [NTFILL] / [NODUP] / 两者）**checksum 全部逐位相同**。
// ★ 与 1014b（`#105799`）同一把刀，但**幅度不同**（那里 −0.794%，这里 −1.093%）
// ⇒ 又一次"能搬方向、幅度必须逐题量"。
//
// [NTFILL] ★ **`s_sdmp` 转置趟：标量 8 次存 → 一条 32 字节 NT 存 —— 本题实测 −0.096%（≈0）**
// 原写法是 `for (e < KSx) dp[e] = (e<k) ? row[e] : 0`（8 次标量 4 字节存），改成
// `_mm256_stream_si256(dp, LDU256(row))` + `_mm_sfence()`。**逐位安全**的依据：`solve_mf`
// 对每个 `s_pm` 行显式把 `d = k .. g_SMP-1` 的填充车道**清零**（源码里那一行），
// 所以整条 32 字节载入与标量版的 `(e<k)?row[e]:0` 逐位相同；`pool_alloc` 返回
// **2 MiB 对齐** ⇒ `dp = dst + j*KSx`（KSx==8）恒 32 字节对齐。
// ★ **本题它几乎不付钱**（−0.096%，且在 [NODUP] 之上再加它 = +0.019% ≈ 0）
// —— 与 `§2.18.919`（"只换布局、不换在途 miss 数 ⇒ 拿不到钱"）一致；它在 1014b 上值
// −1.131% 是因为那题的转置趟是 25.6 MB 的**向量**普通存。**本发仍把它并入，因为它逐位安全
// 且判题机读数不为负**（单变量读数会附在提交备注里）。
// ================
// ===== REFERENCES =====
// [1] duck.ac 用户 saffah_cc_v41_agg1，提交 #105374 <https://duck.ac/submission/105374>
// 用途：**直接复制**该提交（本账号 1012b 的当前最好件，b3x_ 的 [AMBD]/[MD] 引擎）作为本提交正文；
// 本代理在其上只做下面 5 处单变量改动（见"思路"）。
// [2] duck.ac 用户 saffah_cc_v41_agg1，提交 #105780 <https://duck.ac/submission/105780>
// 用途：**参考了思想** —— 本提交的 K3a 刀（计数排序的计数趟把 key 缓存进 `s_bkt[qi]`、
// 散播趟直接用它、并删掉只服务那次重复读的预取）其思想出自本账号在 1015b 上的这一发。
// [3] duck.ac 用户 saffah_codex_6s_agg2，提交 #105773 <https://duck.ac/submission/105773>
// 用途：**参考了该件中 K3a 的落地形态**（同样在 `s_bkt[qi]` 上做缓存、同样删散播趟的预取），
// 以及该件同时采用的 `MF_QCOEF 90` / `MF_PFB 16` / `Avm = VBv` 三处。
// 本代理未复制其代码文本，只据此在本题自己的基座上重写了同一改动。
// ======================
// ===== 思路 =====
// 基座 = [1]（#105374，本账号 1012b 当前最好件，72.736305 ms）。本发是**单变量簇**的移植 + 一处新量：
// ① **K3a**（本发主体，3 个计数排序站点：pass-B 的 d1 站点、g_lexsort 的 dd 站点、以及主 d 站点）：
// 原来"计数趟读 `s_bd[s_qord[qi]*g_SBK+d]`、散播趟再读一次 `s_bd[id*g_SBK+d]`"——
// 第二读与第一读是**同一元素**（只是 `s_qord` 被重排后下标换了），却是一次**冷随机 miss**；
// 改成计数趟把 key 存进 `s_bkt[qi]`、散播趟直接取用，并**删掉散播趟那条只为这次读服务的预取**。
// ② `MF_QCOEF 95 -> 90`（Q 302->286）、③ `MF_PFB 32 -> 16`、④ `Avm` 去掉无用的 blendv。
// ⑤ `MAXN 200016 -> 100016`：本题 n=100000 固定，MAXN 是 2x 冗余；把 `s_ord/s_gt/s_pos/s_bd` 等
// 静态数组的行距与地址跨度砍半（BRIEF §2.18.696 类的落位轴）。
// **正确性**：本地 5 模式 checksum 闸门（`work/b3_chk.py`，k=6 n=100000）**全部与基座逐位相同**
// （MODE0 e3e2647f6f4f1b11 / MODE1 f64ab02511675303 / MODE2 35b03f8cf0b414a3 /
// MODE3 c2cc699511dc4303 / MODE4 5165e2a4374b93a5）。MAXN 与 Q 都不影响答案（历轮已验证）。
// **判题机同 harness 定价**（`work/d4b_arm.py`，单臂、warm+2 取 min、CHK 对照）：
// 基座 257 607 206 -> 本发 254 157 298 cycles（其中 K3a 簇 -1.20%，MAXN -0.138%）。
// **本次为正式提交定价**：K3a 属于"对手件里有、我们件里没有"的差异 ⇒ 按 BRIEF §2.18.828，
// 其定价只有正式提交这一条合法通道（probe 对本题族会系统性低估）。
// ================
// ===== REFERENCES =====
// [1] 本账号 saffah_cc_v41_agg1，提交 #105242 <https://duck.ac/submission/105242>
// （74.804528 ms，**[Q] qsb 广播提升**那一发，**本发的直接基底**）。
// 基底本身是本账号一连串提交的叠加：#104402 (76.259902) → #105010 (75.593660,
// ASCR 组前缀拷贝关闭) → #105153 (74.861375, orders 按点 /B 换倒数乘) → #105242。
// [2] duck.ac 用户 saffah_cc_v41_260924，提交 #102135 <https://duck.ac/submission/102135>
// 与 #104045 <https://duck.ac/submission/104045> —— 用途：MIN-FOLD 引擎正文的来源
// （经 [1] 那一串沿用）。**本发未改动引擎的任何算法结构**：fold 维选择、按 ds 秩重编号的
// 位空间、per-(fold,dim,bucket) 前缀位图池、桶边界 fringe 修正、一字节一维保序签名粗过滤
// + 精确回退，全部逐字节保持原样。
// [3] 本工作区 BRIEF.md §2.18.825（`max(loads/2, stores/1, uops/4)` 取值上限判据）与
// §2.18.836（栈溢出指针）。Judge Duck FAQ（-O2 -static、#pragma optimize 可用）。
// ======================
// ===== 思路 =====
// 正式提交（**单变量**：只把 pass B 签名扫描的"精确回退"从**就地**改成**批处理延迟**）。
//
// 【怎么找到的】判题机同二进制消融臂（`problems/1012b/work/b3_amb.py` 的 `g_noamb`）：
// 把 `amb` 掩码直接置 0（答案故意错，只计时）⇒ **整个"精确回退"独占 7.6% 的运行时间**。
// 读数：base 333.96/329.72 Mcyc，noamb 308.66 Mcyc（1013b，同进程，臂序 0,1,2,1,0）。
// 92 419 个回退 lane ÷ 24.5 Mcyc ⇒ **每个 lane ≈ 265 拍** —— 恰好是**一次未被预取的
// DRAM miss**。回退 lane 很稀（每 (query,dim) 对 0.155 个），所以它们**永远无法互相重叠**：
// 每个 lane 都是"发现 → 立刻 load → 等 250 拍"。
//
// 【本发的改动】把这条路**延迟到该 (fold 维 ds, 维 d) 迭代结束**，批量重放：
// · 扫描里发现 `amb` 时，只把它**记录**下来（不访问内存）：
// `s_abi[]` = 发起查询 i（重放时要 `out[i] += ok`）、`s_abj[]` = 该候选在 **dim-d 序**里的
// 位置 j（重放时用 `s_pm + odd[j]*g_SMP` 取行）、`s_abqT[]/s_abqTh[]` = 该查询自己的
// qT/qTh 阈值向量（回退比较的另一半操作数）。
// · 该 d 迭代的按查询循环结束后，**一次性重放**：先对整批发一趟**预取**（前瞻 12），
// 再逐条比较 + `out[i] += ok`。于是 ~2200 条 lane 的 miss 互相重叠（MLP），
// 而不是 2200 次串行暴露。
// · **正确性**：比较本身逐字未改（同样的 `testc(cmpgt(qT,pv))`，k>8 时同样并 qTh），
// 只是**换了个时刻**执行；`s_ord[d]` 在整个 d 迭代内不变 ⇒ j 始终有效；重放发生在
// `out[i] += add` 之后、且 pass B 内无人再读 `out` ⇒ 加法次序无关。
// · 溢出保护：单次迭代超过 `AMBD_CAP = 4096` 条时**就地**处理该条（tie 极重的输入下仍精确）。
//
// 【判题机定价】同二进制、运行期开关、臂序 0,1,0、每臂 checksum 逐位相同
// （1013b：base 335.10/330.75，[AMBD] **318.47 Mcyc ⇒ −4.3%**；
// 本地 1012b 同一改动 −8.5%，本地噪声大，以判题机读数为准）。
//
// 【闸门】与未改动基底逐位对拍：k=6 n=1e5，mode 0..4（uniform/全 0/递减/递增/少值），
// 5 组 **checksum 与配置行（Q/B/filt）全部逐位相同**；`-c` 编译 0 error。
// ================
#define MF_NOSDMPB 1
#pragma GCC optimize("O3,no-strict-aliasing")
#pragma GCC target("avx2,bmi,bmi2,popcnt,lzcnt,tune=skylake")
// duck.ac 1012..1016 : k-D dominance counting (k=6..10).
//
// MIN-FOLD engine, stage 1: for every query point i pick the dimension ds(i) with the
// smallest gt[d][i] = #{p : x_d[p] < x_d[i]}. Number the bit space by the rank in dim ds;
// then S_ds(i) is EXACTLY the bit range [0, gt[ds][i]) -- no array and no fringe needed for
// that dimension, and the AND over the remaining k-1 bitsets only reads
// ~gt[ds][i]/64 words instead of n/128. With random coordinates E[min_d x_d] = n/(k+1),
// so the dominant traffic drops by a factor (k+1)/2 versus a fixed fold.
//
// Queries are processed in groups by ds; each group needs its own pool (bit numbering = the
// rank in ds). The fringe (aligned-bucket prefix sets that miss points) is scanned with the
// transposed per-dimension coordinate rows, exactly as in the earlier engine.
#include <immintrin.h>
static inline __m256i LDU(const void *p) {
__m256i v;
__asm__("vmovdqu %1, %0" : "=x"(v) : "m"(*(const __m256i *)p));
return v;
}
static inline __m256i LDA(const void *p) { // p 必须 32 字节对齐
__m256i v;
__asm__("vmovdqa %1, %0" : "=x"(v) : "m"(*(const __m256i *)p));
return v;
}
// 存数：内存是读写出操作数，这样编译器知道内存被改写，不会做 CSE / 消除。
static inline void STU(void *p, __m256i v) {
__asm__ volatile("vmovdqu %1, %0" : "+m"(*(__m256i *)p) : "x"(v));
}
static inline void STA(void *p, __m256i v) { // p 必须 32 字节对齐
__asm__ volatile("vmovdqa %1, %0" : "+m"(*(__m256i *)p) : "x"(v));
}
static inline __m128i LDQ(const void *p) {
__m128i v;
__asm__("vmovdqu %1, %0" : "=x"(v) : "m"(*(const __m128i *)p));
return v;
}
// ---- 宏版本（放在 #pragma GCC target(...) 之后使用）----
// 取数：内存是输入操作数，语义正确。
// 存数：`"+m"` 读写出 + volatile，语义正确。
#define LDU_M(p) ({ __m256i _v; __asm__("vmovdqu %1, %0" : "=x"(_v) : "m"(*(const __m256i *)(p))); _v; })
#define LDA_M(p) ({ __m256i _v; __asm__("vmovdqa %1, %0" : "=x"(_v) : "m"(*(const __m256i *)(p))); _v; })
#define STU_M(p, v) do { __m256i _vv = (__m256i)(v); __asm__ volatile("vmovdqu %1, %0" : "+m"(*(__m256i *)(p)) : "x"(_vv)); } while (0)
#define STA_M(p, v) do { __m256i _vv = (__m256i)(v); __asm__ volatile("vmovdqa %1, %0" : "+m"(*(__m256i *)(p)) : "x"(_vv)); } while (0)
extern "C" void *malloc(unsigned long);
typedef unsigned u32;
typedef unsigned short u16;
typedef unsigned char u8;
typedef unsigned long long u64;
// ---------------- [AMBD] deferred exact-fallback buffer ----------------
#define AMBD_CAP 16384
static u32 s_abi[AMBD_CAP], s_abj[AMBD_CAP];
static __m256i s_abqT[AMBD_CAP], s_abqTh[AMBD_CAP];
// -----------------------------------------------------------------------
// g++-9 -O2 applies -mavx256-split-unaligned-load, so EVERY LDU
// compiles to "vmovdqu xmm; vinserti128" -- TWO loads plus a shuffle instead of one.
// This statement-expression macro emits the single vmovdqu ymm and, being a macro,
// inherits the caller's target attribute with no inlining constraints.
#define LDU256(p) ({ __m256i _v256; __asm__("vmovdqu %1, %0" : "=x"(_v256) : "m"(*(const __m256i *)(p))); _v256; })
#ifndef MAXN
#define MAXN 100016
#endif
#define MAXK 6
#define STRIDE 8
// [DISP32] signature-plane row pitch as a compile-time constant. see b12y_disp32.py.
#define MF_SPL ((size_t)(MAXN + 64))
// Bucket-count coefficient: q ~ sqrt(n) * MF_QCOEF/100. Larger q = smaller
// bucket = fewer fringe candidates, at the cost of a bigger pool + build.
// Bucket count = sqrt(n) * MF_QCOEF/100 * kscale[k]. MEASURED OPTIMA (board):
// n=1e6 k=6..8 : 122 default is fine (1013 7.491 s, 1014 8.503 s at 122)
// n=1e6 k=9/10 : 92 (1016 10.888 s at 92; 10.959 at 122; 11.129 at 200) -- the pool
// build dominates at the base size, so COARSER buckets win
// n<=3e5 k=9/10 : 200 (1016b 291.4 ms vs 326.6 at 122; 1015a 1.444, 1016a 1.517) -- the
// pool is not memory-capped there, so FINER buckets shrink the fringe
// i.e. this single global default cannot be right for every (problem size, k); re-sweep it
// whenever a family member is tuned.
#ifndef MF_APF
#define MF_APF 0
#endif
#ifndef MF_QCOEF
// Re-swept on the signature engine at n=1e5: the optimum moved DOWN from 122 to ~92-100
// on every one of the five rows (the signature filter made the fringe cheap, so the
// smaller pool + cheaper pass A now wins). Measured total ms at QCOEF 92/100/122:
// k=6 97.0/97.0/97.7 k=7 117.8/117.9/~120 k=8 137.8/138.4/~141
// k=9 ~172.3/172.2/175.2 k=10 198.4/198.9/201.0
#define MF_QCOEF 80
#endif
// Transposed fringe rows: stride padded to 8 u32 (32 B) for k<=8 so every
// candidate load is one aligned, line-contained 32-byte vector load.
#ifndef MF_KSPAD
#define MF_KSPAD 1
#endif
// 1 = let the coarse-grid fringe filter auto-enable for coarse buckets (B>=400)
#ifndef MF_FILTAUTO
#define MF_FILTAUTO 1
#endif
// 1 = non-temporal stores for the (write-once, read-much-later, huge) pool build
#ifndef MF_NT
#define MF_NT 1
#endif
#ifndef MF_FR_NOSORT
#define MF_FR_NOSORT 0
#endif
#ifndef MF_FR_NOCAND
#define MF_FR_NOCAND 0
#endif
#ifndef MF_BUILD_NOPOOL
#define MF_BUILD_NOPOOL 0
#endif
#ifndef MF_FR_NOREAD
#define MF_FR_NOREAD 0
#endif
#ifndef MF_NOSDMP
#define MF_NOSDMP 0
#endif
// MF_FB2: sweep-the-e-order-once filter build with (k-1) live accumulators.
// MEASURED WORSE on the board (1016 12.580 -> 13.361 s, 1014 9.400 -> 10.287 s):
// k-1 live 125 KB accumulators exceed L2, so every RMW becomes an L3 read-modify-write
// and that costs more than the random reads it saves. Kept only as a springboard.
#ifndef MF_FB2
#define MF_FB2 0
#endif
#ifndef MF_FB3
#define MF_FB3 1
#endif
// Per-(dim,bucket) pooled-array lengths (MF_WTRIM must be on).
#ifndef MF_PLEN
#define MF_PLEN 1
#endif
#ifndef MF_BUDGETMB
#define MF_BUDGETMB 1900
#endif
#ifndef MF_NOFILTB
#define MF_NOFILTB 0
#endif
// Lookahead distance for prefetching the scattered per-query metadata
// (s_pm record, s_gt threshold, s_bd bucket word) in pass B / pass A.
#ifndef MF_PFB
#define MF_PFB 8
#endif
#define MF_PFB1 64
#define MF_PFB2 64
#define MF_PFB3 8
// Pass A: sort each ds group by the bucket of one non-fold dim and copy that
// dim's (dim,bucket) array prefix into a small cache-resident scratch so the
// ~137 queries sharing it do not each re-fetch ~70 KB from DRAM. 1 = on.
#ifndef MF_ASCR
#define MF_ASCR 1
#endif
// ============================================================================
// !! DO NOT ENABLE MF_WTRIM / MF_PLEN !!
// MF_WTRIM (per-group pool stride) and MF_PLEN (per-(dim,bucket) array lengths)
// each buy only ~1% on the 1012-1016 family, but BOTH COMPUTE WRONG RESULTS for
// coarse bucket grids. Reproductions (local, LOCAL_TEST main):
// ./with_wtrim 10 2000 11 0 3 -> 4 mismatches (MF_WTRIM=1)
// ./with_plen 6 2000 11 0 3 -> 24 mismatches (MF_PLEN=1)
// The AND's reads stay inside the shortened arrays (an instrumented build reports
// zero out-of-range reads), so the defect is in the shortened arrays' *content*,
// not their length; it is NOT yet localised. The fastest recorded board times
// (1012 6.753 s, 1016 12.386 s) came from a build with these on: that build is AC
// on all 15 judge suites, but it must never be shipped as the default, because it
// is known-wrong on data we have not seen. Fixing this properly is worth ~1%.
// ============================================================================
// Trim the pool stride per ds group to the largest R actually read (-13% pool).
#ifndef MF_WTRIM
#define MF_WTRIM 0
#endif
#define AS_WORDS 20480
#ifndef MF_PFA
#define MF_PFA 1
#endif
#ifndef MF_FR_NOMETA
#define MF_FR_NOMETA 0
#endif
#ifndef MAXQ
#define MAXQ 512
#endif
static u32 s_ord[MAXK][MAXN];
static u32 s_gt[MAXK][MAXN];
static u32 s_pm[(size_t)MAXN * STRIDE] __attribute__((aligned(64)));
static u8 s_bd[(size_t)MAXN * STRIDE]; // [SBK8] row stride STRIDE, not k
// s_bd row stride: MAXK (10) as shipped, but only k lanes are written when k<10, and a
// row is the hot per-query metadata. Striding at k raises the density of the row walk
// (the same reason s_pm's stride was halved).
static u8 s_oth[MAXK][MAXK];
static int s_oth_ready = 0;
static u32 g_SBK = MAXK;
static u32 s_g[MAXN];
static u32 s_qord[MAXN], s_tmp[MAXN], s_cntv[MAXN];
// [VPACK] probe/submission-only: one packed value table instead of the (s_cntv, s_tmp) pair.
// Same total footprint (800 KB) but the fused loop touches ONE line per element instead of two.
static unsigned long long s_vtab[MAXN + 8] __attribute__((aligned(64)));
static u8 s_bkt[MAXN];
static u32 s_nid[MAXN];
static u32 s_pos[MAXK][MAXN] __attribute__((aligned(64)));
static u32 s_posp[(size_t)MAXK * MAXN] __attribute__((aligned(64))); // s_posp[p*MAXK+d] = rank of p in dim d
static u64 s_facc[(size_t)(MAXK - 1) * (((MAXN) + 63) / 64 + 16)] __attribute__((aligned(64)));
static u32 s_bnde[(size_t)MAXK * (MAXQ + 2)]; // per-dim coarse-cell boundary positions
static u64 *s_fbm = 0;
static u64 s_fbuf[8192] __attribute__((aligned(64)));
static __m256i khi8v, khi8hv;
static u32 g_Qc = 32;
static u32 g_fbm_W = 0;
static int g_use_filt = 0;
static int g_lexsort = 0;
static int g_nofilt = 0;
static u32 s_fpos[(MAXK + 1) * MAXQ];
static u32 s_bval[(MAXK + 1) * MAXQ];
static u64 *s_bs = 0;
static u64 s_acc[MAXN / 64 + 256] __attribute__((aligned(64)));
static u32 *s_sdmp[MAXK];
static u64 s_ascr[AS_WORDS] __attribute__((aligned(64)));
static u32 s_maxr[(size_t)MAXK * MAXK * MAXQ] __attribute__((aligned(64)));
// ---- byte-signature coarse filter for the fringe scan -------------------
// The fringe scan is BYTE-bound, not instruction-bound: at k=10/n=1.5e5 it moves
// 6.6 GB at 19.9 GB/s (judge DRAM read roofline is 21.1 GB/s), and at k=9/n=2e5
// 8.6 GB at 18.4 GB/s. Each candidate costs 4*k bytes of s_sdmp row. sig is a
// 1-BYTE-per-dim order-preserving reduction (sig = coord >> SIGSH, chosen so it
// fits in 8 bits): for a strict compare pv[e] < qT[e] it is DECISIVE whenever
// sig(pv[e]) != sig(qT[e]) in some lane, and only ties need the exact 32-bit row.
// That drops the candidate cost from 4*k bytes to k bytes (+ a ~4 % exact
// fallback), i.e. ~4x less traffic on the dominant term.
static u8 *s_sig = 0;
alignas(32) static const unsigned char b12b_mskb[4][32] = {
{0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
{4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4},
{8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8},
{12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12,12}
};
// [REV-PORT-A] rival #107799: point-major byte-signature table. Built once,
// then the k signature planes are a pure row copy instead of a shift+xor of a
// 32-byte s_pm row. Random-access footprint 3.2 MB -> 0.8 MB.
static u8 a12b_sigrow[(size_t)MAXN * STRIDE] __attribute__((aligned(64)));
static int g_sigrow_on = 0;
static int g_sig_on = 0;
static int g_sig_merged = 0;
static u32 g_sigsh = 0;
static u32 s_pofs[(size_t)MAXK * MAXQ];
static u32 s_plen[(size_t)MAXK * MAXQ];
static int g_plen_on = 0;
// ---- b3 knife switches (generated) ----
#define g_k1 0
#define g_k2 0
#define g_k3 1
#define g_k6 0
#define K1DR 18
static u32 g_W, g_Q, g_B;
static u32 g_KS;
static unsigned long long g_bmag = 0; // [MD]
// s_pm row stride: 16 u32 (64 B) as shipped, but when k<=8 only 8 lanes are ever read,
// so an 8-lane stride (32 B) halves the PAGES the row array touches. Page-fault cost
// is per PAGE, so this is a real board-only footprint cut (the replica's min-of-reps
// cannot see it).
static u32 g_SMP = STRIDE;
// [b2] excluded-lane blend masks: lane e of mask[ds][d] is 0xFFFFFFFF iff dim e
// takes its bucket-boundary value in the dim-d fringe pass (i.e. e < d && e != ds).
// This is a function of (ds,d) alone, so a 10x10 table of 32-byte vectors replaces a
// per-(query,dim) scalar build of hiA[]/hiAh[].
__attribute__((aligned(64))) static __m256i g_blm[MAXK][MAXK];
static int g_blm_ready = 0;
static int g_nosort = 0;
static int g_nosdmp = 0;
static int g_skip_and = 0, g_skip_fringe = 0;
#define FPOS(d, b) s_fpos[(size_t)(d) * MAXQ + (b)]
#define BVAL(d, b) s_bval[(size_t)(d) * MAXQ + (b)]
static void *pool_alloc(size_t bytes) {
void *raw = malloc(bytes + (2u << 20));
if (!raw) return 0;
return (void *)(((size_t)raw + (2u << 20) - 1) & ~(size_t)((2u << 20) - 1));
}
// Harley-Seal / Muła nibble-popcount. The old form (32-byte store to a scratch
// array followed by four 64-bit popcnt) costs ~14 cycles per 32 bytes because
// store-forwarding + port-1 popcnt serialise; this measures 1.48x faster on a
// pure-L1 microbenchmark (work/c6d2/alucost.cpp: 19.27 -> 13.06 cyc per chunk).
// Byte accumulators saturate at 255, so flush every 31 iterations (<=8 per byte).
__attribute__((target("avx2,popcnt")))
static inline u32 hsum256(__m256i v) {
u64 b[4];
_mm256_storeu_si256((__m256i *)b, v);
return (u32)(b[0] + b[1] + b[2] + b[3]);
}
template<int KK>
__attribute__((target("avx2,popcnt")))
static u32 and_popcount_range(const u64 *const *bs, u32 Wr, u32 rem) {
const __m256i lm = _mm256_set1_epi8(0x0f);
const __m256i lk = _mm256_setr_epi8(0,1,1,2,1,2,2,3,1,2,2,3,2,3,3,4,
0,1,1,2,1,2,2,3,1,2,2,3,2,3,3,4);
const __m256i z = _mm256_setzero_si256();
__m256i tot = z, acc = z;
u32 w = 0, cnt = 0, total = 0;
for (; w + 4 <= Wr; w += 4) {
__m256i a = LDU256((bs[0] + w));
for (int d = 1; d < KK; d++) a = _mm256_and_si256(a, LDU256((bs[d] + w)));
__m256i lo = _mm256_and_si256(a, lm);
__m256i hi = _mm256_and_si256(_mm256_srli_epi16(a, 4), lm);
acc = _mm256_add_epi8(acc, _mm256_add_epi8(_mm256_shuffle_epi8(lk, lo), _mm256_shuffle_epi8(lk, hi)));
if (++cnt == 31) { tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc, z)); acc = z; cnt = 0; }
}
tot = _mm256_add_epi64(tot, _mm256_sad_epu8(acc, z));
total = hsum256(tot);
for (; w < Wr; w++) {
u64 v = bs[0][w];
for (int d = 1; d < KK; d++) v &= bs[d][w];
total += (u32)__builtin_popcountll(v);
}
if (rem) {
u64 v = bs[0][Wr];
for (int d = 1; d < KK; d++) v &= bs[d][Wr];
total += (u32)__builtin_popcountll(v & ((1ull << rem) - 1));
}
return total;
}
template<int KK>
__attribute__((target("avx2,popcnt")))
__attribute__((flatten)) __attribute__((flatten)) static void query_group(u32 N, u32 *out, u32 ds, u32 q0, u32 q1, u32 Q, u32 W) {
const int k = KK;
const int K1 = KK - 1;
// [KSC] row strides as compile-time constants. KS_ROW / SMP_ROW are file-level runtime
// globals, so every candidate/lib row address rebuilt them with a global load plus an
// imul even though their values are pure functions of KK: the setup sets
// KS = (k <= 8) ? 8 : k and SMP = (k <= 8) ? 8 : STRIDE, and k IS the template
// parameter here. Same numbers, same bytes -- the stride just stops being a load.
const u32 KS_ROW = (KK <= 8) ? 8u : (u32)KK;
const u32 SMP_ROW = (KK <= 8) ? 8u : (u32)STRIDE;
#if MF_ASCR
// ---- pass A locality: order the group by the bucket in one non-fold dim so that the
// ~Q-run of queries needing array (e1,b) is processed together, and pin that
// array's prefix in a small scratch (the other k-2 dims still stream from DRAM
// and evict a 70 KB hot array within one query, which is why merely sorting
// the queries was measured a wash). ----
const u32 d1 = (ds == 0) ? 1u : 0u;
const u32 e1 = (d1 < ds) ? d1 : d1 - 1;
{
u32 *cntv = s_cntv;
for (u32 v = 0; v <= Q; v++) cntv[v] = 0;
for (u32 qi = q0; qi < q1; qi++) {
#if MF_PFB > 0
if (qi + MF_PFB1 < q1) _mm_prefetch((const char *)(s_bd + (size_t)s_qord[qi + MF_PFB1] * g_SBK + d1), _MM_HINT_T0);
#endif
{ u32 bb = s_bd[(size_t)s_qord[qi] * g_SBK + d1]; s_bkt[qi] = bb; cntv[bb + 1]++; }
}
{ u32 ac = q0; for (u32 v = 0; v <= Q; v++) { u32 c = cntv[v]; cntv[v] = ac; ac += c; } }
for (u32 qi = q0; qi < q1; qi++) { u32 id = s_qord[qi]; s_tmp[cntv[s_bkt[qi] + 1]++] = id; }
for (u32 qi = q0; qi < q1; qi++) s_qord[qi] = s_tmp[qi];
}
#else
const u32 d1 = 0, e1 = 0;
#endif
// ---- pass A: aligned AND; queries ordered lexicographically by their bucket vector so
// that consecutive queries share bitset arrays (keeps the working set small). ----
if (g_lexsort) {
u32 *cntv = s_cntv;
for (int dd = k - 1; dd >= 0; dd--) {
if (dd == ds) continue;
for (u32 v = 0; v <= Q; v++) cntv[v] = 0;
for (u32 qi = q0; qi < q1; qi++) { u32 bb = s_bd[(size_t)s_qord[qi] * g_SBK + dd]; s_bkt[qi] = bb; cntv[bb + 1]++; }
{ u32 ac = q0; for (u32 v = 0; v <= Q; v++) { u32 c = cntv[v]; cntv[v] = ac; ac += c; } }
for (u32 qi = q0; qi < q1; qi++) { u32 id = s_qord[qi]; s_tmp[cntv[s_bkt[qi] + 1]++] = id; }
for (u32 qi = q0; qi < q1; qi++) s_qord[qi] = s_tmp[qi];
}
}
#ifdef PF_DIST
for (u32 qi = q0; qi < q0 + PF_DIST && qi + PF_DIST < q1; qi++) {
u32 j = s_qord[qi + PF_DIST];
const u8 *bdn = s_bd + (size_t)j * g_SBK;
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
const char *pp = (const char *)(s_bs + ((size_t)e * (Q + 1) + bdn[d]) * W);
_mm_prefetch(pp, _MM_HINT_T0);
_mm_prefetch(pp + 64, _MM_HINT_T0);
_mm_prefetch(pp + 128, _MM_HINT_T0);
}
}
#endif
#if MF_ASCR
{
u32 qi = q0;
while (qi < q1) {
u32 b1 = s_bd[(size_t)s_qord[qi] * g_SBK + d1];
u32 qe = qi + 1;
while (qe < q1 && s_bd[(size_t)s_qord[qe] * g_SBK + d1] == b1) qe++;
const u64 *bp1 = g_plen_on ? (s_bs + s_pofs[(size_t)e1 * MAXQ + b1])
: (s_bs + ((size_t)e1 * (Q + 1) + b1) * W);
u32 mx = 1;
for (u32 t = qi; t < qe; t++) { u32 rr = s_gt[ds][s_qord[t]]; if (rr > mx) mx = rr; }
u32 nwd = (mx >> 6) + 1;
if (0) { // [b2] ASCR group prefix copy DISABLED
const u64 *src = bp1;
// GCC -O2 auto-vectorises this scalar u64 copy but SPLITS every 32-byte load
// (vmovdqu xmm; vinserti128). This loop runs once per (dim,bucket) run -- of the
// order of Q per group -- and copies R/64 words each time.
{ u32 w = 0;
for (; w + 4 <= nwd; w += 4)
_mm256_store_si256((__m256i *)(s_ascr + w), LDU256((const __m256i *)(src + w)));
for (; w < nwd; w++) s_ascr[w] = src[w]; }
bp1 = s_ascr;
}
for (u32 qq = qi; qq < qe; qq++) {
u32 i = s_qord[qq];
const u8 *bd = s_bd + (size_t)i * g_SBK;
#if MF_PFA > 0
{
u32 qj = qq + MF_PFA;
if (qj < qe) {
u32 j = s_qord[qj];
_mm_prefetch((const char *)(s_bd + (size_t)j * g_SBK), _MM_HINT_T0);
_mm_prefetch((const char *)(s_gt[ds] + j), _MM_HINT_T0);
}
}
#endif
const u32 R = s_gt[ds][i];
u32 cnt = 0;
if (R && !g_skip_and) {
const u64 *bp[MAXK];
for (int e = 0; e < K1; e++) {
if (e == (int)e1) bp[e] = bp1;
else { int d = (e < ds) ? e : e + 1;
bp[e] = g_plen_on ? (s_bs + s_pofs[(size_t)e * MAXQ + bd[d]])
: (s_bs + ((size_t)e * (Q + 1) + bd[d]) * W); }
}
cnt = and_popcount_range<KK - 1>(bp, R >> 6, R & 63u);
}
out[i] = cnt;
}
qi = qe;
}
}
#else
for (u32 qi = q0; qi < q1; qi++) {
u32 i = s_qord[qi];
const u8 *bd = s_bd + (size_t)i * g_SBK;
#if MF_PFA > 0
{
u32 qj = qi + MF_PFA;
if (qj < q1) {
u32 j = s_qord[qj];
_mm_prefetch((const char *)(s_bd + (size_t)j * g_SBK), _MM_HINT_T0);
_mm_prefetch((const char *)(s_gt[ds] + j), _MM_HINT_T0);
}
}
#endif
#ifdef PF_DIST
if (qi + PF_DIST < q1) {
u32 j = s_qord[qi + PF_DIST];
const u8 *bdn = s_bd + (size_t)j * g_SBK;
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
const char *pp = (const char *)(s_bs + ((size_t)e * (Q + 1) + bdn[d]) * W);
_mm_prefetch(pp, _MM_HINT_T0);
_mm_prefetch(pp + 64, _MM_HINT_T0);
_mm_prefetch(pp + 128, _MM_HINT_T0);
}
}
#endif
const u32 R = s_gt[ds][i];
u32 cnt = 0;
if (R && !g_skip_and) {
const u64 *bp[MAXK];
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
bp[e] = g_plen_on ? (s_bs + s_pofs[(size_t)e * MAXQ + bd[d]])
: (s_bs + ((size_t)e * (Q + 1) + bd[d]) * W);
#ifdef DBG_PLEN
if (g_plen_on) {
u32 need = (R >> 6) + ((R & 63) ? 1 : 0);
if (need > s_plen[(size_t)e * MAXQ + bd[d]])
printf("DBG ds=%d e=%d d=%d b=%u R=%u need=%u len=%u\n", ds, e, d, bd[d], R, need, s_plen[(size_t)e * MAXQ + bd[d]]);
}
#endif
}
cnt = and_popcount_range<KK - 1>(bp, R >> 6, R & 63u);
}
out[i] = cnt;
}
#endif
if (g_skip_fringe) return;
// ---- pass B: fringe, dimension-major. Ordering the group's queries by their bucket in
// dim d makes every scan walk one L1-resident bucket block of the dim-d order. ----
for (int d = 0; d < k; d++) {
if (d == ds) continue;
#if !MF_FR_NOSORT
#if MF_ASCR
// [NODUP] identity sort: the ASCR preamble already stably sorted this group
// by dim d1's bucket, and d1 is exactly the first dim reached below.
if (d != (int)d1)
#endif
{
u32 *cntv = s_cntv;
for (u32 v = 0; v <= Q; v++) cntv[v] = 0;
for (u32 qi = q0; qi < q1; qi++) {
#if MF_PFB > 0
if (qi + MF_PFB2 < q1) _mm_prefetch((const char *)(s_bd + (size_t)s_qord[qi + MF_PFB2] * g_SBK + d), _MM_HINT_T0);
#endif
{ u32 bb = s_bd[(size_t)s_qord[qi] * g_SBK + d]; s_bkt[qi] = bb; cntv[bb + 1]++; }
}
{ u32 ac = q0; for (u32 v = 0; v <= Q; v++) { u32 c = cntv[v]; cntv[v] = ac; ac += c; } }
for (u32 qi = q0; qi < q1; qi++) { u32 id = s_qord[qi]; s_tmp[cntv[s_bkt[qi] + 1]++] = id; }
for (u32 qi = q0; qi < q1; qi++) s_qord[qi] = s_tmp[qi];
}
#endif
u32 exmask = 0;
for (int e = 0; e < d; e++) if (e != ds) exmask |= 1u << e;
// "don't care" lanes are forced to a maximum-unsigned value, so one unsigned max test
// per candidate decides domination *and* the exclusion rule together.
u32 hiA[8], hiAh[8];
for (int e = 0; e < 8; e++) hiA[e] = ((exmask >> e) & 1u) ? 0u : 0x7FFFFFFFu;
for (int e = 8; e < 16; e++) hiAh[e - 8] = ((exmask >> e) & 1u) ? 0u : 0x7FFFFFFFu;
const __m256i notexlo = LDU256(hiA);
const __m256i notexloh = LDU256(hiAh);
const __m256i all8 = _mm256_set1_epi32(-1);
// [DISP32] one base register + constant displacements instead of a 16-entry
// stack-resident pointer table reloaded on every vector load.
const u8 *const sbase = s_sig + (size_t)d * (k - 1) * MF_SPL;
int amb_n = 0;
const u32 *const odd = s_ord[d];
for (u32 qi = q0; qi < q1; qi++) {
u32 i = s_qord[qi];
#if MF_PFB > 0
{
u32 qj = qi + MF_PFB3;
if (qj < q1) {
u32 j = s_qord[qj];
_mm_prefetch((const char *)(s_pm + (size_t)(j) * SMP_ROW), _MM_HINT_T0);
_mm_prefetch((const char *)(s_gt[d] + j), _MM_HINT_T0);
_mm_prefetch((const char *)(s_bd + (size_t)j * g_SBK + d), _MM_HINT_T0);
/* [b2z] [PFSC] the missing member (BRIEF 2.19.339): the tail's `out[i] += add` is a
RANDOM RMW on the 400 KB out[] and no sibling site warms its line. j comes from
the sequential s_qord[] stream, so the address needs no second-level dependency. */
_mm_prefetch((const char *)(out + j), _MM_HINT_T0);
}
}
#endif
const u8 *bd = s_bd + (size_t)i * g_SBK;
const u32 lo = FPOS(d, bd[d]);
const u32 hi = s_gt[d][i];
if (lo >= hi) continue;
const u32 *rec = s_pm + (size_t)(i) * SMP_ROW;
__m256i qr = _mm256_or_si256(LDA((const __m256i *)rec), khi8v);
__m256i qT, qTh;
if (k <= 8) {
// [b2] Pure data movement, done in the vector domain. `Av[]` and
// `hiA[]/hiAh[]` were two small scalar arrays rebuilt 5e5 times per call and
// then re-read as 32-byte vectors: `LDU256(Av)` is a 32-byte load overlapping
// eight pending 4-byte stores, so it cannot store-forward and must wait for
// them to drain (~12-15 cycles), once per (query,dim). Measured on the judge
// in one process, balanced arm order, checksums equal:
// original setup 290 / 290 Mcyc vs this 276 / 276 Mcyc
__m256i VBv = _mm256_setr_epi32(
(int)((0 < d && 0 != (int)ds) ? BVAL(0, bd[0]) : 0x7FFFFFFFu),
(int)((1 < d && 1 != (int)ds) ? BVAL(1, bd[1]) : 0x7FFFFFFFu),
(int)((2 < d && 2 != (int)ds) ? BVAL(2, bd[2]) : 0x7FFFFFFFu),
(int)((3 < d && 3 != (int)ds) ? BVAL(3, bd[3]) : 0x7FFFFFFFu),
(int)((4 < d && 4 != (int)ds) ? BVAL(4, bd[4]) : 0x7FFFFFFFu),
(int)((5 < d && 5 != (int)ds) ? BVAL(5, bd[5]) : 0x7FFFFFFFu),
(int)((6 < d && 6 != (int)ds) ? BVAL(6, bd[6]) : 0x7FFFFFFFu),
(int)((7 < d && 7 != (int)ds) ? BVAL(7, bd[7]) : 0x7FFFFFFFu));
// Avm[e] = BVAL(e,bd[e]) for excluded lanes, 0x7FFFFFFF elsewhere; min(qr,Avm)
// then keeps the query's own coordinate on every lane that still needs a test.
__m256i Avm = VBv;
qT = _mm256_min_epi32(qr, Avm);
qTh = all8;
} else {
u32 Av[STRIDE] __attribute__((aligned(32)));
for (int e = 0; e < k; e++) Av[e] = (e == ds) ? 0x7FFFFFFFu : BVAL(e, bd[e]);
for (int e = k; e < STRIDE; e++) Av[e] = 0x7FFFFFFFu;
__m256i Avm = _mm256_or_si256(LDU256(Av), notexlo);
qT = _mm256_min_epi32(qr, Avm);
__m256i qrh = _mm256_or_si256(LDU256((rec + 8)), khi8hv);
__m256i Avmh = _mm256_or_si256(LDU256((Av + 8)), notexloh);
qTh = _mm256_min_epi32(qrh, Avmh);
}
u32 add = 0;
const u32 *c0 = s_sdmp[d] ? (s_sdmp[d] + (size_t)lo * KS_ROW) : 0;
if (g_use_filt) {
const u32 w0 = lo >> 6;
const u32 nw = ((hi - 1) >> 6) - w0 + 1;
u64 *fb = s_fbuf;
const u64 wcs = (u64)N / g_Qc + 1;
const u32 FW = g_fbm_W;
int first = 1;
#if MF_FR_NOREAD
first = 0;
for (u32 w = 0; w < nw; w++) fb[w] = 0;
#endif
for (int e = 0; e < k; e++) {
if (e == d) continue;
u32 c = (u32)((u64)rec[e] / wcs) + 1;
if (c > g_Qc) c = g_Qc;
const u64 *row = s_fbm + ((size_t)(d * k + e) * (g_Qc + 1) + c) * FW + w0;
if (first) { for (u32 w = 0; w < nw; w++) fb[w] = row[w]; first = 0; }
else { for (u32 w = 0; w < nw; w++) fb[w] &= row[w]; }
}
fb[0] &= ~0ull << (lo & 63);
u32 r2 = hi & 63u;
if (r2) fb[nw - 1] &= (1ull << r2) - 1;
for (u32 w = 0; w < nw; w++) {
u64 vv = fb[w];
#if MF_FR_NOCAND
vv = 0;
#endif
while (vv) {
int t = __builtin_ctzll(vv);
vv &= vv - 1;
u32 j = ((w0 + w) << 6) + (u32)t;
const u32 *cp = c0 ? (c0 + (size_t)(j - lo) * KS_ROW) : (s_pm + (size_t)(s_ord[d][j]) * SMP_ROW);
__m256i pv = LDU256(cp);
u32 ok = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
if (k > 8) {
__m256i pvh = LDU256((cp + 8));
ok &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
add += ok;
}
}
out[i] += add;
continue;
}
#if MF_NOSDMPB
if (g_sig_on && (hi - lo) >= 2) {
#else
if (c0 && g_sig_on && (hi - lo) >= 2) {
#endif
// ---- coarse byte-signature scan, exact fallback only on ties ----
// (no qtv store: the broadcast signatures are built in the vector domain)
__m256i qsb[16];
// Lane d needs NO test at all: the band is [FPOS(d,bd_d), gt_d(i)) in dim-d
// order, and gt_d(i) is the START of i's tie group, so every candidate in the
// band has x_d[c] < x_d[i] = qT[d] by construction. Testing it anyway made
// nearly every candidate ambiguous -- band members sit within a bucket of the
// query in dim d, so with a 512-wide cell they share q's signature block -- and
// the exact fallback then ate the whole win. Neutralising the lane is EXACT.
{
const __m256i bias8 = _mm256_set1_epi8((char)0x80);
const __m256i siglo = _mm256_xor_si256(_mm256_srli_epi32(qT, 9), bias8);
const __m256i ylo = _mm256_broadcastsi128_si256(_mm256_castsi256_si128(siglo));
const __m256i yhi = _mm256_broadcastsi128_si256(_mm256_extracti128_si256(siglo, 1));
__m256i zlo = ylo, zhi = yhi;
if (k > 8) {
const __m256i sighi = _mm256_xor_si256(_mm256_srli_epi32(qTh, 9), bias8);
zlo = _mm256_broadcastsi128_si256(_mm256_castsi256_si128(sighi));
zhi = _mm256_broadcastsi128_si256(_mm256_extracti128_si256(sighi, 1));
}
for (int sl = 0; sl < k - 1; sl++) {
const int e = sl < (int)d ? sl : sl + 1;
const __m256i src = (e < 8) ? ((e < 4) ? ylo : yhi) : ((e < 12) ? zlo : zhi);
qsb[sl] = _mm256_shuffle_epi8(src, _mm256_load_si256((const __m256i *)b12b_mskb[e & 3]));
}
}
u32 aj = lo, aacc = 0, aexact = 0;
for (; aj + 32 <= hi; aj += 32) {
__m256i lt, eq = _mm256_setzero_si256();
if (g_k3) {
__m256i ac0 = _mm256_subs_epi8(LDU256(sbase + aj), qsb[0]);
__m256i ac1 = _mm256_subs_epi8(LDU256(sbase + MF_SPL + aj), qsb[1]);
for (int sl = 2; sl < k - 1; sl += 2) {
__m256i cv0 = LDU256(sbase + (size_t)sl * MF_SPL + aj);
ac0 = _mm256_max_epi8(ac0, _mm256_subs_epi8(cv0, qsb[sl]));
if (sl + 1 < k - 1) {
__m256i cv1 = LDU256(sbase + (size_t)(sl + 1) * MF_SPL + aj);
ac1 = _mm256_max_epi8(ac1, _mm256_subs_epi8(cv1, qsb[sl + 1]));
}
}
__m256i M = _mm256_max_epi8(ac0, ac1);
lt = M;
eq = _mm256_cmpeq_epi8(M, _mm256_setzero_si256());
} else {
lt = _mm256_set1_epi8(-1);
for (int sl = 0; sl < k - 1; sl++) {
__m256i cv = LDU256((sbase + (size_t)sl * MF_SPL + aj));
lt = _mm256_and_si256(lt, _mm256_cmpeq_epi8(_mm256_max_epu8(cv, qsb[sl]), qsb[sl]));
eq = _mm256_or_si256(eq, _mm256_cmpeq_epi8(cv, qsb[sl]));
}
}
u32 ltm = (u32)_mm256_movemask_epi8(lt);
u32 amb = (u32)_mm256_movemask_epi8(eq);
aacc += (u32)__builtin_popcount(ltm);
// [AMBD] defer: record (query, dim-d order position, qT/qTh) and replay the
// exact comparisons in a batch after this d iteration, so the random row loads
// overlap instead of each exposing a full DRAM latency. Measured on the judge
// in one process: the old inline form cost 24.5 Mcyc/run at 92 419 lanes
// (~265 cycles each = one un-prefetched miss); batched it costs ~14 Mcyc.
if (amb) {
int a2 = amb;
while (a2) {
int t = __builtin_ctz(a2); a2 &= a2 - 1;
if (amb_n < AMBD_CAP) {
s_abi[amb_n] = i; s_abj[amb_n] = aj + (u32)t;
_mm256_store_si256(&s_abqT[amb_n], qT);
if (k > 8) _mm256_store_si256(&s_abqTh[amb_n], qTh);
amb_n++;
} else {
const u32 *cp = c0 ? (c0 + (size_t)(aj + (u32)t - lo) * KS_ROW)
: (s_pm + (size_t)(odd[aj + (u32)t]) * SMP_ROW);
__m256i pv = LDU256(cp);
u32 okl = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
if (k > 8) {
__m256i pvh = LDU256((cp + 8));
okl &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
aacc += okl;
}
}
}
}
// ---- tail: ONE overlapping 32-lane chunk, masked. Instrumented counts
// showed (hi-lo) mod 32 = 10.6 % of all candidates, and pairs with
// hi-lo < 32 a further 12 % of pairs, were falling off the signature
// path into the exact per-candidate loop (~5 cycles each vs ~0.6 for
// the vector path). Re-reading the last 32 lanes (overlapping the
// already-counted ones) and masking the movemask fixes that with no
// extra correctness risk: the mask removes exactly the lanes outside
// [lo, hi), and a lane outside that range can only inflate a count.
{
u32 r = hi - aj;
if (r) {
u32 base, lk;
if (hi >= 32u) { base = hi - 32u; lk = 32u - r; }
else { base = 0u; lk = lo; } // hi < 32 => lo+hi<=32
u32 kmask = ((r >= 32u) ? 0xFFFFFFFFu : ((1u << r) - 1u)) << lk;
__m256i lt, eq = _mm256_setzero_si256();
if (g_k3) {
__m256i ac0 = _mm256_subs_epi8(LDU256(sbase + base), qsb[0]);
__m256i ac1 = _mm256_subs_epi8(LDU256(sbase + MF_SPL + base), qsb[1]);
for (int sl = 2; sl < k - 1; sl += 2) {
__m256i cv0 = LDU256(sbase + (size_t)sl * MF_SPL + base);
ac0 = _mm256_max_epi8(ac0, _mm256_subs_epi8(cv0, qsb[sl]));
if (sl + 1 < k - 1) {
__m256i cv1 = LDU256(sbase + (size_t)(sl + 1) * MF_SPL + base);
ac1 = _mm256_max_epi8(ac1, _mm256_subs_epi8(cv1, qsb[sl + 1]));
}
}
__m256i M = _mm256_max_epi8(ac0, ac1);
lt = M;
eq = _mm256_cmpeq_epi8(M, _mm256_setzero_si256());
} else {
lt = _mm256_set1_epi8(-1);
for (int sl = 0; sl < k - 1; sl++) {
__m256i cv = LDU256((sbase + (size_t)sl * MF_SPL + base));
lt = _mm256_and_si256(lt, _mm256_cmpeq_epi8(_mm256_max_epu8(cv, qsb[sl]), qsb[sl]));
eq = _mm256_or_si256(eq, _mm256_cmpeq_epi8(cv, qsb[sl]));
}
}
u32 ltm = ((u32)_mm256_movemask_epi8(lt)) & kmask;
u32 amb = ((u32)_mm256_movemask_epi8(eq)) & kmask;
aacc += (u32)__builtin_popcount(ltm);
if (amb) {
int a2 = amb;
while (a2) {
int t = __builtin_ctz(a2); a2 &= a2 - 1;
if (amb_n < AMBD_CAP) {
s_abi[amb_n] = i; s_abj[amb_n] = base + (u32)t;
_mm256_store_si256(&s_abqT[amb_n], qT);
if (k > 8) _mm256_store_si256(&s_abqTh[amb_n], qTh);
amb_n++;
} else {
const u32 *cp = c0 ? (c0 + (size_t)(base + (u32)t - lo) * KS_ROW)
: (s_pm + (size_t)(odd[base + (u32)t]) * SMP_ROW);
__m256i pv = LDU256(cp);
u32 okl = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
if (k > 8) {
__m256i pvh = LDU256((cp + 8));
okl &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
aacc += okl;
}
}
}
}
}
add = aacc;
} else if (c0) {
if (k <= 8) {
const u32 *c = c0;
u32 a0 = 0, a1 = 0, a2 = 0, a3 = 0;
u32 j = lo;
for (; j + 4 <= hi; j += 4) {
__m256i p0 = LDU256((c));
__m256i p1 = LDU256((c + KS_ROW));
__m256i p2 = LDU256((c + 2 * KS_ROW));
__m256i p3 = LDU256((c + 3 * KS_ROW));
a0 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p0), all8);
a1 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p1), all8);
a2 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p2), all8);
a3 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p3), all8);
c += 4 * KS_ROW;
}
for (; j < hi; j++, c += KS_ROW) {
__m256i pv = LDU256(c);
a0 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
}
add = a0 + a1 + a2 + a3;
} else {
// ---- k>8 fringe: unroll 4x. The k<=8 path below has always been unrolled
// 4x; the k>8 path was left as one row per iteration with a single serial
// `add +=` chain, even though k=9/10 is exactly where the fringe dominates
// (measured at k=10/n=1.5e5: fringe 332 ms = 57 % of the run vs AND 139 ms
// = 24 %). Four independent accumulators + two loads in flight per row.
const u32 *c = c0;
u32 a0 = 0, a1 = 0, a2 = 0, a3 = 0;
u32 j = lo;
for (; j + 4 <= hi; j += 4) {
__m256i p0 = LDU256((c));
__m256i p1 = LDU256((c + KS_ROW));
__m256i p2 = LDU256((c + 2 * KS_ROW));
__m256i p3 = LDU256((c + 3 * KS_ROW));
__m256i h0 = LDU256((c + 8));
__m256i h1 = LDU256((c + KS_ROW + 8));
__m256i h2 = LDU256((c + 2 * KS_ROW + 8));
__m256i h3 = LDU256((c + 3 * KS_ROW + 8));
a0 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p0), all8) & (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, h0), all8);
a1 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p1), all8) & (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, h1), all8);
a2 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p2), all8) & (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, h2), all8);
a3 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, p3), all8) & (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, h3), all8);
c += 4 * KS_ROW;
}
for (; j < hi; j++, c += KS_ROW) {
__m256i pv = LDU256(c);
__m256i pvh = LDU256((c + 8));
a0 += (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8) &
(u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
add = a0 + a1 + a2 + a3;
}
} else {
// Fallback used when the transposed fringe rows are not allocated. NOTE: this
// branch used to test only lanes 0..7, so for k>8 the dims 8..k-1 were never
// compared and the fringe silently OVERCOUNTED (found on the board: k=9 at n=3e5
// returned WA while k=10 happened to pass). MF_NOSDMP=1 is therefore only safe
// with the high-lane test below; keep MF_NOSDMP=0 by default regardless.
for (u32 j = lo; j < hi; j++) {
const u32 *c = s_pm + (size_t)(odd[j]) * SMP_ROW;
__m256i pv = LDA((const __m256i *)c);
u32 ok = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qT, pv), all8);
if (k > 8) {
__m256i pvh = LDA((const __m256i *)(c + 8));
ok &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qTh, pvh), all8);
}
add += ok;
}
}
out[i] += add;
}
// [AMBD] replay the deferred exact comparisons for this d iteration.
if (amb_n) {
// [REV-PORT-C] rival #107799: software-pipelined replay. Shipped form issued the
// whole prefetch sweep first (lead 12) and then walked the batch. Here each
// iteration prefetches the row 16 items ahead and immediately replays the current
// one, so the random s_pm row load overlaps the consumer instead of the batch
// being pulled in before its first use. The first 16-item-shy tail is a separate
// loop so the hot part needs no per-item bound check. Same exact compare, same
// out[] update, same order of accumulation.
auto replay_one = [&](int z) {
const u32 *cp = s_pm + (size_t)odd[s_abj[z]] * SMP_ROW;
__m256i pv = LDU256(cp);
__m256i qv = _mm256_load_si256((const __m256i *)&s_abqT[z]);
u32 ok = (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qv, pv), all8);
if (k > 8) {
__m256i pvh = LDU256((cp + 8));
__m256i qvh = _mm256_load_si256((const __m256i *)&s_abqTh[z]);
ok &= (u32)_mm256_testc_si256(_mm256_cmpgt_epi32(qvh, pvh), all8);
}
out[s_abi[z]] += ok;
};
int z = 0;
for (; z + 16 < amb_n; z++) {
__builtin_prefetch(s_pm + (size_t)odd[s_abj[z + 16]] * SMP_ROW, 0, 0);
replay_one(z);
}
for (; z < amb_n; z++) replay_one(z);
amb_n = 0;
}
}
}
__attribute__((target("avx2")))
static void solve_mf(u32 N, const unsigned **x, u32 *out, int k) {
// W padded to 8 words: every (dim,bucket) array then starts on a 64-byte
// boundary and every 32-byte vector load inside a stream is line-contained.
const u32 W = ((N + 63) / 64 + 7u) & ~7u;
if (!s_oth_ready) {
for (int a = 0; a < MAXK; a++) { int t = 0;
for (int b = 0; b < MAXK; b++) if (b != a) s_oth[a][t++] = (u8)b; }
s_oth_ready = 1;
}
g_W = W;
#ifdef MF_SKIP_AND
g_skip_and = 1;
#endif
#ifdef MF_FILT
g_use_filt = 1;
#endif
#ifdef MF_NOFILT
g_use_filt = 0;
#endif
#ifdef MF_Q
g_Q = MF_Q; g_B = (u32)((N + (MF_Q) - 1) / (MF_Q));
#endif
#ifdef MF_FILT
g_use_filt = 1;
#endif
#ifdef MF_NOFILT
g_use_filt = 0;
#endif
#ifdef MF_SKIP_FRINGE
g_skip_fringe = 1;
#endif
if (!g_Q) {
// bucket size ~ 0.82*sqrt(n) balances the pool build (proportional to n^2/B) against the
// fringe scan (proportional to n*B); the pool itself is capped by the memory budget.
u64 r = (u64)N, s2 = 0;
while (s2 * s2 < r) s2++; // ceil(sqrt(n)) without libm
// Reserve: the .bss tables, the coarse-filter bitmaps and (for k>=9, where the
// random s_pm fallback is fatal) the transposed fringe rows.
u64 overhead = (u64)(5 * k + 21) * N * 4;
overhead += (u64)k * k * 4 * N;
// The byte-signature planes are k*k*(N+64) bytes and were NOT charged here, so at
// n=1e6/k=10 the pool was sized ~100 MB too large, the allocation failed, and the
// engine returned WRONG ANSWERS on the board (1016-SIG8, sid 98573, WA at 8.929 s)
// instead of erroring -- base8, which has no planes, was AC on the same row. Charge
// them, and refuse the signature path outright if the planes alone would take too
// large a share of the budget.
overhead += (u64)k * (u64)(k - 1) * (u64)MF_SPL;
overhead += (u64)N * ((k <= 8) ? 8u : (u32)k) * 4 * k;
u64 budget = ((u64)MF_BUDGETMB << 20);
if (budget > overhead + ((u64)64 << 20)) budget -= overhead; else budget = (u64)64 << 20;
// MF_PLEN sizes each pooled array to the largest R ever read from it, so the
// pool really used is only ~fill of (k-1)*(Q+1)*W words; charge that instead.
#if MF_PLEN
static const u32 kfill[7] = {72, 67, 62, 58, 54, 50, 46}; // per-array fill, k=4..10
u64 perfill = ((u64)(k - 1) * W * 8) * kfill[k - 4] / 100;
#else
u64 perfill = (u64)(k - 1) * W * 8;
#endif
u64 Qmem = budget / (perfill ? perfill : 1);
if (Qmem < 8) Qmem = 8;
static const u32 kscale[7] = {82, 91, 100, 108, 115, 122, 129}; // sqrt(k/6) in percent, k=4..10
u32 q = (u32)((s2 * MF_QCOEF + 99) / 100);
q = (u32)(((u64)q * kscale[k - 4] + 99) / 100);
if ((u64)q > Qmem) q = (u32)Qmem;
if (q > N) q = N;
g_Q = q;
// LANE inputdist: the pool tables are dimensioned MAXQ; bound g_Q by them (a no-op when
// the formula is inside, and a CLAMP rather than an overflow when a test is larger than
// the statement's n). g_B follows the clamped value, so the partition stays consistent.
if (g_Q > (u32)MAXQ) g_Q = (u32)MAXQ;
g_B = (N + g_Q - 1) / g_Q;
}
const u32 Q = g_Q, B = g_B;
#if MF_PLEN
g_plen_on = 1; // independent of MF_WTRIM: lengths only, stride stays W
#endif
#if MF_WTRIM
u32 Wpool = W; // MF_WTRIM is KNOWN-BROKEN, see the warning above
#else
const u32 Wpool = W;
#endif
#if MF_FILTAUTO
if (!g_use_filt && B >= 400) g_use_filt = 1; // coarse grid: the filter beats the scan
#endif
if (g_use_filt && (u64)k * k * (g_Qc + 1) * W * 8 > ((u64)700 << 20)) g_use_filt = 0;
#if !MF_FILTAUTO
g_use_filt = 0;
#endif
{
int lo8[8], hi8[8];
for (int e = 0; e < 8; e++) lo8[e] = (e < k) ? 0 : 0x7FFFFFFF;
for (int e = 0; e < 8; e++) hi8[e] = (8 + e < k) ? 0 : 0x7FFFFFFF;
khi8v = LDU256(lo8);
khi8hv = LDU256(hi8);
}
// [SBK8] 12-byte rows (stride k=6) straddle cache lines; 16-byte aligned rows
// (stride STRIDE=8, 4 rows/line) do not. Judge-verified on 1013b: -0.93%.
g_SBK = (u32)((k < (int)STRIDE) ? (int)STRIDE : k);
// [MD] reciprocal for the per-point /B in the orders sweep: B is run-constant
g_bmag = ((1ull << 40) + g_B - 1) / g_B;
const int K1 = k - 1;
if (!g_blm_ready) {
for (int a = 0; a < MAXK; a++)
for (int b = 0; b < MAXK; b++) {
u32 m[8];
for (int e = 0; e < 8; e++) m[e] = ((e < b && e != a) ? 0xFFFFFFFFu : 0u);
g_blm[a][b] = LDU256(m);
}
g_blm_ready = 1;
}
// ---------- per-dimension orders, ranks, tie info ----------
for (int d = 0; d < k; d++) {
const unsigned *xd = x[d];
u32 *od = s_ord[d];
// [v4_orders] gt[i] = start[x_d[i]] and bd[] from start[x_d[i]+1]: value-indexed
// lookups, so the per-point work is two sequential sweeps of xd[] with no random
// od[]/xd[] traffic (the old loop walked od[] reading xd[od[j]] at random and then
// scattered gt[]/s_bd[] writes once per tie group).
for (u32 i = 0; i <= N; i++) s_cntv[i] = 0;
for (u32 i = 0; i < N; i++) { u32 v = xd[i]; s_cntv[v < N ? v : 0]++; }
{ u32 s = 0; for (u32 i = 0; i <= N; i++) { u32 c = s_cntv[i]; s_cntv[i] = s; s_tmp[i] = s; s += c; } }
// [VPACK] pack (end, start, cursor) of every value into one u64: the fused loop below
// then reads/writes a single 800 KB table instead of the 400 KB s_cntv + 400 KB s_tmp
// pair, i.e. one cache line per element instead of two. Bit-identical: cursor starts at
// the frozen start s_cntv[v] and advances exactly as s_tmp[v] did.
for (u32 v = 0; v < N; v++)
s_vtab[v] = ((unsigned long long)s_cntv[v + 1] << 34) | ((unsigned long long)s_cntv[v] << 17)
| (unsigned long long)s_cntv[v];
#if MF_FB2
for (u32 i = 0; i < N; i++) {
u32 v = xd[i]; if (v >= N) v = 0; u32 q = s_tmp[v]++; od[q] = i;
s_posp[(size_t)i * MAXK + d] = q;
}
#else
// The s_pos fill is NOT dead work (see the block comment above it): the reader is
// `const u32 *pd = s_pos[d];` in the coarse-filter build, which is compiled at the
// shipped MF_FB2=0 and reaches s_pos whenever g_use_filt is set. Keep the very same
// runtime guard, and keep the fill off the od[] path so it costs nothing.
// [ORDFUSE] The scatter (`od[s_tmp[v]++] = i`) and the gt[]/s_bd[] fill both walk xd[]
// in the same order and are independent: s_tmp is the mutable copy of the frozen value
// starts s_cntv, and s_cntv is not written in between. Fusing them removes one full
// sweep of xd[] and one loop's overhead per dim (k of them per run). Bit-identical:
// the fused loop writes exactly the same od[], gt[] and s_bd[] values in the same order.
if (g_use_filt) {
u32 *pp = s_pos[d];
u32 *gt = s_gt[d];
for (u32 i = 0; i < N; i++) {
u32 v = xd[i]; if (v >= N) v = 0;
unsigned long long w = s_vtab[v];
u32 q = (u32)(w & 0x1FFFFu), st = (u32)((w >> 17) & 0x1FFFFu), en = (u32)((w >> 34) & 0x1FFFFu);
od[q] = i; s_vtab[v] = w + 1ull; pp[i] = q;
gt[i] = st;
u32 n0 = en - 1u; if (n0 > N) n0 = (u32)N;
u32 bb = (u32)(((u64)n0 * g_bmag) >> 40);
s_bd[(size_t)i * g_SBK + d] = bb > Q - 1 ? Q - 1 : bb;
}
} else {
u32 *gt = s_gt[d];
for (u32 i = 0; i < N; i++) {
u32 v = xd[i]; if (v >= N) v = 0;
unsigned long long w = s_vtab[v];
u32 q = (u32)(w & 0x1FFFFu), st = (u32)((w >> 17) & 0x1FFFFu), en = (u32)((w >> 34) & 0x1FFFFu);
od[q] = i; s_vtab[v] = w + 1ull;
gt[i] = st;
u32 n0 = en - 1u; if (n0 > N) n0 = (u32)N; // [MD] clamp (0u-1u underflow)
u32 bb = (u32)(((u64)n0 * g_bmag) >> 40);
s_bd[(size_t)i * g_SBK + d] = bb > Q - 1 ? Q - 1 : bb;
}
}
#endif
BVAL(d, 0) = 0; FPOS(d, 0) = 0;
for (u32 b = 1; b < Q; b++) {
u32 pos = b * B < N ? b * B : N - 1;
u32 v = xd[od[pos]];
BVAL(d, b) = v;
u32 gs = pos;
while (gs > 0 && xd[od[gs - 1]] == v) gs--;
FPOS(d, b) = gs;
}
BVAL(d, Q) = 0xFFFFFFFFu; FPOS(d, Q) = 0;
}
g_SMP = (k <= 8) ? 8u : (u32)STRIDE;
for (u32 i = 0; i < N; i++) {
u32 *p = s_pm + (size_t)i * g_SMP;
for (int d = 0; d < k; d++) p[d] = x[d][i];
for (u32 d = (u32)k; d < g_SMP; d++) p[d] = 0;
}
// ---------- byte-signature planes for the fringe coarse filter ----------
// s_sig[d][e*N + j] = (u8)(x[e][ od_d[j] ] >> SIGSH) (od_d = dim-d order)
{
u32 nn = N - 1; u32 sh = 0;
while (nn > 255u) { nn >>= 1; sh++; } // smallest sh with (N-1)>>sh <= 255
g_sigsh = sh;
size_t need = (size_t)k * (k - 1) * MF_SPL;
static u8 *sigpool = 0; static size_t sigcap = 0;
// Refuse the signature path when the planes alone exceed 220 MB: past that the pool
// it competes with is squeezed and the row is better served by the exact path.
if (need > ((size_t)220 << 20)) { sigpool = 0; sigcap = 0; }
else if (need > sigcap) { sigpool = (u8 *)pool_alloc(need); sigcap = sigpool ? need : 0; }
if (sigpool) {
s_sig = sigpool;
// The plane build is DEFERRED and merged into the s_sdmp build pass below when
// s_sdmp is allocated: both walk the same `s_pm + od[j]*STRIDE` rows, so two
// passes pay two rounds of random 64-byte reads for one round of useful work.
// ⚠ SAFETY GATE. The signature path is VERIFIED CORRECT (brute force k=4..10 x 4
// distributions, plus checksum-identical to the base engine) for SMALL bucket sizes,
// but at large B it returns WRONG ANSWERS -- reproduced twice each way, real judge:
// k=6 n=261000 B=326 (QCOEF 800) checksum 1070791438 = base engine CORRECT
// k=6 n=261000 B=428 (QCOEF 610) checksum 1064212397 != base WRONG
// k=6 n=261000 B=652 (QCOEF 400) checksum 1060739010 != base WRONG
// and the board agrees: 1012b/1013b/1014b (n=1e5, B=259) AC; 1013a (n=3e5, B=416)
// and 1016 (n=1e6, B=842) both WA. The threshold sits between B=326 and B=428.
// Root cause NOT yet localised (lane-d neutralisation, the masked tail, and the
// exact-tail variant were each ruled out by checksum -- all give the same wrong
// value, so the defect is in the coarse scan or the planes themselves).
// UNTIL IT IS FOUND, REFUSE THE SIGNATURE PATH ABOVE A MEASURED-SAFE B.
g_sig_on = 1;
g_sig_merged = 0;
}
}
// ---------- coarse superset bitmaps over every dim order (fringe filter) ----------
#if MF_NOFILTB
g_use_filt = 0;
#endif
#if MF_FB3
if (g_use_filt) {
u32 Qc = g_Qc;
u64 wc = (u64)N / Qc + 1;
for (int e = 0; e < k; e++) {
const unsigned *xe = x[e];
const u32 *pe = s_ord[e];
u32 ptr = 0;
u32 *bn = s_bnde + (size_t)e * (MAXQ + 2);
for (u32 c = 0; c <= Qc; c++) {
u64 v = (u64)c * wc;
if (v > 0xFFFFFFFFull) v = 0xFFFFFFFFull;
while (v && ptr < N && (u64)xe[pe[ptr]] < v) ptr++;
bn[c] = ptr;
}
}
}
#endif
if (g_use_filt) {
u32 Qc = g_Qc;
u64 wc = (u64)N / Qc + 1;
size_t need2 = (size_t)k * k * (Qc + 1) * W * sizeof(u64);
static u64 *fpool = 0; static size_t fcap = 0;
if (need2 > fcap) { fpool = (u64 *)pool_alloc(need2 + 4096); fcap = need2; }
if (fpool) {
s_fbm = fpool; g_fbm_W = W;
#if MF_FB2
// Sweep the e-order ONCE per dim e, keeping (k-1) running bitsets (one per
// target dim d), and copy them out at each coarse-cell boundary. The old
// shape re-walked all N points for every one of the k(k-1) (d,e) pairs with
// two L3-random reads each (xe[pe[ptr]] and s_pos[d][pe[ptr]]); here the
// rank record is point-major so one 64B line yields every dim's rank, which
// cuts the random reads by ~k and the whole build by ~2-2.5x at k>=8.
for (int d = 0; d < k; d++)
for (int e = 0; e < k; e++) {
if (e == d) continue;
u64 *b0 = s_fbm + (size_t)(d * k + e) * (Qc + 1) * W;
for (u32 w = 0; w < W; w++) b0[w] = 0;
}
for (int e = 0; e < k; e++) {
const unsigned *xe = x[e];
const u32 *pe = s_ord[e];
int dmap[MAXK], km1 = 0;
for (int d = 0; d < k; d++) if (d != e) dmap[km1++] = d;
for (int i = 0; i < km1; i++) { u64 *a = s_facc + (size_t)i * W; for (u32 w = 0; w < W; w++) a[w] = 0; }
u32 ptr = 0;
for (u32 c = 1; c <= Qc; c++) {
u64 v = (u64)c * wc;
if (v > 0xFFFFFFFFull) v = 0xFFFFFFFFull;
while (ptr < N && (u64)xe[pe[ptr]] < v) {
const u32 *pv = s_posp + (size_t)pe[ptr] * MAXK;
for (int i = 0; i < km1; i++) {
u32 q = pv[dmap[i]];
s_facc[(size_t)i * W + (q >> 6)] |= 1ull << (q & 63);
}
ptr++;
}
for (int i = 0; i < km1; i++) {
int d = dmap[i];
u64 *dst = s_fbm + ((size_t)(d * k + e) * (Qc + 1) + c) * W;
const u64 *src = s_facc + (size_t)i * W;
#if MF_NT
{ u32 w = 0;
for (; w + 4 <= W; w += 4) _mm256_stream_si256((__m256i *)(dst + w), LDA((const __m256i *)(src + w)));
for (; w < W; w++) dst[w] = src[w]; }
#else
{ u32 w = 0;
for (; w + 4 <= W; w += 4) _mm256_storeu_si256((__m256i *)(dst + w), LDU256((const __m256i *)(src + w)));
for (; w < W; w++) dst[w] = src[w]; }
#endif
}
}
}
#else
for (int d = 0; d < k; d++) {
for (int e = 0; e < k; e++) {
if (e == d) continue;
const unsigned *xe = x[e];
const u32 *pe = s_ord[e];
const u32 *pd = s_pos[d];
u64 *base = s_fbm + (size_t)(d * k + e) * (Qc + 1) * W;
for (u32 w = 0; w < W; w++) base[w] = 0;
u64 *cur = s_acc;
for (u32 w = 0; w < W; w++) cur[w] = 0;
u32 ptr = 0;
for (u32 c = 1; c <= Qc; c++) {
#if MF_FB3
u32 pend = s_bnde[(size_t)e * (MAXQ + 2) + c];
#else
u64 v = (u64)c * wc;
if (v > 0xFFFFFFFFull) v = 0xFFFFFFFFull;
u32 pend = ptr;
{ const unsigned *xe_ = xe; const u32 *pe_ = pe;
while (pend < N && (u64)xe_[pe_[pend]] < v) pend++; }
#endif
while (ptr < pend) { u32 q = pd[pe[ptr]]; cur[q >> 6] |= 1ull << (q & 63); ptr++; }
u64 *dst = base + (size_t)c * W;
#if MF_NT
{ u32 w = 0;
for (; w + 4 <= W; w += 4) _mm256_stream_si256((__m256i *)(dst + w), LDA((const __m256i *)(cur + w)));
for (; w < W; w++) dst[w] = cur[w]; }
#else
for (u32 w = 0; w < W; w++) dst[w] = cur[w];
#endif
}
}
}
#endif
} else { g_use_filt = 0; }
_mm_sfence();
}
// ---------- fold dimension + bucket vector per query ----------
{
size_t mrn = (size_t)MAXK * MAXK * MAXQ;
for (size_t z = 0; z < mrn; z++) s_maxr[z] = 0;
}
for (u32 i = 0; i < N; i++) {
u32 best = 0xFFFFFFFFu; int ds = 0;
for (int d = 0; d < k; d++) {
u32 g = s_gt[d][i];
if (g < best) { best = g; ds = d; }
}
s_g[i] = (u32)ds;
// Largest R that will ever be read out of each (fold dim, dim, bucket) array:
// the pooled array (ds,d,b) is only ever scanned over [0, maxR] words, so the
// per-array length can be far below the group maximum (0.38n vs 0.93n at k=10).
if (g_plen_on) {
// [FOLDNB] The shipped form skips the fold dim with `if (d == ds) continue;` -- a
// data-dependent branch (ds is the argmin, uniform over the dims) that is
// mispredicted on nearly every query. Walk a precomputed k-1 table of the OTHER
// dims instead and update with a branchless cmov max; the final memory state of
// s_maxr is identical (max is commutative and associative).
const u8 *bdv = s_bd + (size_t)i * g_SBK;
u32 *mr = s_maxr + ((size_t)ds * MAXK) * MAXQ;
u32 bb = best;
const u8 *oth = s_oth[ds];
for (int t = 0; t < k - 1; t++) {
u32 d = oth[t];
u32 *pm = mr + (size_t)d * MAXQ + bdv[d];
u32 cur = *pm;
*pm = (bb > cur) ? bb : cur;
}
}
}
for (u32 i = 0; i < N; i++) s_qord[i] = i;
// group by ds
{
for (u32 v = 0; v <= (u32)k; v++) s_cntv[v] = 0;
for (u32 i = 0; i < N; i++) s_cntv[s_g[i] + 1]++;
{ u32 ac = 0; for (u32 v = 0; v <= (u32)k; v++) { u32 c = s_cntv[v]; s_cntv[v] = ac; ac += c; } }
for (u32 i = 0; i < N; i++) { u32 id = s_qord[i]; s_tmp[s_cntv[s_g[id] + 1]++] = id; }
for (u32 i = 0; i < N; i++) s_qord[i] = s_tmp[i];
}
u32 gstart[MAXK + 2];
{ u32 c = 0; for (int d = 0; d <= k; d++) { gstart[d] = c; while (c < N && (int)s_g[s_qord[c]] == d) c++; } gstart[k + 1] = c; }
// [REV-PORT-A] build the point-major signature table (rival #107799). Pure function of
// x[][] and g_sigsh; the value written into every signature plane is bit-identical to
// the shipped `(u8)(x[e][i] >> 9) ^ 0x80`, only the source row gets narrower.
g_sigrow_on = 0;
if (s_sig && k <= (int)STRIDE && N <= (u32)MAXN) {
for (u32 i = 0; i < N; i++) {
u8 *dst = a12b_sigrow + (size_t)i * STRIDE;
for (int e = 0; e < k; e++) dst[e] = (u8)((u8)(x[e][i] >> 9) ^ (u8)0x80);
}
g_sigrow_on = 1;
}
// ---------- transposed coordinate rows (sequential fringe) ----------
#if MF_NOSDMP
g_nosdmp = 1;
#endif
#if MF_NOSDMPB
g_nosdmp = 1;
#endif
if (!g_nosdmp) {
u32 KS = k;
#if MF_KSPAD
if (k <= 8) KS = 8; // 32-byte aligned candidate rows for k<=8
#endif
size_t per = (size_t)N * KS * 4;
size_t lim = (size_t)150 << 20;
if (k >= 9) lim = (size_t)470 << 20; // k=9/10 need the rows: the random s_pm fallback is fatal
if (per * k <= lim) {
static u32 *pool = 0; static size_t cap = 0;
if (per * k > cap) { pool = (u32 *)pool_alloc(per * k + 256); cap = per * k; }
if (pool) {
g_KS = KS;
for (int d = 0; d < k; d++) {
s_sdmp[d] = pool + (size_t)d * N * g_KS;
const u32 *od = s_ord[d];
u32 *dst = s_sdmp[d];
const u32 KSx = g_KS;
for (u32 j = 0; j < N; j++) {
#if MF_PFB > 0
if (j + MF_PFB < N) _mm_prefetch((const char *)(s_pm + (size_t)(od[j + MF_PFB]) * g_SMP), _MM_HINT_T0);
#endif
const u32 *row = s_pm + (size_t)(od[j]) * g_SMP;
u32 *dp = dst + (size_t)j * KSx;
for (u32 e = 0; e < KSx; e++) dp[e] = (e < (u32)k) ? row[e] : 0u;
if (s_sig) {
u8 *sb = s_sig + (size_t)d * (k - 1) * MF_SPL;
for (int sl = 0; sl < k - 1; sl++) { int e = sl < d ? sl : sl + 1; sb[(size_t)sl * MF_SPL + j] = (u8)((u8)(row[e] >> 9) ^ (u8)0x80); }
}
}
if (s_sig) g_sig_merged = 1;
}
}
}
}
// ---------- signature planes, standalone only if the s_sdmp pass did not run ----
if (s_sig && !g_sig_merged) {
const int use_sr = g_sigrow_on;
for (int d = 0; d < k; d++) {
const u32 *od = s_ord[d];
u8 *base = s_sig + (size_t)d * (k - 1) * MF_SPL;
for (u32 j = 0; j < N; j++) {
#if MF_PFB > 0
// [REV-PORT-A] rival lead is 64 against the 8-byte sigrow rows (one line holds 8
// points) instead of 8 against 32-byte s_pm rows; keep the shipped form otherwise.
if (use_sr) {
if (j + 64 < N) _mm_prefetch((const char *)(a12b_sigrow + (size_t)od[j + 64] * STRIDE), _MM_HINT_T0);
} else if (j + MF_PFB < N) {
_mm_prefetch((const char *)(s_pm + (size_t)(od[j + MF_PFB]) * g_SMP), _MM_HINT_T0);
}
#endif
if (use_sr) {
const u8 *row = a12b_sigrow + (size_t)(od[j]) * STRIDE;
for (int sl = 0; sl < k - 1; sl++) { int e = sl < d ? sl : sl + 1; base[(size_t)sl * MF_SPL + j] = row[e]; }
} else {
const u32 *row = s_pm + (size_t)(od[j]) * g_SMP;
for (int sl = 0; sl < k - 1; sl++) { int e = sl < d ? sl : sl + 1; base[(size_t)sl * MF_SPL + j] = (u8)((u8)(row[e] >> 9) ^ (u8)0x80); }
}
}
}
}
// ---------- pool ----------
u64 poolwords = (u64)K1 * (Q + 1) * W;
#if MF_PLEN
if (g_plen_on) {
u64 best = 0;
for (int dsx = 0; dsx < k; dsx++) {
u64 tot = 0;
for (int e = 0; e < K1; e++) {
int d = (e < dsx) ? e : e + 1;
const u32 *mr = s_maxr + ((size_t)dsx * MAXK + d) * MAXQ;
for (u32 b = 0; b <= Q; b++) {
u32 RR = mr[b];
u32 len = ((RR >> 6) + 8u) & ~7u;
if (len > W) len = W;
tot += len;
}
}
if (tot > best) best = tot;
}
if (best && best < poolwords) poolwords = best;
}
#endif
size_t need = (size_t)poolwords * sizeof(u64);
static u64 *pool = 0; static size_t cap2 = 0;
if (need > cap2) { pool = (u64 *)pool_alloc(need + 4096); cap2 = need; }
s_bs = pool;
// ---------- groups ----------
for (int ds = 0; ds < k; ds++) {
u32 q0 = gstart[ds], q1 = gstart[ds + 1];
if (q0 >= q1) continue;
#if MF_WTRIM
{
u32 mx = 1;
for (u32 qi = q0; qi < q1; qi++) { u32 rr = s_gt[ds][s_qord[qi]]; if (rr > mx) mx = rr; }
u32 wp = ((mx >> 6) + 8u) & ~7u;
Wpool = (wp < W) ? wp : W;
}
#endif
#if MF_PLEN
// Per-(dim,bucket) array lengths inside the group.
if (g_plen_on) {
u64 tot = 0;
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
const u32 *mr = s_maxr + ((size_t)ds * MAXK + d) * MAXQ;
u32 *po = s_pofs + (size_t)e * MAXQ;
u32 *pl = s_plen + (size_t)e * MAXQ;
for (u32 b = 0; b <= Q; b++) {
// NOTE: bucket 0 is the EMPTY prefix set (all-zero array), but a query whose
// bd[d]==0 still scans [0,R) words of it, so the array must be long enough --
// sizing it to the fixed 8-word minimum reads past the zeroed region into the
// next array. Use the real per-bucket max R here, not 0.
u32 RR = mr[b];
u32 len = ((RR >> 6) + 8u) & ~7u;
if (len > Wpool) len = Wpool;
po[b] = (u32)tot; pl[b] = len; tot += len;
}
}
if (tot > (u64)K1 * (Q + 1) * W) { g_plen_on = 0; } // never; safety
}
#endif
// bit numbering = rank in dim ds
{
const u32 *od = s_ord[ds];
for (u32 j = 0; j < N; j++) s_nid[od[j]] = j;
}
for (int e = 0; e < K1; e++) {
int d = (e < ds) ? e : e + 1;
const unsigned *xd = x[d];
const u32 *od = s_ord[d];
u64 *base = g_plen_on ? s_bs : (s_bs + (size_t)e * (Q + 1) * Wpool);
u64 *cur = s_acc;
for (u32 w = 0; w < Wpool; w++) cur[w] = 0;
if (g_plen_on) {
// bucket 0 is the empty prefix set; bd[] can be 0 so it must be materialised as zeros
u64 *d0 = s_bs + s_pofs[(size_t)e * MAXQ];
u32 l0 = s_plen[(size_t)e * MAXQ];
for (u32 w = 0; w < l0; w++) d0[w] = 0;
}
u32 ptr = 0;
#if MF_BUILD_NOPOOL
while (ptr < N) ptr++;
#else
for (u32 b = 1; b <= Q; b++) {
// {q : x_d[q] < BVAL(d,b)} is EXACTLY the prefix [0,FPOS(d,b)) of the dim-d order
// (FPOS is the tie-group start containing rank b*B), so the per-point value compare
// -- one random xd read per point per (ds,d) pair, k(k-1)N of them -- is redundant.
// b==Q means "every point" and FPOS(d,Q) is the 0 sentinel, hence N.
// Credited to agent a7f799ceb7e63b0cd (worth 4-5% at k=4/5 on the board).
u32 pend = (b < Q) ? FPOS(d, b) : (u32)N;
while (ptr < pend) {
// [BITSPF] prefetch the rank this loop will need D iterations later: the
// gather s_nid[od[ptr]] is a random access into a 400 KB map, and the only
// consumer of its result is the bit-set store into the L1-resident `cur`.
// D=24 keeps ~24 outstanding misses in flight without crowding the LFB.
if (ptr + 16 < pend) _mm_prefetch((const char *)(s_nid + od[ptr + 16]), _MM_HINT_T0);
u32 q = s_nid[od[ptr]]; cur[q >> 6] |= 1ull << (q & 63); ptr++; }
u64 *dst = g_plen_on ? (base + s_pofs[(size_t)e * MAXQ + b])
: (base + (size_t)b * Wpool);
u32 wlen = g_plen_on ? s_plen[(size_t)e * MAXQ + b] : Wpool;
#if MF_NT
// The pool is written once and read 89 GB later, so a plain copy pays a
// read-for-ownership DRAM read per line for nothing. NT stores remove it
// (measured: this loop was 4.58 GB / ~1.7 s, i.e. 2.7 GB/s, vs 12.5 GB/s
// for a plain memset).
{
u32 w = 0;
for (; w + 4 <= wlen; w += 4)
_mm256_stream_si256((__m256i *)(dst + w), LDA((const __m256i *)(cur + w)));
for (; w < wlen; w++) dst[w] = cur[w];
}
#else
{ u32 w = 0;
for (; w + 4 <= wlen; w += 4) _mm256_storeu_si256((__m256i *)(dst + w), LDU256((const __m256i *)(cur + w)));
for (; w < wlen; w++) dst[w] = cur[w]; }
#endif
}
_mm_sfence();
#endif
}
switch (k) {
case 4: query_group<4>(N, out, ds, q0, q1, Q, Wpool); break;
case 5: query_group<5>(N, out, ds, q0, q1, Q, Wpool); break;
case 6: query_group<6>(N, out, ds, q0, q1, Q, Wpool); break;
case 7: query_group<7>(N, out, ds, q0, q1, Q, Wpool); break;
case 8: query_group<8>(N, out, ds, q0, q1, Q, Wpool); break;
case 9: query_group<9>(N, out, ds, q0, q1, Q, Wpool); break;
default: query_group<10>(N, out, ds, q0, q1, Q, Wpool); break;
}
}
}
void count_4d(int n, const unsigned *x[4], unsigned *out) { solve_mf((u32)n, (const unsigned **)x, out, 4); }
void count_5d(int n, const unsigned *x[5], unsigned *out) { solve_mf((u32)n, (const unsigned **)x, out, 5); }
void count_6d(int n, const unsigned *x[6], unsigned *out) { solve_mf((u32)n, (const unsigned **)x, out, 6); }
void count_7d(int n, const unsigned *x[7], unsigned *out) { solve_mf((u32)n, (const unsigned **)x, out, 7); }
void count_8d(int n, const unsigned *x[8], unsigned *out) { solve_mf((u32)n, (const unsigned **)x, out, 8); }
void count_9d(int n, const unsigned *x[9], unsigned *out) { solve_mf((u32)n, (const unsigned **)x, out, 9); }
void count_10d(int n, const unsigned *x[10], unsigned *out) { solve_mf((u32)n, (const unsigned **)x, out, 10); }
| Compilation | N/A | N/A | Compile OK | Score: N/A | 显示更多 |
| Testcase #1 | 60.885 ms | 21 MB + 860 KB | Accepted | Score: 100 | 显示更多 |