/* References:
- saffah_cc_v41_agg1, https://duck.ac/submission/123022: copied latest1024 byte engine, retaining complete attribution and license notices.
- saffah_codex_6a_agg3, https://duck.ac/submission/125867: our six-bit Gray table builder and paired query traversal from4096 byte multiplication.
Idea: Use twelve five-bit groups plus one four-bit group, adapting our https://duck.ac/submission/125931 Gray builder and paired queries to the padded parity input. This shrinks the active table to12.5KiB so it can coexist with the8KiB A masks and8KiB output in L1; the additional two lookups replace256 table entries. Transfer the register-only correction to this engine's65-row padded parity input. Replace the all-word four-bit lookup slab with one20.5KiB slab per input word and paired two-chain queries. The prior local six-bit experiment built tables by reading old entries; this version uses Gray XORs and no table reads during construction. Purpose: reduce table query count and build latency while preserving the accepted folded multiply kernel and every parity bit.
*/
// ===== REFERENCES (mc2i_ 席, 2026-10-03) =====
// [1] duck.ac 本账号 saffah_cc_v41_agg1，提交 #122880
// <https://duck.ac/submission/122880>（6.837403 ms，本账号当前最好件 = mine）
// 用途：**本发基座**（正文一字未改，只加 4 条 `prefetchw`）。件 =
// problems/mmmc1k/work/mc2g__tmp/sub_p384cm.cpp（= work/mc2i__tmp/base.cpp）。
// [2] 全队已兑现律「**`prefetchw` 必须用内联汇编**」——出处：x11z_ 席 1011e6 封条：
// gcc9 把 `__builtin_prefetch(p, 1, ·)`（rw=1）编成 `prefetchnta`、**不取独占行**
// ⇒ 随后的 RMW/普通 store 仍付整次 RFO；改内联汇编 `prefetchw` 后 1011e6 单刀
// **−16.129 ms**（相位口径 −4.4%）。用途：**刀形的出处**（本席按本题实测重新定价）。
// [3] BRIEF 在册「按 DRAM 字节计价（37 ms/GB）· 每元素每省 1 B ≈ 6.9 ms 板面 ·
// 写 1 B 付 2 B 流量（RFO）⇒ 省写 1 B 的价 = 省读 2 B」(n08h_ 1008e8) +
// 「板面 = 探针 + 页错税 + stdout + 静态表」(e7d_ 1004e7)。
// 用途：**机理依据** —— 打包器面板写是「首次写未驻留行」= 每 32B store 付一次
// 整行 RFO（读为写，37 ms/GB 量级）。
// ======================
// ===== 思路 =====
// **单变量：给两个打包器的面板**写目的地**加 4 条内联汇编 `prefetchw`（数值、布局、
// 循环序、访存宽度一字未改；只提前把目的行取成独占态）。**
// (a) `buildA4` 内层 k 循环（k += 64，每步写 `dst+2k` 共 128 B = 2 行）：
// `prefetchw dst+2k+256` / `dst+2k+320`（= 前 2 个 k 步的 2 行，128 B 全覆盖）。
// (b) `buildB4` 内层 t 循环（t += 32，每步写 `o+2t` 共 64 B = 1 行）：
// `prefetchw o+2t+256` / `o+2t+320`（= 前 2 个 g 步的 2 行）。
// 机理：`A0b` 128 KB/叶、`BI3b` 64 KB/叶，每次都是**整块重写**且上一叶之后已被
// mk_s 的 8 MB 级 L1 面板读逐出 ⇒ 每个 32B store 都打在**未驻留行**上 = 付整行 RFO
// （把 32 B 写放大成 64 B 流量 + L3 往返延迟）。`prefetchw` 把该行的独占取得提前到
// store 之前几百字节，store 命中独占行 ⇒ RFO 延迟对 store 缓冲不再可见。
// **判题机实测（`tools/tryrun.py` 通道，同进程 rdtscp 包夹整次 `matrix_multiply`，
// 轮转臂序、同批锚 base；本批同臂重复散布 ≤0.005%）**：
// base 6.488578 ms（r0 6.538462）
// +A4+A4 两条 prefetchw（vPFA4）**−0.384%** · +B4 两条（vPFB4）**−0.747%** ·
// 两条都加（**本发**）**−1.116% ~ −1.133%**（r0 亦同号：6.5385 → 6.4645）
// · 距离扫描：前 1 步 −1.157%（r0 6.4763）· 前 2 步 −1.133%（r0 **6.4645**，
// 本发取此：板面是**单次冷跑**，以 r0 为准）· 前 3 步 −1.071% · 每步只 1 行 −0.637%
// （⇒ 覆盖不足，2 行/步是必需的）。
// 另加 BP2/A1b/epilogue/XT4 四处目的地的 prefetchw（vPFAB2/3/4）**均为 0**（饱和）。
// 逐位等价：只加 hint、不动任何数据通路 ⇒ 与 #122880 输出逐字节相同（闸门 6 形态
// `gate_main.inc` 6/6 IDENTICAL + NONDEG_ALL=1；4 个变体各自通过）。
// **本轮判负留档（勿重走，均已逐位闸门）**：题面指定的「`mk_s` K 维分块 + 部分和、
// 摆脱 16 KB 切片 L1 常驻」方向**经实测不可行**：(1) 判题机消融臂 `dNoBc`（删 16 条
// `vpbroadcastd`/轮）= **+0.039% = 零** ⇒ **载入端口压力不是本题的税**，任何「减载入」
// 或「改 B 面板分块粒度」的重写都无收益；(2) `vExpl`（每 32B 块先 `vmovdqa` 再两条
// 寄存器源 madd，载入 80→48/轮）= **+0.531%**（fused uop 148→180 撞 4 宽 rename）；
// (3) `dNoAdd`（删 64 条 `vpaddw`，诊断臂）= **−7.44%** ⇒ 本题的**唯一货币是 ALU 向量
// op 数**，而 `mk_s` 的 (64 madd + 64 add)/迭代 = **每 32 MAC 付 2 ALU op = AVX2
// `vpmaddubsw`+`vpaddw` 硬地板**（3 端口 × 16 MAC/op = 48 MAC/拍），实测 92.4% 已达；
// K 维分块/部分和只会**增加** ALU op（归约/部分和搬运）；其目的（解开生产者行分带）
// 已被 mc2h_ 直接实验否掉（源读 9.2→3 MB 实测 **+0.196%**，流量守恒）。(4) 其余判负：
// G16 展开（环体 ×2）**+0.329%** · 去 mk_s 预取 **+0.188%**（预取承重）· 64B 对齐 +0.006%。
// ======================
// ===== REFERENCES (mc2g_ 席, 2026-10-03) =====
// [1] duck.ac 本账号 saffah_cc_v41_agg1，提交 #122256
// <https://duck.ac/submission/122256>（6.852133 ms，本账号当前最好件 / mine）
// 用途：**本发基座**（正文一字未改，只改一处数据布局）。件 =
// problems/mmmc1k/work/mc2a__tmp/mc2a_sub1.cpp（= work/mc2c__tmp/base.cpp，
// md5 87b7ea4798368c0ae2135d1a598382ca，本席复跑逐位闸门通过）。
// [2] 本队三题已兑现的「7 条同相流行距填充/页对齐」布局刀：
// * mmmc4k 席 `mc4a_`（封条 (1830)）：`M[7][体积+64]` → `+24576`，判题机 tryrun
// 实测 −3.79 ms（−0.86%）；同批扫描 0 → +1.94（比 +64 更差）· 4096 −2.32 ·
// 8192 −3.19 · 16384 −3.77 · 24576 −3.79 · 32768 −2.88。
// * mmmc2k 席 `mc2f_`（封条 (1908)）：`S2M[7][1024*1024+64]` → `+24576`。
// * 1014b 席 `b14am_`（封条）：`s_sig` 逐维行距 → 4096 B 页对齐，板面 −0.698 ms。
// 用途：**刀形的出处**（本席按本题实测重新定价，不跨题外推）。
// [3] 本账号 `mc2c_`/`mc2d_`/`mc2e_` 席（problems/mmmc1k/notes.md 末三节）在册结论：
// `mk_s` 真下界 = (64 madd + 64 add)/3 = 683 拍/调用（`vpbroadcastd m32` 是载入
// uop、不占 ALU 端口）⇒ 板面 = mk_s 地板 4.760 ms（69.5%）+ 非 mk_s 2.092 ms；
// 非 mk_s 纯延迟绑定 ⇒ 「减 uop/减调用数」整族零收成。用途：本席立项依据 ——
// 缺口只能从**非 mk_s 的内存形态**里拿，故取最便宜的布局刀先试。
// ======================
// ===== 思路 =====
// **单变量：`M2[7][256*256]` → `M2[7][256*256+P]`（只在数组第二维加常数 P，代码一字未改）。**
// 机理（本题独有、与姊妹题同形）：n=1024 走 `strassen1` ⇒ 7 × `strassen2_512x`，每次
// (a) `submm512_4` 把 7 个 256×256 叶乘积写进 `M2[0..6]`（每个叶只写自己那一块）；
// (b) `qwall` **同时读 M2[0..6] 七条流、索引完全相同**（`M2[k]+ i*h`，h=256）
// —— 七条流的地址差恰为 65536 = 2^16：
// L1 组号 = (addr>>6)&63 ⇒ 2^16 是 64 的整数倍 ⇒ **七条流全落在同一个 L1 组**
// （再加 qwall 的 4 条 C 目的流 e1..e4，其偏移 0/256/h*n/h*n+h 也都是 4096 的
// 整数倍 ⇒ 同一组共 9~11 条线 > 8 路）；页内偏移亦完全相同 ⇒ DRAM 通道/体
// 哈希位不分摊。⇒ 原版 `P=0` 正是「七条同相流挤同一组」的最坏形态。
// 扫描（判题机「立即测试」通道，`tools/tryrun.py` 轮转臂序、每臂 2~3 发、臂内取
// r1..r5 最小值；**同一进程内 rdtscp 包夹整次 matrix_multiply**，故只看引擎、不含驱动）：
// 臂内重复噪声 0.004%（r1..r5 散布 0.03%）⇒ ≥0.05% 效应可判。
// P=0（原版）6.536779 ms · 64 −0.387% · 128 −0.073% · 192 −0.37% · **256 −0.650%** ·
// 288 −0.32% · 320 −0.425% · 352 −0.29% · **384 −0.66%** · 416 −0.29% · 448 −0.36% ·
// 512 −0.634% · 768 −0.33% · 1024 −0.263% · 2048 −0.321% · 4032 −0.417% ·
// 4096 **+0.132%**（4096 的整数倍不分离 L1 组 ⇒ 与机理一致）· 8192/16384/24576/32768 −0.18~−0.24%。
// **对照（证明差距不是「BSS 长大」）**：同批加 `dummy[86016]` / `dummy[172032]` 两个
// 无关 BSS 尾数组（BSS 增量分别 ≈ P=12288/P=24576 的 7 倍）⇒ −0.006% / +0.019%
// （噪声底）⇒ **效应全部来自 M2 自身行距**，与堆/输入数组的物理位置无关。
// 逐位等价：每个 M2 块内部的行距（h=256）与全部索引一字未动，`M2[k]` 仍是行指针；
// 七块只被整体平移 P —— 输出与 #122256 逐字节相同（下方闸门）。
// **第二处（同一发内）：`corr4_256` 四趟合并成一趟。** 原式 = 对 4 个 kb 组
// （kb=0,64,128,192）各建一次 16×16 的 `XT4_256` nibble 表 + 各扫一遍 256 行、
// 每行对 `Pbb` 做一次 `load^xor^store`（另有开头 `memset(Pbb,0,8KB)`）。本发把
// 四张表同时建好（`XT4_256[4][16][16]` = 32 KB），行循环只跑**一遍**、四个累加器
// 在寄存器里 XOR 成最终值后**一次 store**（memset 与全部 RMW 消失）。GF(2) 加法
// 可交换可结合、且 Pbb 每调用被全量（n*W*8 B）重写 ⇒ **逐位相同**（闸门同 6/6）。
// 判题机实测：单这一处 −0.015%（噪声级）；叠在布局刀上 −0.04%（3 发同号）。
// ======================
// ===== REFERENCES (mc2a_ 席, 2026-10-03) =====
// [1] duck.ac 本账号 saffah_cc_v41_agg1，提交 #122138
// <https://duck.ac/submission/122138>（6.883430 ms，本账号当前最好件）
// 用途：**本发基座**。件 = problems/mmmc1k/work/mc1y__tmp/f1.cpp。
// [2] duck.ac 本账号 saffah_cc_v41_agg1 席 `mc1z_`（problems/mmmc1k/notes.md 末节）自产候选
// `problems/mmmc1k/work/mc1z__tmp/cand_qw.cpp`（= #122138 正文 + qwall 融合）。
// 用途：本发采纳其 **qwall 融合**（4 趟象限写合并成 1 趟按行扫；`qw<>` 4 次调用/
// k1 改成 `qwall<M1,M2>` 1 次；#122138 里 12 块 M2 读/k1 → 7 块）。逐字节闸门由该席
// 6/6 通过、并由本席复跑通过（见下）。其余一字未动。
// ======================
// ===== 思路（本发单变量 2 处，全部逐位等价） =====
// 本席（mc2a_）承接交棒「打包器的 L3 延迟」轴，先做分件与相位实测（本机 gcc9 + rdtsc 相位仪）：
// a4=5.6% · b4=5.6% · corr=4.7% · mk_s=74.5% · epi≈0 · qwall=4.2~5.1%（整题）。
// 打包器热/冷对照：buildA4 25.0K 拍/调用（与 mk_s 交错）↔ 17~18K 拍/调用（源数据热）
// ⇒ 交错的惩罚 = 28%，但源象限按 k2 的并集最大 1 MB > 判题机 L2（256 KB）⇒ 行条带化
// 形态不可得（7 份 A 面板 896 KB / 7 份 BI 面板 448 KB 都放不进 L2），本席实测的
// `(i 外, g 内)` 换序更把 D1mr 从 6.8M 打到 43M（BI 切片被逐出 L1）⇒ **现有循环序是最优**。
// 单变量 A：`qwall` 加软件预取（7 条 M2 源流各 +2 行、8 条目的行 +2 行）。
// 单变量 B：`corr4_256` 的 16 项 nibble 表**构表**改用"4 条连续 BP2 行的基向量 + 11 条 xor"
// （原式对 15 个条目逐个重算 `(k>>6)*(65*W)+(k&63)*W` 地址，~10 条标量 op/条目；
// 4g+3<=63 ⇒ 4 个基 k 恒在同一 64 行组内 ⇒ 地址是 b0+j*W 的等差，可一次算出）。
// 两处**均不改数值、不改内存布局**（结果逐位相同）。Ir 85,201,866 → 85,070,840（−0.154%）。
// ======================
// ===== 闸门（逐字节，本席复跑） =====
// `mc1z__tmp/gate_main.inc`（splitmix64 六形态 × n=1024，两个哈希 + 四象限和 + 非退化断言）
// ⇒ 本发件与 cand_qw/#122138 **6/6 IDENTICAL**；阳性对照（KDEST[3][1].M 1→−1）**DIFFER** ✓。
// ======================
// ===== 判题机定价（tools/probe.py --reps 3 --json，同窗口同驱动 r0/热调用） =====
// base #122138: r0=23,832,542 · 热均值 23,628,446 · 整进程 29.9745 ms
// 本发件 : r0=23,637,652 · 热均值 23,483,255 · 整进程 29.8318 ms
// （cand_qw 单臂：r0=23,649,504 · 热均值 23,498,883，与本发件差别 = 单变量 A+B）
// ======================
// ===== REFERENCES =====
// [1] duck.ac 本账号 saffah_cc_v41_agg1，提交 #120323
// <https://duck.ac/submission/120323>（6.945315 ms，本账号当前最好件）
// 用途：**本发基座**。件 = problems/mmmc1k/work/mc1x_tmp/mc1x_sub_fold.cpp
// （= #106545 + fold 单变量，逐字节闸门通过）。本发正文 = 该件正文 + 下述单变量改动。
// [2] duck.ac 用户 saffah_codex_6a_agg3，提交 #120247
// <https://duck.ac/submission/120247>（6.976894 ms，本题当前 T）
// 用途：**本发所叠加改动的出处** —— `mk_s` 内层循环每轮 A 值 **4→8**：
// (a) `long cnt=n>>3` → `n>>4`；
// (b) 循环体由 4 组（广播偏移 0/4/8/12，`addq $16 %[x0]`、`addq $16 %[x1]`、
// `addq $512 %[bi]`）扩为 8 组（偏移 0/4/8/12/16/20/24/28，
// `addq $32`、`addq $32`、`addq $1024`）；其余逐字全同。
// 取证：去注释逐行 `diff` 后，该件与**我方 #106545** 只差这**一个内容块**（2 个 hunk，
// 其余 92 行全是注释/头），即"我方旧件正文 + 对手单点改动"。该改动在**未折叠**基座上的
// 判题机当量 = #106545 7.058140 → #120247 6.976894 = **−1.152%**。
// [3] 本账号自产：`c1w_` 席（problems/mmmc1k/notes.md，2026-09-28 16:2x 节）的 **fold 刀**
// （= #120323 相对 #106545 的唯一改动：16 条 `vmovdqa N(%[bi]),%%ymm4` 折进紧随两条
// `vpmaddubsw` 的 src2 内存操作数）。判题机实测 −1.598%（112.825 µs）。
// [4] Intel SDM / Agner Fog 端口表 + 本席计数：`mk_s` 折叠后每轮 = 32 `vpmaddubsw`(p01) +
// 32 `vpaddw`(p015) + 8 `vpbroadcastd`(p23) + 32 madd 内存操作数(载入 p23) + 5 控制。
// ======================
// ===== 思路 =====
// **试验性提交（目的：给"抄对手 + 复刻自家刀"的组合在判题机上定价）**。
// 单变量：把对手 #120247 的 4→8 展开**叠在我方 fold 刀之上**（其余一字不动）。
// 1. 逐位等价：展开只是把同一批 (广播, k 块, madd, add) 拷贝一遍、把指针步长与计数同步翻倍，
// madd/add 的**相对次序、每次的地址与操作数**一字未动 ⇒ 与 #120323 输出逐字节相同
// （本席闸门：6 种数据形态 × n=256/512/1024，全部 IDENTICAL，含阳性对照）。
// 2. 为什么还要发（机理与判题机定价的冲突）：我方 fold 已把融合域条数 93→77、前端下界
// 23.25→19.25 拍/轮，**伏在端口下界 21.33 之下** ⇒ 按端口模型，展开省下的那 2.5 条/轮
// 应当是**零剩余**（本机两臂配对实测正是 ≈0：−0.40% / −0.15% / −0.07%，符号不稳定）。
// 但 #120247 证明**判题机对"更少的前端条数"的定价高于本机**（同一改动本机 −0.83% ↔
// 判题机 −1.15%，比值 1.39），所以"展开在判题机上仍有残值"是一个**只能由判题机裁决**的
// 假设（§2.18.349「读数属于某一代/某一台机器」）。
// 3. 若假设成立，本发应比 #120323 再快一个展开当量的一部分；若假设不成立，本发与 #120323
// 同级（不劣化 mine：判题机取历史最优）。
// ======================
// ===== mc1y_ 本发单变量：把顶层 Strassen 折进 C（删掉 MMS[7][512*512] = 1792 KiB） =====
// 起因（本席实测，先证伪上一手交棒的前提）：`mc1y_diag.py` 逐类删指令（本机 gcc9，min-of-3）：
// 删 16 条广播 = **+0.0%** ；删 64 条 madd 内存操作数（改 register）= **+0.0%** ；
// 删 64 条 `vpaddw` = **−8.5%**。⇒ `mk_s` 循环是 **ALU/端口受限**，不是载入/前端受限。
// 配套：`mc1y_split.py` 桩测 = 打包器 11.3% · `cmb` 9.0% · 其余（含 `mk_s`）≈80%。
// ⇒ mc1x_ 交棒的「B 面板 4-k 布局」其所依赖的两项（广播 8→4、载入 40→34）**都是免费的**，
// 而该布局本身要把 A 面板的 dword 周期变成 4 个不同 A 值 ⇒ 每 madd 只覆盖 8 列
// （偶/奇 16 位道分别是 k0k1 / k2k3 两个半和，必须回并），或需补一条 p5 洗牌 ⇒
// **在 ALU/p5 受限的循环里净亏**（旧 `vpbroadcastw` 刀的 −3.72% 正是 p5 之争）。
// ⇒ 改攻内存轴：判题机对**被写页**收费（本席实测标定 0.34 µs/页，与 BRIEF「0.244 µs/页 +
// 省 1 MB = 省 61 µs」同族），而 MMS 占本引擎 3516 KiB 足迹的 **51%**。
// 单变量：`strassen2_512(MMS[k1], ...)` + 末尾 4 条 `cmb_nt` ⇒ 改成模板 `strassen2_512x<K1>`
// 把该 512×512 积块的 4 个象限**直接写进 C 的对应目的块**（表 `KDEST[7][2]`：R/Cc/模式，
// 模式 0=store 1=add −1=sub 即 P-form 里该块的符号），一个有 2 个目的的 k1 由 `qw<>` 的
// 一次源读取同时喂两个目的 ⇒ 源读量不变；`MMS` 与 `cmb_nt` 通道整体消失。
// 闸门（逐字节，`c1w_dump.cpp` 驱动，6 数据形态 × n=256/512/1024）: 18/18 IDENTICAL
// （n=1024 覆盖被改路径；阳性对照把 KDEST[3][1].M 由 1 改 −1 ⇒ n=1024 **DIFFER** ⇒ 闸门活）。
// 定价（**判题机**，`tools/probe.py` 同文件同臂、单调用 rdtsc 内计）：
// 冷调用 24.211 Mticks → **23.788 Mticks（−1.75%）** · 热调用 23.455 → 23.575 Mticks（+0.51%）
// · 引擎自付页费（冷−热）217 µs → **61 µs** · `mem_kb` 5588 → **3796（−1792 KiB）**。
// 本机 gcc9 min-of-N：5.856 → 5.612（−4.2%，本机只作相对量）。
// 预期判题机当量 ≈ **−120 µs**（单次冷调用口径，页费已含在冷调用里，不与页模型重复计）。
// 否决记录（本席实测定负、已弃）：目的块首写改 `_mm256_stream_si256` = **+41 µs**（NT 版热调用
// 23.718 vs 23.575 Mticks）⇒ 保留普通 store。
// ======================================================================================
/* References:
* saffah_cc_v41_agg1, https://duck.ac/submission/106545 : copied split-byte
* Strassen engine, original microkernel and original reference notices below.
* No independent license notice in the original public source.
* Idea: double the hot kernel's K unrolling from4 to8 pairs, halving pointer-update
* and loop-control instructions per product. Official generated-input profiling
* found the leaf multiply stage dominates79 percent of cycles and is front-end bound.
* Experiment: reduce loop overhead with identical arithmetic and output layout.
*/
// ===== REFERENCES =====
// [1] 本席 c1w_（mmmc1k，2026-09-28）。基线 = 本账号 #106371（work/l1x_dt_sub.cpp）。
// 出处（机制）：判题机多臂诊断臂（前席 nn_，notes「诊断臂：把那 8 条 vpbroadcastw 换成
// vmovdqu ⇒ −13%~−19%」）+ 本机端口标定（notes nn_：「vpbroadcastw 1/cycle，只 p5」）。
// [2] Intel SDM / Agner：Haswell..Coffee Lake 的载入端口只能原生广播 32/64/128 位；
// `vpbroadcastw m16` = 2 uop（载入 p23 + 洗牌 p5），`vpbroadcastd m32` = 1 uop（p23）。
// ======================
// ===== 思路 =====
// **把 mk_s 每轮的 8 条 `vpbroadcastw m16` 换成 `vpbroadcastd m32`**（单变量：只改广播指令与
// A 面板的存放形态，其它一字不动）。
// A 面板（A0b）把每个 16 位 A 对**存两遍**（pair P 落在字节 4P），于是
// `vpbroadcastd 4p(%[x0])` == 旧的 `vpbroadcastw 2p(%[x0])` **逐位相同**（32B 操作数一致）。
// 目的：把 8 条 p5 uop 从内层循环里拿掉。内层是 p5 饱和的（8 条广播 + vpaddw 在 p5 上的份额
// ⇒ 端口下界 max((32+A1)/2, 8+A5), A1+A5=32 ⇒ 24 拍/轮）；拿掉广播后下界降到 max(23.25 前端,
// 21.33 端口) ⇒ 理论上把那个 24 拍的墙拆掉。
// 代价：A 面板 2×（buildA4 的写流量 64KB→128KB/次，读不变 256KB/次）⇒ 预计 +60 µs 内，远小于收益。
// ===== REFERENCES =====
// [1] duck.ac 用户 saffah_codex_6s_agg2，提交 #106361 <https://duck.ac/submission/106361>
// 用途：本发唯一改动的**出处**（其头部自述 "dual_temp_revadd"，官方 Try It 12.029930 -> 11.867511 ms，
// 输出 hash 相同）。取证：`diff`（去注释）后它与本账号 #105098 只差 **三类**改动 ——
// (A) `mk_s` 内联 asm 里 16 处 `vpmaddubsw ... %ymm5` 与随后的 `vpaddw %ymm5, %ymm1[2345], ...`
// **同时**改成 `%ymm6`（= 其自述的 dual-temp 排程）；← **本发只取这一处**
// (B) `corr4_256` 的 16 路 xor 展开**退回**成 rolled 循环；
// (C) `corr4_256` 的 memset 删除 + `kb==0` 直存（本队史测：在**我方基座**上"去 memset" = **+59.9 µs**）。
// ⇒ (B)(C) 是**他的旧形态**、对**我方基座**是负向 ⇒ 他的 7.324042 不能直接归因给 (A)。
// [2] duck.ac 用户 saffah_cc_v41_agg1（本账号），提交 #105098 <https://duck.ac/submission/105098>
// 用途：**本发的基座**（7.490983 ms）。本发正文逐字节 = 该件正文 + (A) 的 16 处替换，其余一律还原。
// ======================
// ===== 思路 =====
// **本发 = 只取对手 #106361 的命名刀 (A)，把它的两处旧形态还原成我方 #105098 的形态**（单变量，逐位闸门）。
//
// 机理：`mk_s` 每读一个 32B b 分块，要为**两列**各做一次 `vpmaddubsw`，两次都把结果写进**同一个** `%ymm5`
// 再各自 `vpaddw` 累加。改后：**累加进 `%ymm12..%ymm15`（第二列）的那一路改用 `%ymm6`** ⇒
// 两列乘积不再共用同一个临时寄存器（他自述 "dual temporary before adds"）。
//
// 取证定位（§2.19.018 第一动作，先做再动手）：
// `diff <(grep -v '^\s*//' 我方#105098) <(grep -v '^\s*//' 对手#106361)` = 64 行，**全部**是这 16 对的
// `ymm5 -> ymm6`；此外 3 处 corr 侧差异**是他在退回我方的既有优化**（见 [1]）。
//
// 本发的单变量证明：`diff <(去注释 我方#105098) <(去注释 本发)` = **恰好 64 行 / 32 处替换**，无其它差异。
//
// 闸门（本地，逐字节）：
// `g++ -O2 -static -std=c++17`（引擎在前、dump driver 在后**同一 TU**）
// n=1024 LCG 数据 ⇒ 两件 stdout 的 1,048,576 字节 **`cmp` 完全一致**；
// 同进程 2 臂探针（rdtsc，R=6）：FNV = `34fd0a86984a4fc3` **两臂相同**（= 本 lane 历史常数 ✓）。
// ★ 本地探针读数在噪声内（min −2.93% / med +1.1% 反号）⇒ **本机 gcc15 不作为定价依据**；
// 定价依据 = 对手 #106361 的**判题机**实测：7.578116 → 7.324042，且 (B)(C) 对我方基座为负
// ⇒ **(A) 对我方的保守下界 = 7.490983 − 7.324042 = 0.166941 ms（−2.23%）** ✓
// ================
/* mmmc_k1: int8 GEMM C = A*B mod 256, EXACT, AVX2.
Split A = A0 + 128*A1 (A0 = A & 0x7F unsigned, A1 = top bit of each byte).
term1 = sum_k A0[i][k]*B[k][j] -> vpmaddubsw, no saturation (127*128*2 < 32768)
term2 = 128*(A1 . Bpar) mod 256 -> GF(2) matmul, 64 bits/word (Bpar = B&1)
Layout (lanes = columns): accum lane m is column j0+m.
x = vpbroadcastw of the 2-byte pair A0[i][2p..2p+1]
y = 32B load from BI[p] at 2*j0 ; BI[p][2j]=B[2p][j], BI[p][2j+1]=B[2p+1][j]
*/
#define MC4_KNIVES_DE 1
#include <immintrin.h>
#pragma GCC optimize("O3,unroll-loops,no-strict-aliasing,rename-registers")
#pragma GCC target("avx2,bmi,bmi2,popcnt,lzcnt")
#include <stdint.h>
#include <string.h>
#include <stdio.h>
#define AVX2 __attribute__((target("avx2")))
#define NMAX 4096
static int8_t A0b[2*NMAX*NMAX] __attribute__((aligned(64)));
static uint64_t A1b[NMAX*NMAX/8] __attribute__((aligned(64)));
static uint64_t Bpb[NMAX*NMAX/8] __attribute__((aligned(64)));
static uint64_t Pbb[NMAX*NMAX/8] __attribute__((aligned(64)));
static uint8_t PBb[NMAX*NMAX] __attribute__((aligned(64)));
static void generic_mm(int n, const int8_t*A, const int8_t*B, int8_t*C){
for (int i=0;i<n;i++) for (int j=0;j<n;j++){
int s=0; for(int k=0;k<n;k++) s += (int)A[(size_t)i*n+k]*(int)B[(size_t)k*n+j];
C[(size_t)i*n+j]=(int8_t)s;
}
}
AVX2 static inline __m256i bcast16(const int8_t*p){
__m256i v; __asm__("vpbroadcastw %1, %0" : "=x"(v) : "m"(*(const short*)p)); return v;
}
AVX2 static inline __m256i packlo8(__m256i a, __m256i b){
const __m256i m=_mm256_set1_epi16(0x00FF);
return _mm256_permute4x64_epi64(
_mm256_packus_epi16(_mm256_and_si256(a,m), _mm256_and_si256(b,m)), 0xD8);
}
AVX2 static void build_A0(int n, const int8_t*A){
const __m256i m=_mm256_set1_epi8(0x7F);
size_t tot=(size_t)n*n;
for (size_t i=0;i<tot;i+=32)
_mm256_store_si256((__m256i*)(A0b+i),
_mm256_and_si256(_mm256_load_si256((const __m256i*)(A+i)), m));
}
AVX2 static void build_A1bits(int n, const int8_t*A){
int W=n/64;
for (int i=0;i<n;i++){
const int8_t*row=A+(size_t)i*n; uint64_t*w=A1b+(size_t)i*W;
for (int k=0;k<n;k+=64){
uint32_t lo=(uint32_t)_mm256_movemask_epi8(_mm256_load_si256((const __m256i*)(row+k)));
uint32_t hi=(uint32_t)_mm256_movemask_epi8(_mm256_load_si256((const __m256i*)(row+k+32)));
w[k>>6]=(uint64_t)lo|((uint64_t)hi<<32);
}
}
}
AVX2 static void build_Bpbits(int n, const int8_t*B){
const __m256i one=_mm256_set1_epi8(1);
int W=n/64;
for (int k=0;k<n;k++){
const int8_t*row=B+(size_t)k*n; uint64_t*w=Bpb+(size_t)k*W;
for (int j=0;j<n;j+=64){
__m256i a=_mm256_slli_epi16(_mm256_and_si256(_mm256_load_si256((const __m256i*)(row+j)),one),7);
__m256i b=_mm256_slli_epi16(_mm256_and_si256(_mm256_load_si256((const __m256i*)(row+j+32)),one),7);
uint32_t lo=(uint32_t)_mm256_movemask_epi8(a);
uint32_t hi=(uint32_t)_mm256_movemask_epi8(b);
w[j>>6]=(uint64_t)lo|((uint64_t)hi<<32);
}
}
}
/* P[i] = XOR of Bp[k] over k with A1[i][k] set. NY = n/256 ymm per row. */
AVX2 static void corr_bits(int n){
int W=n/64, NY=W>>2;
for (int i=0;i<n;i++){
__m256i acc[16];
for(int t=0;t<NY;t++) acc[t]=_mm256_setzero_si256();
const uint64_t*row=A1b+(size_t)i*W;
for(int w=0;w<W;w++){
uint64_t b=row[w];
while(b){
int t=__builtin_ctzll(b); b&=b-1;
const uint64_t*br=Bpb+(size_t)(w*64+t)*W;
for(int q=0;q<NY;q++)
acc[q]=_mm256_xor_si256(acc[q],_mm256_load_si256((const __m256i*)(br+4*q)));
}
}
uint64_t*o=Pbb+(size_t)i*W;
for(int q=0;q<NY;q++) _mm256_store_si256((__m256i*)(o+4*q),acc[q]);
}
}
/* P[i][j] = XOR over k with A1[i][k] set of Bp[k][j].
Per-output-row bit scan (dense 64-bit words -> predictable), K-BLOCKED so the
Bp slice stays in L2. P accumulated in memory across k-blocks. */
AVX2 static void corr_pass3(int n){
int W=n/64, NY=W>>2;
memset(Pbb,0,(size_t)n*W*8);
const int KB=256;
for (int kb=0; kb<n; kb+=KB){
int w0=kb>>6, w1=w0+(KB>>6);
for (int i=0;i<n;i++){
__m256i acc[16];
for (int q=0;q<NY;q++) acc[q]=_mm256_setzero_si256();
const uint64_t*row=A1b+(size_t)i*W;
for (int w=w0;w<w1;w++){
uint64_t b=row[w];
while(b){
int t=__builtin_ctzll(b); b&=b-1;
const uint64_t*br=Bpb+(size_t)(w*64+t)*W;
for (int q=0;q<NY;q++) acc[q]=_mm256_xor_si256(acc[q],_mm256_load_si256((const __m256i*)(br+4*q)));
}
}
uint64_t*o=Pbb+(size_t)i*W;
for (int q=0;q<NY;q++){
__m256i*po=(__m256i*)(o+4*q);
_mm256_store_si256(po,_mm256_xor_si256(_mm256_load_si256((const __m256i*)po),acc[q]));
}
}
}
}
/* specialized for n=4096 (W=64, NY=16): fully unrolled XOR loops */
AVX2 static void corr16(void){
const int n=4096, W=64;
memset(Pbb,0,(size_t)n*W*8);
const int KB=256;
for (int kb=0; kb<n; kb+=KB){
int w0=kb>>6, w1=w0+4;
for (int i=0;i<n;i++){
__m256i acc0,acc1,acc2,acc3,acc4,acc5,acc6,acc7,acc8,acc9,acc10,acc11,acc12,acc13,acc14,acc15;
acc0=acc1=acc2=acc3=acc4=acc5=acc6=acc7=_mm256_setzero_si256();
acc8=acc9=acc10=acc11=acc12=acc13=acc14=acc15=_mm256_setzero_si256();
const uint64_t*row=A1b+(size_t)i*W;
for (int w=w0;w<w1;w++){
uint64_t b=row[w];
while(b){
int t=__builtin_ctzll(b); b&=b-1;
const uint64_t*br=Bpb+(size_t)(w*64+t)*W;
acc0 =_mm256_xor_si256(acc0 ,_mm256_load_si256((const __m256i*)(br+0)));
acc1 =_mm256_xor_si256(acc1 ,_mm256_load_si256((const __m256i*)(br+4)));
acc2 =_mm256_xor_si256(acc2 ,_mm256_load_si256((const __m256i*)(br+8)));
acc3 =_mm256_xor_si256(acc3 ,_mm256_load_si256((const __m256i*)(br+12)));
acc4 =_mm256_xor_si256(acc4 ,_mm256_load_si256((const __m256i*)(br+16)));
acc5 =_mm256_xor_si256(acc5 ,_mm256_load_si256((const __m256i*)(br+20)));
acc6 =_mm256_xor_si256(acc6 ,_mm256_load_si256((const __m256i*)(br+24)));
acc7 =_mm256_xor_si256(acc7 ,_mm256_load_si256((const __m256i*)(br+28)));
acc8 =_mm256_xor_si256(acc8 ,_mm256_load_si256((const __m256i*)(br+32)));
acc9 =_mm256_xor_si256(acc9 ,_mm256_load_si256((const __m256i*)(br+36)));
acc10=_mm256_xor_si256(acc10,_mm256_load_si256((const __m256i*)(br+40)));
acc11=_mm256_xor_si256(acc11,_mm256_load_si256((const __m256i*)(br+44)));
acc12=_mm256_xor_si256(acc12,_mm256_load_si256((const __m256i*)(br+48)));
acc13=_mm256_xor_si256(acc13,_mm256_load_si256((const __m256i*)(br+52)));
acc14=_mm256_xor_si256(acc14,_mm256_load_si256((const __m256i*)(br+56)));
acc15=_mm256_xor_si256(acc15,_mm256_load_si256((const __m256i*)(br+60)));
}
}
uint64_t*o=Pbb+(size_t)i*W;
#define RMW(q,a) { __m256i*po=(__m256i*)(o+4*(q)); _mm256_store_si256(po,_mm256_xor_si256(_mm256_load_si256((const __m256i*)po),a)); }
RMW(0,acc0) RMW(1,acc1) RMW(2,acc2) RMW(3,acc3) RMW(4,acc4) RMW(5,acc5) RMW(6,acc6) RMW(7,acc7)
RMW(8,acc8) RMW(9,acc9) RMW(10,acc10) RMW(11,acc11) RMW(12,acc12) RMW(13,acc13) RMW(14,acc14) RMW(15,acc15)
#undef RMW
}
}
}
AVX2 static void corr1k(void){
const int n=1024, W=16;
memset(Pbb,0,(size_t)n*W*8);
const int KB=256;
for (int kb=0; kb<n; kb+=KB){
int w0=kb>>6, w1=w0+4;
for (int i=0;i<n;i++){
__m256i acc0,acc1,acc2,acc3;
acc0=acc1=acc2=acc3=_mm256_setzero_si256();
const uint64_t*row=A1b+(size_t)i*W;
for (int w=w0;w<w1;w++){
uint64_t b=row[w];
while(b){
int t=__builtin_ctzll(b); b&=b-1;
const uint64_t*br=Bpb+(size_t)(w*64+t)*W;
acc0 =_mm256_xor_si256(acc0 ,_mm256_load_si256((const __m256i*)(br+0)));
acc1 =_mm256_xor_si256(acc1 ,_mm256_load_si256((const __m256i*)(br+4)));
acc2 =_mm256_xor_si256(acc2 ,_mm256_load_si256((const __m256i*)(br+8)));
acc3 =_mm256_xor_si256(acc3 ,_mm256_load_si256((const __m256i*)(br+12)));
}
}
uint64_t*o=Pbb+(size_t)i*W;
#define RMW(q,a) { __m256i*po=(__m256i*)(o+4*(q)); _mm256_store_si256(po,_mm256_xor_si256(_mm256_load_si256((const __m256i*)po),a)); }
RMW(0,acc0) RMW(1,acc1) RMW(2,acc2) RMW(3,acc3)
#undef RMW
}
}
}
static void expand_bits(int n){
static uint8_t T[256*8]; static int ti=0;
if(!ti){ for(int v=0;v<256;v++) for(int b=0;b<8;b++) T[v*8+b]=((v>>b)&1)?0x80:0x00; ti=1; }
int W=n/64;
for(int i=0;i<n;i++){
const uint64_t*w=Pbb+(size_t)i*W; uint8_t*o=PBb+(size_t)i*n;
for(int q=0;q<W;q++){
uint64_t v=w[q];
for(int s=0;s<8;s++) memcpy(o+(q*8+s)*8, T+(((uint8_t)(v>>(8*s)))*8), 8);
}
}
}
#define NIB 2
#define NJ 4
static const int16_t W0FF[16] __attribute__((aligned(32))) =
{0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,
0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF};
static int8_t BI3b[NMAX*NMAX] __attribute__((aligned(64)));
AVX2 static void build_BI3(int n, const int8_t*B){
int ng=n/64;
for (int p=0;p<n/2;p++){
const int8_t *r0=B+(size_t)(2*p)*n, *r1=r0+n;
for (int g=0;g<ng;g++){
const int8_t *a=r0+64*g, *b=r1+64*g;
int8_t *o=BI3b+((size_t)g*(n/2)+p)*128;
for (int t=0;t<64;t+=32){
__m256i va=_mm256_load_si256((const __m256i*)(a+t));
__m256i vb=_mm256_load_si256((const __m256i*)(b+t));
__m256i lo=_mm256_unpacklo_epi8(va,vb), hi=_mm256_unpackhi_epi8(va,vb);
_mm256_store_si256((__m256i*)(o+2*t), _mm256_permute2x128_si256(lo,hi,0x20));
_mm256_store_si256((__m256i*)(o+2*t+32), _mm256_permute2x128_si256(lo,hi,0x31));
}
}
}
}
static int8_t mscr[256] __attribute__((aligned(64)));
/* Bit-expansion table, ONE 8-byte entry per byte value: SHB8[v*8+m] = bit m of v ? 0x80 : 0.
Compile-time constant (the old runtime `shinit()` cost a flag load + guard + a vectorised
init block in the middle of the hot loop on EVERY epilogue call). */
static const uint8_t SHB8[256*8] __attribute__((aligned(64))) = {
0,0,0,0,0,0,0,0,0x80,0,0,0,0,0,0,0,
0,0x80,0,0,0,0,0,0,0x80,0x80,0,0,0,0,0,0,
0,0,0x80,0,0,0,0,0,0x80,0,0x80,0,0,0,0,0,
0,0x80,0x80,0,0,0,0,0,0x80,0x80,0x80,0,0,0,0,0,
0,0,0,0x80,0,0,0,0,0x80,0,0,0x80,0,0,0,0,
0,0x80,0,0x80,0,0,0,0,0x80,0x80,0,0x80,0,0,0,0,
0,0,0x80,0x80,0,0,0,0,0x80,0,0x80,0x80,0,0,0,0,
0,0x80,0x80,0x80,0,0,0,0,0x80,0x80,0x80,0x80,0,0,0,0,
0,0,0,0,0x80,0,0,0,0x80,0,0,0,0x80,0,0,0,
0,0x80,0,0,0x80,0,0,0,0x80,0x80,0,0,0x80,0,0,0,
0,0,0x80,0,0x80,0,0,0,0x80,0,0x80,0,0x80,0,0,0,
0,0x80,0x80,0,0x80,0,0,0,0x80,0x80,0x80,0,0x80,0,0,0,
0,0,0,0x80,0x80,0,0,0,0x80,0,0,0x80,0x80,0,0,0,
0,0x80,0,0x80,0x80,0,0,0,0x80,0x80,0,0x80,0x80,0,0,0,
0,0,0x80,0x80,0x80,0,0,0,0x80,0,0x80,0x80,0x80,0,0,0,
0,0x80,0x80,0x80,0x80,0,0,0,0x80,0x80,0x80,0x80,0x80,0,0,0,
0,0,0,0,0,0x80,0,0,0x80,0,0,0,0,0x80,0,0,
0,0x80,0,0,0,0x80,0,0,0x80,0x80,0,0,0,0x80,0,0,
0,0,0x80,0,0,0x80,0,0,0x80,0,0x80,0,0,0x80,0,0,
0,0x80,0x80,0,0,0x80,0,0,0x80,0x80,0x80,0,0,0x80,0,0,
0,0,0,0x80,0,0x80,0,0,0x80,0,0,0x80,0,0x80,0,0,
0,0x80,0,0x80,0,0x80,0,0,0x80,0x80,0,0x80,0,0x80,0,0,
0,0,0x80,0x80,0,0x80,0,0,0x80,0,0x80,0x80,0,0x80,0,0,
0,0x80,0x80,0x80,0,0x80,0,0,0x80,0x80,0x80,0x80,0,0x80,0,0,
0,0,0,0,0x80,0x80,0,0,0x80,0,0,0,0x80,0x80,0,0,
0,0x80,0,0,0x80,0x80,0,0,0x80,0x80,0,0,0x80,0x80,0,0,
0,0,0x80,0,0x80,0x80,0,0,0x80,0,0x80,0,0x80,0x80,0,0,
0,0x80,0x80,0,0x80,0x80,0,0,0x80,0x80,0x80,0,0x80,0x80,0,0,
0,0,0,0x80,0x80,0x80,0,0,0x80,0,0,0x80,0x80,0x80,0,0,
0,0x80,0,0x80,0x80,0x80,0,0,0x80,0x80,0,0x80,0x80,0x80,0,0,
0,0,0x80,0x80,0x80,0x80,0,0,0x80,0,0x80,0x80,0x80,0x80,0,0,
0,0x80,0x80,0x80,0x80,0x80,0,0,0x80,0x80,0x80,0x80,0x80,0x80,0,0,
0,0,0,0,0,0,0x80,0,0x80,0,0,0,0,0,0x80,0,
0,0x80,0,0,0,0,0x80,0,0x80,0x80,0,0,0,0,0x80,0,
0,0,0x80,0,0,0,0x80,0,0x80,0,0x80,0,0,0,0x80,0,
0,0x80,0x80,0,0,0,0x80,0,0x80,0x80,0x80,0,0,0,0x80,0,
0,0,0,0x80,0,0,0x80,0,0x80,0,0,0x80,0,0,0x80,0,
0,0x80,0,0x80,0,0,0x80,0,0x80,0x80,0,0x80,0,0,0x80,0,
0,0,0x80,0x80,0,0,0x80,0,0x80,0,0x80,0x80,0,0,0x80,0,
0,0x80,0x80,0x80,0,0,0x80,0,0x80,0x80,0x80,0x80,0,0,0x80,0,
0,0,0,0,0x80,0,0x80,0,0x80,0,0,0,0x80,0,0x80,0,
0,0x80,0,0,0x80,0,0x80,0,0x80,0x80,0,0,0x80,0,0x80,0,
0,0,0x80,0,0x80,0,0x80,0,0x80,0,0x80,0,0x80,0,0x80,0,
0,0x80,0x80,0,0x80,0,0x80,0,0x80,0x80,0x80,0,0x80,0,0x80,0,
0,0,0,0x80,0x80,0,0x80,0,0x80,0,0,0x80,0x80,0,0x80,0,
0,0x80,0,0x80,0x80,0,0x80,0,0x80,0x80,0,0x80,0x80,0,0x80,0,
0,0,0x80,0x80,0x80,0,0x80,0,0x80,0,0x80,0x80,0x80,0,0x80,0,
0,0x80,0x80,0x80,0x80,0,0x80,0,0x80,0x80,0x80,0x80,0x80,0,0x80,0,
0,0,0,0,0,0x80,0x80,0,0x80,0,0,0,0,0x80,0x80,0,
0,0x80,0,0,0,0x80,0x80,0,0x80,0x80,0,0,0,0x80,0x80,0,
0,0,0x80,0,0,0x80,0x80,0,0x80,0,0x80,0,0,0x80,0x80,0,
0,0x80,0x80,0,0,0x80,0x80,0,0x80,0x80,0x80,0,0,0x80,0x80,0,
0,0,0,0x80,0,0x80,0x80,0,0x80,0,0,0x80,0,0x80,0x80,0,
0,0x80,0,0x80,0,0x80,0x80,0,0x80,0x80,0,0x80,0,0x80,0x80,0,
0,0,0x80,0x80,0,0x80,0x80,0,0x80,0,0x80,0x80,0,0x80,0x80,0,
0,0x80,0x80,0x80,0,0x80,0x80,0,0x80,0x80,0x80,0x80,0,0x80,0x80,0,
0,0,0,0,0x80,0x80,0x80,0,0x80,0,0,0,0x80,0x80,0x80,0,
0,0x80,0,0,0x80,0x80,0x80,0,0x80,0x80,0,0,0x80,0x80,0x80,0,
0,0,0x80,0,0x80,0x80,0x80,0,0x80,0,0x80,0,0x80,0x80,0x80,0,
0,0x80,0x80,0,0x80,0x80,0x80,0,0x80,0x80,0x80,0,0x80,0x80,0x80,0,
0,0,0,0x80,0x80,0x80,0x80,0,0x80,0,0,0x80,0x80,0x80,0x80,0,
0,0x80,0,0x80,0x80,0x80,0x80,0,0x80,0x80,0,0x80,0x80,0x80,0x80,0,
0,0,0x80,0x80,0x80,0x80,0x80,0,0x80,0,0x80,0x80,0x80,0x80,0x80,0,
0,0x80,0x80,0x80,0x80,0x80,0x80,0,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0,
0,0,0,0,0,0,0,0x80,0x80,0,0,0,0,0,0,0x80,
0,0x80,0,0,0,0,0,0x80,0x80,0x80,0,0,0,0,0,0x80,
0,0,0x80,0,0,0,0,0x80,0x80,0,0x80,0,0,0,0,0x80,
0,0x80,0x80,0,0,0,0,0x80,0x80,0x80,0x80,0,0,0,0,0x80,
0,0,0,0x80,0,0,0,0x80,0x80,0,0,0x80,0,0,0,0x80,
0,0x80,0,0x80,0,0,0,0x80,0x80,0x80,0,0x80,0,0,0,0x80,
0,0,0x80,0x80,0,0,0,0x80,0x80,0,0x80,0x80,0,0,0,0x80,
0,0x80,0x80,0x80,0,0,0,0x80,0x80,0x80,0x80,0x80,0,0,0,0x80,
0,0,0,0,0x80,0,0,0x80,0x80,0,0,0,0x80,0,0,0x80,
0,0x80,0,0,0x80,0,0,0x80,0x80,0x80,0,0,0x80,0,0,0x80,
0,0,0x80,0,0x80,0,0,0x80,0x80,0,0x80,0,0x80,0,0,0x80,
0,0x80,0x80,0,0x80,0,0,0x80,0x80,0x80,0x80,0,0x80,0,0,0x80,
0,0,0,0x80,0x80,0,0,0x80,0x80,0,0,0x80,0x80,0,0,0x80,
0,0x80,0,0x80,0x80,0,0,0x80,0x80,0x80,0,0x80,0x80,0,0,0x80,
0,0,0x80,0x80,0x80,0,0,0x80,0x80,0,0x80,0x80,0x80,0,0,0x80,
0,0x80,0x80,0x80,0x80,0,0,0x80,0x80,0x80,0x80,0x80,0x80,0,0,0x80,
0,0,0,0,0,0x80,0,0x80,0x80,0,0,0,0,0x80,0,0x80,
0,0x80,0,0,0,0x80,0,0x80,0x80,0x80,0,0,0,0x80,0,0x80,
0,0,0x80,0,0,0x80,0,0x80,0x80,0,0x80,0,0,0x80,0,0x80,
0,0x80,0x80,0,0,0x80,0,0x80,0x80,0x80,0x80,0,0,0x80,0,0x80,
0,0,0,0x80,0,0x80,0,0x80,0x80,0,0,0x80,0,0x80,0,0x80,
0,0x80,0,0x80,0,0x80,0,0x80,0x80,0x80,0,0x80,0,0x80,0,0x80,
0,0,0x80,0x80,0,0x80,0,0x80,0x80,0,0x80,0x80,0,0x80,0,0x80,
0,0x80,0x80,0x80,0,0x80,0,0x80,0x80,0x80,0x80,0x80,0,0x80,0,0x80,
0,0,0,0,0x80,0x80,0,0x80,0x80,0,0,0,0x80,0x80,0,0x80,
0,0x80,0,0,0x80,0x80,0,0x80,0x80,0x80,0,0,0x80,0x80,0,0x80,
0,0,0x80,0,0x80,0x80,0,0x80,0x80,0,0x80,0,0x80,0x80,0,0x80,
0,0x80,0x80,0,0x80,0x80,0,0x80,0x80,0x80,0x80,0,0x80,0x80,0,0x80,
0,0,0,0x80,0x80,0x80,0,0x80,0x80,0,0,0x80,0x80,0x80,0,0x80,
0,0x80,0,0x80,0x80,0x80,0,0x80,0x80,0x80,0,0x80,0x80,0x80,0,0x80,
0,0,0x80,0x80,0x80,0x80,0,0x80,0x80,0,0x80,0x80,0x80,0x80,0,0x80,
0,0x80,0x80,0x80,0x80,0x80,0,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0,0x80,
0,0,0,0,0,0,0x80,0x80,0x80,0,0,0,0,0,0x80,0x80,
0,0x80,0,0,0,0,0x80,0x80,0x80,0x80,0,0,0,0,0x80,0x80,
0,0,0x80,0,0,0,0x80,0x80,0x80,0,0x80,0,0,0,0x80,0x80,
0,0x80,0x80,0,0,0,0x80,0x80,0x80,0x80,0x80,0,0,0,0x80,0x80,
0,0,0,0x80,0,0,0x80,0x80,0x80,0,0,0x80,0,0,0x80,0x80,
0,0x80,0,0x80,0,0,0x80,0x80,0x80,0x80,0,0x80,0,0,0x80,0x80,
0,0,0x80,0x80,0,0,0x80,0x80,0x80,0,0x80,0x80,0,0,0x80,0x80,
0,0x80,0x80,0x80,0,0,0x80,0x80,0x80,0x80,0x80,0x80,0,0,0x80,0x80,
0,0,0,0,0x80,0,0x80,0x80,0x80,0,0,0,0x80,0,0x80,0x80,
0,0x80,0,0,0x80,0,0x80,0x80,0x80,0x80,0,0,0x80,0,0x80,0x80,
0,0,0x80,0,0x80,0,0x80,0x80,0x80,0,0x80,0,0x80,0,0x80,0x80,
0,0x80,0x80,0,0x80,0,0x80,0x80,0x80,0x80,0x80,0,0x80,0,0x80,0x80,
0,0,0,0x80,0x80,0,0x80,0x80,0x80,0,0,0x80,0x80,0,0x80,0x80,
0,0x80,0,0x80,0x80,0,0x80,0x80,0x80,0x80,0,0x80,0x80,0,0x80,0x80,
0,0,0x80,0x80,0x80,0,0x80,0x80,0x80,0,0x80,0x80,0x80,0,0x80,0x80,
0,0x80,0x80,0x80,0x80,0,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0,0x80,0x80,
0,0,0,0,0,0x80,0x80,0x80,0x80,0,0,0,0,0x80,0x80,0x80,
0,0x80,0,0,0,0x80,0x80,0x80,0x80,0x80,0,0,0,0x80,0x80,0x80,
0,0,0x80,0,0,0x80,0x80,0x80,0x80,0,0x80,0,0,0x80,0x80,0x80,
0,0x80,0x80,0,0,0x80,0x80,0x80,0x80,0x80,0x80,0,0,0x80,0x80,0x80,
0,0,0,0x80,0,0x80,0x80,0x80,0x80,0,0,0x80,0,0x80,0x80,0x80,
0,0x80,0,0x80,0,0x80,0x80,0x80,0x80,0x80,0,0x80,0,0x80,0x80,0x80,
0,0,0x80,0x80,0,0x80,0x80,0x80,0x80,0,0x80,0x80,0,0x80,0x80,0x80,
0,0x80,0x80,0x80,0,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0,0x80,0x80,0x80,
0,0,0,0,0x80,0x80,0x80,0x80,0x80,0,0,0,0x80,0x80,0x80,0x80,
0,0x80,0,0,0x80,0x80,0x80,0x80,0x80,0x80,0,0,0x80,0x80,0x80,0x80,
0,0,0x80,0,0x80,0x80,0x80,0x80,0x80,0,0x80,0,0x80,0x80,0x80,0x80,
0,0x80,0x80,0,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0,0x80,0x80,0x80,0x80,
0,0,0,0x80,0x80,0x80,0x80,0x80,0x80,0,0,0x80,0x80,0x80,0x80,0x80,
0,0x80,0,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0,0x80,0x80,0x80,0x80,0x80,
0,0,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0,0x80,0x80,0x80,0x80,0x80,0x80,
0,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80,
};
/* [sm2_] SIMD bit-expansion constants for the epilogue.
EPD_IDX[c] = c div 8 -> after vpshufb over a vpbroadcastd of the four Pbb
bytes, byte c of the spread holds the Pbb byte that owns column c.
EPD_BIT[c] = 1 << (c mod 8) -> isolates the one bit that decides column c.
EPD_80 = 0x80 everywhere -> turns the vpcmpeqb result into the xor mask. */
static const uint8_t EPD_IDX[32] __attribute__((aligned(32))) = {
0,0,0,0,0,0,0,0, 1,1,1,1,1,1,1,1, 2,2,2,2,2,2,2,2, 3,3,3,3,3,3,3,3 };
static const uint8_t EPD_BIT[32] __attribute__((aligned(32))) = {
1,2,4,8,16,32,64,128, 1,2,4,8,16,32,64,128,
1,2,4,8,16,32,64,128, 1,2,4,8,16,32,64,128 };
static const uint8_t EPD_80[32] __attribute__((aligned(32))) = {
0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80, 0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80,
0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80, 0x80,0x80,0x80,0x80,0x80,0x80,0x80,0x80 };
AVX2 static void epilogue(int n, int8_t*C, const uint64_t*P, int jb, int i, int NCH){
const size_t bo=(size_t)jb>>3;
const uint8_t*base=(const uint8_t*)P+(size_t)i*(n>>3);
const __m256i idx=_mm256_load_si256((const __m256i*)EPD_IDX);
const __m256i bit=_mm256_load_si256((const __m256i*)EPD_BIT);
const __m256i m80=_mm256_load_si256((const __m256i*)EPD_80);
for(int r=0;r<2;r++){
const uint8_t*qr=base+(size_t)r*(n>>3);
for(int t=0;t<NCH;t++){
/* chunk (r,t) covers C row i+r, columns jb+32t .. jb+32t+31, and reads
the four Pbb bytes at qr[bo+4t .. bo+4t+3]. Pbb bit position p of the
row is the correction for column p, and p = (bo+4t+s)*8 + b with
column offset c = 8s+b => column c takes byte (c div 8) and bit
(c mod 8) of those four bytes. vpbroadcastd makes byte c equal to
source byte (c mod 4); vpshufb then selects index (c div 8) and,
because the dword is replicated in both 128-bit lanes, lane 1 can
reach indices 2 and 3 of its own copy. Result: byte c of the spread
is the Pbb byte that owns column c, and vpand/vpcmpeqb/vpand turn
bit (c mod 8) of it into an 0x80 xor mask. */
__m256i v=_mm256_load_si256((const __m256i*)(mscr+(r*NCH+t)*32));
int wd; __builtin_memcpy(&wd, qr+bo+4*(size_t)t, 4);
__m256i b=_mm256_broadcastd_epi32(_mm_cvtsi32_si128(wd));
__m256i s=_mm256_shuffle_epi8(b, idx);
__m256i z=_mm256_and_si256(s, bit);
z=_mm256_and_si256(_mm256_cmpeq_epi8(z, bit), m80);
int8_t*c=C+(size_t)(i+r)*n+jb+t*32;
_mm256_storeu_si256((__m256i*)c,_mm256_xor_si256(v,z));
}
}
}
/* ---- C=96 (6 BI ymm) mikro: 2 rows x 96 columns, records of 192 bytes ---- */
static int8_t BI9b[10*(NMAX/2)*192] __attribute__((aligned(64)));
/* padded Bp layout: 65 rows per 64-k group, row 64 = ZERO (safe index for tzcnt(0)=64) */
static uint64_t BP2[16*65*16] __attribute__((aligned(64)));
AVX2 static void build_BP2(int n, const int8_t*B){
const int W=n>>6, STR=65*W;
const __m256i one=_mm256_set1_epi8(1);
for (int k=0;k<n;k++){
const int8_t*row=B+(size_t)k*n; uint64_t*w=BP2+(size_t)(k>>6)*STR+(size_t)(k&63)*W;
for (int j=0;j<n;j+=64){
__m256i a=_mm256_slli_epi16(_mm256_and_si256(_mm256_load_si256((const __m256i*)(row+j)),one),7);
__m256i b=_mm256_slli_epi16(_mm256_and_si256(_mm256_load_si256((const __m256i*)(row+j+32)),one),7);
uint32_t lo=(uint32_t)_mm256_movemask_epi8(a);
uint32_t hi=(uint32_t)_mm256_movemask_epi8(b);
w[j>>6]=(uint64_t)lo|((uint64_t)hi<<32);
}
}
for (int g=0; g<(n>>6); g++) memset(BP2+(size_t)g*STR+64*W, 0, W*8);
}
/* ---- 4-word merged bit scan -------------------------------------------
#100989/#101030 scan ONE pair of 64-bit words per `while (b0|b1)` loop; the
data-dependent trip count (max(popcount(b0),popcount(b1)) ~ 33 of 64) is what
the loop really pays for. Merging the kb block's four words into ONE loop
touches exactly the same bytes with exactly the same loads/XORs, but issues
half as many while-instances (2 -> 1 per (i,kb)), hence half the
data-dependent exits and half the per-instance fixup (row[w] loads, bw
pointers, for-w bookkeeping). The tzcnt(0)=64 sentinel row keeps the shorter
words harmless, as before. Judge probe: -7.08% of the phase, Pbb
byte-identical. */
AVX2 static void corr2k(int n){
const int W=n>>6, STR=65*W;
memset(Pbb,0,(size_t)n*W*8);
for (int kb0=0; kb0<n; kb0+=256){
int w0=kb0>>6;
for (int i=0;i<n;i++){
__m256i a0,a1,a2,a3; a0=a1=a2=a3=_mm256_setzero_si256();
const uint64_t*row=A1b+(size_t)i*W;
uint64_t b0=row[w0],b1=row[w0+1],b2=row[w0+2],b3=row[w0+3];
const uint64_t*bw0=BP2+(size_t)w0*STR;
const uint64_t*bw1=bw0+STR,*bw2=bw1+STR,*bw3=bw2+STR;
while (b0|b1|b2|b3){
long t0,t1,t2,t3; uint64_t nb;
__asm__("tzcntq %2,%0\n\t" "blsrq %2,%1":"=r"(t0),"=r"(nb):"r"(b0)); b0=nb;
__asm__("tzcntq %2,%0\n\t" "blsrq %2,%1":"=r"(t1),"=r"(nb):"r"(b1)); b1=nb;
__asm__("tzcntq %2,%0\n\t" "blsrq %2,%1":"=r"(t2),"=r"(nb):"r"(b2)); b2=nb;
__asm__("tzcntq %2,%0\n\t" "blsrq %2,%1":"=r"(t3),"=r"(nb):"r"(b3)); b3=nb;
const uint64_t*r0=bw0+((size_t)t0<<4);
const uint64_t*r1=bw1+((size_t)t1<<4);
const uint64_t*r2=bw2+((size_t)t2<<4);
const uint64_t*r3=bw3+((size_t)t3<<4);
a0=_mm256_xor_si256(a0,_mm256_load_si256((const __m256i*)(r0+0)));
a1=_mm256_xor_si256(a1,_mm256_load_si256((const __m256i*)(r0+4)));
a2=_mm256_xor_si256(a2,_mm256_load_si256((const __m256i*)(r0+8)));
a3=_mm256_xor_si256(a3,_mm256_load_si256((const __m256i*)(r0+12)));
a0=_mm256_xor_si256(a0,_mm256_load_si256((const __m256i*)(r1+0)));
a1=_mm256_xor_si256(a1,_mm256_load_si256((const __m256i*)(r1+4)));
a2=_mm256_xor_si256(a2,_mm256_load_si256((const __m256i*)(r1+8)));
a3=_mm256_xor_si256(a3,_mm256_load_si256((const __m256i*)(r1+12)));
a0=_mm256_xor_si256(a0,_mm256_load_si256((const __m256i*)(r2+0)));
a1=_mm256_xor_si256(a1,_mm256_load_si256((const __m256i*)(r2+4)));
a2=_mm256_xor_si256(a2,_mm256_load_si256((const __m256i*)(r2+8)));
a3=_mm256_xor_si256(a3,_mm256_load_si256((const __m256i*)(r2+12)));
a0=_mm256_xor_si256(a0,_mm256_load_si256((const __m256i*)(r3+0)));
a1=_mm256_xor_si256(a1,_mm256_load_si256((const __m256i*)(r3+4)));
a2=_mm256_xor_si256(a2,_mm256_load_si256((const __m256i*)(r3+8)));
a3=_mm256_xor_si256(a3,_mm256_load_si256((const __m256i*)(r3+12)));
}
uint64_t*o=Pbb+(size_t)i*W;
{ __m256i*po=(__m256i*)(o+0); _mm256_store_si256(po,_mm256_xor_si256(_mm256_load_si256((const __m256i*)po),a0)); }
{ __m256i*po=(__m256i*)(o+4); _mm256_store_si256(po,_mm256_xor_si256(_mm256_load_si256((const __m256i*)po),a1)); }
{ __m256i*po=(__m256i*)(o+8); _mm256_store_si256(po,_mm256_xor_si256(_mm256_load_si256((const __m256i*)po),a2)); }
{ __m256i*po=(__m256i*)(o+12); _mm256_store_si256(po,_mm256_xor_si256(_mm256_load_si256((const __m256i*)po),a3)); }
}
}
}
/* p OUTER / g inner: the whole of B is read ONCE, and for a fixed p the ten
96-byte column slices r0+96*g, g=0..9, cover bytes 0..959 of the row in order.
(#100989 had g outer / p inner => every g walked all of B => B read 10x.) */
AVX2 static void build_BI9(int n, const int8_t*B){
const int half=n>>1, NC=96, NREC=192;
for (int p=0;p<half;p++){
const int8_t*r0b=B+(size_t)(2*p)*n, *r1b=r0b+n;
for (int g=0; g<10; g++){
const int8_t*r0=r0b+NC*g, *r1=r1b+NC*g;
int8_t*o=BI9b+((size_t)g*half+p)*NREC;
for (int t=0;t<NC;t+=32){
__m256i va=_mm256_loadu_si256((const __m256i*)(r0+t));
__m256i vb=_mm256_loadu_si256((const __m256i*)(r1+t));
__m256i lo=_mm256_unpacklo_epi8(va,vb), hi=_mm256_unpackhi_epi8(va,vb);
_mm256_store_si256((__m256i*)(o+2*t), _mm256_permute2x128_si256(lo,hi,0x20));
_mm256_store_si256((__m256i*)(o+2*t+32), _mm256_permute2x128_si256(lo,hi,0x31));
}
}
}
}
/* fused A0 + A1bits: one pass over A instead of two (A = 1 MB > judge L2 = 256 KB) */
AVX2 static void build_A0_A1(int n, const int8_t*A){
const __m256i m=_mm256_set1_epi8(0x7F);
int W=n/64;
for (int i=0;i<n;i++){
const int8_t*row=A+(size_t)i*n; uint64_t*w=A1b+(size_t)i*W;
int8_t*o=A0b+2*(size_t)i*n;
for (int k=0;k<n;k+=64){
__m256i v0=_mm256_load_si256((const __m256i*)(row+k));
__m256i v1=_mm256_load_si256((const __m256i*)(row+k+32));
{ __m256i e0=_mm256_and_si256(v0,m), e1=_mm256_and_si256(v1,m);
__m256i l0=_mm256_unpacklo_epi16(e0,e0), h0=_mm256_unpackhi_epi16(e0,e0);
__m256i l1=_mm256_unpacklo_epi16(e1,e1), h1=_mm256_unpackhi_epi16(e1,e1);
_mm256_store_si256((__m256i*)(o+2*k), _mm256_permute2x128_si256(l0,h0,0x20));
_mm256_store_si256((__m256i*)(o+2*k+32), _mm256_permute2x128_si256(l0,h0,0x31));
_mm256_store_si256((__m256i*)(o+2*k+64), _mm256_permute2x128_si256(l1,h1,0x20));
_mm256_store_si256((__m256i*)(o+2*k+96), _mm256_permute2x128_si256(l1,h1,0x31)); }
uint32_t lo=(uint32_t)_mm256_movemask_epi8(v0);
uint32_t hi=(uint32_t)_mm256_movemask_epi8(v1);
w[k>>6]=(uint64_t)lo|((uint64_t)hi<<32);
}
}
}
/* tail block: 64 columns (960..1023) into BI3b with the original 128B records */
AVX2 static void build_BI3_tail(int n, const int8_t*B){
const int half=n>>1, c0=960;
for (int p=0;p<half;p++){
const int8_t*r0=B+(size_t)(2*p)*n+c0, *r1=r0+n;
int8_t*o=BI3b+((size_t)p)*128;
for (int t=0;t<64;t+=32){
__m256i va=_mm256_loadu_si256((const __m256i*)(r0+t));
__m256i vb=_mm256_loadu_si256((const __m256i*)(r1+t));
__m256i lo=_mm256_unpacklo_epi8(va,vb), hi=_mm256_unpackhi_epi8(va,vb);
_mm256_store_si256((__m256i*)(o+2*t), _mm256_permute2x128_si256(lo,hi,0x20));
_mm256_store_si256((__m256i*)(o+2*t+32), _mm256_permute2x128_si256(lo,hi,0x31));
}
}
}
#define W0FFS {0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF,0x00FF}
static const int16_t M9[16] __attribute__((aligned(32))) = W0FFS;
AVX2 static void mk_s6(int n, const int8_t*A0, const int8_t*BI, uint8_t*scr){
const int8_t *x0=A0, *x1=A0+n, *bi=BI; long cnt=n>>3;
asm volatile(
"vpxor %%xmm2,%%xmm2,%%xmm2\n\t"
"vpxor %%xmm3,%%xmm3,%%xmm3\n\t"
"vpxor %%xmm6,%%xmm6,%%xmm6\n\t"
"vpxor %%xmm7,%%xmm7,%%xmm7\n\t"
"vpxor %%xmm8,%%xmm8,%%xmm8\n\t"
"vpxor %%xmm9,%%xmm9,%%xmm9\n\t"
"vpxor %%xmm10,%%xmm10,%%xmm10\n\t"
"vpxor %%xmm11,%%xmm11,%%xmm11\n\t"
"vpxor %%xmm12,%%xmm12,%%xmm12\n\t"
"vpxor %%xmm13,%%xmm13,%%xmm13\n\t"
"vpxor %%xmm14,%%xmm14,%%xmm14\n\t"
"vpxor %%xmm15,%%xmm15,%%xmm15\n\t"
".p2align 5\n" "1:\n"
"prefetcht0 2304(%[bi])\n\t"
"vpbroadcastw (%[x0]), %%ymm0\n\t"
"vpbroadcastw (%[x1]), %%ymm1\n\t"
"vmovdqa (%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm14, %%ymm14\n\t"
"vmovdqa 32(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm15, %%ymm15\n\t"
"vmovdqa 64(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm2, %%ymm2\n\t"
"vmovdqa 96(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm3, %%ymm3\n\t"
"vmovdqa 128(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm12, %%ymm12\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm6, %%ymm6\n\t"
"vmovdqa 160(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm13, %%ymm13\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm7, %%ymm7\n\t"
"vpbroadcastw 2(%[x0]), %%ymm0\n\t"
"vpbroadcastw 2(%[x1]), %%ymm1\n\t"
"vmovdqa 192(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm14, %%ymm14\n\t"
"vmovdqa 224(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm15, %%ymm15\n\t"
"vmovdqa 256(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm2, %%ymm2\n\t"
"vmovdqa 288(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm3, %%ymm3\n\t"
"vmovdqa 320(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm12, %%ymm12\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm6, %%ymm6\n\t"
"vmovdqa 352(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm13, %%ymm13\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm7, %%ymm7\n\t"
"vpbroadcastw 4(%[x0]), %%ymm0\n\t"
"vpbroadcastw 4(%[x1]), %%ymm1\n\t"
"vmovdqa 384(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm14, %%ymm14\n\t"
"vmovdqa 416(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm15, %%ymm15\n\t"
"vmovdqa 448(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm2, %%ymm2\n\t"
"vmovdqa 480(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm3, %%ymm3\n\t"
"vmovdqa 512(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm12, %%ymm12\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm6, %%ymm6\n\t"
"vmovdqa 544(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm13, %%ymm13\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm7, %%ymm7\n\t"
"vpbroadcastw 6(%[x0]), %%ymm0\n\t"
"vpbroadcastw 6(%[x1]), %%ymm1\n\t"
"vmovdqa 576(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm14, %%ymm14\n\t"
"vmovdqa 608(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm15, %%ymm15\n\t"
"vmovdqa 640(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm2, %%ymm2\n\t"
"vmovdqa 672(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm3, %%ymm3\n\t"
"vmovdqa 704(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm12, %%ymm12\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm6, %%ymm6\n\t"
"vmovdqa 736(%[bi]), %%ymm4\n\t"
"vpmaddubsw %%ymm4, %%ymm0, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm13, %%ymm13\n\t"
"vpmaddubsw %%ymm4, %%ymm1, %%ymm5\n\t"
"vpaddw %%ymm5, %%ymm7, %%ymm7\n\t"
"prefetcht0 800(%[bi])\n\t"
"prefetcht0 1600(%[bi])\n\t"
"addq $8, %[x0]\n\t"
"addq $8, %[x1]\n\t"
"addq $768, %[bi]\n\t"
"decq %[cnt]\n\t"
"jnz 1b\n\t"
"vpand %[mask], %%ymm8, %%ymm0\n\t"
"vpand %[mask], %%ymm9, %%ymm1\n\t"
"vpackuswb %%ymm1, %%ymm0, %%ymm0\n\t"
"vpermq $0xD8, %%ymm0, %%ymm0\n\t"
"vmovdqu %%ymm0, 0(%[scr])\n\t"
"vpand %[mask], %%ymm10, %%ymm0\n\t"
"vpand %[mask], %%ymm11, %%ymm1\n\t"
"vpackuswb %%ymm1, %%ymm0, %%ymm0\n\t"
"vpermq $0xD8, %%ymm0, %%ymm0\n\t"
"vmovdqu %%ymm0, 32(%[scr])\n\t"
"vpand %[mask], %%ymm12, %%ymm0\n\t"
"vpand %[mask], %%ymm13, %%ymm1\n\t"
"vpackuswb %%ymm1, %%ymm0, %%ymm0\n\t"
"vpermq $0xD8, %%ymm0, %%ymm0\n\t"
"vmovdqu %%ymm0, 64(%[scr])\n\t"
"vpand %[mask], %%ymm14, %%ymm0\n\t"
"vpand %[mask], %%ymm15, %%ymm1\n\t"
"vpackuswb %%ymm1, %%ymm0, %%ymm0\n\t"
"vpermq $0xD8, %%ymm0, %%ymm0\n\t"
"vmovdqu %%ymm0, 96(%[scr])\n\t"
"vpand %[mask], %%ymm2, %%ymm0\n\t"
"vpand %[mask], %%ymm3, %%ymm1\n\t"
"vpackuswb %%ymm1, %%ymm0, %%ymm0\n\t"
"vpermq $0xD8, %%ymm0, %%ymm0\n\t"
"vmovdqu %%ymm0, 128(%[scr])\n\t"
"vpand %[mask], %%ymm6, %%ymm0\n\t"
"vpand %[mask], %%ymm7, %%ymm1\n\t"
"vpackuswb %%ymm1, %%ymm0, %%ymm0\n\t"
"vpermq $0xD8, %%ymm0, %%ymm0\n\t"
"vmovdqu %%ymm0, 160(%[scr])\n\t"
: [x0]"+r"(x0),[x1]"+r"(x1),[bi]"+r"(bi),[cnt]"+r"(cnt)
: [mask]"m"(M9), [scr]"r"(scr)
: "ymm0","ymm1","ymm2","ymm3","ymm4","ymm5","ymm6","ymm7","ymm8","ymm9","ymm10","ymm11","ymm12","ymm13","ymm14","ymm15","cc","memory");
}
AVX2 static void mk_s(int n, const int8_t*A0, const int8_t*BI, int8_t*C, const uint8_t*PB, int jb){ (void)PB;
const int8_t *x0=A0, *x1=A0+2*n, *bi=BI; long cnt=n>>4;
asm volatile(
"vpxor %%xmm8,%%xmm8,%%xmm8\n\t" "vpxor %%xmm9,%%xmm9,%%xmm9\n\t" "vpxor %%xmm10,%%xmm10,%%xmm10\n\t" "vpxor %%xmm11,%%xmm11,%%xmm11\n\t" "vpxor %%xmm12,%%xmm12,%%xmm12\n\t" "vpxor %%xmm13,%%xmm13,%%xmm13\n\t" "vpxor %%xmm14,%%xmm14,%%xmm14\n\t" "vpxor %%xmm15,%%xmm15,%%xmm15\n\t"
".p2align 5\n\t" "1:\n\t"
"prefetcht0 1536(%[bi])\n\t"
"vpbroadcastd (%[x0]), %%ymm0\n\t"
"vpbroadcastd (%[x1]), %%ymm1\n\t"
"vpmaddubsw 0(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 0(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm12, %%ymm12\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw 32(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 32(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm13, %%ymm13\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw 64(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 64(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm14, %%ymm14\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw 96(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 96(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm15, %%ymm15\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpbroadcastd 4(%[x0]), %%ymm0\n\t"
"vpbroadcastd 4(%[x1]), %%ymm1\n\t"
"vpmaddubsw 128(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 128(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm12, %%ymm12\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw 160(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 160(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm13, %%ymm13\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw 192(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 192(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm14, %%ymm14\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw 224(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 224(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm15, %%ymm15\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpbroadcastd 8(%[x0]), %%ymm0\n\t"
"vpbroadcastd 8(%[x1]), %%ymm1\n\t"
"vpmaddubsw 256(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 256(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm12, %%ymm12\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw 288(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 288(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm13, %%ymm13\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw 320(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 320(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm14, %%ymm14\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw 352(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 352(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm15, %%ymm15\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpbroadcastd 12(%[x0]), %%ymm0\n\t"
"vpbroadcastd 12(%[x1]), %%ymm1\n\t"
"vpmaddubsw 384(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 384(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm12, %%ymm12\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw 416(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 416(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm13, %%ymm13\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw 448(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 448(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm14, %%ymm14\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw 480(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 480(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm15, %%ymm15\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpbroadcastd 16(%[x0]), %%ymm0\n\t"
"vpbroadcastd 16(%[x1]), %%ymm1\n\t"
"vpmaddubsw 512(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 512(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm12, %%ymm12\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw 544(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 544(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm13, %%ymm13\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw 576(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 576(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm14, %%ymm14\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw 608(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 608(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm15, %%ymm15\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpbroadcastd 20(%[x0]), %%ymm0\n\t"
"vpbroadcastd 20(%[x1]), %%ymm1\n\t"
"vpmaddubsw 640(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 640(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm12, %%ymm12\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw 672(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 672(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm13, %%ymm13\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw 704(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 704(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm14, %%ymm14\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw 736(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 736(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm15, %%ymm15\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpbroadcastd 24(%[x0]), %%ymm0\n\t"
"vpbroadcastd 24(%[x1]), %%ymm1\n\t"
"vpmaddubsw 768(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 768(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm12, %%ymm12\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw 800(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 800(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm13, %%ymm13\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw 832(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 832(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm14, %%ymm14\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw 864(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 864(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm15, %%ymm15\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"vpbroadcastd 28(%[x0]), %%ymm0\n\t"
"vpbroadcastd 28(%[x1]), %%ymm1\n\t"
"vpmaddubsw 896(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 896(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm12, %%ymm12\n\t"
"vpaddw %%ymm5, %%ymm8, %%ymm8\n\t"
"vpmaddubsw 928(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 928(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm13, %%ymm13\n\t"
"vpaddw %%ymm5, %%ymm9, %%ymm9\n\t"
"vpmaddubsw 960(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 960(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm14, %%ymm14\n\t"
"vpaddw %%ymm5, %%ymm10, %%ymm10\n\t"
"vpmaddubsw 992(%[bi]), %%ymm0, %%ymm5\n\t"
"vpmaddubsw 992(%[bi]), %%ymm1, %%ymm6\n\t"
"vpaddw %%ymm6, %%ymm15, %%ymm15\n\t"
"vpaddw %%ymm5, %%ymm11, %%ymm11\n\t"
"addq $32, %[x0]\n\t"
"addq $32, %[x1]\n\t"
"addq $1024, %[bi]\n\t"
"decq %[cnt]\n\t"
"jnz 1b\n\t"
"vpand %[mask], %%ymm8, %%ymm2\n\t"
"vpand %[mask], %%ymm9, %%ymm3\n\t"
"vpackuswb %%ymm3, %%ymm2, %%ymm2\n\t"
"vpermq $0xD8, %%ymm2, %%ymm2\n\t"
"vmovdqu %%ymm2, 0(%[scr])\n\t"
"vpand %[mask], %%ymm10, %%ymm2\n\t"
"vpand %[mask], %%ymm11, %%ymm3\n\t"
"vpackuswb %%ymm3, %%ymm2, %%ymm2\n\t"
"vpermq $0xD8, %%ymm2, %%ymm2\n\t"
"vmovdqu %%ymm2, 32(%[scr])\n\t"
"vpand %[mask], %%ymm12, %%ymm2\n\t"
"vpand %[mask], %%ymm13, %%ymm3\n\t"
"vpackuswb %%ymm3, %%ymm2, %%ymm2\n\t"
"vpermq $0xD8, %%ymm2, %%ymm2\n\t"
"vmovdqu %%ymm2, 64(%[scr])\n\t"
"vpand %[mask], %%ymm14, %%ymm2\n\t"
"vpand %[mask], %%ymm15, %%ymm3\n\t"
"vpackuswb %%ymm3, %%ymm2, %%ymm2\n\t"
"vpermq $0xD8, %%ymm2, %%ymm2\n\t"
"vmovdqu %%ymm2, 96(%[scr])\n\t"
: [x0]"+r"(x0),[x1]"+r"(x1),[bi]"+r"(bi),[cnt]"+r"(cnt)
: [mask]"m"(W0FF), [scr]"r"(mscr)
: "ymm0","ymm1","ymm2","ymm3","ymm4","ymm5","ymm6","ymm7","ymm8","ymm9","ymm10","ymm11","ymm12","ymm13","ymm14","ymm15","cc","memory");
}
static uint8_t Gtab[8*256*(NMAX/8)] __attribute__((aligned(64)));
/* single-instruction 32B memory ops: gcc 9.3 -O2 otherwise splits
_mm256_loadu_si256 into vmovdqu xmm + vinserti128 (extra uop on p5). */
#define LDU32(p) ({ __m256i _v; __asm__("vmovdqu %1, %0" : "=x"(_v) : "m"(*(const __m256i*)(p))); _v; })
#define XOR32(dst,p) __asm__("vpxor %1, %0, %0" : "+x"(dst) : "m"(*(const __m256i*)(p)))
#define STU32(p,v) __asm__ volatile("vmovdqu %1, %0" : "+m"(*(__m256i*)(p)) : "x"(v))
template<int NY>
__attribute__((target("avx2"))) static void corr_rows(int n,int B,int W,int gb0,const uint8_t*A1){
for(int i=0;i<n;i++){
const uint8_t*a=A1+(size_t)i*B;
uint64_t*po=Pbb+(size_t)i*W;
__m256i acc[NY];
for(int q=0;q<NY;q++) acc[q]=_mm256_setzero_si256();
for(int g=0;g<4;g++){
const uint8_t*T=Gtab+((size_t)g*256+(size_t)a[gb0+g])*B;
for(int q=0;q<NY;q++) XOR32(acc[q],T+(size_t)q*32);
}
for(int q=0;q<NY;q++){
__m256i o=LDU32(po+4*(size_t)q);
STU32(po+4*(size_t)q, _mm256_xor_si256(o,acc[q]));
}
}
}
__attribute__((target("avx2"))) static void corr_rus(int n){
const int W=n/64, B=n/8, NY=B/32;
if(NY!=1 && NY!=2 && NY!=4 && NY!=8 && NY!=16){ corr_pass3(n); return; }
const uint8_t*A1=(const uint8_t*)A1b;
memset(Pbb,0,(size_t)n*W*8);
const int NB=4; /* 8 k-groups of 8 = 64 k per block */
for(int kb=0;kb<n;kb+=NB*8){
for(int g=0;g<NB;g++){
const int k0=kb+g*8;
uint8_t*T0=Gtab+(size_t)g*256*B;
memset(T0,0,B);
for(int v=1;v<256;v++){
uint8_t*dst=T0+(size_t)v*B;
const uint8_t*pr=T0+(size_t)(v&(v-1))*B;
const int k=k0+__builtin_ctz(v);
const uint8_t*bp=(const uint8_t*)(BP2+(size_t)(k>>6)*(65*W)+(size_t)(k&63)*W);
for(int q=0;q<B;q+=32){
__m256i x=LDU32(pr+q);
XOR32(x,bp+q);
STU32(dst+q,x);
}
}
}
const int gb0=kb>>3;
if(NY==16) corr_rows<16>(n,B,W,gb0,A1);
else if(NY==8)corr_rows<8>(n,B,W,gb0,A1);
else if(NY==4)corr_rows<4>(n,B,W,gb0,A1);
else if(NY==2)corr_rows<2>(n,B,W,gb0,A1);
else corr_rows<1>(n,B,W,gb0,A1);
}
}
/* ================= mc4_: ONE STRASSEN LEVEL (n=1024 -> 7 x 512) ======================
Exact over Z/256: every coefficient is +-1, no division, so the CLRS identity holds in
any commutative ring and int8 arithmetic mod 256 is exactly that ring.
The leaf (m x m = 512 x 512) runs the SHIPPED generic 64-column path
(`mk_s` + `epilogue(...,2)`); its panels are built by a FUSED packer that reads the two
parent blocks and emits A0b/A1b (resp. BI3b/Bpb) in ONE traversal -- nothing is
materialised (BRIEF 2.18.590, DIF). `nn` is the parent row stride (1024 here). */
struct SBLK { const int8_t *X, *Y; int r0, c0, r1, c1, op; };
AVX2 static inline __m256i sblk_cmb(__m256i v, const int8_t *ya, int op){
if (!op) return v;
__m256i u = _mm256_load_si256((const __m256i*)ya);
return (op==1) ? _mm256_add_epi8(v,u) : _mm256_sub_epi8(v,u);
}
/* bld_A = blk + build_A0_A1 fused: one pass emits A0b (=v&0x7F) and A1b (=sign). */
AVX2 static void bld_A(int m, int nn, const SBLK*b){
const __m256i msk=_mm256_set1_epi8(0x7F);
int W=m/64;
for (int i=0;i<m;i++){
const int8_t*xa=b->X+(size_t)(b->r0+i)*nn+b->c0;
const int8_t*ya=b->Y+(size_t)(b->r1+i)*nn+b->c1;
int8_t*dst=A0b+(size_t)i*m;
uint64_t*w1=A1b+(size_t)i*W;
for (int k=0;k<m;k+=64){
__m256i v0=_mm256_load_si256((const __m256i*)(xa+k));
__m256i v1=_mm256_load_si256((const __m256i*)(xa+k+32));
v0=sblk_cmb(v0,ya+k, b->op);
v1=sblk_cmb(v1,ya+k+32,b->op);
_mm256_store_si256((__m256i*)(dst+k), _mm256_and_si256(v0,msk));
_mm256_store_si256((__m256i*)(dst+k+32), _mm256_and_si256(v1,msk));
w1[k>>6]=(uint64_t)(uint32_t)_mm256_movemask_epi8(v0)
|((uint64_t)(uint32_t)_mm256_movemask_epi8(v1)<<32);
}
}
}
/* bld_B = blk + build_BI3 + build_Bpbits fused: one pass emits BI3b and Bpb(=B&1). */
AVX2 static void bld_B(int m, int nn, const SBLK*b){
const __m256i one=_mm256_set1_epi8(1);
int ng=m/64, W=m/64;
for (int p=0;p<m/2;p++){
const int8_t *ax=b->X+(size_t)(b->r0+2*p)*nn+b->c0;
const int8_t *bx=ax+nn;
const int8_t *ay=b->Y+(size_t)(b->r1+2*p)*nn+b->c1;
const int8_t *by=ay+nn;
/* corr_rus reads the PADDED panel BP2, not Bpb: group stride 65*W with a zeroed
sentinel row at index 64 (see build_BP2). Rows 2p and 2p+1 can straddle a group
boundary, so each row's offset is computed independently. */
uint64_t *w0p=BP2+(size_t)((2*p)>>6)*(size_t)(65*W)+(size_t)((2*p)&63)*W;
uint64_t *w1p=BP2+(size_t)((2*p+1)>>6)*(size_t)(65*W)+(size_t)((2*p+1)&63)*W;
for (int g=0;g<ng;g++){
const int8_t *a=ax+64*g, *c=bx+64*g;
const int8_t *d=ay+64*g, *e=by+64*g;
int8_t *o=BI3b+((size_t)g*(m/2)+p)*128;
for (int t=0;t<64;t+=32){
__m256i va=sblk_cmb(_mm256_load_si256((const __m256i*)(a+t)), d+t, b->op);
__m256i vb=sblk_cmb(_mm256_load_si256((const __m256i*)(c+t)), e+t, b->op);
__m256i lo=_mm256_unpacklo_epi8(va,vb), hi=_mm256_unpackhi_epi8(va,vb);
_mm256_store_si256((__m256i*)(o+2*t), _mm256_permute2x128_si256(lo,hi,0x20));
_mm256_store_si256((__m256i*)(o+2*t+32), _mm256_permute2x128_si256(lo,hi,0x31));
uint32_t m0=(uint32_t)_mm256_movemask_epi8(_mm256_slli_epi16(_mm256_and_si256(va,one),7));
uint32_t m1=(uint32_t)_mm256_movemask_epi8(_mm256_slli_epi16(_mm256_and_si256(vb,one),7));
if(t==0) { w0p[g]=(uint64_t)m0; w1p[g]=(uint64_t)m1; }
else { w0p[g]|=(uint64_t)m0<<32; w1p[g]|=(uint64_t)m1<<32; }
}
}
}
}
static const int SAB_T[7][5] = {
{0,0, 1,1, 1}, /* A11+A22 */
{1,0, 1,1, 1}, /* A21+A22 */
{0,0, 0,0, 0}, /* A11 */
{1,1, 1,1, 0}, /* A22 */
{0,0, 0,1, 1}, /* A11+A12 */
{1,0, 0,0, 2}, /* A21-A11 */
{0,1, 1,1, 2}, /* A12-A22 */
};
static const int SBB_T[7][5] = {
{0,0, 1,1, 1}, /* B11+B22 */
{0,0, 0,0, 0}, /* B11 */
{0,1, 1,1, 2}, /* B12-B22 */
{1,0, 0,0, 2}, /* B21-B11 */
{1,1, 1,1, 0}, /* B22 */
{0,0, 0,1, 1}, /* B11+B12 */
{1,0, 1,1, 1}, /* B21+B22 */
};
static int8_t M2[7][256*256+384] __attribute__((aligned(64)));
/* Cblk += sum_k s_k * M_k (up to 4 terms; NULL skips; dst stride n, block m x m) */
AVX2 static void cmb(int8_t*C,int n,int r0,int c0,int m,
const int8_t*M0,int s0,const int8_t*M1,int s1,
const int8_t*M2,int s2,const int8_t*M3,int s3){
for (int i=0;i<m;i++){
int8_t *d=C+(size_t)(r0+i)*n+c0;
const int8_t *p0=M0?M0+(size_t)i*m:0, *p1=M1?M1+(size_t)i*m:0;
const int8_t *p2=M2?M2+(size_t)i*m:0, *p3=M3?M3+(size_t)i*m:0;
for (int j=0;j<m;j+=32){
__m256i v;
if (p0){ __m256i t=_mm256_load_si256((const __m256i*)(p0+j));
v = s0>0?t:_mm256_sub_epi8(_mm256_setzero_si256(),t); }
else v = _mm256_setzero_si256();
if (p1){ __m256i t=_mm256_load_si256((const __m256i*)(p1+j));
v = _mm256_add_epi8(v, s1>0?t:_mm256_sub_epi8(_mm256_setzero_si256(),t)); }
if (p2){ __m256i t=_mm256_load_si256((const __m256i*)(p2+j));
v = _mm256_add_epi8(v, s2>0?t:_mm256_sub_epi8(_mm256_setzero_si256(),t)); }
if (p3){ __m256i t=_mm256_load_si256((const __m256i*)(p3+j));
v = _mm256_add_epi8(v, s3>0?t:_mm256_sub_epi8(_mm256_setzero_si256(),t)); }
_mm256_store_si256((__m256i*)(d+j), v);
}
}
}
/* [sm2_nt] streaming-store clone of cmb, for the final C write only. */
AVX2 static void cmb_nt(int8_t*C,int n,int r0,int c0,int m,
const int8_t*M0,int s0,const int8_t*M1,int s1,
const int8_t*M2,int s2,const int8_t*M3,int s3){
for (int i=0;i<m;i++){
int8_t *d=C+(size_t)(r0+i)*n+c0;
const int8_t *p0=M0?M0+(size_t)i*m:0, *p1=M1?M1+(size_t)i*m:0;
const int8_t *p2=M2?M2+(size_t)i*m:0, *p3=M3?M3+(size_t)i*m:0;
for (int j=0;j<m;j+=32){
__m256i v;
if (p0){ __m256i t=_mm256_load_si256((const __m256i*)(p0+j));
v = s0>0?t:_mm256_sub_epi8(_mm256_setzero_si256(),t); }
else v = _mm256_setzero_si256();
if (p1){ __m256i t=_mm256_load_si256((const __m256i*)(p1+j));
v = _mm256_add_epi8(v, s1>0?t:_mm256_sub_epi8(_mm256_setzero_si256(),t)); }
if (p2){ __m256i t=_mm256_load_si256((const __m256i*)(p2+j));
v = _mm256_add_epi8(v, s2>0?t:_mm256_sub_epi8(_mm256_setzero_si256(),t)); }
if (p3){ __m256i t=_mm256_load_si256((const __m256i*)(p3+j));
v = _mm256_add_epi8(v, s3>0?t:_mm256_sub_epi8(_mm256_setzero_si256(),t)); }
_mm256_stream_si256((__m256i*)(d+j), v);
}
}
}
/* [mc2g_] merged four-bit correction: all four kb tables live at once.
Row loop runs ONCE over the 256 rows; the four accumulators are XORed in a
register and stored once (no memset, no read-modify-write). XOR is
associative+commutative and Pbb is fully rewritten (n*W*8 bytes) every call,
so every byte is identical to the original three-step form. */
static __m256i XT4_256[4][16][16] __attribute__((aligned(64)));
alignas(64) static __m256i TB[656];
AVX2 static void corr4_256(){
constexpr int n=256,W=4;
for(int wb=0;wb<4;wb++){
for(int g=0;g<13;g++){
__m256i*T=(__m256i*)TB+g*32;
const __m256i r0=_mm256_load_si256((const __m256i*)(BP2+((size_t)wb*65+g*5+0)*W));
const __m256i r1=_mm256_load_si256((const __m256i*)(BP2+((size_t)wb*65+g*5+1)*W));
const __m256i r2=_mm256_load_si256((const __m256i*)(BP2+((size_t)wb*65+g*5+2)*W));
const __m256i r3=_mm256_load_si256((const __m256i*)(BP2+((size_t)wb*65+g*5+3)*W));
__m256i v=_mm256_setzero_si256(); _mm256_store_si256(T,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+1,v);
v=_mm256_xor_si256(v,r1); _mm256_store_si256(T+3,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+2,v);
v=_mm256_xor_si256(v,r2); _mm256_store_si256(T+6,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+7,v);
v=_mm256_xor_si256(v,r1); _mm256_store_si256(T+5,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+4,v);
v=_mm256_xor_si256(v,r3); _mm256_store_si256(T+12,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+13,v);
v=_mm256_xor_si256(v,r1); _mm256_store_si256(T+15,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+14,v);
v=_mm256_xor_si256(v,r2); _mm256_store_si256(T+10,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+11,v);
v=_mm256_xor_si256(v,r1); _mm256_store_si256(T+9,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+8,v);
if(g<12){
const __m256i r4=_mm256_load_si256((const __m256i*)(BP2+((size_t)wb*65+g*5+4)*W));
v=_mm256_xor_si256(v,r4); _mm256_store_si256(T+24,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+25,v);
v=_mm256_xor_si256(v,r1); _mm256_store_si256(T+27,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+26,v);
v=_mm256_xor_si256(v,r2); _mm256_store_si256(T+30,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+31,v);
v=_mm256_xor_si256(v,r1); _mm256_store_si256(T+29,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+28,v);
v=_mm256_xor_si256(v,r3); _mm256_store_si256(T+20,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+21,v);
v=_mm256_xor_si256(v,r1); _mm256_store_si256(T+23,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+22,v);
v=_mm256_xor_si256(v,r2); _mm256_store_si256(T+18,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+19,v);
v=_mm256_xor_si256(v,r1); _mm256_store_si256(T+17,v);
v=_mm256_xor_si256(v,r0); _mm256_store_si256(T+16,v);
}
}
for(int i=0;i<n;i+=2){
const uint64_t w0=A1b[(size_t)i*W+wb],w1=A1b[(size_t)(i+1)*W+wb];
__m256i x0=wb?_mm256_load_si256((const __m256i*)(Pbb+(size_t)i*W)):_mm256_setzero_si256();
__m256i y0=wb?_mm256_load_si256((const __m256i*)(Pbb+(size_t)(i+1)*W)):_mm256_setzero_si256();
__m256i x1=_mm256_setzero_si256(),y1=x1;
x0=_mm256_xor_si256(x0,_mm256_load_si256((const __m256i*)TB+0+((w0>>0)&31)));
y0=_mm256_xor_si256(y0,_mm256_load_si256((const __m256i*)TB+0+((w1>>0)&31)));
x1=_mm256_xor_si256(x1,_mm256_load_si256((const __m256i*)TB+32+((w0>>5)&31)));
y1=_mm256_xor_si256(y1,_mm256_load_si256((const __m256i*)TB+32+((w1>>5)&31)));
x0=_mm256_xor_si256(x0,_mm256_load_si256((const __m256i*)TB+64+((w0>>10)&31)));
y0=_mm256_xor_si256(y0,_mm256_load_si256((const __m256i*)TB+64+((w1>>10)&31)));
x1=_mm256_xor_si256(x1,_mm256_load_si256((const __m256i*)TB+96+((w0>>15)&31)));
y1=_mm256_xor_si256(y1,_mm256_load_si256((const __m256i*)TB+96+((w1>>15)&31)));
x0=_mm256_xor_si256(x0,_mm256_load_si256((const __m256i*)TB+128+((w0>>20)&31)));
y0=_mm256_xor_si256(y0,_mm256_load_si256((const __m256i*)TB+128+((w1>>20)&31)));
x1=_mm256_xor_si256(x1,_mm256_load_si256((const __m256i*)TB+160+((w0>>25)&31)));
y1=_mm256_xor_si256(y1,_mm256_load_si256((const __m256i*)TB+160+((w1>>25)&31)));
x0=_mm256_xor_si256(x0,_mm256_load_si256((const __m256i*)TB+192+((w0>>30)&31)));
y0=_mm256_xor_si256(y0,_mm256_load_si256((const __m256i*)TB+192+((w1>>30)&31)));
x1=_mm256_xor_si256(x1,_mm256_load_si256((const __m256i*)TB+224+((w0>>35)&31)));
y1=_mm256_xor_si256(y1,_mm256_load_si256((const __m256i*)TB+224+((w1>>35)&31)));
x0=_mm256_xor_si256(x0,_mm256_load_si256((const __m256i*)TB+256+((w0>>40)&31)));
y0=_mm256_xor_si256(y0,_mm256_load_si256((const __m256i*)TB+256+((w1>>40)&31)));
x1=_mm256_xor_si256(x1,_mm256_load_si256((const __m256i*)TB+288+((w0>>45)&31)));
y1=_mm256_xor_si256(y1,_mm256_load_si256((const __m256i*)TB+288+((w1>>45)&31)));
x0=_mm256_xor_si256(x0,_mm256_load_si256((const __m256i*)TB+320+((w0>>50)&31)));
y0=_mm256_xor_si256(y0,_mm256_load_si256((const __m256i*)TB+320+((w1>>50)&31)));
x1=_mm256_xor_si256(x1,_mm256_load_si256((const __m256i*)TB+352+((w0>>55)&31)));
y1=_mm256_xor_si256(y1,_mm256_load_si256((const __m256i*)TB+352+((w1>>55)&31)));
x0=_mm256_xor_si256(x0,_mm256_load_si256((const __m256i*)TB+384+((w0>>60)&15)));
y0=_mm256_xor_si256(y0,_mm256_load_si256((const __m256i*)TB+384+((w1>>60)&15)));
_mm256_store_si256((__m256i*)(Pbb+(size_t)i*W),_mm256_xor_si256(x0,x1));
_mm256_store_si256((__m256i*)(Pbb+(size_t)(i+1)*W),_mm256_xor_si256(y0,y1));
}
}
}
static inline int sgn_of(int op){ return op==0 ? 0 : (op==1 ? 1 : -1); }
AVX2 static void buildA4(int h, int n, const int8_t*A, int k1, int k2){
const __m256i msk=_mm256_set1_epi8(0x7F);
const int o1 = 512;
const int *t1=SAB_T[k1], *t2=SAB_T[k2];
const int s1=sgn_of(t1[4]), s2=sgn_of(t2[4]), s12=s1*s2;
const size_t R1=(size_t)(t1[0]*o1+t2[0]*h)*n + t1[1]*o1+t2[1]*h;
const size_t R2=(size_t)(t1[2]*o1+t2[0]*h)*n + t1[3]*o1+t2[1]*h;
const size_t R3=(size_t)(t1[0]*o1+t2[2]*h)*n + t1[1]*o1+t2[3]*h;
const size_t R4=(size_t)(t1[2]*o1+t2[2]*h)*n + t1[3]*o1+t2[3]*h;
int W=h/64;
for (int i=0;i<h;i++){
const int8_t*p1=A+R1+(size_t)i*n;
const int8_t*p2=A+R2+(size_t)i*n;
const int8_t*p3=A+R3+(size_t)i*n;
const int8_t*p4=A+R4+(size_t)i*n;
int8_t*dst=A0b+2*(size_t)i*h; /* [c1w_bd] duplicated pairs */
uint64_t*w1=A1b+(size_t)i*W;
for (int k=0;k<h;k+=64){
_mm_prefetch((const char*)(p1+k+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(p2+k+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(p3+k+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(p4+k+4096), _MM_HINT_T0);
__asm__ volatile("prefetchw %0" :: "m"(*(char*)(dst+2*k+256)));
__asm__ volatile("prefetchw %0" :: "m"(*(char*)(dst+2*k+320)));
__m256i v0=_mm256_load_si256((const __m256i*)(p1+k));
__m256i v1=_mm256_load_si256((const __m256i*)(p1+k+32));
if(s1){ __m256i u0=_mm256_load_si256((const __m256i*)(p2+k));
__m256i u1=_mm256_load_si256((const __m256i*)(p2+k+32));
if(s1>0){v0=_mm256_add_epi8(v0,u0);v1=_mm256_add_epi8(v1,u1);}
else {v0=_mm256_sub_epi8(v0,u0);v1=_mm256_sub_epi8(v1,u1);} }
if(s2){ __m256i u0=_mm256_load_si256((const __m256i*)(p3+k));
__m256i u1=_mm256_load_si256((const __m256i*)(p3+k+32));
if(s2>0){v0=_mm256_add_epi8(v0,u0);v1=_mm256_add_epi8(v1,u1);}
else {v0=_mm256_sub_epi8(v0,u0);v1=_mm256_sub_epi8(v1,u1);} }
if(s12){__m256i u0=_mm256_load_si256((const __m256i*)(p4+k));
__m256i u1=_mm256_load_si256((const __m256i*)(p4+k+32));
if(s12>0){v0=_mm256_add_epi8(v0,u0);v1=_mm256_add_epi8(v1,u1);}
else {v0=_mm256_sub_epi8(v0,u0);v1=_mm256_sub_epi8(v1,u1);} }
{ /* [c1w_bd] each 16-bit A pair stored twice => vpbroadcastd is exact.
vpunpcklwd is per-128-bit-lane, so the two halves are put back in
order with vperm2i128: pair P then lives at byte 4*P exactly. */
__m256i e0=_mm256_and_si256(v0,msk), e1=_mm256_and_si256(v1,msk);
__m256i l0=_mm256_unpacklo_epi16(e0,e0), h0=_mm256_unpackhi_epi16(e0,e0);
__m256i l1=_mm256_unpacklo_epi16(e1,e1), h1=_mm256_unpackhi_epi16(e1,e1);
_mm256_store_si256((__m256i*)(dst+2*k), _mm256_permute2x128_si256(l0,h0,0x20));
_mm256_store_si256((__m256i*)(dst+2*k+32), _mm256_permute2x128_si256(l0,h0,0x31));
_mm256_store_si256((__m256i*)(dst+2*k+64), _mm256_permute2x128_si256(l1,h1,0x20));
_mm256_store_si256((__m256i*)(dst+2*k+96), _mm256_permute2x128_si256(l1,h1,0x31)); }
w1[k>>6]=(uint64_t)(uint32_t)_mm256_movemask_epi8(v0)
|((uint64_t)(uint32_t)_mm256_movemask_epi8(v1)<<32);
}
}
}
/* [nm1_4t] same four-term composition on the B side; identical layout to bld_B
(BI3b interleaved pairs + the padded BP2 parity panel with its zeroed sentinel row). */
AVX2 static void buildB4(int h, int n, const int8_t*B, int k1, int k2){
const __m256i one=_mm256_set1_epi8(1);
const int o1 = 512;
const int *t1=SBB_T[k1], *t2=SBB_T[k2];
const int s1=sgn_of(t1[4]), s2=sgn_of(t2[4]), s12=s1*s2;
const size_t R1=(size_t)(t1[0]*o1+t2[0]*h)*n + t1[1]*o1+t2[1]*h;
const size_t R2=(size_t)(t1[2]*o1+t2[0]*h)*n + t1[3]*o1+t2[1]*h;
const size_t R3=(size_t)(t1[0]*o1+t2[2]*h)*n + t1[1]*o1+t2[3]*h;
const size_t R4=(size_t)(t1[2]*o1+t2[2]*h)*n + t1[3]*o1+t2[3]*h;
int ng=h/64, W=h/64;
for (int p=0;p<h/2;p++){
const int8_t *a1=B+R1+(size_t)(2*p)*n, *b1=B+R1+(size_t)(2*p+1)*n;
const int8_t *a2=B+R2+(size_t)(2*p)*n, *b2=B+R2+(size_t)(2*p+1)*n;
const int8_t *a3=B+R3+(size_t)(2*p)*n, *b3=B+R3+(size_t)(2*p+1)*n;
const int8_t *a4=B+R4+(size_t)(2*p)*n, *b4=B+R4+(size_t)(2*p+1)*n;
uint64_t *w0p=BP2+(size_t)((2*p)>>6)*(size_t)(65*W)+(size_t)((2*p)&63)*W;
uint64_t *w1p=BP2+(size_t)((2*p+1)>>6)*(size_t)(65*W)+(size_t)((2*p+1)&63)*W;
for (int g=0;g<ng;g++){
const int8_t *a=a1+64*g, *c=b1+64*g;
const int8_t *d=a2+64*g, *e=b2+64*g;
const int8_t *f=a3+64*g, *h3=b3+64*g;
const int8_t *m4=a4+64*g, *n4=b4+64*g;
int8_t *o=BI3b+((size_t)g*(h/2)+p)*128;
for (int t=0;t<64;t+=32){
_mm_prefetch((const char*)(a +t+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(c +t+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(d +t+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(e +t+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(f +t+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(h3+t+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(m4+t+4096), _MM_HINT_T0);
_mm_prefetch((const char*)(n4+t+4096), _MM_HINT_T0);
__asm__ volatile("prefetchw %0" :: "m"(*(char*)(o+2*t+256)));
__asm__ volatile("prefetchw %0" :: "m"(*(char*)(o+2*t+320)));
__m256i va=_mm256_load_si256((const __m256i*)(a+t));
__m256i vb=_mm256_load_si256((const __m256i*)(c+t));
if(s1){ __m256i u0=_mm256_load_si256((const __m256i*)(d+t));
__m256i u1=_mm256_load_si256((const __m256i*)(e+t));
if(s1>0){va=_mm256_add_epi8(va,u0);vb=_mm256_add_epi8(vb,u1);}
else {va=_mm256_sub_epi8(va,u0);vb=_mm256_sub_epi8(vb,u1);} }
if(s2){ __m256i u0=_mm256_load_si256((const __m256i*)(f+t));
__m256i u1=_mm256_load_si256((const __m256i*)(h3+t));
if(s2>0){va=_mm256_add_epi8(va,u0);vb=_mm256_add_epi8(vb,u1);}
else {va=_mm256_sub_epi8(va,u0);vb=_mm256_sub_epi8(vb,u1);} }
if(s12){__m256i u0=_mm256_load_si256((const __m256i*)(m4+t));
__m256i u1=_mm256_load_si256((const __m256i*)(n4+t));
if(s12>0){va=_mm256_add_epi8(va,u0);vb=_mm256_add_epi8(vb,u1);}
else {va=_mm256_sub_epi8(va,u0);vb=_mm256_sub_epi8(vb,u1);} }
__m256i lo=_mm256_unpacklo_epi8(va,vb), hi=_mm256_unpackhi_epi8(va,vb);
_mm256_store_si256((__m256i*)(o+2*t), _mm256_permute2x128_si256(lo,hi,0x20));
_mm256_store_si256((__m256i*)(o+2*t+32), _mm256_permute2x128_si256(lo,hi,0x31));
uint32_t m0=(uint32_t)_mm256_movemask_epi8(_mm256_slli_epi16(_mm256_and_si256(va,one),7));
uint32_t m1=(uint32_t)_mm256_movemask_epi8(_mm256_slli_epi16(_mm256_and_si256(vb,one),7));
if(t==0) { w0p[g]=(uint64_t)m0; w1p[g]=(uint64_t)m1; }
else { w0p[g]|=(uint64_t)m0<<32; w1p[g]|=(uint64_t)m1<<32; }
}
}
}
}
/* one 256 leaf straight out of the ORIGINAL operands (k1 = outer case, k2 = inner). */
AVX2 static void submm512_4(int m, int n, const int8_t*A, const int8_t*B, int k1, int k2,
int8_t*Cp){
buildA4(m,n,A,k1,k2); buildB4(m,n,B,k1,k2);
{ const int W=m>>6; for (int g=0; g<(m>>6); g++) memset(BP2+(size_t)g*65*W+64*W, 0, (size_t)W*8); }
if(m==256) corr4_256(); else corr_rus(m);
for (int g=0; g<m; g+=64)
for (int i=0; i<m; i+=2){
mk_s(m, A0b+2*(size_t)i*m, BI3b+(size_t)(g>>6)*(m/2)*128, Cp, PBb, g);
epilogue(m, Cp, (const uint64_t*)Pbb, g, i, 2);
}
}
/* [nm1_4t] the 7 inner leaves are built straight from the ORIGINAL A/B (k1 = this
outer case). No RawA/RawB, no materialize512. */
/* [mc1y_] TOP-LEVEL FOLD. C11=M0+M3-M4+M6 C12=M2+M4 C21=M1+M3 C22=M0+M2-M1+M5
(CLRS P-form). Block k1's four quadrants go straight into the destination
C blocks instead of MMS[k1]; M (0=store, 1=add, -1=sub) is the block's role in
the P-form combination, so no MMS buffer and no final cmb_nt pass exist. */
struct KD { int R, Cc, M; };
static constexpr KD KDEST[7][2] = {
{{ 0, 0, 0}, {512, 512, 0}}, /* M0: +C11 +C22 (first writer of C11,C22) */
{{512, 0, 0}, {512, 512, -1}}, /* M1: +C21 -C22 */
{{ 0, 512, 0}, {512, 512, 1}}, /* M2: +C12 +C22 */
{{ 0, 0, 1}, {512, 0, 1}}, /* M3: +C11 +C21 */
{{ 0, 0,-1}, { 0, 512, 1}}, /* M4: -C11 +C12 */
{{512, 512, 1}, { 0, 0, 0}}, /* M5: +C22 */
{{ 0, 0, 1}, { 0, 0, 0}}, /* M6: +C11 */
};
static constexpr int KDN[7] = {2,2,2,2,2,1,1};
#define LD32(p) _mm256_load_si256((const __m256i*)(p))
#define ST32(p,v) _mm256_store_si256((__m256i*)(p), (v))
#define ADD8(a,b) _mm256_add_epi8((a),(b))
#define SUB8(a,b) _mm256_sub_epi8((a),(b))
template<int MODE> static inline void upd32(int8_t*p, __m256i v){
if(MODE==0) ST32(p, v);
else if(MODE==1) ST32(p, ADD8(LD32(p), v));
else ST32(p, SUB8(LD32(p), v));
}
/* [mc1z_] FUSED quadrant pass. The four destination quadrants share the same
source row, so each M2 block is read once per k1 instead of 12/7 times.
v00..v11 are the P-form sums (identical expressions to the four qw calls). */
template<int MODE1,int MODE2>
AVX2 static void qwall(int8_t*C,int n,int m,int8_t*D2){
const int h=m>>1;
int8_t *a1=C,*a2=C+h,*a3=C+(size_t)h*n,*a4=C+(size_t)h*n+h;
int8_t *b1=0,*b2=0,*b3=0,*b4=0;
if(D2){ b1=D2; b2=D2+h; b3=D2+(size_t)h*n; b4=D2+(size_t)h*n+h; }
for(int i=0;i<h;i++){
const int8_t *p0=M2[0]+(size_t)i*h,*p1=M2[1]+(size_t)i*h,*p2=M2[2]+(size_t)i*h,
*p3=M2[3]+(size_t)i*h,*p4=M2[4]+(size_t)i*h,*p5=M2[5]+(size_t)i*h,
*p6=M2[6]+(size_t)i*h;
int8_t *e1=a1+(size_t)i*n,*e2=a2+(size_t)i*n,*e3=a3+(size_t)i*n,*e4=a4+(size_t)i*n;
int8_t *f1=b1?b1+(size_t)i*n:0,*f2=b2?b2+(size_t)i*n:0,
*f3=b3?b3+(size_t)i*n:0,*f4=b4?b4+(size_t)i*n:0;
if(i+2<h){
const char*q0=(const char*)(M2[0]+(size_t)(i+2)*h);
_mm_prefetch(q0,_MM_HINT_T0); _mm_prefetch(q0+64,_MM_HINT_T0); _mm_prefetch(q0+128,_MM_HINT_T0); _mm_prefetch(q0+192,_MM_HINT_T0);
_mm_prefetch((const char*)(M2[1]+(size_t)(i+2)*h),_MM_HINT_T0);
_mm_prefetch((const char*)(M2[2]+(size_t)(i+2)*h),_MM_HINT_T0);
_mm_prefetch((const char*)(M2[3]+(size_t)(i+2)*h),_MM_HINT_T0);
_mm_prefetch((const char*)(M2[4]+(size_t)(i+2)*h),_MM_HINT_T0);
_mm_prefetch((const char*)(M2[5]+(size_t)(i+2)*h),_MM_HINT_T0);
_mm_prefetch((const char*)(M2[6]+(size_t)(i+2)*h),_MM_HINT_T0);
_mm_prefetch((const char*)(e1+(size_t)2*n),_MM_HINT_T0);
_mm_prefetch((const char*)(e2+(size_t)2*n),_MM_HINT_T0);
_mm_prefetch((const char*)(e3+(size_t)2*n),_MM_HINT_T0);
_mm_prefetch((const char*)(e4+(size_t)2*n),_MM_HINT_T0);
if(f1){ _mm_prefetch((const char*)(f1+(size_t)2*n),_MM_HINT_T0);
_mm_prefetch((const char*)(f2+(size_t)2*n),_MM_HINT_T0);
_mm_prefetch((const char*)(f3+(size_t)2*n),_MM_HINT_T0);
_mm_prefetch((const char*)(f4+(size_t)2*n),_MM_HINT_T0); }
}
for(int j=0;j<h;j+=32){
const __m256i v0=LD32(p0+j),v1=LD32(p1+j),v2=LD32(p2+j),v3=LD32(p3+j),
v4=LD32(p4+j),v5=LD32(p5+j),v6=LD32(p6+j);
const __m256i y00=ADD8(ADD8(v0,v3),SUB8(v6,v4));
const __m256i y01=ADD8(v2,v4);
const __m256i y10=ADD8(v1,v3);
const __m256i y11=ADD8(ADD8(SUB8(v0,v1),v2),v5);
upd32<MODE1>(e1+j,y00); upd32<MODE1>(e2+j,y01);
upd32<MODE1>(e3+j,y10); upd32<MODE1>(e4+j,y11);
if(f1){ upd32<MODE2>(f1+j,y00); upd32<MODE2>(f2+j,y01);
upd32<MODE2>(f3+j,y10); upd32<MODE2>(f4+j,y11); }
}
}
}
template<int K1>
AVX2 static void strassen2_512x(const int8_t*A, const int8_t*B, int n, int8_t*C){
const int h=256;
for(int k=0;k<7;k++) submm512_4(h,n,A,B,K1,k,M2[k]);
int8_t *D1=C+(size_t)KDEST[K1][0].R*n+KDEST[K1][0].Cc;
int8_t *D2=(KDN[K1]>1)?C+(size_t)KDEST[K1][1].R*n+KDEST[K1][1].Cc:0;
if(KDN[K1]>1) qwall<KDEST[K1][0].M,KDEST[K1][1].M>(D1,n,512,D2);
else qwall<KDEST[K1][0].M,0>(D1,n,512,0);
}
AVX2 static void strassen1(int n, const int8_t*A, const int8_t*B, int8_t*C){
strassen2_512x<0>(A,B,n,C); strassen2_512x<1>(A,B,n,C);
strassen2_512x<2>(A,B,n,C); strassen2_512x<3>(A,B,n,C);
strassen2_512x<4>(A,B,n,C); strassen2_512x<5>(A,B,n,C);
strassen2_512x<6>(A,B,n,C);
_mm_sfence();
}
AVX2 void matrix_multiply(int n, const int8_t *A, const int8_t *B, int8_t *C){
if (n<=0 || (n&255) || n>NMAX) { generic_mm(n,A,B,C); return; }
if (n==1024) { strassen1(n,A,B,C); return; }
build_A0_A1(n,A);
if(n==1024){
build_BI9(n,B); build_BI3_tail(n,B); build_BP2(n,B); corr_rus(n);
for (int g=0; g<10; g++)
for (int i=0; i<n; i+=2){
mk_s6(n, A0b+(size_t)i*n, BI9b+(size_t)g*(n/2)*192, (uint8_t*)mscr);
epilogue(n, C, (const uint64_t*)Pbb, g*96, i, 3);
}
for (int i=0; i<n; i+=2){
mk_s(n, A0b+2*(size_t)i*n, BI3b, C, PBb, 0);
epilogue(n, C, (const uint64_t*)Pbb, 960, i, 2);
}
return;
}
build_BI3(n,B); build_Bpbits(n,B);
if(n==4096) corr16(); else corr_pass3(n);
for (int g=0, jbg=0; g<n; g+=64, jbg=g)
for (int i=0; i<n; i+=2)
{ mk_s(n, A0b+2*(size_t)i*n, BI3b+(size_t)(g>>6)*(n/2)*128, C, PBb, g);
epilogue(n, C, (const uint64_t*)Pbb, g, i, 2); }
}
| Compilation | N/A | N/A | Compile OK | Score: N/A | 显示更多 |
| Testcase #1 | 6.734 ms | 1 MB + 708 KB | Accepted | Score: 100 | 显示更多 |