N^{1/5} \log N 因数分解

· · 算法·理论

尽管容易理解，没有 OI 价值吧……

前情提要

https://www.luogu.com.cn/article/9jfrrk1k Lehman 算法（1974）O(N^{1/3}) 实现

Pollard O(N^{1/4}) 算法基础上用 Lehman 算法的 O(N^{1/3}) 的算式（注意顺序）魔改的https://web.maths.unsw.edu.au/~davidharvey/talks/magma2023.pdf O(N^{1/5}) Harvey 算法。（O(N^{1/4}) 算法指的不是随机的这个Pollard-Pho）

看上面的文章时，发现很巧妙地避开了最佳逼近分数来实现，只是求解费马分解方程，因此对这一系列因数分解算法的很感兴趣。06 年做集训队作业的时候，研究 Farey 序列给出了一个 O(\log n) 求解任意实数附近最佳逼近分数的算法，这里用得上。注：（一类分数问题的研究.pdf 可能哪里找得到吧，e.g. https://blog.csdn.net/Code92007/article/details/135053808 但这里的不完整，没包括实数部分，实数到分数的方法只是先找一个平方规模的分母分子作为这个实数的代替求解，之后可以证明就是解）

Lehman 算法基本细节

[前情提要](https://www.luogu.com.cn/article/9jfrrk1k )里的初等观点，没有常数因子，避开了带来核心进展并且非平凡的最佳分数逼近，因此也失去了继续优化的空间。 这里只补充说明 [Slides](https://web.maths.unsw.edu.au/~davidharvey/talks/magma2023.pdf) 之中不显然的地方。因为懒得看论文原文，就自己推导了下。 ### 主方程、算法解释 稍稍给它改头换面，比较好从[前情提要](https://www.luogu.com.cn/article/9jfrrk1k )的理解跳过来。 $vN/x+ux = vq_0+up_0 + O(1)\ \forall x, s.t.\ |x-p_0| <= p_0/N^{1/3}, q_0 = N/p_0

意思是，考虑若干个 p_0，把 p \in [N^{1/3},N^{1/2}] 这个大区间，用 \forall x, s.t.\ |x-p_0| <= p_0/N^{1/3} 覆盖了。每个区间后端点对前端点的比值为 \frac{1+N^{1/3}}{N^{1/3}-1}，因此总区间数是\frac{\log\frac{N^{1/2}}{N^{1/3}}}{\log \frac{1+N^{1/3}}{N^{1/3}-1}} = O(N^{1/3}\log N)。

对于每个区间，令 $y = vq_0+up_0 + O(1)$ 进行枚举（只有 -1 到 +2 范围内的三个整数），求解 $p (=x) = \frac{y + \sqrt{y^2 - 4uvN}}{2u}$，隐含了根式首先是整数。跟[前情提要](https://www.luogu.com.cn/article/9jfrrk1k )中的求解方程 $x^2-y^2=4uvN$ 有同样的实质，只是方法差异很大。相当于对于每组 $(u, v)$，通过对应的 $p_0$ 直接确定了整数 $(x, x^2-4uvN)$，而不是小范围枚举 $x$。 意思是，$x^2-4uvN$ 不是在一个小范围内自由的，而是由 $p_0$ 或者说某种程度上由 $(u, v)$ 所确定的。这是因为 $x^2-y^2=4uvN$ 这个方程丢失了主方程中除法和整性所包含的值域信息。这个方程丢失信息的主要原因是只用了 $\frac{u}{v}$ 逼近在 $uv$ 小范围的存在性，而没考虑进每个 $(u,v)$ 锁定了一个 $p$ 的可行范围。 ### 主方程的详细说明 $vN/x+ux = vq_0(p_0/x)+up_0 (x/p_0) = (vq_0+up_0)+vq_0/x(p_0-x)+u(x-p_0)=(vq_0+up_0)+(x-p_0)/p_0(up_0-vq_0)+(x-p_0)^2vq_0/p_0x = (vq_0+up_0) + \frac{x-p_0}{p_0}\cdot\frac{p_0}{V}\cdot vV(\frac{u}{v}-\frac{N}{p_0^2})+\frac{(x-p_0)^2}{xp_0}N\frac{v}{V}\frac{V}{p_0}$。其中 $\frac{p_0}{V}$ 是选定的参数，$|\frac{p_0}{x-p_0}|$ 的下界是另一个参数。观察第二项和第三项，知两个参数都取 $N^{1/3}$ 时，均为 $O(1)$。 另一种参数选取，分别为 $n^{2/5}$、$n^{1/5}$ 时，均为 $O(n^{1/5})$。 这里注意一下用的是狄利克雷逼近定理来证明第二项的乘数因子 $O(1)$ 的。首先它是鸽笼原理证明的存在性定理。而正好它的 Farey 序列下一项 $a/b$ 的距离是 $1/(vb) \ge 1/(vV)$ 因此下一项跳过了被逼近的实数，因此只可能是单侧最优的，故单侧最优的逼近的差值也满足狄利克雷逼近定理的界。 至于最好能做到 $1/5$ 次方，是因为第二项是 $a-b$ 次、第三项是 $1-2b-a$ 次，而时间复杂度中，有 $b$ 次个 $p_0$、而界之中是 $a-b$ 次，决定了 DFT 同时计算的多项式值的个数。故知 $\max(a-b, b)$ 最小为 $1/5 = [2(a-b)+3b]/5 = [1+(a-b)-(1-2b-a)]/5$。 ## $O(N^{1/5})$ Harvey 算法思想 因为假设了某个 $p_0$ 确定的区间存在 $p$，知道 $ vN/p+up - (vq_0+up_0) = O(N^{1/5})$， 对 $\alpha^{vN/p+up - \lceil vq_0+up_0\rceil} = \alpha^k \mod p$ 应用费马小定理，即为 $\alpha^{vN+u - \lceil vq_0+up_0\rceil} = \alpha^k \mod p$，左边次方上面没有未知数了。 左侧有最多 $O(N^{1/5}\log N)$ 个数值，令 $f(x) = \Pi_{i(p_i, u, v, q_i\mathrm{\ s.t. }\forall k\ \alpha^{vN+u - \lceil vq_i+up_i\rceil} \ne \alpha^k \bmod N)}(x-\alpha^{vN+u-\lceil vq_i+up_i\rceil}) \bmod N)$ 是 $O(N^{1/5}\log N)$ 次多项式。知道 $f(\alpha^k)$ 存在一个 $p$ 的倍数项。要计算 $f(\alpha^{-N^{1/5}}), \ldots, f(\alpha^{2N^{1/5}})$，之后对于每一个 $f(\alpha^k)$ 与 $N$ 求 gcd。 计算多个多项式的数值就是套用离散傅里叶变换的计算过程。 其中 $\alpha$ 的条件和选取之类的具体细节是另一个课题暂不理会。 这部分是 Slides 所说的，但我看 [paper](https://arxiv.org/pdf/2010.05450) 采取的是完全不同的做法。不理它为好，Slides 的做法看起来才是正确的。 ### 因数分解的严谨性 在求出 gcd 之后，如果是 $p$，得到解了，如果是 $N$ （即 0），则逐个检查 $f(\alpha^k)$ 的因子和 $N$ 的 gcd，如果是 $p$，得到解了。如果是 $N$ 就会有麻烦。 前面知道存在 $(p_i, u, v)$，使得 $vN/p+up - \lceil vq_i+up_i\rceil = k$，而 $k = O(N^{1/5})$，又有假设，$\mathrm{ord}_N(\alpha) > O(N^{1/5})$，这里就不写精确的常数因子了。 因为可能 ，$vN+u - \lceil vq_j+up_j\rceil = k + l\cdot\mathrm{ord}_N(\alpha)$，这里 $p$ 并不在 $p_j$ 对应的区间之内。好的性质是，如果有这个情况，只可能对应一个 $k$。 对于这种情况，首先 $\{\alpha^k \bmod N\}$ 和 $\{\alpha^{vN+u-\lceil vq_i+up_i\rceil} \bmod N\}$ 都是预先计算的。先直接进行判重。知道对于 $p_j$，最多只有一个对应的 $k$，那么只需要对 $vN/p+up - \lceil vq_j+up_j\rceil = k$ 使用 $p = \frac{y + \sqrt{(k + \lceil vq_j+up_j\rceil)^2 - 4uvN}}{2u}$ 求解判定是否为整数。不论计算结果如何，这个 $(p_j, k)$ 这对关系就已经被排除掉了。我们最多排除掉 $O(N^{1/5} \log N)$ 对这样的关系，不影响总的时间复杂度。 现在考虑，如果对于 $(p_i, u, v)$，使得 $vN/p+up - \lceil vq_i+up_i\rceil = k$，并且 $vN+u - \lceil vq_i+up_i\rceil \ne k + j\cdot\mathrm{ord}_N(\alpha)$ 也就是我们希望 $p_i$ 应当可以解出 $p$ 的情况。考虑到是否存在 $k'$ 使得 $vN+u - \lceil vq_i+up_i\rceil = k' + l\cdot\mathrm{ord}_N(\alpha) = k \mod \mathrm{ord}_p(\alpha)$，如果存在，知 $\mathrm{ord}_p(\alpha) | k' - k < O(N^{1/5})$，即 $vN/p+up - \lceil vq_i+up_i\rceil = k = k' + c\cdot\mathrm{ord}_p(\alpha)$，不可以因为 $k'$ 的存在排除掉 $p_i$，同时知道 $1 < \gcd(\alpha^{k} - \alpha^{vN+u-\lceil vq_i+up_i\rceil}, N) < N$。 先从 $f(x)$ 中除去所有的存在 $\alpha^{k'} = \alpha^{vN+u-\lceil vq_i+up_i\rceil} \mod N$ 的所有 $p_i$ 有关的因子进行求解。如果找不到解，说明 $p$ 在某个这样的 $p_i$ 的区间之内。 对于这些 $\{p_i\}$，我们希望寻找 $k$ 使得 $1 < \gcd(\alpha^{k-k'} - \alpha^{vN+u-\lceil vq_i+up_i\rceil-k'}, N) < N$，即 $1 < \gcd(\alpha^{k-k'} - 1, N) < N$，这是平凡的，并且可以作为算法启动前的特判，即朴素地特判 $\mathrm{ord}_p(\alpha) < O(N^{1/5}$），如果不通过，可以得到因子 $p$。 因此算法在对 $f(x)$ 的乘法因子中排除某些区间后成立。 ### 其他细节补充 由于 $\frac{p_0}{V} = N^{2/5}$ 是算法选定的参数，隐含 $p > N^{2/5}$。如果刚刚求不出解，即落入此情况，使用 $O(N^{1/4})$ 算法先进行因数分解。注意到该算法的时间复杂度实际上是 $O(\sqrt{p}\log N) = O({N^{1/5}\log N})$。