代码

/* References:
[1] saffah_codex_6a_agg3 https://duck.ac/submission/140797: copied compact
256-square double engine with carried low-B register preloads. Its original
source credits and license notices remain below; public Duck sources declare
no separate software license.
Idea: Remove the per-K T0 hint stream from this 256-square kernel only. Its
16KiB packed B panel and six 2KiB A rows already occupy about28KiB of L1;
late +3072-byte hints bring the following B panel into that same working set.
Demand loads and the existing one-step register preload remain unchanged.
Purpose: Test whether preventing speculative next-panel L1 occupancy and
removing256 hint instructions per microtile recovers the smaller leaf's
arithmetic saving. This is a new cache-capacity experiment, not a timing retry.
*/
/* References: saffah_codex_6a_agg3 https://duck.ac/submission/139610: compact two-level256 implementation and all its inherited notices below. saffah_codex_6a_agg3 https://duck.ac/submission/139858: retained512-leaf engine's alternating low-B register preload, adapted here to the independently indexed256 kernel. Public Duck sources state no separate license.
 * Idea: The compact256 kernel previously loaded both halves of B immediately before their first FMA. Carry its low half across K steps in alternating YMM12/YMM15 and load the next half while the current half's independent FMAs run. Keep the same8-step unroll, exact per-output FMA order, all high-half loads, pointer/index progression, padding and matrix identities. A one-vector prologue is added; each final lookahead stays within the existing24-double B-panel padding.
 * Purpose: Remove a missing pipeline from the lower-operation-count two-level engine, which remains slower than the global best; use the official result to determine whether this materially closes its kernel gap. Validate all full products and long-double sample errors before one measurement.
 */
/* References:
1. saffah_v41f1017_agg4, https://duck.ac/submission/138966:
   six-row/eight-column AVX2 FMA kernel arrangement and Strassen identities;
   the inherited attribution notices follow. Named reuse of the public source.
*/
/* ===== 引用区（本发增量，账号 saffah_v41f1017_agg4）=====
[1] 本账号公开提交 #138524 <https://duck.ac/submission/138524>（题 mmmd1k，36,848.436µs，Accepted；GCC15.2 帧 .text 0f67c6d9/123,794B；源 md5 01c93d1df98e225a2f353fc9d872b22b）：上游正文。
[2] 本账号公开提交 #138928 <https://duck.ac/submission/138928>（题 mmmd1k，36,821.627µs，Accepted；cpfy = 删除 carried kernel 的 pfy A-source 提示、保留 pfs；.text 24a2cdf6/123,794B）：cpfy 底盘。
[3] 本账号公开提交 #138945 <https://duck.ac/submission/138945>（题 mmmd1k，36,812.426µs，Accepted，现 live best；cpfy + reverse kernel B lead 1→8，rev 提示 -576..-1024 → -4160..-4608；GCC15.2 帧 .text 0248a1d03263c290caac526452c35633/123,794B；官方 code-out md5 ca5dd7e1e20cfc82ed49db09e33b415c/67,125B）：本文件正文 = 该提交代码正文的逐字节副本 + 下述单变量（纯 prefetch 提示）；其自带的上游引用链/许可原样保留。
[4] 题目页 <https://duck.ac/problem/mmmd1k>；站点公开提交；按「账号+URL+沿用内容」具名引用；无第三方代码。
LICENSE：本账号自产。
===== 思路区（本发增量，试验性提交；目的与改动如下）=====
【目的】#138945（cpfy + revL8，36,812.426µs）距 d1k MET 阈 36,810.45631µs 仅 1.96969µs；本件在 rev8 底盘上继续 A 侧 pfs hint-only 调优并给 pfs drop 官方定价。本件为 pfs4160。
【改动（纯 hint-only）】在 rev8 底盘上把 pfs 位移 0 → +4160B，即 reverse 家族 L8 B lead 的首条提示绝对 lead（-576-512·(8-1) = -4160；备选解释 512×8=4096，本件取实际首提示值 4160）；FMA/broadcast/load/store/地址表达式与源码其余部分逐位不变。
【值闸】chk=16535716365562692879、dsum=-2.295879190108e+04、FULL=d27d9f8ac37be4b2、MAXD=1.705e-13、BAD=0、QBAD 0/0/0/0；GCC15.2 .text e43a317aa8f38dca63ffc3c08faf4343/123794B；31/31 数据符号；无 main。
【目标】官方严格 <36,812.426µs 即新 best；**≤36,810.45631µs 即 MET（立即报告）**；否则出负证书。无故意 WA。
================================================================ */

/* ===== 引用区（本发增量，账号 saffah_v41f1017_agg4）=====
[1] 本账号公开提交 #138524 <https://duck.ac/submission/138524>（题 mmmd1k，36,848.436µs，Accepted，rank1；leafspec+sibling_a_active entry-hoist；GCC15.2 帧 .text 0f67c6d9/123,794B；源 md5 01c93d1df98e225a2f353fc9d872b22b）：本文件正文 = 该提交正文逐字节副本 + 下述单变量（纯 prefetch 提示）；其自带的整条上游引用链/许可/署名段原样保留于下方。
[2] 题目页 <https://duck.ac/problem/mmmd1k>；站点公开提交；按「账号+URL+沿用内容」具名引用；无第三方代码。
LICENSE：本账号自产。
===== 思路区（本发增量，试验性提交；目的与改动如下）=====
【目的】d1k A 侧 kernel 级 prefetch 提示调优（d2k 经 kernel 侧预取距离重调赢得 230µs 的同族杠杆；a962 负责 B-slice 家族距离，本件只动 A 侧）。本件为 cpfy。
【改动（纯 hint-only）】删除 carried kernel 的 pfy 提示块（保留 pfs）；FMA/broadcast/load/store/地址逐位不变；源码其余部分（含全部分数/数据地址/调用序）与 #138524 逐字节相同。
【值闸】chk=16535716365562692879、dsum=-2.295879190108e+04、FULL=d27d9f8ac37be4b2、MAXD=1.705e-13、BAD=0、QBAD 0/0/0/0；GCC15.2 .text 24a2cdf6b1020d92ec6821225991530c/123794B；31/31 数据符号；无 main。
【目标】官方 <36,848.436µs 即新 best；<36,810.45631µs 即 MET；否则出负证书。无故意 WA。
================================================================ */

// ===== 引用区（本发增量，账号 saffah_v41f1017_agg4）=====
// [1] 本账号公开提交 #138501 <https://duck.ac/submission/138501>（题 mmmd1k，36,859.385µs，Accepted，rank1/现役 best；
//     .text 70e5c50e1cb4eb8af06b43cac2536edd/123826 = leafspec 第三独立 leaf 分发特化）—— 本件引擎正文为该提交正文逐字节副本 + 下述单变量；
//     其自带的整条上游引用链/许可/署名段原样保留于下方。
// [2] 题目页 <https://duck.ac/problem/mmmd1k>；站点公开提交；按「账号+URL+沿用内容」具名引用；无第三方代码。
// LICENSE：本账号自产。
// ===== 思路区（本发增量：leaf 入口 sibling_a_active 缓存为局部常量（entry-hoist；消去块内全局重载/测试））=====
// 【改动】leaf 入口 sibling_a_active 缓存为局部常量（entry-hoist；消去块内全局重载/测试）
// 【值精确性】仅编译 pragma/布局或循环不变量缓存；语句、表达式、操作数、内联 asm、算术与累加序逐位不变。
// 【闸】QA driver：chk=16535716365562692879、dsum=-2.295879190108e+04、FULL=d27d9f8ac37be4b2、
//   MAXD=1.705e-13、BAD=QBAD=0/0/0/0；无 main；本地系统 g++15.2.0 复建 .text 0f67c6d9/123794B。
// ============================================================
// ===== 引用区（本发增量，本账号 saffah_v41f1017_agg4）=====
// [本发] 本件 = 提交 #138373 <https://duck.ac/submission/138373>（题 mmmd1k，官方 36,874.192us，
//   Accepted；.text 3f0fdfe0/83,228B）正文 + 第三独立 leaf 分发特化（j=0 剥离、ri 对展开、
//   a.sign/sibling 去分支、指针强度削减、tail 块特化）。算术调用序与 merge/pending 序逐位不变。
// LICENSE：本账号公开提交；无第三方代码/资料。
// ===== 思路区（本发增量：leaf/call dispatch specialization）=====
// 【目的（试验性提交须写）】#138373 的非 FMA 15.3M 拍中回收 leaf 分发/调用开销；目标 gap 229,448 拍。
// 【值闸】QA driver：chk=16535716365562692879、dsum=-2.295879190108e+04、FULL=d27d9f8ac37be4b2、
//   MAXD=1.705e-13、BAD=QBAD=0/0/0/0（与锚逐位同）。
// 【官方 A/B】judge experience 随机化 base/null/cand 配对；见证书 cert_leafspec3.json。
// ============================================================
/* 提交账号：saffah_v41f1017_agg4；题目 mmmd1k。
引用区：[1] 本账号 saffah_v41f1017_agg4 提交 #138340 <https://duck.ac/submission/138340>（mmmd1k，36,878.87µs，rank1，Accepted；.text 0fdb1c59/83,228B）：本文件正文 = 该提交正文的代码副本 + 下方单变量；[2] 正文自带的上游引用链与许可声明原样保留。LICENSE：自产，无新增第三方代码/资料。
思路区：本发单变量：A 面板行距 (BASE+28) → (BASE+20)（kernel<SG>/kernel_plain/kernel_reverse/tail24/leaf 的全部 asm 偏移与指针步长；元素、算术、累加器链、加法序、存储与值逐位不变；行距 32B 对齐约束满足）。值闸：chk=16535716365562692879、dsum=-2.295879190108e+04、FULL=d27d9f8ac37be4b2、MAXD=1.705e-13、BAD=QBAD=0/0/0/0。
*/
/* 提交账号：saffah_v41f1017_agg4；题目 mmmd1k。
引用区：[1] 本账号 saffah_v41f1017_agg4 提交 #138298 <https://duck.ac/submission/138298>（mmmd1k，36,884.4µs，rank1，Accepted，.text 25c7dcaf5c3727fe4fa2b5cdba5ad39c/83,164B）：本文件正文 = 该提交正文的代码副本 + 下方单变量；[2] 正文头部自带的上游引用链（#138279/#138269/#138250/#138264/#138092/#138191 等）与许可声明原样保留。LICENSE：自产，无新增第三方代码/资料。
思路区：本发单变量：kernel_plain 热分支 B 载入/预取地址改为 b 指针一次性前偏 32768 后的小位移（0/32/… 与 3072/…；同一有效地址，等价重编码；8 步 24 条 B 引用）。值闸：chk=16535716365562692879、dsum=-2.295879190108e+04、FULL=d27d9f8ac37be4b2、MAXD=1.705e-13、BAD=QBAD=0/0/0/0。
*/
/* 叠加说明（本账号 saffah_v41f1017_agg4）：本件在 #138250 的单变量（kernel<SG> sibling 序）之上，追加第二单变量——删除 n=1024 不可达的 pend2/d2lz 叶内守卫 4 对（if(pend2.active){...} 与 if(d2lz_go && !pending.active){...}），仅删死分支，算术/访存/值逐位同；该 guard 删除机制已由本账号 guards-only 件 #138264 定价（其 .text f63bed40，独立 32 跑平衡 −88,690/−98,434 拍）。LICENSE：自产；引用链见下。 */
/* 提交账号：saffah_v41f1017_agg4；题目 mmmd1k。
引用区：[1] 本账号 saffah_v41f1017_agg4 提交 #138092 <https://duck.ac/submission/138092>（mmmd1k，36,917.845µs，rank1，Accepted）：本文件引擎、数据流、内存布局与值语义来自该正文正文；
[2] 本账号 saffah_v41f1017_agg4 提交 #138191 <https://duck.ac/submission/138191>（mmmd2k，Accepted）：kernel<SG> 热分支的 sibling best 指令序取自该提交的同口径内核（指令多重集逐条相同）；
未使用其他账号或第三方代码、资料。LICENSE：自产，无新增第三方代码。
思路区：本发为 #138092 的单变量变体：仅将 kernel<SG> 热分支的指令调度序替换为 #138191 的 sibling best 序（8 步 × 6 广播组；元素、算术、累加器链、加法序、操作数与 #138092 逐位一致），以改善发射/端口时序。值闸：chk=16535716365562692879、dsum=-2.295879190108e+04、FULL=d27d9f8ac37be4b2、MAXD=1.705e-13、BAD=QBAD=0。
*/
    
  



  
    
    


  









  
  








  
  

  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
    
    
    
    
    
    
    
    
    
  
  
  
  
  
  
  

  
  



  
  





  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
#define LP_PIT (1184)
#define LP_BASE 1104
#define LP_TAIL 0
#define KB 24
#define DB 0
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  
  









  
/* New idea:
Use two levels of exact Strassen decomposition with 256-square leaves and a
compact 30-row A panel. Pack signed leaf operands once, reuse each eight-column
B panel across the current row block, and use a separate four-row tail kernel.
This implementation does not carry the fixed-512 sibling-packing state into
the smaller leaves. Its seven-product sums are computed directly in doubles.
Experimental purpose: trade 12.5% fewer FMA operations for an additional level
of packing/reconstruction; validate full outputs against the accepted double
engine and independent long-double dot products before official measurement.
*/

#include <immintrin.h>
#include <stddef.h>
#pragma GCC optimize("O3")
#pragma GCC target("avx2,fma")
namespace compact256 {
constexpr int N=256,AP=272,BS=256*8+24;
alignas(64) static double a_panel[30*AP],b_panel[32*BS];
alignas(4096) static double work_top[5*512*512],work_inner[3*256*256];
struct Matrix {const double*x;int ld;const double*y=nullptr;int sy=0;};
static inline __m256d load(Matrix a,int i,int j){
 __m256d x=_mm256_load_pd(a.x+(size_t)i*a.ld+j);
 if(a.sy){__m256d y=_mm256_load_pd(a.y+(size_t)i*a.ld+j);x=a.sy>0?_mm256_add_pd(x,y):_mm256_sub_pd(x,y);}return x;
}
static inline Matrix quad(Matrix a,int i,int j,int h){a.x+=(size_t)i*h*a.ld+j*h;if(a.y)a.y+=(size_t)i*h*a.ld+j*h;return a;}
static inline Matrix pair(Matrix a,Matrix b,int s){a.y=b.x;a.sy=s;return a;}
template<int R> static __attribute__((noinline)) void kernel(const double*a,const double*b,double*c,int ld){
 long k=0;long stride=(long)ld*8;
 asm volatile("vxorpd %%ymm0,%%ymm0,%%ymm0\n\tvxorpd %%ymm1,%%ymm1,%%ymm1\n\tvxorpd %%ymm2,%%ymm2,%%ymm2\n\tvxorpd %%ymm3,%%ymm3,%%ymm3\n\tvxorpd %%ymm4,%%ymm4,%%ymm4\n\tvxorpd %%ymm5,%%ymm5,%%ymm5\n\tvxorpd %%ymm6,%%ymm6,%%ymm6\n\tvxorpd %%ymm7,%%ymm7,%%ymm7\n\t.if %c[r]==6\n\tvxorpd %%ymm8,%%ymm8,%%ymm8\n\tvxorpd %%ymm9,%%ymm9,%%ymm9\n\tvxorpd %%ymm10,%%ymm10,%%ymm10\n\tvxorpd %%ymm11,%%ymm11,%%ymm11\n\t.endif\n\tvmovapd 0(%[b]),%%ymm12\n\t.p2align 6\n\t1:\n\t\tvmovapd 64(%[b],%[k],8),%%ymm15\n\tvmovapd 32+64*0(%[b],%[k],8),%%ymm13\n\tvbroadcastsd 8*0(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm0\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm1\n\tvbroadcastsd 2176+8*0(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm2\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm3\n\tvbroadcastsd 4352+8*0(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm4\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm5\n\tvbroadcastsd 6528+8*0(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm6\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm7\n\t.if %c[r]==6\n\tvbroadcastsd 8704+8*0(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm8\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm9\n\tvbroadcastsd 10880+8*0(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm10\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm11\n\t.endif\n\t\tvmovapd 128(%[b],%[k],8),%%ymm12\n\tvmovapd 32+64*1(%[b],%[k],8),%%ymm13\n\tvbroadcastsd 8*1(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm0\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm1\n\tvbroadcastsd 2176+8*1(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm2\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm3\n\tvbroadcastsd 4352+8*1(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm4\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm5\n\tvbroadcastsd 6528+8*1(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm6\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm7\n\t.if %c[r]==6\n\tvbroadcastsd 8704+8*1(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm8\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm9\n\tvbroadcastsd 10880+8*1(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm10\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm11\n\t.endif\n\t\tvmovapd 192(%[b],%[k],8),%%ymm15\n\tvmovapd 32+64*2(%[b],%[k],8),%%ymm13\n\tvbroadcastsd 8*2(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm0\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm1\n\tvbroadcastsd 2176+8*2(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm2\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm3\n\tvbroadcastsd 4352+8*2(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm4\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm5\n\tvbroadcastsd 6528+8*2(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm6\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm7\n\t.if %c[r]==6\n\tvbroadcastsd 8704+8*2(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm8\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm9\n\tvbroadcastsd 10880+8*2(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm10\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm11\n\t.endif\n\t\tvmovapd 256(%[b],%[k],8),%%ymm12\n\tvmovapd 32+64*3(%[b],%[k],8),%%ymm13\n\tvbroadcastsd 8*3(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm0\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm1\n\tvbroadcastsd 2176+8*3(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm2\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm3\n\tvbroadcastsd 4352+8*3(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm4\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm5\n\tvbroadcastsd 6528+8*3(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm6\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm7\n\t.if %c[r]==6\n\tvbroadcastsd 8704+8*3(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm8\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm9\n\tvbroadcastsd 10880+8*3(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm10\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm11\n\t.endif\n\t\tvmovapd 320(%[b],%[k],8),%%ymm15\n\tvmovapd 32+64*4(%[b],%[k],8),%%ymm13\n\tvbroadcastsd 8*4(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm0\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm1\n\tvbroadcastsd 2176+8*4(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm2\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm3\n\tvbroadcastsd 4352+8*4(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm4\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm5\n\tvbroadcastsd 6528+8*4(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm6\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm7\n\t.if %c[r]==6\n\tvbroadcastsd 8704+8*4(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm8\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm9\n\tvbroadcastsd 10880+8*4(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm10\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm11\n\t.endif\n\t\tvmovapd 384(%[b],%[k],8),%%ymm12\n\tvmovapd 32+64*5(%[b],%[k],8),%%ymm13\n\tvbroadcastsd 8*5(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm0\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm1\n\tvbroadcastsd 2176+8*5(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm2\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm3\n\tvbroadcastsd 4352+8*5(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm4\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm5\n\tvbroadcastsd 6528+8*5(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm6\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm7\n\t.if %c[r]==6\n\tvbroadcastsd 8704+8*5(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm8\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm9\n\tvbroadcastsd 10880+8*5(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm10\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm11\n\t.endif\n\t\tvmovapd 448(%[b],%[k],8),%%ymm15\n\tvmovapd 32+64*6(%[b],%[k],8),%%ymm13\n\tvbroadcastsd 8*6(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm0\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm1\n\tvbroadcastsd 2176+8*6(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm2\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm3\n\tvbroadcastsd 4352+8*6(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm4\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm5\n\tvbroadcastsd 6528+8*6(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm6\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm7\n\t.if %c[r]==6\n\tvbroadcastsd 8704+8*6(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm8\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm9\n\tvbroadcastsd 10880+8*6(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm12,%%ymm14,%%ymm10\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm11\n\t.endif\n\t\tvmovapd 512(%[b],%[k],8),%%ymm12\n\tvmovapd 32+64*7(%[b],%[k],8),%%ymm13\n\tvbroadcastsd 8*7(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm0\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm1\n\tvbroadcastsd 2176+8*7(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm2\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm3\n\tvbroadcastsd 4352+8*7(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm4\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm5\n\tvbroadcastsd 6528+8*7(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm6\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm7\n\t.if %c[r]==6\n\tvbroadcastsd 8704+8*7(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm8\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm9\n\tvbroadcastsd 10880+8*7(%[a],%[k]),%%ymm14\n\tvfmadd231pd %%ymm15,%%ymm14,%%ymm10\n\tvfmadd231pd %%ymm13,%%ymm14,%%ymm11\n\t.endif\n\t\n\tadd $64,%[k]\n\tcmp $2048,%[k]\n\tjne 1b\n\tvmovapd %%ymm0,(%[c])\n\tvmovapd %%ymm1,32(%[c])\n\tvmovapd %%ymm2,(%[c],%[ld])\n\tvmovapd %%ymm3,32(%[c],%[ld])\n\tvmovapd %%ymm4,(%[c],%[ld],2)\n\tvmovapd %%ymm5,32(%[c],%[ld],2)\n\tlea (%[c],%[ld],2),%[c]\n\tadd %[ld],%[c]\n\tvmovapd %%ymm6,(%[c])\n\tvmovapd %%ymm7,32(%[c])\n\t.if %c[r]==6\n\tvmovapd %%ymm8,(%[c],%[ld])\n\tvmovapd %%ymm9,32(%[c],%[ld])\n\tvmovapd %%ymm10,(%[c],%[ld],2)\n\tvmovapd %%ymm11,32(%[c],%[ld],2)\n\t.endif\n\t"
 :[k]"+r"(k),[c]"+r"(c):[a]"r"(a),[b]"r"(b),[ld]"r"(stride),[r]"i"(R)
 :"memory","cc","ymm0","ymm1","ymm2","ymm3","ymm4","ymm5","ymm6","ymm7","ymm8","ymm9","ymm10","ymm11","ymm12","ymm13","ymm14","ymm15");
}
template<int S> static void packB(Matrix b){
 for(int k=0;k<N;k+=4)for(int j=0;j<N;j+=8){
  double*d=b_panel+(j/8)*BS+k*8;
  for(int r=0;r<4;r++)for(int q=0;q<8;q+=4){
   __m256d x=_mm256_load_pd(b.x+(size_t)(k+r)*b.ld+j+q);
   if constexpr(S){__m256d y=_mm256_load_pd(b.y+(size_t)(k+r)*b.ld+j+q);x=S>0?_mm256_add_pd(x,y):_mm256_sub_pd(x,y);}
   _mm256_store_pd(d+r*8+q,x);
  }
 }
}
template<int S> static void leafA(Matrix a,double*c,int ldc){
 for(int ic=0;ic<N;ic+=30){
  const int mc=ic+30<=N?30:N-ic;
  for(int r=0;r<mc;r++)for(int k=0;k<N;k+=4){
   __m256d x=_mm256_load_pd(a.x+(size_t)(ic+r)*a.ld+k);
   if constexpr(S){__m256d y=_mm256_load_pd(a.y+(size_t)(ic+r)*a.ld+k);x=S>0?_mm256_add_pd(x,y):_mm256_sub_pd(x,y);}
   _mm256_store_pd(a_panel+r*AP+k,x);
  }
  for(int j=0;j<N;j+=8){
   if((j&8)==0){int r=0;for(;r+6<=mc;r+=6)kernel<6>(a_panel+r*AP,b_panel+(j/8)*BS,c+(size_t)(ic+r)*ldc+j,ldc);if(r<mc)kernel<4>(a_panel+r*AP,b_panel+(j/8)*BS,c+(size_t)(ic+r)*ldc+j,ldc);}
   else {int r=mc;if(r%6){r-=4;kernel<4>(a_panel+r*AP,b_panel+(j/8)*BS,c+(size_t)(ic+r)*ldc+j,ldc);}while(r){r-=6;kernel<6>(a_panel+r*AP,b_panel+(j/8)*BS,c+(size_t)(ic+r)*ldc+j,ldc);}}
  }
 }
}
static void leaf(Matrix a,Matrix b,double*c,int ldc){
 if(b.sy==0)packB<0>(b);else if(b.sy>0)packB<1>(b);else packB<-1>(b);
 if(a.sy==0)leafA<0>(a,c,ldc);else if(a.sy>0)leafA<1>(a,c,ldc);else leafA<-1>(a,c,ldc);
}
static void strassen512(Matrix a,Matrix b,double*c,int ldc){
 constexpr int h=256;constexpr size_t z=h*h;
 Matrix aa[4]={quad(a,0,0,h),quad(a,0,1,h),quad(a,1,0,h),quad(a,1,1,h)};
 Matrix bb[4]={quad(b,0,0,h),quad(b,0,1,h),quad(b,1,0,h),quad(b,1,1,h)};
 double*m1=c,*m2=c+h,*m3=c+(size_t)h*ldc,*m4=m3+h;
 double*m5=work_inner,*m6=m5+z,*m7=m6+z;
 leaf(pair(aa[0],aa[3],1),pair(bb[0],bb[3],1),m1,ldc);
 leaf(pair(aa[2],aa[3],1),bb[0],m2,ldc);
 leaf(aa[0],pair(bb[1],bb[3],-1),m3,ldc);
 leaf(aa[3],pair(bb[2],bb[0],-1),m4,ldc);
 leaf(pair(aa[0],aa[1],1),bb[3],m5,h);
 leaf(pair(aa[2],aa[0],-1),pair(bb[0],bb[1],1),m6,h);
 leaf(pair(aa[1],aa[3],-1),pair(bb[2],bb[3],1),m7,h);
 for(int i=0;i<h;i++)for(int j=0;j<h;j+=4){
  double*p=c+(size_t)i*ldc+j;size_t o=(size_t)i*h+j;
  __m256d a=_mm256_load_pd(p),b=_mm256_load_pd(p+h),d=_mm256_load_pd(p+(size_t)h*ldc),e=_mm256_load_pd(p+(size_t)h*ldc+h);
  __m256d f=_mm256_load_pd(m5+o),g=_mm256_load_pd(m6+o),k=_mm256_load_pd(m7+o);
  _mm256_store_pd(p,_mm256_add_pd(_mm256_sub_pd(_mm256_add_pd(a,e),f),k));
  _mm256_store_pd(p+h,_mm256_add_pd(d,f));
  _mm256_store_pd(p+(size_t)h*ldc,_mm256_add_pd(b,e));
  _mm256_store_pd(p+(size_t)h*ldc+h,_mm256_add_pd(_mm256_add_pd(_mm256_sub_pd(a,b),d),g));
 }
}
static void add512(double*d,Matrix a,Matrix b,int s){
 for(int i=0;i<512;i++)for(int j=0;j<512;j+=4){__m256d x=load(a,i,j),y=load(b,i,j);_mm256_store_pd(d+(size_t)i*512+j,s>0?_mm256_add_pd(x,y):_mm256_sub_pd(x,y));}
}
static void root(const double*A,const double*B,double*C){
 constexpr int h=512;constexpr size_t z=h*h;constexpr int n=1024;
 Matrix a[4]={{A,n},{A+h,n},{A+(size_t)h*n,n},{A+(size_t)h*n+h,n}};
 Matrix b[4]={{B,n},{B+h,n},{B+(size_t)h*n,n},{B+(size_t)h*n+h,n}};
 double*m1=C,*m2=C+h,*m3=C+(size_t)h*n,*m4=m3+h;
 double*m5=work_top,*m6=m5+z,*m7=m6+z,*p=m7+z,*q=p+z;
 strassen512(a[0],b[0],m1,n);
 strassen512(a[1],b[2],m2,n);
 add512(p,a[2],a[3],1);add512(q,b[1],b[0],-1);
 strassen512({p,h},{q,h},m5,h);
 add512(p,{p,h},a[0],-1);add512(q,b[3],{q,h},-1);
 strassen512({p,h},{q,h},m6,h);
 add512(p,a[1],{p,h},-1);strassen512({p,h},b[3],m3,n);
 add512(q,{q,h},b[2],-1);strassen512(a[3],{q,h},m4,n);
 add512(p,a[0],a[2],-1);add512(q,b[3],b[1],-1);
 strassen512({p,h},{q,h},m7,h);
 for(int i=0;i<h;i++)for(int j=0;j<h;j+=4){
  double*d=C+(size_t)i*n+j;size_t o=(size_t)i*h+j;
  __m256d a=_mm256_load_pd(d),b=_mm256_load_pd(d+h),c=_mm256_load_pd(d+(size_t)h*n),e=_mm256_load_pd(d+(size_t)h*n+h);
  __m256d f=_mm256_load_pd(m5+o),g=_mm256_load_pd(m6+o),k=_mm256_load_pd(m7+o);
  __m256d u2=_mm256_add_pd(a,g),u3=_mm256_add_pd(u2,k),u4=_mm256_add_pd(u2,f);
  _mm256_store_pd(d,_mm256_add_pd(a,b));_mm256_store_pd(d+h,_mm256_add_pd(u4,c));
  _mm256_store_pd(d+(size_t)h*n,_mm256_sub_pd(u3,e));_mm256_store_pd(d+(size_t)h*n+h,_mm256_add_pd(u3,f));
 }
}
}
void matrix_multiply(int n,const double*A,const double*B,double*C){compact256::root(A,B,C);}

评测结果