解決33問題──將33寫成3個整數的立方和

這篇文章內容翻譯自論文 Cracking the problem with 33，論文研究了方程 $x^3+y^3+z^3=k$ 在一些小的 $k$ 值的解，並首次將33寫成了3個整數的立方和。完成中文能夠查看項目 qiwihui/cracking-the-problem-with-33。截止到目前，100之內的天然數就剩下42尚未找到關於立方和的整數解了！

Answer to the Ultimate Question of Life, the Universe, and Everything. -- 42

<!--more-->

如下是論文正文翻譯：

---

解決33問題

做者：ANDREW R. BOOKER

摘要 受到Tim Browning和Brady Haran的 Numberphile 視頻"未解決的33問題"的啓發，咱們研究了方程 $x^3+y^3+z^3=k$ 在一些小的 $k$ 值的解。 咱們找到了 $k=33$ 的第一個已知解。

1. 簡介

令 $k$ 爲正整數，其中 $k \equiv ±4(\mod 9)$。 而後Heath-Brown [HB92] 推測 有無限多的三元組 $(x，y，z) \in \mathbb{Z}^3$ 知足

$$k = x^3 + y^3 + z^3. \quad \text{(1)}$$

早在1954年就開始對（1）進行各類數值研究 [MW55]；請參閱 [BPTYJ07]，瞭解截至2000年的這些研究的歷史。 自那時起進行的計算因爲Elkies [Elk00] 而被算法所主導。咱們所知道的最新內容是Huisman [Hui16] 的論文， 該論文肯定了（1）的全部解，其中 $k \le 1000$ 且 $\max\{|x|,|y|,|z|\}\le 10^15$。 特別是，Huisman報告說除了13個 $k \le 1000$ 的值之外的全部解決方案都是已知的：

$$33, 42, 114, 165, 390, 579, 627, 633, 732, 795, 906, 921, 975. \quad \text{(2)}$$

Elkies的算法經過使用格基減小（lattice basis reduction）在Fermat曲線 $X^3+Y^3=1$ 附近尋找有理點來工做；它很是適合同時找到許多 $k$ 值的解。 在本文中，咱們描述了一種在k值肯定時更有效的不一樣方法。 它的優勢是能夠找到全部具備 最小 座標界限的解，而不是Elkies算法中的最大座標。 這老是產生搜索範圍的非平凡的擴張（nontrivial expansion），由於除了能夠單獨考慮的有限多個例外以外，還有

$$\max \{|x|,|y|,|z|\} > \sqrt[3]{2} \min \{|x|,|y|,|z|\}$$

此外，根據經驗，一般狀況是其中一個變量比其餘變量小得多，所以咱們但願實際上增益更大。

咱們的策略相似於一些早期的方法（特別參見 [HBLtR93]，[Bre95]，[KTS97] 和 [BPTYJ07]）， 而且基於觀察：$k-z^3=x^3+y^3$ 的任何解都具備 $x+y$ 做爲一個因子。 相對於早期研究，咱們的主要貢獻是注意到，經過一些時間空間權衡，運行時間在高度邊界內很是接近線性， 而且在現代64位計算機上實現時很是實用。

更詳細地說，假設 $(x，y，z)$ 是（1）的解，而且不失通常性，假設 $|x| \ge |y| \ge |z|$。 而後咱們有

$$k-z^{3}=x^{3}+y^{3}=(x+y)(x^{2}-x y+y^{2})$$

若是 $k-z^3=0$ 則 $y=-x$，而且 $x$ 的每一個值都產生一個解。 不然，設 $d=|x+y|=|x|+y \operatorname{sgn} x$， 咱們看到 $d$ 能夠除 $|k-z^3|$ 而且

$$begin{aligned}
begin{aligned}
frac{left|k-z^{3}right|}{d} &=x^{2}-x y+y^{2}=x(2 x-(x+y))+y^{2} \
&=|x|(2|x|-d)+(d-|x|)^{2}=3 x^{2}-3 d|x|+d^{2}
end{aligned}
end{aligned}$$

獲得

$$\{x, y\}=\left\{\frac{1}{2} \operatorname{sgn}\left(k-z^{3}\right)\left(d \pm \sqrt{\frac{4|k-z^{3}|-d^{3}}{3 d}}\right)\right\}$$

所以，給定 $z$ 的候選值，經過遍歷 $|k-z^3|$ 的全部除數，有一個有效的程序來查找 $x$ 和 $y$ 的全部相應值。這個基本算法在假設整數分解的時間複雜度的標準啓發式（standard heuristics）下，已經能在 時間 $O(B^{1+\varepsilon})$ 內找到知足 $\min\{|x|,|y|,|z|\}\ge B$ 的全部解。 在下一節中，咱們將解釋如何避免因子分解並更有效地實現相同目的。

感謝 感謝Roger Heath-Brown提供了有用的意見和建議。

2. 方法

爲了便於表示，咱們假設 $k \equiv ±3(\mod 9)$；請注意，這適用於（2）中的全部 $k$。 因爲上述基本算法對於尋找小解是合理的，所以咱們將假設 $|z|>\sqrt{k}$。 此外，若是咱們將（1）專門用於 $y=z$ 的解，那麼咱們獲得Thue方程 $x^3+2y^3=k$，這是有效可解的。 使用 PARI/GP [The18] 中的Thue求解器，咱們驗證了（2）中的 $k$ 不存在這樣的解。 所以，咱們能夠進一步假設 $y \ne z$。

因爲 $|z|>\sqrt{k} \ge \sqrt[3]{k}$，咱們有

$$\operatorname{sgn} z=-\operatorname{sgn}(k-z^{3})=-\operatorname{sgn}(x^{3}+y^{3})=-\operatorname{sgn} x.$$

一樣，由於 $x^3 + z^3 = k-y^3$ 和 $|y|\ge |z|$， 咱們有 $\operatorname{sgn} y=-\operatorname{sgn} x=\operatorname{sgn} z$。將（1）的兩邊乘以 $-\operatorname{sgn} z$，咱們獲得

$$|x|^{3}-|y|^{3}-|z|^{3}=-k \operatorname{sgn} z \quad \text{(4)}$$

令 $\alpha=\sqrt[3]{2}-1$，而且 $d=|x+y|=|x|-|y|$。 若是 $d \ge \alpha |z|$ 則

$$begin{aligned}
begin{aligned}
-k operatorname{sgn} z &=|x|^{3}-|y|^{3}-|z|^{3} geq(|y|+alpha|z|)^{3}-|y|^{3}-|z|^{3} \
&=3 alpha(alpha+2)(|y|-|z|) z^{2}+3 alpha(|y|-|z|)^{2}|z| \
& geq 3 alpha(alpha+2)|y-z| z^{2}
end{aligned}
end{aligned}$$

因爲 $3 \alpha(\alpha+2)>1$， 這與咱們的假設不相容，即 $y \ne z$ 和 $|z|>\sqrt{k}$。 所以咱們必然有 $0<d<\alpha|z|$。

接下來，減小（4）模3並回想咱們的假設 $k \equiv ±3(\mod 9)$，咱們有

$$d=|x|-|y| \equiv|z| \quad(\mod 3).$$

設 $\epsilon\in\{±1\}$ 使得 $k \equiv 3 \epsilon(\mod 9)$。 而後，因爲每一個立方數都與 $0$ 或 $±1(mod 9)$ 相等， 咱們必然有 $x \equiv y \equiv z \equiv \epsilon(\mod 3)$， 所以 $\operatorname{sgn} z=\epsilon(\frac{|z|}{3})=\epsilon(\frac{d}{3})$。 基於（3），當且僅當 $d | z^{3}-k$ 以及 $3d(4|z^{3}-k|-d^3) = 3d(4\epsilon(\frac{d}{3})(z^{3}-k)-d^{3})$ 是平方數時， 咱們獲得（1）的解。

總之，找到（1）的全部解而且知足 $|x| \ge |y| \ge |z|>\sqrt{k}$，$y \ne z$ 和 $|z|\le B$，對於每一個與3互質的 $d\in\mathbb{Z}\cap(0,\alpha B)$，解決如下系統就足夠了：

$$begin{aligned}
begin{aligned}
&{frac{d}{sqrt[3]{2}-1}<|z| le B, quad operatorname{sgn} z=epsilonleft(frac{d}{3}right), quad z^{3} equiv k quad(mod d)} \
&{3 dleft(4 epsilonleft(frac{d}{3}right)(z^{3}-k)-d^{3}right)=square} & text{(5)}
end{aligned}
end{aligned}$$

咱們解決這個問題的方法很簡單：咱們經過它們的主要因子分解遞歸地計算 $d$ 的值， 並應用中國剩餘定理來將 $z^{3} \equiv k(\mod d)$ 的解減小到素數模冪的狀況下， 其中標準算法能夠適用。設 $r_{d}(k)=\# \left\{z(\mod d):z^{3} \equiv k(\mod d)\right\}$ 表示 $k$ 模 $d$ 的立方根數。經過標準分析估計，因爲 $k$ 不是立方數，咱們有

$$\sum_{d \le \alpha B} r_{d}(k) \ll_{k} B$$

啓發式地，計算對全部素數 $p\le \alpha B$ 的 $z^{3} \equiv k(\mod p)$ 的解 能夠用 $[0, \alpha B]$ 上的整數在 $O(B)$ 算術運算來完成； 見例如 [[NZM91]，§2.9，練習8]中描述的算法。假設這一點，能夠看出， 使用Montgomery的批量反轉技巧[[Mon87]，§10.3.1]，計算對全部正整數 $p\le \alpha B$ 的 $z^{3} \equiv k(\mod p)$ 的根的剩餘工做能夠再次用 $O(B)$ 算術運算完成。

所以，咱們能夠在線性時間內計算知足（5）的第一行的全部 $z$， 做爲算術進展（arithmetic progressions）的並集。爲了檢測最後一行的解，有一個快速的方法來肯定 $\Delta :=3d\left(4\epsilon(\frac{d}{3})(z^{3}-k)-d^{3}\right)$ 是一個平方數 相當重要。咱們首先注意到對於固定 $d$，這種狀況減小到在橢圓曲線上找到積分點； 特別是，令 $X=12d|z|$ 和 $Y=(6d^2|x-y|$，從（3）中咱們看到（X，Y）位於Mordell曲線上

$$Y^{2}=X^{3}-2(6 d)^{3}\left(d^{3}+4 \epsilon\left(\frac{d}{3}\right) k\right). \quad \text{(6)}$$

所以，對於固定 $d$，存在至多有限多個解，而且它們能夠被有效地約束。 對於 $d$ 的一些小值，找到（6）上的全部積分點並檢查是否產生任何知足（1）的解是切實可行的。 例如，使用Magma[[BCFS18]，§128.2.8]中的積分點函數（functionality），咱們驗證瞭如（2）中的 $k$ 和 $d \le 40$ 狀況下沒有解， 除了 $(k, d)\in\{(579,29),(579,34),(975,22)\}$。

接下來咱們天然注意到一些同餘和可分性約束：

引理 設 $z$ 爲（5）的解，設 $p$ 爲素數， 設 $s=ord_p d$，$t=ord_p(z^3-k)$。則

(i) $z \equiv \frac{4}{3} k\left(2-d^{2}\right)+9(k+d)(\mod 18)$；
(ii) 若是 $p \equiv 2 (\mod 3)$ 則 $t \le 3s$；
(iii) 若是 $t \le 3s$ 則 $s \equiv t (\mod 2)$；
(iv) 若是 $ord_p k \in \{1,2\}$ 則 $s \in \{0,ord_p k\}$。

證實 令 $\Delta=3d\left(4\epsilon(\frac{d}{3})(z^3-k)-d^3\right)$， 令 $\delta=(\frac{d}{3})$，咱們有 $|z| \equiv d \equiv \delta(\mod 3)$， 觀察到 $(\delta+3 n)^{3} \equiv \delta+9 n(\mod 27)$，模27，咱們有

$$begin{aligned}
begin{aligned}
frac{Delta}{3 d} &=4 epsilon deltaleft(z^{3}-kright)-d^{3}=4|z|^{3}-d^{3}-4 epsilon delta k \
& equiv 4[delta+3(|z|-delta)]-[delta+3(d-delta)]-4 epsilon delta k=3(4|z|-d)-delta[18+4(epsilon k-3)] \
& equiv 3(4|z|-d)-d[18+4(epsilon k-3)]=12|z|-9 d-4 epsilon d k \
& equiv 3|z|-4 epsilon d k
end{aligned}
end{aligned}$$

這消失了模9，因此爲了使 $\Delta$ 成爲平方數，它也必須消除mod 27。 因而

$$z=\epsilon \delta|z| \equiv \frac{4 \delta d k}{3} \equiv \frac{4(2-d^{2}) k}{3} \quad(\mod 9)$$

減小（1）模2咱們獲得 $z \equiv k+d(\mod 2)$，這獲得（i）。

接下來設 $u=p^{-s} d$ 和 $v=p^{-t} \epsilon \delta(z^{3}-k)$，這樣就有

$$\Delta=3\left(4 p^{s+t} u v-p^{4 s} u^{4}\right)$$

若是 $3s<t$ 則 $p^{-4 s} \Delta \equiv-3 u^{4}(\mod 4 p)$， 可是當 $p \equiv 2(\mod 3)$ 時這是不可能的，由於 $-3$ 不是 $4p$ 的平方模。所以，在這種狀況下咱們必須 $t<3s$。

接下來假設 $t<3s$。 咱們考慮如下狀況，涵蓋全部可能性：

• 若 $p = 3$ 則 $s = t = 0$，那麼 $s \equiv t(\bmod 2)$。
• 若 $p \ne 3$ 且 $3s > t+2 \operatorname{ord}_{p} 2$， 則 $\operatorname{ord}_{p} \Delta=s+t+2 \operatorname{ord}_{p} 2$，那麼 $s \equiv t(\mod 2)$。
• 若 $3s\in\{t, t+2\}$ 則 $s \equiv t(\bmod 2)$。
• 若是 $p=2$ 且 $3s = t + 1$ 則 $2^{-4 s} \Delta=3(2 u v-u^{4}) \equiv 3(\bmod 4)$，這是不可能的。

所以，在任何狀況咱們得出結論 $s \equiv t(\mod 2)$。

最後，假設 $p|k$ 和 $p \not | 3k$。若是 $s=0$ 則無需證實的，因此假設否則。 因爲 $d | z^{3}-k$，咱們必須有 $d | k$，由於

$$0 < s \le t=\operatorname{ord}_{p}(z^{3}-k)=\operatorname{ord}_{p} k<3 s$$

經過部分（iii）得出 $s \equiv \operatorname{ord}_{p} k(\mod 2)$， 所以 $s=\operatorname{ord}_{p} k$。

所以，一旦 $z(\mod d)$ 的殘差類（residue class）固定， 則其殘差模 $lcm(d,18)$ 是肯定的。還要注意，條件（ii）和（iii）對於測試 $p=2$ 是有效的。

然而，即便有這些優化，也有 $\ll B\log B$ 對 $d, z$ 知足（5）的第一行和引理的結論（i）和（iv）。 所以，爲了實現比 $O(B\log B)$ 更好的運行時間，須要從一開始就消除一些 $z$ 值。 咱們經過標準的時間空間交換來實現這一目標。確切地說，設置 $P=3(\log \log B)(\log \log \log B)$， 而且讓 $M=\prod_{5 \le p \le P} p$ 是區間 $[5, P]$ 之間的素數的乘積。 根據素數定理，咱們獲得 $\log M=(1+o(1)) P$。若是 $\Delta$ 是平方數， 那麼對於任意素數 $p|M$ 咱們有

$$\left(\frac{\Delta}{p}\right)=\left(\frac{3 d}{p}\right)\left(\frac{|z|^{3}-c}{p}\right) \in\{0,1\} \quad \text{(7)}$$

其中 $c \equiv \epsilon\left(\frac{d}{3}\right) k+\frac{d^{3}}{4}$。 當 $\operatorname{lcm}(d, 18) \le \alpha B / M$ 時， 咱們首先爲每一個殘差類 $|z|(\bmod M)$ 計算該函數， 而且僅選擇對於每一個 $p|M$ 知足（7）的那些殘基。 由Hasse約束，容許的殘差的數量最多爲

$$\frac{M}{2^{\omega(M /(M, d))}} \prod_{p | \frac{M}{(M, d)}}\left(1+O\left(\frac{1}{\sqrt{p}}\right)\right)=\frac{M}{2^{\omega(M /(M, d))}} e^{O(\sqrt{P} / \log P)}$$

所以，要考慮的 $z$ 值的總數最多爲

$$begin{aligned}
begin{array}{l}{
sum_{operatorname{lcm}(d, 18) le frac{alpha B}{M}} r_{d}(k)left[M+frac{e^{O(sqrt{P} / log P)}}{2^{omega(M /(M, d))}} frac{alpha B}{d}right] +sum_{d le alpha B, {lcm}(d, 18) le frac{alpha B}{M}} frac{r_{d}(k) alpha B}{d}} \
{ll_{k} B log M+frac{e^{O(sqrt{P} / log P)}}{2^{omega(M)}} sum_{g | M} frac{2^{omega(g)} r_{g}(k)}{g} sum_{d^{prime} le frac{alpha B}{9 g M}} frac{r_{d^{prime}}(k) alpha B}{d^{prime}}} \
{ll_{k} B log M+B log B frac{e^{O(sqrt{P} / log P)}}{2^{omega(M)}} prod_{p | M}left(1+frac{2 r_{p}(k)}{p}right)} \
{ll B P+frac{B log B}{2^{(1+o(1)) P / log P}} ll B(log log B)(log log log B)
}end{array}
end{aligned}$$

對於沒有以這種方式消除的 $z$，咱們遵循相似的策略， 其中一些其餘輔助模 $M^{\prime}$ 由較大的素數組成，以加速平方測試。 咱們預先計算模爲 $M^{\prime}$ 的立方數表和Legendre符號模 $p|M^{\prime}$， 所以將測試（7）簡化爲了表查找。只有當全部這些測試都經過時， 咱們才能在多精度算術中計算 $\Delta$ 並應用通常的平方檢驗，這種狀況對於一小部分候選值來講都是如此。 事實上，咱們指望Legendre測試的數量平均有限，因此總的來講， 找到全部解決方案的 $|z| \le B$ 應該要求不超過 $O_k(B(\log \log B)(\log \log \log B))$ 次表查找和對 $[0, B]$ 中整數的算術運算。

所以，當 $B$ 符合機器字大小時，咱們預計運行時間幾乎是線性的，這就是咱們在實踐中觀察到的 $B<2^{64}$。

3. 實現

咱們在C中實現了上述算法，其中有一些內聯彙編程序來源於由Ben Buhrow [Buh19] 編寫的Montgomery算法 [Mon85]， 以及Kim Walisch的用於枚舉素數的 primesieve 庫 [Wal19]。

該算法天然地在具備超過 $\sqrt{\alpha B}$ 的素因子和 具備 $\sqrt{\alpha B}$ -平滑的素數的 $d$ 的值之間分配。 前一組 $d$ 消耗超過運行時間的三分之二，但更容易並行化。 咱們在布里斯托大學高級計算研究中心的大規模並行集羣Bluecrystal Phase 3上運行了這一部分。 對於平滑的 $d$，咱們使用了一個單獨的32核和64核節點的小集羣。

咱們搜索了知足 $k \in \{33,42\}$ 和 $\min\{|x|, |y|, |z|\} \le 10^16$ 的（1）的解，找到了如下結果：

$$33 = 8 866 128 975 287 528^3 +（-8 778 405 442 862 239)^3 +（-2 736 111 468 807 040)^3$$

總計算在三個星期的實際時間中大約使用了15個核年。

參考文獻

（略）

School of Mathematics, University of Bristol, University Walk, Bristol, BS8 1TW, United Kingdom

E-mail address: <andrew.booker@bristol.ac.uk>

---

博客參考：

• 人類第一次將33寫成了3個整數的立方和