利用现代优化与 AlphaEvolve 改进矩阵乘法指数

Improving the matrix multiplication exponent with modern optimization and AlphaEvolve

arXiv: 2608.16884v1

论文信息

标题: Improving the matrix multiplication exponent with modern optimization and AlphaEvolve

作者: Emilien Dupont, Marvin Eisenberger, Borislav Kozlovskii, et al.

发布日期: 2026-08-17

arXiv ID: 2608.16884v1

PDF 链接: 下载 PDF

3 分钟速览

  • 研究问题：本文研究矩阵乘法渐近复杂度指数 ω\omega 的上界，目标是改进此前由组合损失分析得到的 ω<2.371339\omega<2.371339。
  • 核心方法：把组合损失分析中的非凸优化问题重构为可微、可并行的形式，用 JAX 实现梯度下降，并用 Sinkhorn-Knopp 算法处理最大熵分布；随后用 AlphaEvolve 自动改进优化程序。
  • 关键结果：在 q=5q=5、最大递归层数 ℓ∗=4\ell^*=4 的设置下，经过严格有理算术验证，得到 ω<2.371177\omega<2.371177，比此前上界改进约 1.62×10−41.62\times10^{-4}。
  • 主要局限：该改进幅度较小，作者明确表示更大的改进可能需要新的数学思想；方法仍依赖数值优化和事后严格验证。
  • 适合读者：理论计算机科学、代数复杂度、矩阵乘法算法、机器学习优化以及自动算法发现方向的研究者。

论文背景和研究动机

矩阵乘法是计算机科学中最基础的运算之一。其渐近复杂度用指数 ω\omega 刻画：两个 n×nn\times n 矩阵相乘可以在 O(nω+o(1))\mathcal{O}(n^{\omega+o(1)}) 次算术运算内完成。Strassen 在 1969 年首次证明 ω<2.81\omega<2.81，此后数十年的大量工作不断压低 ω\omega 的上界。

过去约 40 年中，所有主要改进都依赖激光方法（laser method）。该方法并不直接构造快速矩阵乘法算法，而是通过张量分解与组合分析间接推出复杂度上界。当前最好结果来自激光方法的一种精化版本——组合损失分析（combination loss analysis）。该技术把一个计算机辅助证明转化为一个非凸优化问题：只要找到该问题的可行解，就能得到 ω\omega 的上界。

此前最好上界 ω<2.371339\omega<2.371339 来自 Alman、Duan、Williams、Xu、Xu、Zhou 的工作《More asymmetry yields faster matrix multiplication》。该工作使用最大递归层数 ℓ∗=3\ell^*=3，并用 SNOPT 软件包中的序列二次规划（SQP）求解。本文作者认为，该优化问题本身仍有改进空间，尤其是现代机器学习中的可微优化、自动微分和 GPU 并行技术尚未被充分引入。于是，他们尝试用梯度下降方法处理更大规模的问题实例，并借助 AlphaEvolve 进一步优化求解算法。

核心方法和技术细节

组合损失分析的核心优化问题可以用一棵带参数的树来描述。树中每个非根节点对应一个形状 (sX,sY,sZ)(s_X,s_Y,s_Z) 和一个区域 r∈[6]r\in[6]。优化变量主要是概率分布：根节点和正形状节点上有区域分布 AA 和形状分裂分布 α\alpha；零形状节点上有完整分裂分布 β\beta；第 2 层特殊节点上还有标量 μ∈[0,1/2]\mu\in[0,1/2]。

这些分布通过自顶向下的质量传播、完整分裂分布的递归构造、保留指数（retained exponent）和局部矩阵尺寸等派生量，最终汇总为目标函数与约束。具体地，优化问题为：

minimizeΩsubject toEtotal+Mtotal⋅Ω≥2ℓ∗−1log⁡(q+2),\begin{array}{ll} \text{minimize} & \Omega \\ \text{subject to} & E_{\text{total}} + M_{\text{total}}\cdot \Omega \ge 2^{\ell^*-1}\log(q+2), \end{array}

其中 EtotalE_{\text{total}} 是总保留指数，MtotalM_{\text{total}} 是总矩阵尺寸。根据论文引用的定理，任何满足该约束的可行解都意味着 ω≤Ω\omega\le\Omega。

本文没有沿用 SQP 求解器，而是将问题改造为适合 GPU 并行和自动微分的可微计算图。主要技术包括：

  1. 概率分布的无约束参数化：所有概率分布用 logits 表示，再通过 softmax 得到概率向量。这样原来的单纯形约束被自动满足，问题变为无约束优化。

  2. 最大熵分布的可微求解：派生量中出现 HDmax⁡(ρ)H_D^{\max}(\rho)，即在边缘分布固定时最大化熵的上确界。该量没有解析表达式。作者将其视为最优传输中的熵正则化问题，用 Sinkhorn-Knopp 算法迭代求解，并用隐式微分进行反向传播，以获得稳定梯度。

  3. 张量化表示与硬件并行：与先前在图上逐节点消息传递的实现不同，作者引入幻影节点（phantom nodes）和掩码，将高度不均匀的图结构嵌入到多维张量中；再将节点聚类为若干 “阶段”。这样整个计算可以利用 JAX 在 GPU 上并行执行。

  4. 提高递归层数：此前工作使用 ℓ∗=3\ell^*=3，可优化参数约 2.5 万；本文提升到 ℓ∗=4\ell^*=4，可优化参数增长到约 700 万。该提升是以双指数级复杂度增长为代价的，但梯度方法和大规模并行使其可行。

  5. AlphaEvolve 改进优化算法：AlphaEvolve 是一个编码智能体，能够修改并进化优化程序。本文让 AlphaEvolve 改写优化算法，每次运行约 5 小时（单 GPU）以输出 ω\omega 上界，并通过进化逐步降低该上界。作者使用其 “进化构造” 特性，让每代优化从父算法找到的最优解附近继续搜索。

创新点和贡献

本文的主要贡献不在于提出新的矩阵乘法理论，而在于将现代数值优化和自动算法发现引入到一个此前由传统 SQP 主导的计算机辅助证明流程中。

第一个创新是问题重构与可微求解。以前的最大熵分布和对应拉格朗日乘子被当作自由变量与其余参数联合优化；本文则把它们改为通过 Sinkhorn-Knopp 算法动态求解，并用隐式微分传递梯度。这种思路来自最优传输和机器学习，在此问题中带来了稳定性和可扩展性。

第二个创新是并行化表示。将图结构映射为高维张量，并引入幻影节点和阶段聚类，使原本难以并行的非均匀图计算能够在 GPU 上高效执行。这使得 ℓ∗=4\ell^*=4 的约 700 万参数优化成为可能。

第三个创新是 AlphaEvolve 的应用。AlphaEvolve 不仅被用来寻找更好的参数，还被用来改写求解算法本身。摘要中给出的数据表明：单独使用梯度优化，相较此前最佳上界改进约 0.97×10−40.97\times10^{-4}；加入 AlphaEvolve 后，总改进提升到约 1.62×10−41.62\times10^{-4}。

第四个贡献是严格验证。数值优化得到的浮点解不能直接作为数学上界。作者将解舍入为有理数，确保最大熵证书仍有效，再用精确有理算术重新计算所有派生量，并给每个对数计算有理数界，从而消除浮点误差。论文指出，验证代码和发现的最优解将随仓库发布。

实验结果分析

论文给出的最终结果是：

ω<2.371177.\omega < 2.371177.

这一结果来自 q=5q=5、ℓ∗=4\ell^*=4 的优化问题实例。此前最好上界为 2.3713392.371339，因此改进绝对量约为 1.62×10−41.62\times10^{-4}。

从历史角度看，作者在讨论中指出，该改进幅度与过去 40 年间多数改进大致相当。表 1 显示，从 Coppersmith–Winograd 的 2.3754772.375477 到 Williams 的 2.3728732.372873，再到后续若干工作，每次改进通常也是以 10−410^{-4} 甚至更小的量级推进。因此，2.371339→2.3711772.371339\to2.371177 虽然绝对值小，但在这个高度精细化的竞赛中仍属于可报告进展。

需要注意的是，这个数字不是某个具体算法在有限矩阵尺寸下的运行时间测量，而是渐近复杂度指数的一个严格上界。它的有效性来自两件事：一是优化问题本身的定理保证，即任意可行解都对应 ω\omega 的上界；二是本文的事后有理算术验证，确保证书没有因浮点误差而失效。

在方法层面，单用梯度优化的结果是改进约 0.97×10−40.97\times10^{-4}，而 AlphaEvolve 将改进幅度扩大到约 1.62×10−41.62\times10^{-4}。这说明自动算法搜索在这一问题上提供了额外的、不可忽略的收益。

局限与待解决问题

本文作者在讨论中明确承认，这类方法可能只会带来 “适度” 的进一步改进。更大的进步可能需要新的数学思想，而不是更精细的优化或更高的 ℓ∗\ell^*。

从方法本身看，优化问题是非凸的。梯度下降和 AlphaEvolve 都无法保证找到全局最优解；最终结果只是一个经过验证的可行解，不代表该优化问题的最小 Ω\Omega。同时，ℓ∗\ell^* 增大时，参数数量和计算量呈双指数增长。本文成功处理了 ℓ∗=4\ell^*=4，但论文未说明进一步增大 ℓ∗\ell^* 是否可行。

此外，最大熵分布的计算依赖 Sinkhorn-Knopp 迭代，并且验证步骤需要额外构造有效证书。论文中的 Lemma 1 给出了证书误差界：若近似误差为 ε\varepsilon，则真实最大熵在 [H(y),H(y)+2ε][H(y), H(y)+2\varepsilon] 区间内。这保证了严谨性，但也增加了工程复杂度。

最后，AlphaEvolve 的使用虽然带来收益，但其运行成本较高：每次执行约需单 GPU 5 小时，多次进化意味着大量计算资源。论文没有给出 AlphaEvolve 过程的整体资源消耗，因此该方法在类似问题上的可复现性和成本效益仍有待后续实践检验。