Number Theory
202404
0 前言
离散数学是本书的重点，而整数又是离散数学的中心议题。数论 是讨论整数性质的重要数学分支，因此我们要来探索它。
——《具体数学》第四章
标有 * 的为拓展内容，或者说比较难的题目，但它们都非常有趣。
部分题目的代码是洛谷的提交记录，阅读这篇文章的大多数同学都该能看（若有需要可以找我要）。
没看懂的地方直接问，不排除我会存在手误。
1 基本概念
1.1 整除
两个整数相除，我们自然想去讨论商是否是一个整数。小学生常常期待这个，因为整数比分数写起来好看。
这引出了数论中的基础概念：整除。
如果 \(m \gt 0\) 且 \(n /m\) 是一个整数，我们就说 \(m\) 整除 \(n\)（或者 \(n\) 被 \(m\) 整除）。这个性质奠定了整个数论的基础，所以我们赋予它一个特殊的记号会更方便，记为
同样的，如果 \(m\) 不整除 \(n\)，我们就写成 \(m \nmid n\)。
1.2 整值函数
整数是离散数学的支柱，我们常常需要将分数或者任意的实数转换到整数。 ——《具体数学》第三章
我们首先来讨论 底（floor，最大整数）函数和 顶（ceiling，最小整数）函数，对于所有实数 \(x\)，其定义如下：
俗称 下取整 和 上取整，底和顶也会在数论中有许多应用，在后续的部分内容中我们也会用到。
\(x\) 和 \(\lfloor x \rfloor\) 之间的差称之为 \(x\) 的 分数部分，它在应用中经常出现，所以也值得拥有自己的记号：
我们有时称 \(\lfloor x \rfloor\) 是 \(x\) 的 整数部分，因为 \(x= \lfloor x \rfloor + \{ x\}\)，这也是一个实数的表示方法。
尝试在这里就证明一个非常有趣且有用的事实：如果 \(m\) 和 \(n\) 是整数且分母 \(n\) 为正，则
设 \(f(x)\) 是任意一个具有如下性质且在一个实数区间连续的单调递增函数
于是只要 \(f(x),f(\left \lceil x \right \rceil),f(\left \lfloor x \right \rfloor)\) 有定义，我们就有
画出图来容易发现这是显然的，这里不做证明。
我们发现其实 \(f(x)=\frac {x+m} n\) 也具有上面要求的 \(f(x)\) 的性质，于是我们就证得了上面的结论。
1.3 带余除法
当 \(m\) 和 \(n\) 是正整数时，\(n\) 被 \(m\) 除的商是 \(\lfloor n/m \rfloor\)。而大多数情况我们都不能做到 \(m\) 整除 \(n\)，所以这个除法的余数也有一个简单的记号，很方便，我们称它是 \(n \bmod m\)。基本公式
告诉我们，可以将 \(n \bmod m\) 表示成 \(n - m \lfloor n/m \rfloor\)，我们可以将它推广到任意实数：
\(y=0\) 呢？为了避免用 \(0\) 作除数，也为了完整起见，我们可以定义：
这里的 \(\text{mod}\) 就被定义成了一个二元运算，也就是 取模。容易发现 \(x \bmod y\) 是在 \([0,y-1]\) 中的一个值。
我们约定 \(\text{mod}\) 比加法或者减法的优先级高。同时，分配律 是 \(\text{mod}\) 最重要的代数性质，对于所有实数 \(c,x,y\)，我们有
容易从定义证明这个法则，因为如果 \(cy \neq 0\)，则有
且模数为零时该式显然也为真（两边都是 \(cx\)）。
关于 \(\text{mod}\)，还有一些非常显然的性质。
- \((a+b) \bmod p = ((a \bmod p) + (b \bmod p)) \bmod p\)。
- \((ab) \bmod p= ((a \bmod p)(b \bmod p)) \bmod p\)。
- \(a \bmod bk \bmod b = a \bmod b\)。
前两条定义了 \(\text{mod}\) 的加法和乘法，第三条则是关于模数的一个转化。
这些性质都是耳熟能详的，这里也不做证明。
1.4 同余关系
模运算是数论提供的一种重要工具，我们在上面把它当成二元运算使用，而在这里我们也可以把 \(\text{mod}\) 运用到整个方程上面，为此，我们使用稍不同的记号会更加方便：
由于 \(x \bmod m\) 与 \(x\) 相差 \(m\) 的倍数，因而我们可以用另一种方式来解读同余式：
同余符号 \(\equiv\) 看起来很像 \(=\)，因为同余式与方程非常相像。例如，我们将同余式元素相加减，仍保持同余关系：
乘法同样有效，只要处理的对象是整数：
证明直接作差：\(ac-bd=(a-b)c+b(c-d)\)。反复利用乘法性质，我们可以取幂：
这样一来，我们对方程所习惯做的大多数代数运算对同余式都可以运用，但并不是所有运算。除法运算显然不在其中。
以下部分可能会用到之后要讲的内容。
如果 \(ad \equiv bd \pmod m\)，我们不能总断言 \(a \equiv b \pmod m\)。例如 \(3 \times 2 \equiv 5 \times 2 \pmod 4\)，但 \(3 \not \equiv 5 \pmod 4\)。
然而在 \(d\) 与 \(m\) 互素的情形中，我们可以挽救这一消元的性质：
它的证明是直接在式子两边乘上 \(d\) 的逆元，进而，我们可以把这个东西推广到更为一般的法则，它尽可能小地改变模：
进一步观察改变模的想法，我们可以得到另外一些式子：
反过来，如果我们知道对于两个小的模数有 \(a \equiv b\)，是否能断定对于一个更大的模数有 \(a \equiv b\) 呢？这个规则是
因为如果 \(a-b\) 是 \(m\) 和 \(n\) 的倍数，就一定是 \(\text{lcm}(n,m)\) 的倍数。
注意关注这个式子的特殊情形 \(m \bot n\)，这在之后的 中国剩余定理 颇有应用。
2 欧几里得
有着我们熟悉的 \(\gcd(a,b)\) 和 \(\operatorname{lcm}(a,b)\)。
2.1 欧几里得算法（gcd）
欧几里得算法，又称辗转相除法，常用递归得到最大公因数。
而这是容易简单证明的：
设 \(a=kb+r\)，\(a,b\) 的最大公因数为 \(d\)，那么一定满足 \(d|a,d|b\)。
于是 \(r= a-kb\)，我们将两边同时 \(\div d\) 就可以得到 \(\frac r d =\frac a d - \frac{kb} d\)。
由于 \(d|a,d|b\)，所以上面式子的右边一定是整数，于是 \(\frac r d\) 也是一个整数。
这样一来，\(d\) 也是 \(a \bmod b = r\) 的因数，得证。
CF1806F2 GCD Master (hard version)
给定 \(n, m\) 和一个长度为 \(n\) 的序列 \(\{a_i\} (a_i\leq m)\)。
定义一次对一个长度为 \(m\) 的序列的操作为，选择序列中两个下标 \(1\leq i < j \leq m\)，删去 \(a_i\) 与 \(a_j\)，然后在序列末端加入 \(\gcd(a_i, a_j)\)。
例如，对于 \([7, 6, 2]\)，一次操作可以选择下标 \(2\) 与 \(3\)，这样操作后，序列变成 \([7, 2]\)。
给定 \(k\)，求对序列 \(\{a_i\}\) 执行 \(k\) 次操作后得到序列中的数的和的最大值。
\(1 \le K \lt n \le 10^6,1 \le m \le 10^{18}\)。
感谢🐨的供题。
首先我们考虑有相同的数如何处理？
发现我们对相同的数操作，相当于删去一个，进而在后面的讨论中我们发现其实相同的数并无意义，所以我们最后枚举删了几个相同的数就可以了，于是我们就只需要对不相同的数进行操作即可。
考虑每一次操作，删除两个数 \(x,y\) 并加入 \(\gcd(x,y)\)，而当我们再把 \(\gcd(x,y)\) 和 \(z\) 操作时就等价于删除 \(x,y,z\) 并加入 \(\gcd(x,y,z)\)。
那么这个问题就被转化成把数分成若干组，删去这一组并加上它们的 \(\gcd\)，对于一组 \(k\) 个元素，它们需要的操作次数是 \(k-1\) 次。
我们如何让分组最优？不难猜到把所有的分成一组。
假设现在我们存在两组 \(S\) 和 \(T\)，它们的 \(\gcd\) 分别是 \(x,y\) 满足 \(x \ge y\)，那么此时的贡献是
我们尝试合并两组，也就是把 \(S\) 中最大的元素拿走，\(\gcd\) 变成 \(\gcd(x,y)\)，由于数都是不同的，所以 \(S\) 中的 \(mx \ge 2x\)，那么现在的贡献就是
由于 \(mx \ge 2x\)，所以
一定是成立的，所以合并两组一定更优，于是我们就得到了这一结论。
这样我们就可以取枚举 \(\gcd\)，选出其和最小的集合，再枚举一下删了多少个相同的数，直接减去它们的和就可以在 \(\mathcal O(m \log m)\) 的时间复杂度解决问题。（我们发现把相同的数加到删除的集合和直接删除是本质相同的，所以若删去 \(p\) 个相同的数，我们只需要把前 \(p\) 小的提出来删掉，具体可以见代码的 \(b\) 数组）
但在这道题中，依次枚举 \(\gcd\) 很明显是不能行的，那么我们希望有一个方法去调优，假设当前选了 \(b_1,b_2,\cdots ,b_t\) 这样的 \(t\) 个数，按照升序排序，当前的 \(\gcd\) 为 \(d\)，我们考虑替换掉一个数，是否可以使得答案更优？
假设我们现在加入一个没有选的最小的数 \(a_x\)，用它来换掉 \(b_t\)，那么贡献就是
因为元素互不相同，所以 \(b_t \ge b_{t-1} + d\)，如果 \(a_x \lt b_{t-1}\)，那么我们就有
于是这样的调整一定是更优的。
这就意味这什么呢？
假设我们选了 \(t\) 个数，那么有原数组的前 \(t-1\) 小值。
也就是说我们只需要枚举选了多少个数和最大的数选的是什么就可以了。但这样还是会 T。
考虑把 \(a\) 排序之后做一个前缀的 \(\gcd\)，容易知道它一定是不断递减且成阶梯状的，段数是 \(\mathcal O(\log m)\) 级别的，所以你只需要枚举当前的前缀 \(\gcd\) 是那个段，每一段暴力处理一遍，时间复杂度就是 \(\mathcal O(n \log m)\) 的。代码。
双倍经验：CF1806F1 GCD Master (easy version)。
2.2 扩展欧几里得算法（exgcd）
扩展欧几里得算法，即 exgcd，旨在求不定方程 \(ax+by=\gcd(a,b)\) 的一组整数解 \(x,y\)。
关于一个不定方程 \(ax+by=c\) 是否有解，我们需要证明一个引理：裴蜀定理。
如果 \(a,b\) 均为整数，则有整数 \(x,y\) 满足 \(ax+by=\gcd(a,b)\)。
证明是简单的，我们设 \(d=ax+by\) 容易发现 \(d\) 一定是 \(\gcd(a,b)\) 的倍数，那么最小的 \(d\) 就是 \(\gcd(a,b)\)。
容易发现，如果 \(ax+by=c\) 中 \(c\) 不是 \(\gcd(a,b)\) 的倍数，这个不定方程是无解的，反之一定有解。
同时，根据 裴蜀定理 我们可以得到一个推论：
两个整数 \(a,b\) 互质，当且仅当 \(ax+by=1\) 存在解 \(x,y\)。
那如何求这组方程 \(ax+by=\gcd(a,b)\) 的解呢？
我们考虑欧几里得算法的递归过程 \(\gcd(a,b)=\gcd(b,a \bmod b)\)。
当递归结束时 \(b=0\)，我们返回 \(a\)，此时有 \(ax+by=\gcd(a,b) \to x=1,y=0\)。
所以我们可以尝试先递归去求 \(bx + (a \bmod b)y=\gcd(b,a\bmod b)\)。
那么
于是这就可以在欧几里得算法中解决了。
int exgcd(int a,int b,int &x,int &y){
if(b==0) return x=1,y=0,a;
int d=exgcd(b,a%b,y,x);
y-=a/b*x;return d;
}
那么我们如何从一组解推广到更多组解呢？（这样可以方便找到一组最小正整数解）
我们用 exgcd 得到了一组解 \(x,y\)，满足 \(ax+by=\gcd(a,b)\)。
对于 \(x,y\) 的两边我们分别 \(-kb,+ka\)，容易发现 \(x-kb,y+ka\) 也是一组合法解，进而，我们将 \(ka,kb\) 分别 \(\div \gcd(a,b)\) 得到的 \(-kb/\gcd(a,b),+ka/\gcd(a,b)\) 依然是一组合法解。
所以 \(x\) 的最小正整数解是 \((x\%m+m)\%m\) 其中 \(m=b/\gcd(a,b)\)。
通过这样的调整，我们就可以找到想要的解了。
而回到最初，如果 \(\gcd(a,b)|c\)，我们把最后的答案 \(x,y\) 直接乘上 \(\frac c {\gcd(a,b)}\) 即可。
P1082 [NOIP2012 提高组] 同余方程
求关于 \(x\) 的同余方程 \(ax \equiv 1 \pmod b\) 的最小正整数解。
整理一下方程，即将同余拆开，我们可以得到
那么下面这个东西是直接可以用 exgcd 解决的，于是你就做完啦！
那如何找到最小整数解呢？
根据上面的转化，其实就是找到最小的 \(x-kb\)，于是直接用 \(x\) 对 \(b\) 取模即可。代码。
然后你就会求逆元了！
ABC340F S = 1
给你整数 \(X\) 和 \(Y\) ，它们至少满足 \(X \neq 0\) 和 \(Y \neq 0\) 中的一个。
请找出一对满足以下所有条件的整数 \((A, B)\) 。如果不存在这样的一对，请报告。
- \(-10^{18} \leq A, B \leq 10^{18}\)
- 在 \(xy\) 平面上，顶点为 \((0, 0), (X, Y), (A, B)\) 的三角形的面积为 \(1\) 。
\(-10^{17} \le X,Y \le 10^{17}\)。
机翻的，还能看。一场非常近的 ABC，有些人在场上。
考虑如何计算这三个点组成的三角形面积呢？
用四边形面积直接减去三个角上的三角形，于是可以得到
而由于我们只考虑了一种象限的情况，所以 \(S= |\frac {ay-bx} 2|\)。
题目要求 \(|\frac {ay-bx} 2 | =1\)，则 \(|ay-bx| = 2\)，于是可以直接用 exgcd 得到一组特解了。Hanghang 的代码。代码。
2.3 类欧几里得算法
设
其中 \(a,b,c,n\) 是常数，我们需要一个 \(O(\log n)\) 的算法。
这个式子和《具体数学》3.5 中的例三相当像啊，也就是说当 \(n=c\) 时，这个式子是具有封闭形式的。
但是先不急，我们先把类欧几里得算法推完再来讨论这个问题。
如果 \(a \ge c\) 或者 \(b \ge c\)，那么我们可以通过取模以减小范围。
我们就把 \(a \ge c\) 或者 \(b \ge c\) 的情况转化成了 \(a \lt c,b \lt c\) 的情况了！
接下来，用一些处理和式和处理底和顶的技巧
设 \(m= \left \lfloor \frac {an+b} c \right \rfloor\)，则
这就是一个递归的式子了，交换了 \(a,c\) 的位置，于是又可以重复上次的过程：先取模，再递归。
这也就是欧几里得算法的辗转相除过程，于是类欧几里得名字就是这样来的。
接下来我们对于类欧几里得算法进行一个推广，尝试去求另外两个和式：
首先来推导 \(g\)，用上面类似的方法，先来考虑 \(a \ge c\) 或者 \(b \ge c\) 的情况。
接下来考虑 \(a \lt c, b \lt c\) 的情况，和上面 \(f\) 较为类似.
设 \(t=\left \lfloor \frac {jc+c-b-1} a \right \rfloor,m=\left \lfloor \frac {an+b} c \right \rfloor\)，那么
接着，我们推导 \(h\)，方法也是类似的。
考虑 \(a \ge c\) 或者 \(b \ge c\) 的情况。
考虑 \(a\lt c,b \lt c\)，但是平方并不好处理，于是我们对其进行一个转化
那么
还是设 \(m=\left \lfloor \frac {an+b} c\right \rfloor,t=\left \lfloor \frac {jc+c-b-1} a \right \rfloor\)，则对前面的进行求和
三个函数一起进行递归，于是你就可以通过 P5170 【模板】类欧几里得算法 啦！代码。
const ll mod=998244353,i2=499122177,i6=166374059;
ll T,a,b,c,n;
struct node{
node () {f=g=h=0;}
ll f,g,h;
}ans;
node calc(ll a,ll b,ll c,ll n){
ll ac=a/c,bc=b/c,m=(a*n+b)/c,in=n*(n+1)%mod*(2*n+1)%mod*i6%mod;
node res;
if(a==0){
res.f=bc*(n+1)%mod;
res.g=bc*n%mod*(n+1)%mod*i2%mod;
res.h=bc*bc%mod*(n+1)%mod;
return res;
}
if(a>=c||b>=c){
res.f=n*(n+1)%mod*i2%mod*ac%mod+(n+1)*bc%mod;res.f%=mod;
res.g=in*ac%mod+n*(n+1)%mod*i2%mod*bc%mod;res.g%=mod;
res.h=ac*ac%mod*in%mod+bc*bc%mod*(n+1)%mod+ac*bc%mod*n%mod*(n+1)%mod;res.h%=mod;
node e=calc(a%c,b%c,c,n);
res.h+=e.h+2*bc%mod*e.f%mod+2*ac%mod*e.g%mod;
res.g+=e.g;res.f+=e.f;
res.h%=mod;res.f%=mod;res.g%=mod;
return res;
}
node e=calc(c,c-b-1,a,m-1);
res.f=(n*m%mod-e.f+mod)%mod;
res.g=(m*n%mod*(n+1)%mod-e.h-e.f+mod+mod)%mod*i2%mod;
res.h=(n*m%mod*(m+1)%mod-2*e.g%mod-2*e.f%mod-res.f%mod+3*mod)%mod;
return res;
}
P5179 Fraction
给你四个正整数 \(a,b,c,d\)，求一个最简分数 \(\frac p q\) 满足 \(\frac a b \lt \frac p q \lt \frac c d\)。
若有多组解，输出 \(q\) 最小的一组，若仍有多组解，输出 \(p\) 最小的一组。
\(a,b,c,d \le 10^9\)。
一道用到类欧类似思想的题目。
什么时候我们可以不加思考地直接算出这个目标分数？
发现当 \(\frac a b \lt 1,\frac c d \gt 1\)，我们直接返回 \(\frac 1 1\) 就可以了，那么现在我们考虑如何进行转化。
举一个例子：
对其取倒数，可以得到
把两边 \(\ge 1\) 的部分换掉
我们把这两边的分数变成更小的真分数，继续循环：
从而一步一步带回即可。代码。
这道题中我们把假分数转成真分数的操作和类欧类似。
ABC283Ex Popcount Sum
\(T\) 组询问，每组求 \(1 \sim n\) 内， \(\bmod m\) 余 \(r\) 的数在二进制下 \(1\) 的个数的和。
\(1 \le T \le 10^5,1 \le m \le n \le 10^9\)。
对 \(x\) 在二进制下 \(1\) 的个数进行简单转化，代码中常用 (x>>i)&1 那么转化成数学表达式就是
而由于我们只统计 \(mj+r\)，于是整个式子其实就是
这个东西就是容易用类欧几里得算法在 \(\mathcal O(T \log^2 n)\) 的时间复杂度内求出了。代码。
*当 \(n=c-1\) 时的封闭形式
我们来尝试求
的封闭形式，也就是 \(f\) 中 \(n=c\) 的情况。这是有趣的，可以教给我们一些其他的东西。
前置知识：来自《具体数学》3.25 和 3.26。
\[n= \left \lfloor \frac n m \right \rfloor + \left \lfloor \frac {n+1} m \right \rfloor + \cdots + \left \lfloor \frac{n+m-1} m \right \rfloor \]如何证明这个东西？
你可以看成把 \(n\) 个物品分成 \(m\) 组，那么一定有一些组是 \(\left \lfloor \frac {n} m\right \rfloor\) 个物品，另外一些就是 \(\left \lceil \frac n m \right \rceil\) 个。
而从小到大排好序后，第 \(i\) 组刚好有 \(\left \lfloor \frac {n+i-1} m \right \rfloor\) 个元素。
同时，我们还可以对这个东西进行推广，用 \(\left \lfloor mx \right \rfloor\) 替换 \(n\)，就可以得到一个对所有实数 \(x\) 都成立的恒等式：
\[\left \lfloor mx \right \rfloor=\left \lfloor x \right \rfloor + \left \lfloor x + \frac 1 m \right \rfloor + \cdots + \left \lfloor x + \frac{m-1}m \right \rfloor \]
首先可以考虑打表，或许你可以发现一些东西，具体可以去看《具体数学》3.5。
而推导时我们可以用刚才推导类欧几里得类似的方法：对 \(nk\) 取模。
把它拆成整数部分和小数部分就是
考虑把三项分开来求和，对于第一项和第三项，我们都涉及到了 \(kn \bmod m\)。
那么证明处理这个神秘东西呢？
容易发现这是简单的，设 \(d=\gcd(n,m)\)，根据上面同余关系中的一些证明，则
于是它所得到的就是把按照某种次序排列的数 \(0,d,2d,\cdots ,m-d\) 重复 \(d\) 次，在后面我们会更为细致地讨论它。
所以这个式子的第一部分的和就是
最后一步就用到了上面前置知识中的东西。
而对于第三部分，我们发现它其实就是一个等差数列重复了 \(d\) 次，所以它的和就是
中间的第二部分，也是非常简单的，同样是一个等差数列，它的和就是
于是我们要求的封闭形式就找出来了
如果对封闭形式稍作处理，就可以让它关于 \(m,n\) 对称，变成
就得到这样的互反律
类欧几里得算法，完。
3 同余关系
在前置知识中提到了同余关系的相关法则，那么在这里我们希望对它的理解更为深刻一些。
3.1 费马小定理
费马小定理：
我们知道，\(p-1\) 个数 \(n \bmod p,2n \bmod p,\cdots,(p-1)n \bmod p\) 就是数 \(1,2,\cdots,p-1\)（按照某种次序排列），于是我们把他们乘在一起就是
这个式子是在 \(\bmod p\) 的意义下成立的，所以我们可以得到
由于 \((p-1)!\) 和 \(p\) 是互素的，所以两边可以同时 \(\div (p-1)!\)，证毕。
同时，费马小定理还有一个更通用的表达方式
3.2 逆元
如果一个线性同余方程 \(ax \equiv 1 \pmod b\) 有解，则乘 \(x\) 是 \(a \bmod b\) 的逆元，记作 \(a^{-1}\)。
容易发现，这个东西可以用上面我们讨论过的 exgcd 直接求解，同时当 \(b\) 是质数时，我们也可以直接用 费马小定理 直接快速幂解决。
这样两种方法都可以做到单次求逆元 \(O(\log n)\) 的时间复杂度。
而根据前面我们对 exgcd 是否有解的讨论，容易发现 \(a\) 存在逆元的充要条件是 \(a \bot b\)。
可是有些时候每次都单独求逆元并不能满足我们所需要的时间复杂度，所以这里就有了另外两种求多个数逆元的方法：
-
线性求逆元：
首先我们知道 \(1^{-1}=1 \pmod p\)，我们尝试线性推出 \(i^{-1}\)。
那么令 \(k=\left \lfloor \frac p i \right \rfloor,j=p \bmod i\)，于是 \(p=ki+j\) 有 \(ki+j \equiv 0 \pmod p\)。
我们将式子两边同时乘上 \(i^{-1}j^{-1}\) 可以得到：
\[\begin{aligned} kj^{-1} + i^{-1} & \equiv 0 & \pmod p \\ i^{-1} & \equiv -kj^{-1} & \pmod p \\ i^{-1} & \equiv -\left \lfloor \frac p i \right \rfloor (p \bmod i)^{-1} & \pmod p \end{aligned} \]这样就可以做到线性递推求逆元了。时间复杂度 \(O(n)\)，就可以通过 P3811 【模板】模意义下的乘法逆元 啦！
-
线性求 \(n\) 个数 \(a_i\) 的逆元：利用前缀积求逆元
设 \(s\) 是 \(a\) 的前缀积，那么如果我们知道了 \(s_n\) 的逆元，于是可以从 \(n\) 倒推回去：\(s_i^{-1} = {s_{i+1}}^{-1} \times a_i\)。
这样每个数的逆元就是 \({a_i}^{-1} = s_i^{-1} \times {s_{i-1}}\)，时间复杂度 \(\mathcal O(n + \log p)\)，这种方法在求阶乘的逆元中常常用到。
3.3 二次探测定理
二次探测定理：
如果 \(p\) 是一个奇素数，则方程 \(x^2 \equiv 1 \pmod p\) 仅有两个解：\(x=1\) 和 \(x=p-1\)。
我们先不急着讨论这个东西，先来讨论 \(x^2 \equiv 1 \pmod m\) 有多少个解？
首先我们应该考虑 \(m\) 是素数幂 \(p^k\) 的情形，这里 \(k \gt 0\)，此时同余式 \(x^2 \equiv 1\) 可以写成
所以 \(p\) 必定整除 \(x-1\) 或者 \(x+1\)，或者整除它们两者。但是 \(p\) 不可能同时整除 \(x-1\) 和 \(x+1\)，除非 \(p=2\)，如果 \(p \gt 2\)，那么 \(p^k \mid (x-1)(x+1) \Longleftrightarrow p^k \mid (x-1) 或者 p^k \mid (x+1)\)，所以恰好有两个解 \(x \equiv -1\) 和 \(x \equiv 1\)。
\(p=2\) 的情形稍有不同。如果 \(2^k \mid (x-1)(x+1)\)，那么 \(x-1\) 和 \(x+1\) 中的一个能被 \(2\) 整除，但不能被 \(4\) 整除，所以另外一个必定能被 \(2^{k-1}\) 整除。
这就意味着当 \(k \ge 3\) 时我们有四个解，即 \(x \equiv \pm 1\) 和 \(x \equiv 2^{k-1} \pm 1\)。
现在，\(x^2 \equiv 1\pmod m\) 当且仅当对 \(m\) 的完全分解式中所有满足 \(m_p \gt 0\) 的素数 \(p\) 都有 \(x^2 \equiv 1 \pmod {p^{m_p}}\)。
每个素数都是独立于其他素数的，对于 \(x \bmod {p^{m_p}}\)，除了 \(p=2\) 的情形，皆有两种可能性。
于是如果 \(m\) 有 \(r\) 个不同的素因子，则 \(x^2 \equiv 1\) 的解的总数是 \(2^r\)，除了当 \(m\) 是偶数时要做修正之外，一般来说，精确的解数是
在这个推导过程中，我们已经证明了 二次探测定理 了！
3.4 威尔逊定理
接下来引出一个逆元的应用——威尔逊定理，他会在之后有关素数的判定中有所应用。
这个定理的一半是平凡的：如果 \(n \gt 1\) 不是素数，它就有一个素因子 \(p\)，\(p\) 必然也是 \((n-1)!\) 的一个因子，所以 \((n-1)!\) 不可能与 \(-1\) 同余。
具体来说，在 1.4 我们引入了一个改变模数的式子，也就是如果 \((n-1)!\) 对模 \(n\) 与 \(-1\) 同余，那么就应该与模 \(p\) 与 \(-1\) 同余，但是事实并非这样。
威尔逊定理的另一半说的是，\((p-1)! \equiv -1 \pmod p\)，我们可以用 逆元配对 来证明这一半的结论。
因为 \(p\) 是素数，如果 \(n \bot p\)，那么就一定存在一个 \(n'\) 满足
这里 \(n\) 是 \(n'\) 的逆元，\(n'\) 也是 \(n\) 的逆元。所以如果 \(n' \neq n\)，那么它们两个就可以进行配对，从而直接抵消掉了。
那哪些数的逆元就是它本身呢？
容易发现其实就是问 \(x^2 \equiv 1 \pmod p\) 的解，这就是上面我们证明的二次探测定理。
当 \(p\) 是奇素数时，只会有 \(x=1\) 和 \(x=p-1\)，那么整个式子就是 \((p-1)! \equiv (p-1) \equiv -1 \pmod p\)，证毕。
而 \(p=2\) 时，直接带进去也是非常显然的，这样我们就证明出来了。
3.5 有趣的遗留问题
在上面推类欧几里得算法的一个拓展时，我们希望证明 \(m\) 个数
按照某种次序恰好组成了 \(\frac m d\) 个数
的 \(d\) 分复制，其中 \(d=\gcd(n,m)\)。
或许上面已经给出了一些不太严谨的简要证明，我们在这里希望把它写得更完整一些。
首先证明得到前面 \(\frac m d\) 个值的 \(d\) 份复制是简单的，根据 1.4 同余关系的一些等式，我们可以知道
从而当 \(0 \le k \lt \frac md\) 时，我们得到这些值的 \(d\) 份复制.
那么现在我们必须指出，那 \(\frac md\) 个数就是 \(\{ 0,d,2d,\cdots ,m-d \}\) 按照某种次数排列。
我们记 \(m = m'd,n=n'd\)，那么根据 1.4 中的分配律就有 \(kn \bmod m = d(kn' \bmod m')\)，所以当 \(0 \le k \lt m'\) 时出现的那些之是 \(d\) 去乘
而这种情形就是 \(m' \bot n'\) 的情形。
根据 鸽巢原理，我们就可以轻松证明出这 \(m'\) 个数就是 \(0,1,\cdots,m'-1\)，这个遗留问题就解决了！
它可以被理解成在一个环上面每次跳 \(n\) 步能跳到的点，在 OI 中也会有许多应用。
4 线性同余方程组
求解形如
的 线性同余方程组，在 OI 中有着广泛的应用。
4.1 中国剩余定理（CRT）
中国剩余定理全称 Chinese Remainder Theorem，简称 CRT。它用于求解 模数两两互质 的线性同余方程组，即对于任意 \(i \neq j\)，都有 \(n_i \bot n_j\)。
首先我们来看一个例子：
在《孙子算经》中有这样一个问题：“今有物不知其数，三三数之剩二（除以 \(3\) 余 \(2\)），五五数之剩三（除以 \(5\) 余 \(3\)），七七数之剩二（除以 \(7\) 余 \(2\)），问物几何？
CRT 的步骤是这样的：
-
找出三个数：
从 \(3\) 和 \(5\) 的公倍数中找出被 \(7\) 除余 \(1\) 的最小数 \(15\)；
从 \(3\) 和 \(7\) 的公倍数中找出被 \(5\) 除余 \(1\) 的最小数 \(21\)；
从 \(7\) 和 \(5\) 的公倍数中找出被 \(3\) 除余 \(1\) 的最小数 \(70\)。
-
用 \(15 \times 2\) （\(2\) 为最终除以 \(7\) 的余数），用 \(21 \times 3\)（\(3\) 为最终除以 \(5\) 的余数），用 \(70 \times 2\)。
-
把以上三个数相加得到 \(233\)，再用它除以 \(3,5,7\) 的最小公倍数 \(105\) ，得到余数 \(23\) 即为最小满足条件的数。
它为什么对的呢？
我们假设 \(m_1\) 是一个 \(\equiv 2 \pmod 3\) 的数， \(m_2\) 是 \(\equiv 3 \pmod 5\) 的数， \(m_3\) 是 \(\equiv 2 \pmod 7\) 的数。
从 \(m_1\) 的角度出发，已知 \(m_1 \equiv 2\pmod 3\)，如果需要 \((m_1 +m_2 +m_3) \equiv 2 \pmod 3\) ，那么 \(m_2,m_3\) 都是 $3 $ 的倍数。
同理我们可以得到：
- \(m_1\) 除以 \(3\) 余 \(2\)，且是 \(5\) 和 \(7\) 的公倍数。
- \(m_2\) 除以 \(5\) 余 \(3\)，且是 \(3\) 和 \(7\) 的公倍数。
- \(m_3\) 除以 \(7\) 余 \(2\)，且是 \(3\) 和 \(5\) 的公倍数。
所以 CRT 用数学表示出来就是：
设 \(N=n_1 \times n_2 \times \dots \times n_k\)，\(N_i\) 表示 \(N \div n_i\)。
\(x\equiv (a_1 N_1 N_1^{-1}+a_2 N_2 N_2^{-1}+ \dots + a_kN_kN_K^{-1}) \pmod N\)。
因为 \(N_i {N_i}^{-1}\)（\({N_i}^{-1}\) 是关于 \(\bmod n_i\) 的逆元） 就满足了是其他数的公倍数，且 \(\equiv 1 \pmod {n_i}\) ，于是再用 \(a_i\) 乘上它即可。
#include<bits/stdc++.h>
using namespace std;
#define ll long long
const int N=15;
int n;
ll p,a[N],b[N],m,x,y,ans=0;
void exgcd(ll a,ll b ,ll &x,ll &y){
if(b==0){x=1,y=0;return;}
exgcd(b,a%b,y,x);y-=a/b*x;
}
int main(){
scanf("%d",&n);m=1ll;
for(int i=1;i<=n;i++){
scanf("%lld%lld",&a[i],&b[i]);
m*=a[i];
}
for(int i=1;i<=n;i++){
p=m/a[i];
exgcd(p,a[i],x,y);
ans+=b[i]*p*(x<0?x+a[i]:x);
}
printf("%lld\n",ans%m);
return 0;
}
4.2 扩展中国剩余定理（exCRT）
如果就是单纯的 线性同余方程组 怎么解呢？也就是 \(n\) 并不两两互质。
这就是 exCRT 所解决的问题，它能解决的东西更为普遍，故使用范围也会更广。
先考虑只有两个方程的简单情况。
容易发现 \(x=a_1+n_1 k_1=a_2 + n_2k_2\)，我们把右边的两个式子提出来并进行移项
那么这就是一个可以用 exgcd 求解的方程了。
如果 \(\gcd(n_1,n_2) \nmid a_2-a_1\)，那么线性同余方程组无解。
否则，我们可以得到新的方程，也就是合并了两个方程：
推广到 \(n\) 个方程的情况，这样依次合并下去，最后的答案就是合并出来的 \(a\) 了。
我们就成功求解了 线性同余方程组！
P4777 【模板】扩展中国剩余定理（EXCRT） 的代码。
ll mul(ll a,ll b,ll m){
ll res=0ll;
while(b){if(b&1) res=(res+a)%m;a=(a+a)%m;b>>=1;}
return res;
}
ll exgcd(ll a,ll b, ll &x,ll &y){
if(b==0){x=1,y=0;return a;}
ll d=exgcd(b,a%b,y,x);y-=(a/b)*x;return d;
}
ll wrk(){
ll ans=a[1],M=b[1];
for(int i=2;i<=n;i++){
ll A=M,B=b[i],C=(a[i]-ans%B+B)%B;
ll gcd=exgcd(A,B,x,y),bg=B/gcd;
if(C%gcd!=0) return -1;
x=mul(x,C/gcd,bg);
ans+=x*M;M*=bg;
ans=(ans%M+M)%M;
}
return (ans%M+M)%M;
}
5 素数
接下来我们来讨论素数。
5.1 素数的例子
如果一个正整数 \(p\) 恰好只有两个因子，即 \(1\) 和 \(p\) ，那么这个数就称为 素数。
而其他那些还有非平凡因子的数都称为 合数，每一个大于 \(1\) 的整数要么是素数，要么是合数。
任何正整数 \(n\) 都可以表示成素数的乘积：
而 算数基本定理 表示：仅有一种方法将 \(n\) 按照素数非减的次序写成素数的乘积。（证明略）
关于素数的其他定理：
-
第 \(n\) 个素数 \(P_n\) 大约是 \(n\) 的自然对数的 \(n\) 倍：
\[P_n \sim n \ln n \] -
对于不超过 \(x\) 的素数个数 \(\pi(x)\):
\[\pi(x) \sim \frac{x}{\ln x} \]
证明过于复杂，这里不给出证明，大家可以自行了解。
还有一种与素数相关的数：梅森数。
形如
的数，其中 \(p\) 总表示素数。
如果 \(n\) 是合数，那么 \(2^n -1\) 不可能是素数，因为 \(2^{km}-1\) 以 \(2^m-1\) 作为一个因子：
但是当 \(p\) 是素数时，\(2^p-1\) 并不总是素数，\(2^{11}-1=2047=23 \times 89\) 是最小的这类非素数。
5.2 互素关系
当 \(\gcd(n,m)= 1\) 时，整数 \(m\) 和 \(n\) 没有公共的素因子，我们就称它们是 互素的。
我们写成
而一个分数是最简分数，当且仅当 \(m \bot n\)，通常我们通过消去分子分母的最大公约数来简化成最简分数，于是我们有
那我们有没有什么好的方法来构造满足 \(m \bot n\) 的全部非负的分数 \(\frac m n\) 组成的集合呢？
答案是肯定的，他被成为 Stern-Brocot 树，其思想是从两个分数 \(\left( \frac 0 1 ,\frac 1 0 \right)\) 出发，然后依照你希望的次数重复下面的操作
- 在两个相邻的分数 \(\frac m n\) 和 \(\frac {m'} {n'}\) 之间插入 \(\frac {m + m'} {n + n'}\)。
新的分数就是 \(\frac {m+m'}{n+n'}\) 称为 \(\frac m n\) 和 \(\frac {m'} {n'}\) 的 中位分数，例如
接下来就是
于是我们把这个分数阵列可以看成是一棵无线的二叉树构造，就是这样
每一个分数都是 \(\frac {m+m'} {n+n'}\)，其中 \(\frac m n\) 是位于左上方且离它最近的祖先，而 \(\frac{m'} {n'}\) 则是右上方离它最近的祖先。
为什么这种构造是对的呢？为什么每次的 \(\frac {m+m'} {n+n'}\) 是最简分数？而又只会出现一次？
因为我们发现，在这个构造中任意一个阶段的相邻分数都满足
首先这个关系在 \(1 \times 1 - 0 \times 0=1\) 是成立的，而每一次插入一个新的中位分数 \(\frac {m+m'} {n+n'}\) 时，我们需要检验的就是
括号展开就一样了，而这也让我们发现
一个中位分数不一定在原先两个的中间，但它的确位于它们两个中间的某个地方，于是我们就不可能得到相同的分数。
但是我们现在仍然存在一个问题——会遗漏么？
这个问题时简单的，由于它是无限的，所以我们会不断用上面的单调性去逼近它，这样就不会漏。
这里我们再引入一个概念，阶为 \(N\) 的法里级数 \(\mathscr F_n\)，它是介于 \(0\) 和 \(1\) 之间的分母不超过 \(N\) 的所有最简分数组成的集合，且按照递增的次序排列，那么当 \(N=6\) 时就是
容易发现它是 S-B 树 的子树，好像也没什么好讲的。
事实上，我们可以把 S-B 树 看成一个表示有理数的 数系，因为每一个正的最简分数都恰好出现一次。
我们用字母 \(L,R\) 表示它在这一层往哪边走，就可以唯一确定一个数。
我们记当前的字符串为 \(S\)，那么就有
例如 \(f(LRRL)= \frac 5 7\)。
这个东西和矩阵乘法的构造是比较像的，也就是我们尝试去构造矩阵
那么我们知道
然后我们就发现其实
于是我们就可以通过一个 \(LR\) 字符串定位出它所定义的值了。
还有一种方法就是改变 \(m,n\)，而不是保持状态 \(S\)，我们发现
因为 \(RS\) 像 \(S\)，只是将上面一行加到了下面一行上，也就是
类似的性质对 \(L\) 也有，于是我们就有
这样的定位就可以不用矩阵乘法。
关于互素的讨论到此结束，更多可以参见《具体数学》！
P8058 [BalkanOI2003] Farey 序列
把所有分子和分母都 \(\le n\) 的 最简真分数 从小到大排成一行，形成的序列称为 Farey 序列。
求出 \(n\) 所对应的 Farey 序列中第 \(k\) 小的数。
\(2 \le n \le 4 \times 10^4,1 \le k \le\) 符合条件分数的个数。
sto smb orz
给出一个 smb 的亚线性做法。
有了上面的一些简单铺垫，我们不难想到在 S-B 树 上面进行二分。
那么对于一个最简分数 \(\frac x y\)，有多少个数比它小呢？
容易列出式子
对这个式子进行莫反，可以得到
我们发现前面是可以用杜教筛处理的（其实这道题中不用），而后面这一块是类欧几里得的板子，于是这个值我们是可以在 \(\mathcal O(n^{\frac 2 3} + \sqrt n \log n)\) 的时间内完成的。
而我们如何去二分到这个值呢？
不难想到直接在 S-B 树 上面走，每次枚举走到左子树还是走到右子树，但是在 \(k=1\) 的情况我们会走 \(n\) 步，很明显是不满足条件的。
我们考虑把往一个方向走的连续段缩起来，也就是只枚举在哪些地方拐（eg：原先一直往左子树走，现在变成往右子树走），设 \(F(n)\) 表示拐 \(n\) 次能到达的分母最小值，那么容易得到 \(F(n) = F(n-1)+F(n-2)\)，也就是一直拐 \(n\) 次，那么左边就是 \(F(n-1)\)，右边是 \(F(n-2)\)。
那么这样我们只会拐 \(\log n\) 次，每次走多长可以用倍增实现，于是我们就做完了。
总时间复杂度 \(\mathcal O(n^{\frac 2 3} + \sqrt n \log^3 n)\)。代码。
ABC333G Nearest Fraction
给定一个小于 \(1\) 的正实数 \(r\) 和一个正整数 \(n\)。 要求在满足 \(0 \le p \le q \le n\) 和 \(\gcd(p,q)=1\) 的前提下，找到使 \(|r-\frac{q}{p}|\) 最小的二元组 \((p,q)\) 。 如果存在多个这样的二元组 \((p,q)\)，输出 \(\frac{q}{p}\) 值最小的那个。
数据范围：\(1\le n \le 10^{10}\)，\(r\in (0,1)\) 且最多包含 18 位有效数字。
首先 \(r\) 容易被转化成分数形式，把它设为 \(m\)，我们考虑用 S-B 树 来解决这个问题。
对于当前的区间 \([l,r]\)，我们考虑在树上二分来找到最近的位置，那么每一次一定是形如先从 \(r\) 一直往左跳到第一个 \(\lt m\) 的地方，再把这个点设成左端点，又不断往右跳。
这样我们就可以不断减小这个区间，从而不断逼近 \(m\)。
每一次搜索都用两个答案 \(ans1,ans2\) 分别记录在 \(m\) 两侧得到的最近点，这样就做完了。
用上一道题类似的倍增优化可以做到 \(\mathcal O(\log^2 n)\)。代码。
双倍经验：P1298 最接近的分数
SP26073 DIVCNT1 - Counting Divisors
\(\sigma_0(i)\) 表示 \(i\) 的约数个数。
求
\[S(n) = \sum_{i=1}^n \sigma_0(i) \]多测，\(T \le 10^5,n \le 2^{63}\)。
比较难的应用。
首先我们可以对式子进行转化
那么问题就转化成求
我们尝试把这个反比例函数画出图来，那么我们要求的其实就是下图中红色部分中的点数（from x义x）：
不妨将其称作 R 区域。
我们希望用一些向量去不断拟合这个反比例函数，于是就想到了这个有理数系：S-B 树。
于是我们需要一个单调栈，里面维护的向量斜率单调递减，进行这样的操作：
- 每次从栈顶取出一个元素，不断加到当前点 \((x,y)\) 上面，直到它差一步就走进 R 区域。
由于我们保证了斜率递减，且向量的斜率表示出来已经是最简分数，所以每一步都是对答案有贡献的，也就是图中棕色区域中的点（from x义x）：
容易发现这样的点数一定不会多算，可以做到不重不漏。（在斜线下的点一定在反比例函数中，并且不会多）
如何计算点数？（我们假设每次上都要计算最下面的一行，从而不能计算上面的一行）
画出图来就是这个样子：
把它分成蓝色和橙色的两个部分，而橙色部分一定是上下对称的，所以假设当前向量为 \(L\)，当前起点在 \((x,y)\)，终点在 \((x+L.x,y-L.y)\)，则贡献是
- 不断弹出栈顶，使得攻栈顶刚好走不进 R 区域，称它为 \(r\)。设最近一次弹的是 \(l\)（有可能是上一步的栈顶向量）。
现在就是，加上 \(l\) 刚好走进 R，加上 \(r\) 刚好走不进 R，就是这个样子（from x义x）：
- 开始二分，用 S-B 树 得到中间向量 \((l_x+r_x,l_y+r_y)\)，我们称之为 \(mid\)。
-
如果 \(mid\) 走一步不会进 R，那么 \(mid\) 进栈并把 \(r = mid\) 继续二分。
-
否则
- 如果 \(r \le f'(x+x_{mid})\)，那么继续二分肯定都走不出 R 了，直接停止二分，就是这个样子（from x义x）：
可以发现如果 \(r \le f'(x+x_mid)\)，\(mid,r\) 的中间向量 \(midmid\) 显然也一定会进 R，而且因为 \(f'\) 单调增，\(r \le f'(x+x_{midmid})\) 也一定满足。
- 否则 $l=mid$ 继续二分。
当 \(y \le \sqrt[3] n\) 时可以直接暴力。
这样的时间复杂度是 \(\mathcal O(\sqrt[3] n \log n)\)。（笔者不会证）
注意到其他的一些东西：
-
最开始的点我们选择 \(\left(\left \lfloor \frac n B \right \rfloor,B+1\right)\)，因为这样才能做到不重不漏（细节会少很多）。
-
初始栈中的向量是 \((1,0)\) 和 \((1,1)\)。
这样我们就做完了。代码。
5.3 素数筛法
最重要的知识点，几乎所有题目都需要筛出 \(1 \sim n\) 的所有素数。
首先我们存在一个非常显然的做法——埃氏筛法。
埃氏筛用所有已经筛出的素数倍数标记合数：从 \(2\) 到 \(n\) 枚举所有 \(i\)，若 \(i\) 未被标记，则 \(i\) 是素数，将 \(i\) 除 \(i\) 以外的倍数标记为合数。
埃氏筛的精髓在于复杂度的证明：不超过 \(n\) 的指数的倒数之和为 \(\mathcal O(\ln \ln n)\) 级别。
也就是
这说明埃氏筛的时间复杂度是 \(\mathcal O(n \ln \ln n)\)。
这里是来自 djq 的证明：
因为每个数只会被其素因子筛到，所以
根据 \(d(i)\) 的计算式
根据 \(2^x\) 的凸性和琴生不等式得
两边同时取对数
因此
埃氏筛已经非常接近与线性了，但是它也还不是最优的，我们存在一个 线性筛素数 的做法。
它的思想和埃氏筛类似，埃氏筛的优秀之处在于仅用质数的倍数筛去合数，但合数会被多个质因子筛到。让每个合数仅被晒一次就能做到线性。
考虑用每一个合数 最小质因数 筛它：从 \(2\) 到 \(n\) 枚举所有数 \(i\)。对于每个 \(i\)，令其最小的质因子为 \(p\)，则对于不大于 \(p\) 的质数 \(q\)，\(iq\) 的最小质因子为 \(q\)。将所有 \(iq\) 标记为合数，则每个合数 \(c\) 仅在 \(i = \frac c m\) 时以 \(im\) 的形式晒去，其中 \(m\) 是 \(c\) 的最小质因子。
综上，有如下步骤：
-
从小到大遍历当前筛出的所有素数 \(pr_j\)，将 \(i \times pr_j\) 标记为合数。
-
若 \(pr_j \mid i\)，退出循环，因为 \(pr_j \mid i \times pr_k(k \gt j)\)，所以 \(i \times pr_k\) 的最小质因子为 \(pr_j\) 不是 \(pr_k\)，再筛就重复了。
时间复杂度线性，于是就可以通过 P3383 【模板】线性筛素数。
for(int i=2;i<N;i++){
if(!vis[i]) prm[++cnt]=i;
for(int j=1;j<=cnt&&prm[j]*i<N;j++){
vis[prm[j]*i]=true;
if(i%prm[j]==0) break;
}
}
5.4 素数判定（Miller Rabin 素性测试）
如何判断一个数是否是素数？
小素数的判定是容易的：直接试除法即可，时间复杂度 \(O(\sqrt{n})\)。
至于大素数判定……
费马曾提出一种方法：费马素性测试。
它基于费马小定理：设 \(n\) 是素数， \(a\) 是正整数且与 \(n\) 互素，那么有 \(a^{n-1} \equiv 1 \pmod n\)。
为了测试 \(n\) 是否是素数，我们在 \(1 \sim n\) 中随机选择一个基数 \(a\) ，而 \(a\) 不需要与 \(n\) 互素。
-
如果 \(a^{n-1} \equiv 1 \pmod n\) 不成立，那么 \(n\) 肯定不是素数。
-
如果 \(a^{n-1} \equiv 1 \pmod n\) 成立，那么 \(n\) 很大概率是素数，尝试的 \(a\) 越多， \(n\) 是素数的概率就越大。
但是从第 \(2\) 种情况看出来费马素数测试不是完全正确的。
对于某个 \(a\) 值，总有一些合数被误判通过了测试；不同的 \(a\) 值，被误判的合数不太一样。
特别地，有一些合数，不管选什么 \(a\) 值，都能通过测试，这种数叫 Carmicheal 数 \((561,1105,1729,……)\)。
不过这种数很少，\(1\) 亿以内只有 \(255\) 个，而越到大的越稀疏，所以费马素性测试几乎不会出错，编码也很简单。
由于上面 费马素性测试 的不完美，我们就引入了 Miller-Rabin 素性测试。
把费马素性测试改进一下，它是已知最快的随机素数测试算法。
原理也非常简单，就是费马素性测试排除 Carmicheael 数。
因为并且 Carmicheael 数在进行次数较多的二次探测定理逆否命题判定后，都会被判定出来。
如果 \(n > 2\)，且 \(n\) 是奇数，测试它是否为素数。
根据费马测试，如果 \(a^{n-1}\equiv 1 \pmod n\) 不成立，一定不是素数。
而由于 \(n-1\) 过大，我们考虑一个技巧：
把 \(n-1\) 表示成幂的形式，令 \(n-1=2^tu\) 其中 \(u\) 是奇数，\(t\) 是正整数。
于是有：
所以计算出 \(a^u\) 依次乘二，这样每次平方同时使用二次探测定理进行判定。
容易发现，对于多个随机的基值 \(a\)，假设随机了 \(s\) 个，做了 \(s\) 次测试，出错的概率是 \(2^{-s}\) 。
当 \(s=50\) 时，出错的概率已经可以忽略不计了。
考虑做了 \(s\) 次运算，每次做 \(t\) 次快速幂，在 \(n\) 很大的时候速度极快！
ll mul(ll a,ll b,ll m){ll z=(long double)a/m*b,res=(a*b-z*m)%m;return (res+m)%m;}
ll add(ll a,ll b,ll n){return (a+b>=n)?a+b-n:a+b;}
ll qpow(ll a,ll b,ll m){
ll res=1ll;
while(b){
if(b&1) res=mul(res,a,m);
a=mul(a,a,m);
b>>=1;
}
return res;
}
const ll prm[12]={2,3,5,7,11,13,17,19,23,29,31};
bool miller_rabin(ll n){
if(n<2) return false;
for(int i=0;i<=10;i++){
if(n==prm[i]) return true;
if(n%prm[i]==0) return false;
}
ll u=n-1,t=0;
while(!(u&1)) u>>=1,t++;
for(int i=0;i<=10;i++){
ll a=prm[i],x1=qpow(a,u,n),x2;
for(int j=1;j<=t;j++){
x2=mul(x1,x1,n);
if(x2==1&&x1!=1&&x1!=n-1) return false;//二次探测定理
x1=x2;
}
if(x1!=1) return false;
}
return true;
}
5.5 分解质因数（Pollard-rho）
判断完了素数，我们来尝试分解质因数。
众所周知，可以用试除法做到 \(O(\sqrt{n})\)，但是同样 \(n\) 一旦大了就不行了。
所以我们考虑一种更高消的方法——Pollard-rho 启发式方法 分解质因数。
事实上，我们还有一种奇妙的方法来寻找一个合数的因子（from TQX）:
scanf("%d",&n);
vector<int> fac;
for(int i=1;i<=???;++i){
int a=rand()%(n-1)+1;
if(n%a==0) fac.push_back(a);
}
这个代码试图通过随即一些数并判断它们是否是 \(n\) 的因数，很明显这是一个非常蠢的方法，因为单次操作找到 \(n\) 的因数的概率实在是太小了。
这个算法非常的劣质，但我们的 Pollard_rho 正是基于这个算法而来的，它的核心思想就是随机。
而我们知道 生日悖论：
如果现在有 \(1000\) 个数，在里面随机到 \(37\) 的概率非常小，但如果我们找两个数 \(i,j\) 使 \(|i-j|=37\) 的概率就要大一倍；
而当 \(k=30\) 时，概率已经超过 \(50\%\) ，当 \(k=100\) 时，概率已经到了 \(99.99\%\)。
这个悖论告诉我们：组合随机采样比单个个体满足条件要高得多，这样可以提高正确率。
于是，Pollard 就提出了一个随机数的生成方法：
其中 \(c\) 可以随机生成，在随机生成了一个 \(x_1\) 后， \(x2=f(x1),x3=f(x2),\cdots\)。
这个函数的随机性非常好，但是进行了数次生成之后就会出现一个环，在坐标系中表示出来就变成了这个样子：
发现这个图像很像 \(\rho\)，于是这就是它的名字由来。
那么考虑如何判环？
可以用 Floyd 的方法：
假设我们有 \(2\) 个函数生成的值 \(a,b\)，我们每次让 \(a\) 走一步， \(b\) 走两步，当 \(a,b\) 相等时就出现了环。
于是我们可以用这样随机出来的数相邻两个的差与 \(n\) 取 gcd 来分解 \(n\)，即如果得到的 gcd \(\neq 1\) ，那么就直接返回。
但是 gcd 的时间复杂度时 \(O(\log n)\) ，每一次都去调用一定会很慢。
在计算时，我们将每次产生的 \(abc(a-b)\) 相乘并积累下来，最后直接判断这个乘积与 \(n\) 的 \(\gcd\)。
但是，如果某一时刻累积下来的样本的乘积为 \(0\) 了，为例不让样本丢失，我们直接退出循环进行判断即可。
每次计算的阈值可以倍增，并且加上一个上限，使用 \(128\) 是一种不错的选择（玄学）。
至此，你就可以通过 P4718 【模板】Pollard-Rho 了。代码。
ll mul(ll a,ll b,ll m){ll z=(long double)a/m*b,res=(a*b-z*m)%m;return (res+m)%m;}
ll add(ll a,ll b,ll n){return (a+b>=n)?a+b-n:a+b;}
ll qpow(ll a,ll b,ll m){
ll res=1ll;
while(b){
if(b&1) res=mul(res,a,m);
a=mul(a,a,m);
b>>=1;
}
return res;
}
const ll prm[12]={2,3,5,7,11,13,17,19,23,29,31};
bool miller_rabin(ll n){
if(n<2) return false;
for(int i=0;i<=10;i++){
if(n==prm[i]) return true;
if(n%prm[i]==0) return false;
}
ll u=n-1,t=0;
while(!(u&1)) u>>=1,t++;
for(int i=0;i<=10;i++){
ll a=prm[i],x1=qpow(a,u,n),x2;
for(int j=1;j<=t;j++){
x2=mul(x1,x1,n);
if(x2==1&&x1!=1&&x1!=n-1) return false;
x1=x2;
}
if(x1!=1) return false;
}
return true;
}
ll gcd(ll a,ll b){return (b==0)?a:gcd(b,a%b);}
ll f(ll a,ll c,ll n){return add(mul(a,a,n),c,n);}
ll rad(ll n){return 1ll*rand()*rand()%n*rand()%n+1ll;}
ll pollard_rho(ll n){
if(n==4) return 2;
ll x=rad(n-1),c=rad(n-1),y=x;
x=f(x,c,n);y=f(x,c,n);
for(int lim=1;x!=y;lim=min(128,lim<<1)){
ll cur=1ll,nxt=0;
for(int i=1;i<lim;i++,cur=nxt){
nxt=mul(cur,abs(add(x,n-y,n)),n);
if(!nxt) break;//x=y
x=f(x,c,n);y=f(f(y,c,n),c,n);
}
ll d=gcd(cur,n);
if(d!=1) return d;
}
return n;
}
void sol(ll n){
ll d=pollard_rho(n);
while(d==n) d=pollard_rho(n);
if(miller_rabin(d)) ans=max(d,ans);
else sol(d);
if(miller_rabin(n/d)) ans=max(ans,n/d);
else sol(n/d);
}
void wrk(){
n=read();
if(miller_rabin(n)){puts("Prime");return;}
ans=0;sol(n);print(ans);
}
P4358 [CQOI2016] 密钥破解
一种非对称加密算法的密钥生成过程如下：
1.任选两个不同的质数\(p,q\)
2.计算\(N=p \times q\)，\(r=(p-1)(q-1)\)
3.选取小于\(r\)，且与\(r\)互质的整数\(e\)
4.计算整数\(d\)，使得\(ed≡1(mod r)\)
5.二元组\((N,e)\)称为公钥，二元组\((N,d)\)称为私钥。
当需要加密消息\(n\)时，（假设\(n\)是一个小于\(N\)的整数，因为任何格式的消息都可转为整数表示），使用公钥\((N,e)\)，按照
\[n^e≡c(mod N) \]运算，可得到密文\(c\)
对密文\(c\)解密时，用私钥\((N,d)\)，按照
\[c^d≡n(mod N) \]运算，可得到原文 \(n\)。算法正确性证明省略。
由于用公钥加密的密文仅能用对应的私钥解密，而不能用公钥解密，因此称为非对称加密算法。通常情况下，公钥由消息的接收方公开，而私钥由消息的接收方自己持有。这样任何发送消息的人都可以用公钥对消息加密，而只有消息的接收方自己能够解密消息。
现在，你的任务是寻找一种可行的方法来破解这种加密算法，即根据公钥破解出私钥，并据此解密密文。
\(e,N,c \le 2^{62},c<N\)。
题意写得非常长，但是并不难。
发现只要我们找到 \(p,q\)，就做完了，然而这个 \(p,q\) 可以直接交给 Pollard_rho 即可。代码。
P4714 「数学」约数个数和
给你一个正整数 \(N\)，请计算 \(N\) 的 \((\)所有约数的\()\times K\) 约数个数和。
答案可能很大，请输出对 \(998244353\) 取模的结果。
$1 \leq N \leq 10 ^ {18}, 0 \leq K \leq 10^{18} $。
首先，对于 \(K=0\) 的情况，\(f_0(x) = \sigma(x)\) 是积性函数。
而当 \(K \gt 0\) 时，我们有
这相当于一个狄利克雷卷积，而根据狄利克雷卷积卷积的性质，我们可以知道 \(f_i(n)\) 也是积性函数。
于是现在问题就变成了如何求 \(f_K(p^q)\)，因为我们把 \(n\) 分解之后就可以得到若干个 \(f_K(p^q)\) 的乘积。
发现这是简单的，等价于组合数学中在 \(q+1\) 个空位中插入 \(K+1\) 个板，可以重复。这时经典问题，你把每两个板之间就加一个球就可以得到
根据对称公式，相当于
而 \(q \le 62\)，所以可以直接用下降幂暴力计算。
至于分解质因数，交给 Pollard_Rho 就可以了。代码。
6 离散对数问题
离散对数问题就是解决模 \(p\) 意义下的 \(\log_ab\)，也就等价于解离散对数方程
这样的 \(x\) 就是 \(\log_a b\)，\(x\) 可能存在多个，也可能一个都不存在，我们希望找到这样的一个 \(x\)。
当 \(a \bot p\) 的时候，我们可以用大步小步算法解决，反之则可以用扩展大步小步算法解决。
6.1 大步小步算法（BSGS）
大步小步算法全称 Baby Step Giant Step，简称 BSGS，适用于 \(a \bot p\) 的情形。
整体思路用到了 根号平衡 的思想，设块长为 \(M\)，那么最终的 \(x\) 一定可以被表示成 \(x=AM-B,0 \le B \lt M\)。
于是这个方程就可以被表示成
只要我们直接暴力枚举 \(A\) 和 \(B\)，用哈希表（pb_ds）存一下是否出现过，就可以用 \(\mathcal O(\max(A,B))\) 的时间复杂度之内计算出答案。
根据欧拉定理（稍后会讲），有解就一定在 \([0,p-1]\) 中有至少一个解，所以这里 \(0 \le B \lt M,0 \le A \le \left \lceil \frac p M \right \rceil\)。
那么设块长为 \(\left \lceil \sqrt p \right \rceil\) 时，复杂度最优，为 \(\mathcal O(\left \lfloor \sqrt p \right \rfloor)\)。
于是你就可以通过 P3846 [TJOI2007] 可爱的质数/【模板】BSGS 啦！代码。
int BSGS(int a,int b,int H){
int A=1,B=sqrt(H)+1;a%=H,b%=H;
if(H==1||b==1) return 0;
if(!a) return !b?1:-1;
gp_hash_table<int,int> mp;
for(int i=1;i<=B;i++) mp[1ll*(A=1ll*A*a%H)*b%H]=i;
for(int i=1,AA=A;i<=B;i++,AA=1ll*AA*A%H)
if(mp.find(AA)!=mp.end()) return i*B-mp[AA];
return -1;
}
容易发现你从小到大枚举时就可以得到最小的非负整数解，而关于 \(x\) 的循环性质我们将会在后面阶与原根中讨论（循环的长度为 \(\delta_a(p)\)）。
在这里还有一些板题：
-
P2485 [SDOI2011] 计算器 ：综合一点的板题。代码。
P3306 [SDOI2013] 随机数生成器
最近小 W 准备读一本新书，这本书一共有 \(p\) 页，页码范围为 \(0 \sim p-1\)。
小 W 很忙，所以每天只能读一页书。为了使事情有趣一些，他打算使用 NOI2012 上学习的线性同余法生成一个序列，来决定每天具体读哪一页。
我们用 \(x_i\) 来表示通过这种方法生成出来的第 \(i\) 个数，也即小 W 第 \(i\) 天会读哪一页。这个方法需要设置 \(3\) 个参数 \(a,b,x_1\)，满足 \(0\leq a,b,x_1\lt p\)，且 \(a,b,x_1\) 都是整数。按照下面的公式生成出来一系列的整数：
\[x_{i+1} \equiv a \times x_i+b \pmod p \]其中 \(\bmod\) 表示取余操作。
但是这种方法可能导致某两天读的页码一样。
小 W 要读这本书的第 \(t\) 页，所以他想知道最早在哪一天能读到第 \(t\) 页，或者指出他永远不会读到第 \(t\) 页。
\(1 \le T \le 50,0 \le a,b,x_1,t \lt p,2 \le p \le 10^9,p\) 是质数
考虑先全部左移一维，那么得到的数其实就是 \(x_0\)，我们令 \(s_i = \frac {x_i} {a^i}\)，那么根据递推式，我们可以得到 \(a^i s_i = a^i s_{i-1} + b\)，所以 \(s_i = s_{i-1} + \frac b {a^i}\) 且 \(s_0 = x_0\)，则 \(s_i = x_0 + \sum_{j=1}^i \frac b {a^j}\)。
于是我们就可以将 \(x_i\) 表示出来：
这就可以用 BSGS 直接完成啦——注意特判 \(a=0,a=1\) 的情况。
总时间复杂度 \(O(T\sqrt p)\)。代码。
P4884 多少个 1？
给定整数 \(K\) 和质数 \(m\)，求最小的正整数 \(N\)，使得 \(11\cdots1\)(\(N\) 个 \(1\))\(\equiv K \pmod m\)。
说人话：就是 \(111\cdots 1111 \bmod m = K\)。
\(6\leq m\leq 10^{11}\)，\(0< K< m\)，保证 \(m\) 是质数。
a 如果我们直接把 $11111 \cdots$ 分解出来是 $10^0+10^1+10^2 + \cdots$，这是非常不好处理的。
我们考虑将这个数表示成一个数，也就是它乘 \(9\) 之后就会变成 \(\frac{10^n - 1} 9\)。
那么此时这个方程就变成了
于是这个东西就可以直接用 BSGS 完成了。
注意可能在乘的过程中爆 long long。代码。
CF1106F Lunar New Year and a Recursive Sequence
有一串 \(n\) 个数的数列，给你 \(b_1\sim b_k\)。当 \(i>k\) 时：
\[f_i=(\prod\limits_{j=1}^kf_{i-j}^{b_j})\bmod{998244353} \]已知 \(f_1=f_2=\cdots=f_{k-1}=1,f_n=m\)，问 \(f_k\) 可能是多少。
\(k \leq 100\)，\(n \leq 10^9\)。
会用到原根相关知识，若不会原根建议先看后面原根部分。
有关于 \(\prod\) 的式子是非常不好求的，我们考虑把它转化成有关加法的式子从而可以通过矩阵快速幂来完成。
而这种需要让我们想到了 原根，众所周知，\(3\) 是 \(998244353\) 的一个原根，我们记它为 \(G\)，并记 \(f_n = G^{g_n}\)，于是我们就把式子转化成了
把指数提下来，同时模数变成了 \(\phi (p)\)，我们就可以得到
那么只要我们能求得 \(g_i\) 那么求得 \(f_i\) 也就是非常容易的了。
发现这个东西刚好是矩阵乘法的形式，也就是我们可以构造 \(k \times k\) 矩阵
从而得到
于是当我们知道 \(g_k\) 时，可以快速算出 \(g_n\) 的值。
而由于前面都是 \(0\)，所以实际上我们有
其中 \(A^{n-k}[k][k]\) 是一个可以用矩阵快速幂处理出来的常数。
由于 \(998244352\) 并不是质数，所以当我们知道 \(g_n\) 的时候 \(g_k\) 可以直接用 exGCD 求出。
现在问题就转化成如何求出 \(g_n\)，由于 \(G^{g_n} \equiv f_n \pmod {998244353}\)，容易发现这时直接可以用 BSGS 求出的。
于是这道题就做完了，实现时注意时模 \(998244352\) 还是 \(998244353\) 即可。代码。
6.2 扩展大步小步算法（exBSGS）
接下来我们来讨论更为普遍的问题 \(a^x \equiv b \pmod p\)，如果没有 \(a \bot p\) 怎么处理？
从已知到未知： 由于 \(a \bot p\) 的情形我们已经可以解决了，那么我们现在考虑把 \(a,p\) 从不互质凑到互质。
我们运用同余方程的运算法则，我们把两边同时除以 \(d= \gcd(a,p)\) 就可以得到
由于 \(\frac a d\) 和 \(\frac p d\) 一定是互质的，所以在模 \(\frac p d\) 的意义下 \(\frac a d\) 是又逆元的，于是我们可以把 \(\frac a d\) 除过去，就是
然而这时的 \(a\) 和 \(\frac pd\) 并不一定互质，比如说 \(a=2,p=4\) 的情况，那么我们就重复上面的操作直到 \(a \bot p\)，而这个操作只会最多进行 \(\log p\) 次，这时容易证明的，在 \(a=2,p=2^k\) 取到极值。
而注意到每除依次，\(x\) 就会 \(-1\)。于是当操作 \(k\) 次后，我们需要计算其实是 \(x-k\)。此时如果 \(b=0\)，直接返回 \(k\) 就可以了，这时需要特判的。
容易发现这样的时间复杂度是 \(\mathcal O( \sqrt p + \log^2 p )\)，后者有每一次 exgcd 求逆元的时间复杂度。
这样就可以通过 P4195 【模板】扩展 BSGS/exBSGS 和 SP3105 MOD - Power Modulo Inverted 了。代码。
int gcd(int a,int b){return b==0?a:gcd(b,a%b);}
void exgcd(int a,int b,int &x,int &y){
if(b==0) return x=1,y=0,void();
exgcd(b,a%b,y,x),y-=x*(a/b);
}
int inv(int a,int p){exgcd(a,p,x,y);return (x%p+p)%p;}
int BSGS(int a,int b,int H){
int A=1,B=sqrt(H)+1;a%=H,b%=H;
if(H==1||b==1) return 0;
if(!a) return !b?1:-1;
gp_hash_table<int,int> mp;
for(int i=1;i<=B;i++) mp[1ll*(A=1ll*A*a%H)*b%H]=i;
for(int i=1,AA=A;i<=B;i++,AA=1ll*AA*A%H)
if(mp.find(AA)!=mp.end()) return i*B-mp[AA];
return -1;
}
int exBSGS(int a,int b,int H){
int d=gcd(a,H),cnt=0;a%=H,b%=H;
while(d>1){
if(b==1) return cnt;
if(b%d) return -1;
H/=d,b=1ll*(b/d)*inv(a/d,H)%H;
d=gcd(H,a%=H),++cnt;
}
int res=BSGS(a,b,H);
return ~res?res+cnt:-1;
}
P5345 【XR-1】快乐肥宅
给定 \(n\) 组 \(k_i,r_i,g_i\)，求关于 \(x\) 的方程组
\[\begin{cases} {k_1}^x & \equiv r_1 & \pmod {g_1} \\ {k_2}^x & \equiv r_2 & \pmod {g_2} \\ \cdots \\ {k_n}^x & \equiv r_n & \pmod {g_n} \end{cases} \]（定义 \(0^0=1\)）在 \([0,10^9]\) 范围内的最小整数解，或判断在这个范围内无解。
\(1 \le n \le 10^3,1 \le k_i,r_i \le g_i \le 10^7\)。
难绷。
容易发现，我们的解题过程分成两部分：求解每一个方程的解 和 求线性同余方程组的解。
首先考虑第一部分，如何求解每一个方程的解呢？
根据 exBSGS 的求解过程可以知道，解的情况一定分成了三种：无解、唯一解和无穷解。
画出图来，解一定是一个 \(\rho\) 形的，也就是这个样子（from 小粉兔）
其中 \(\rho\) 的尾巴是我们依次 \(\div \gcd\) 的过程，这种是只有一个解；
而进入环上就变成了无穷解，而这个解是形如线性同余方程；
反之，没有出现过的就是无解。
尾巴的长度是可以直接求 \(a^x \equiv 1 \pmod b\) 的最小正整数解（这里和上面的板子会略有区别），这里的 a 是 exBSGS 化出来的值。那么我们就可以用一次 exBSGS 的过程求出第一个解的位置 \(a\) 和循环节的长度 \(b\)，这就转化成了方程 \(x \equiv a \pmod b\)。
反之对于只有一个解的情况，我们直接记它为答案，最后特判一下就可以了。
接下来考虑第二部分，求 \(n\) 个线性同余方程组的解。
由于不保证模数互质，所以这个需要用 exCRT 解决，而 \(10^3\) 个 \(10^7\) 是容易爆 long long 的，如何判断呢？
考虑 exCRT 的过程，中间每一次我们都是 ans+=x*M，那么只要这次的 \(M \gt 10^9\)，就不需要继续算后面的部分了，因为下一次一定就暴了。
所以这个时候我们只需要判断当前的解对于后面的方程是否都成立即可。
这样就做完了，时间复杂度为 \(n\) 次 exBSGS 的复杂度。具体实现中有细节。代码。
7 数论函数基础
数论函数是数论中相当重要的一环，我们先来将一些基本的函数。
7.1 相关定义
-
数论函数：定义域为正整数的函数称为 数论函数。因其在所有正整数处均有定义，故可视作数列。OI 中常见的数论函数的陪域（即可能的取值范围）为整数。
-
加性函数：若对于任意 \(a,b \in \mathbb N_+\) 且 \(a \bot b\) 均有 \(f(ab) = f(a) + f(b)\)，则称 \(f\) 为 加性函数。注意区分代数中的加性函数。
-
积性函数：若对于任意 \(a,b \in \mathbb N_+\) 且 \(a \bot b\) 均有 \(f(ab) = f(a)f(b)\)，则称 \(f\) 为 积性函数。容易发现 \(f(1)=1\) 是必要条件。
-
完全积性函数：若对于任意 \(a,b \in \mathbb N_+\) 均有 \(f(ab)=f(a)f(b)\)，则称 \(f\) 为 完全积性函数。完全积性函数一定是积性函数。
-
数论函数的 加法：对于数论函数 \(f,g\)，\(f+g\) 表示 \(f\) 和 \(g\) 对应位置相加，即 \((f+g)(x) = f(x) + g(x)\)。
-
数论函数的 数乘：对于数 \(c\) 和数论函数 \(f\)，\(c \cdot f\) 表示 \(f\) 的各个位置乘 \(c\)，即 \((c \cdot f)(x) = c \cdot f(x)\)。一般简记作 \(cf\)。
-
数论函数的 点乘：对于数论函数 \(f,g\)，\(f \cdot g\) 表示 \(f\) 和 \(g\) 对应位置相乘，即 \((f \cdot g) (x) = f(x)g(x)\)。为与狄利克雷卷积符号 \(*\) 作区分，点乘符号通常不省略。
7.2 常见函数
设 \(n\) 的唯一分解为 \(\prod_{i=1}^m {p_i}^{c_i}\)，以下是一些常见的数论函数：
-
单位函数：\(\epsilon (n) = [n=1]\)。完全积性函数。
-
常数函数：\(1(n) = 1\)。完全积性函数。
-
恒等函数：\(\text{id}_k(n) = n^k\)。\(\text{id}_1(n)\) 记作 \(\text{id}(n)\)。完全积性函数。
-
除数函数：\(\sigma_k(n)= \sum_{d \mid n} d^k\)。\(\sigma_0(n)\) 表示 \(n\) 的约数个数，记作 \(\tau (n)\) 或 \(d(n)\)。\(\sigma_1(n)\) 表示 \(n\) 的约约数和，记作 \(\sigma (n)\)。\(\sigma_k(n)\) 有计算式
根据乘法分配律，\(n\) 的所有约数的 \(k\) 次方之和可以写成
等比数列求和后即可。
-
欧拉函数：\(\varphi(n) = \sum_{i=1}^n [i \bot n]\)，表示 \(n\) 以内与 \(n\) 互质的数的个数。后面我们会有所介绍。
-
本质不同质因子个数函数：\(\omega (n) = \sum_{p \in \mathbb P} [p \mid n]\)，表示 \(n\) 的本质不同质因子个数。
-
莫比乌斯函数：
以上所有函数除 \(\omega\) 是加性函数外，其余都是积性函数。根据计算式及积性函数的定义易证。
7.3 线性筛积性函数
在前面的素数部分，我们分析了线性筛素数，它提供了在线性时间内筛出具有特殊性质的积性函数在 \(1 \sim n\) 处所有取值的基本框架。
只要积性函数 \(f\) 可在 \(\mathcal O(1)\) 时间内计算任意指数幂处的取值 \(f(p^k)\)，就适用线性筛。（这只是 \(f\) 可线性筛的 充分但不必要 条件。存在更弱的条件使得 \(f\) 可线性筛但并不常见，如 \(\mathcal O(k)\) 计算 \(f(p^k)\)，这会在下一章介绍）。
根据积性函数的性质，只要预先求出 \(low_i\) 表示 \(i\) 的最小质因子 \(p\) 的最高次幂 \(p^{v_p(i)}\)，对于 \(i \neq p^k\)，即可使用 \(f(low_i)f(\frac {i} {low_i})\) 计算 \(f(i)\)。其中 \(v_p(i)\) 表示质因数次数符号。
下面放一个框架，from Alex_Wei：
for(int i = 2; i < N; i++) {
if(!vis[i]) pr[++pcnt] = i, f[i] = ..., low[i] = i; // 单独算 f(p)
for(int j = 1; j <= pcnt && i * pr[j] < N; j++) {
vis[i * pr[j]] = 1;
if(i % pr[j] == 0) { // i 与 p 不互质
low[i * pr[j]] = low[i] * pr[j];
if(i == low[i]) f[i * pr[j]] = ...; // i = p ^ k，单独算 f(p ^ {k + 1})
else f[i * pr[j]] = f[i / low[i]] * f[low[i * pr[j]]];
break;
}
low[i * pr[j]] = pr[j];
f[i * pr[j]] = f[i] * f[pr[j]]; // i 与 p 互质，f(ip) = f(i)f(p)
}
}
8 狄利克雷卷积
钢筋混凝土。
8.1 定义与性质
两个数论函数 \(f,g\) 的狄利克雷卷积为一个新的数论函数，记作 \(f * g\)。定义
上面简记成 \(h = f *g\)。按照定义是计算狄利克雷卷积，复杂度为 \(\mathcal O(n \ln n)\)。
狄利克雷卷积还有一些有趣的性质。
交换律，分配律，结合律。
这三个东西都可以用定义式稍微推一下就证明出来了。
单位元： \(\epsilon * f = f\)。
逆元：若函数 \(f,g\) 满足 \(f * g =\epsilon\)，则称 \(f,g\) 互为逆元。
积性函数总是有且仅有一个逆元。
\(f\) 可逆当且仅当 \(f(1) \neq 0\)。
下面写出证明和逆元的计算式：
设 \(g= f^{-1}\)。当 \(f(1)=0\) 时，\(f(1)g(1) = 1\) 与 \(f(1)g(1) = \epsilon (1) = 1\) 矛盾。当 \(f(1) \neq 0\) 时，\(g(1) = \frac 1 {f(1)}\)，对于 \(n \gt 1\) 的 \(g(n)\)，我们可以用已知的 \(g(1) \sim g(n-1)\) 推出
这同时也说明了逆元的 唯一性。
\(f= g\) 的充要条件是 \(f*h = g*h\)，其中 \(h(1) \neq 0\)。
直接乘上 \(h\) 的逆元即可。
积性函数的狄利克雷卷积也是积性函数。
设有两个积性函数 \(f,g\)，令 \(h =f * g\)，设 \(a,b\) 满足 \(a \bot b\)。
第二行中，由于 \(a \bot b\)，所以满足 \(d \mid a,t \mid b,dt = T\) 的 \(d,t\) 只有一对，可以从 \(w\) 中分别提取属于 \(a,b\) 的质因子得到。
积性函数的逆元也是积性函数。
设有积性函数 \(f\)，它的逆是 \(g\)。设 \(a,b\) 满足 \(a \bot b\)，证 \(g(ab) = g(a)g(b)\)。
使用归纳法。当 \(a\) 或 \(b\) 为 \(1\) 时，根据计算式显然成立。假设我们对 \(a' \lt a,b' \lt b\) 的 \((a',b')\) 都已经完成了证明，根据求你计算式可得
综合后面的两个性质，两个积性函数的积与商都是积性函数，但和与差不是。
8.2 线性筛狄利克雷卷积
给出 \(f,g\) 如何快速计算 \(h = f*g\) 呢？
我们直接用定义式计算，把枚举因数变成枚举倍数，根据调和级数，时间复杂度为 \(O(n \ln n)\)。
如果 \(f,g\) 是积性函数，我们可以尝试以下的优化，写出 \(h\) 的表达式，有
对于第一项和第三项，我们都可以在线性筛中总代价 \(\mathcal O(n)\) 求出，而现在问题就来到了第二项。
如果 \(f,g\) 在质数幂处的取值已经求出，则需要 \(\mathcal O(k)\) 的时间计算，我们尝试估算它的复杂度。考虑所有 \(\le \sqrt[k] n\) 的质数对复杂度产生 \(O(k)\) 的贡献，因此
\(x^2 \sqrt[x] n\) 随着 \(n\) 增大，而 \(x^2\) 的增大速度，远小于 \(\sqrt[x] n\) 的减小速度，于是从这一点入手，我们考虑证明 \(x^2 \sqrt [x] n \le O(n)\)。
当 \(x=1\) 是显然成立，于是考虑 \(x \in [2,\log_2 n]\) 时 \(x^2\) 的最大值 \(\log_2^2 n\) 与 \(\sqrt [x] n\) 的最大值 \(\sqrt n\) 之积，因为当 \(x \to + \infty\) 时 \(\log_2 x\) 是 \(\sqrt x\) 的高阶无穷小，所以 \(\mathcal O(\sqrt n \log^2 n) \le \mathcal O(n)\)。
因此，
综上，使用线性筛求出两个在质数幂处取值已知的积性函数的狄利克雷卷积在 \(1 \sim n\) 处的取值的时间复杂度是 \(\mathcal O(n)\)。
这样我们就得到了积性函数线性筛更弱的条件：可以 \(\mathcal O(k)\) 时间计算质数幂处的取值。
P6222 「P6156 简单题」加强版
给出 \(k\)，多次给出 \(n\)，求
\[\sum_{i=1}^n \sum_{j=1}^n (i+j)^K \mu^2(\gcd(i,j)) \gcd(i,j) \]\(1 \le n \le 10^7,K \le 2^{31},T \le 10^4\)。
会用到稍后要讲的莫比乌斯反演。
枚举 gcd，再进行莫比乌斯反演
我们需要线性筛预处理出 \(f=(d \times \mu^2) * \mu\)，并对 \(T^k\) 做前缀和，整除分块求解，时间复杂度 \(\mathcal O(n \frac {\log k} {\log n})\)。
这就需要用到我们上面讨论的线性筛狄利克雷卷积。代码。
简单版本：P6156 简单题。
8.3 狄利克雷前缀和
对于任意数论函数 \(f\) 卷常数函数 \(1\) 等价于对 \(f\) 做 狄利克雷前缀和：令 \(g = f * 1\)，则 \(g(n) = \sum_{d \mid n} f(d)\)。
对每个 \(n\) 计算给定数论函数在其所有因数处的取值之和有很好的实际含义，因此狄利克雷前缀和是相当重要的。
首先容易想到一种暴力做法，直接枚举每一个数的倍数，这样就可以做到 \(\mathcal O(n \ln n)\)。
但是并不是最优的，于是我们有了下面的算法。
将每个数 \(n\) 写成一个无穷的序列 \(a_n = \{ c_1,c_2,\cdots \}\)，表示 \(n = \prod p_i^{c_i}\)，其中 \(p_i\) 表示第 \(i\) 个质数。因为 \(x \mid y\) 的充要条件是 \(a_x (c_i) \le a_y(c_i)\)，所以 \(f * 1\) 可以看成对下标做关于其无穷序列的高维前缀和，即
根据高位前缀和的求法，枚举每一维并将所有下标关于该维做前缀和，可以得到狄利克雷前缀和的实现方法：
初始令 \(x_i = f(i)\)，从小到大枚举每一个质数 \(p_i\)，枚举 \(k\)，将 \(x_{p_ik}\) 加上 \(x_k\)，相当于 \(k\) 的贡献给到 \(a_k(i)\) 加上 \(1\) 之后的下标。
最终得到的 \(x\) 即为 \(g\)。根据小于 \(n\) 的素数倒数之和是 \(\ln \ln n\) 这一结论（埃氏筛法），狄利克雷卷积的时间复杂度是 \(\mathcal O(n \ln \ln n)\)。
于是你就可以通过 P5495 【模板】Dirichlet 前缀和 了。代码。
#include <bits/stdc++.h>
using namespace std;
#define uint unsigned int
const int N=2e7+5;
uint seed,ans,a[N];
int n,p[N>>3],cnt=0;
bool vis[N];
uint rd(){seed^=seed<<13,seed^=seed>>17,seed^=seed<<5;return seed;}
void init(){
for(int i=2;i<N;i++){
if(!vis[i]) p[++cnt]=i;
for(int j=1;j<=cnt&&i*p[j]<N;j++){
vis[i*p[j]]=true;
if(i%p[j]==0) break;
}
}
}
int main(){
cin>>n>>seed;init();
for(int i=1;i<=n;i++) a[i]=rd();
for(int i=1;i<=cnt;i++) for(int j=1;p[i]*j<=n;j++) a[p[i]*j]+=a[j];
for(int i=1;i<=n;i++) ans^=a[i];
cout<<ans;return 0;
}
9 整除分块
整除分块又称为数论分块，因其解决的问题与整除密切相关而得名。数论分块是用于求解形如
的合适，前提为 \(f\) 的前缀和可以快速计算。
9.1 介绍
如果 \(\left \lfloor \frac n i \right \rfloor\) 的数量不多，那么我们可以考虑转化成贡献形式，把原式写成若干个 \(g\left( \left \lfloor \frac n i \right \rfloor \right)\) 乘上一段 \(f\) 的和。那么不同的 \(\left \lfloor \frac n i \right \rfloor\) 到底又多少个呢？
\(\forall i \in [1,n],n \in \mathbb N_+\)，不同的 \(\left \lfloor \frac n i \right \rfloor\) 至多只有 \(2 \sqrt n\) 个。
证明是简单的：
-
当 \(i \le \sqrt n\) 时，\(\left \lfloor \frac n i \right \rfloor\) 只有 \(\sqrt n\) 个。
-
当 \(i \gt \sqrt n\) 时，\(\left \lfloor \frac n i \right \rfloor \le \sqrt n\)，也只有 \(\sqrt n\) 个。
根据这条结论，我们枚举 \(\mathcal O(\sqrt n)\) 种整除值 \(d\)，求出最小和最大的 \(i\) 使得 \(\left \lfloor \frac n i \right \rfloor = d\)，分别记作 \(l,r\)，那么原式可以写成
只要我们能对 \(f\) 进行快速的前缀和预处理，在 \(g\) 处能快速得到值，那么整个和式就可以在 \(\mathcal O(\sqrt n)\) 的时间复杂度内求得。
那么问题就转化成求使得 \(\left \lfloor \frac n i \right \rfloor = d\) 的最大和最小的 \(i\)。
转化一下，发现 \(\left \lfloor \frac n i \right \rfloor = d\) 对 \(i\) 有两条限制：
和
两边把 \(d\) 除过去就可以得到
那么我们不重不漏地枚举所有整除值时，可以对于每一个 \(i\) 计算 \(\left \lfloor \frac n i \right \rfloor\) 的 \(l\) 和 \(r\)，借助上面的式子跳过这个极长的连续段。
具体地，我们令当前枚举到的 \(i\) 为 \(l\)，此时整除值为 \(d = \left \lfloor \frac n i \right \rfloor\)。因为使得 \(\left \lfloor \frac n i \right \rfloor =d\) 的最大的 \(i\) 等于 \(\left \lfloor \frac n d \right \rfloor\)，所以
将 \(g(d)(s(r)-s(l-1))\) 算入答案，下一次 \(l=r+1\) 表示跳过了 \(l \sim r\) 这一段，如果 \(l \gt n\) 就退出即可。（\(s\) 是前缀和）
这样每一次的整除值都可以被算一次，时间复杂度 \(\mathcal O(\sqrt n)\)。
-
当 \(i\) 的上界不为 \(n\) 时，要特殊处理：若 \(m \le n\)，每次要把 \(r\) 和 \(m\) 取最小值；反之要特殊处理 \(\left \lfloor \frac n i \right \rfloor =0\) 的情况。
-
写的时候不要写成了 \(r \le n\) 了！（有些人因为这个常常调很久）
9.2 拓展
拓展也是较为简单的。
我们尝试将 向下取整 转换成 向上取值。
对于左边界 \(l\)，求出使得 \(\left \lceil \frac n l \right \rceil = \left \lceil \frac n r \right \rceil\) 的最大的 \(r\)。
设 \(k=\left \lceil \frac n l \right \rceil\)，那么
最后一步是因为 \(n,k\) 都是正整数，于是我们每次有
注意到 \(\left \lceil \frac n l \right \rceil = 1\) 时上界为 \(+ \infty\)，所以要去实际上界。
另外一个拓展就是我们尝试进行 高维数论分块。
当和式中出现若干个下取整，形如
时，我们只需要每次令
即可。最后要对 \(n\) 取 \(\min\)。
时间复杂度 \(\sum \sqrt {n_j}\)。
因为我们把存在 \(n_j\) 满足 \(\left \lfloor \frac {n_j} i \right \rfloor \neq \left \lfloor \frac {n_j} {i+1} \right \rfloor\) 的位置 \(i\) 视作断点，总断点数量为每个下取整式的断点数量相加，于是就有该时间复杂度。
9.3 例题
几道简单题。
CF1603C Extreme Extension
对于一个数列，定义一次操作为选择数列中任何一个元素 \(a_i\) ，将这个位置替换成两个正整数 \(x,y\) 满足 \(x+y=a_i\)
定义一个数列的极端值为将这个数列变成 不降 数列的最小操作次数
给定一个长度为 \(n\) 的数列 \(a\) （ \(1\leq a_i\leq 10^5\) ），让你求它的所有非空子段的极端值之和，对 \(\text{998244353}\) 取模
一个测试点中有多组数据，保证所有数据 \(n\) 的规模总和不超过 \(10^5\)
P2260 [清华集训2012] 模积和
求
\[\sum_{i=1}^{n} \sum_{j=1}^{m} (n \bmod i) \times (m \bmod j), i \neq j \]mod 19940417 的值
\(1 \le n,m \le 10^9\)。
P3579 [POI2014] PAN-Solar Panels
对于 \(n\) 组询问，每组询问给定四个整数 \(a,b,c,d\)，从区间 \([a,b]\) 和 \([c,d]\) 中任意选取两个整数 \(x\) 和 \(y\)，求 \(\gcd(x,y)\) 的最大值是多少。
\(1 \le n \le 1000\)，\(1 \le a \le b \le 10^9\)，\(1 \le c\le d \le 10^9\)。
10 欧拉函数
欧拉函数是一类非常重要的数论函数。
10.1 定义与性质
欧拉函数的定义为在 \([1,n]\) 中与 \(n\) 互质的整数的个数，记作 \(\varphi(n)\)。
显然，当 \(n\) 是质数时，有 \(\varphi(n) = n-1\)。
若 \(p\) 为质数，则 \(\varphi(p^k) = (p-1)\times p^{k-1}\)。
在 \([1,p^k]\) 中不是 \(p\) 的倍数的数都和 \(p^k\) 互质。因此 \(\varphi(p^k) = p^k - p^{k-1} = (p-1)p^{k-1}\)。
\(\varphi\) 是 积性函数。即若 \(\varphi(ab) = \varphi(a) \times \varphi(b)\)。
设与 \(a\) 互质的数为 \(\{ a_1,a_2,\cdots,a_{\varphi(a)} \}\)，那么在 \([1,ab]\) 内与 \(a\) 互质的整数可以表示为 \(i \times a + a_j(0 \le i \lt b,1 \le j \le \varphi(a))\)。
因为 \(a \bot b\)，所以 \(ia (0 \le i \lt b)\) 在模 \(b\) 意义下互不相同，即对于每一个 \(a_j\)，\(\{ i \times a + a_j \}(0 \le i \lt b)\) 均构成 \(b\) 的完全剩余系 \(\{ 0,1,\cdots,b-1 \}\)，即随着 \(i\) 从 \(0\) 变到 \(b-1\)，\(i \times a + a_j\) 模 \(b\) 的余数取遍 \(0 \sim b-1\)。
因此，对于每一个 \(a_j\)，满足 \(b \bot i \times a + a_j\) 的 \(i\) 的个数为 \(\varphi(b)\)，所以与 \(ab\) 互质的数的个数为 \(\varphi(ab) = \varphi(a) + \varphi(b)\)。得证。
根据算术基本定理，设 \(n\) 被唯一分解为 \(\prod_{i=1}^m {p_i}^{c_i}\)，则 \(\varphi(n) = n \times \prod_{i=1}^m \frac {p_i-1} {p_i}\)。
根据前面两个性质
得证。
当然，我们还可以用 容斥原理 来证明。
考虑去掉所有被至少一个 \(n\) 的质因子整除的数，再加上被至少两个 \(n\) 的质因子整除的数，以此类推，我们可以得到
其中 \(S\) 是 \(n\) 的质因子集合的子集。容易发现，这个式子就是
这样我们也可以反推到积性函数上面。设 \(a = \prod_{i=1}^{m_1} {p_i}^{c_i},b= \prod_{i=1}^{m_2} {q_i}^{d_i}\)。由于 \(a \bot b\)，所以 \(p_i,q_i\) 互不相同。因此
于是有
若 \(a \mid b\)，则 \(\varphi(ab) = a \times \varphi(b)\)。
根据前面 \(\varphi\) 的计算式，容易发现这时比较显然的。
于是我们也容易得到这样的结论：
若 \(a \mid b\)，则 \(\varphi(a) \mid \varphi(b)\)。
欧拉反演：
\[\sum_{d \mid n} \varphi(d) = n \]
考虑 \(\varphi(n)\) 的定义是与 \(n\) 互质的数的个数，那么容易发现 \(\varphi(n)\) 就是分母为 \(n\) 的最简基本分数的个数，而法里级数 \(\mathscr F_n\) 既包含分母不超过 \(n\) 的 \([0,1]\) 范围内的所有最简分数，我们尝试从这一点来论证欧拉反演。
我们把法里级数中非零的部分提出来，例如当 \(n=12\) 时，在约分成最简基本分数前，分母为 \(12\) 的所有基本分数集合是
约分得到
我们可以根据它们的分母将这些分数分组，就可以得到
对此，我们如何就是呢？
发现 \(12\) 的每一个因子都出现在分母上，与其一起出现的还有 \(\varphi(d)\) 个分子，从而
也就是
这就是 欧拉反演。
同时它可以被表示成狄利克雷卷积的形式，也就是 \(\varphi * 1 =\text{id}\)，这在杜教筛中颇有应用。
\(2 \mid \varphi(n) (n \gt 2)\)。
假设 \(x \bot n\)，那么一定有 \((n-x) \bot n\)。于是与 \(n\) 互质的数是两两配对的，唯一的特殊情况就是 \(x = n-x\) 的时候。
容易发现这种情况显然是满足条件的，即 \(x=\frac n 2\)，那么 \(\gcd(n,x) = x \neq 1\)。得证。这个性质在莫反的某道题中有所应用。
使得 \(\gcd(n,x) = d\) 且 \(1 \le x \le n\) 的 \(x\) 的个数为 \(\varphi \left( \frac nd \right)\)。
使得这个条件成立的条件可以被表示成 \(\frac n d \bot \frac x d\)，同时 \(1 \le \frac x d \le \frac n d\)，从而答案就是 \(\varphi \left( \frac n d \right)\)。得证
如果我们把条件转移到 \(0 \le x \lt n\)，也是成立的，因为我们之前定义了 \(\gcd(0,n) = \gcd(n,n) = n\)。
10.2 线性筛欧拉函数
根据上面欧拉函数的计算式，我们可以得到它的线性筛法，若 \(p \mid n\)，则
于是代码就是简单的
phi[1]=1;
for(int i=2;i<N;i++){
if(!vis[i]) p[++cnt]=i,phi[i]=i-1;
for(int j=1;j<=cnt&&i*p[j]<N;j++){
vis[i*p[j]]=true;
if(i%p[j]==0){
phi[i*p[j]]=phi[i]*p[j];
break;
}
else phi[i*p[j]]=phi[i]*(p[j]-1);
}
}
10.3 欧拉定理
欧拉发现费马的定理可以用下面的方式推广到非素数的模：
这和费马小定理是相当类似的，于是我们可以尝试用类似于费马小定理的方法来证明它。
\(\varphi(m)\) 个数 \(\{ kn \mod m \mid k \bot m 且 0 \le k \lt m \}\) 就是 \(\{ k \mid k \bot m 且 0 \le k \lt m \}\) 按照某种次序排列。
将它们相乘并用 \(\prod_{0 \le k \lt m,k \mid m} k\) 来除即可。
而当 \(m\) 是素数时，带入就可以得到费马小定理。
10.4 扩展欧拉定理
再把欧拉定理进行一个扩展，可以得到
感性理解：
当 \(\gcd(a,m) \neq 1\) 时，不断乘 \(a\) 会让 \(a^i\) 和 \(m\) 的 \(\gcd\) 不断变大，这到某个特定的 \(a^r\) 之后 \(\gcd\) 就不再变化，因为此时 \(\gcd(a^r,m)\) 受到 \(m\) 的每个与 \(a\) 相同的质因子的次数限制。
令这个 \(\gcd = d\)，考虑 \(a^r \bmod m,a^{r+1} \bmod m,\cdots\)，显然它们均是 \(d\) 的倍数，且除以 \(d\) 之后会取遍所有与 \(\frac m d\) 互质的数，这时因为 \(a \bot \frac m d\)。
根据欧拉定理，\(a^k \bmod \frac m d\) 有循环节 \(\varphi \left( \frac m d \right)\)，因此 \(\left( a^k \bmod {\frac m d} \right) \times d\) 即 \(a^k \bmod m\)，也有循环节 \(\varphi \left( \frac m d \right)\)，故从 \(a^r\) 开始，\(a^{r+k} \bmod m\) 有循环节 \(\varphi(p)\)。
比较抽象，非常感性。
严谨证明 摘自 OI-wiki：
-
\(a\) 的从 \(0\) 次，\(1\) 次到 \(b\) 次幂模 \(m\) 构成的序列中，存在 \(r,s\)，使得前 \(r\) 个数（即从 \(a^0 \bmod m\) 到 \(a^{r-1} \bmod m\)）互不相同，从 \(r\) 开始，每 \(s\) 个数就循环依次。
在 3.5 中，我们对这个遗留问题进行了讨论。
-
\(a\) 为素数的情况，该是成立。
证明：
-
若模 \(m\) 不能被 \(a\) 整除，那么有 \(\gcd(a,m) = 1\) 成立，根据欧拉定理容易证明。
-
若模 \(m\) 能被 \(a\) 整除，那么存在 \(r\) 和 \(m'\) 使得 \(m = a^r m'\)，且 \(\gcd(a,m')= 1\) 成立。所以根据欧拉定理有 \(a^{\varphi(m')} \equiv 1 \pmod {m'}\)。
又由于 \(\gcd(a^r,m') = 1\)，所以根据欧拉函数的计算式，容易得到 \(\varphi(m) = \varphi(m') \times (a-1)a^{r-1}\)，即我们有 \(\varphi(m') \mid \varphi(m)\)。
所以 \(a^{\varphi(m')} \equiv 1 \pmod {m'},\varphi(m') \mid \varphi(m) \Longrightarrow a^{\varphi(m)} \equiv 1 \pmod {m'}\)，即 \(a^{\varphi(m)} = km'+1\)，两边同时乘 \(a^r\)，得 \(a^{r+\varphi(m)} = km+a^r\)，因为 \(m = a^r m'\)。
所以 \(m\) 中素因子 \(a\) 的次数 \(r\) 满足：\(a^r \equiv a^{r+\varphi(m)} \pmod {m}\)。我们可以简单变换形式，得到 推论：
\[ b \gt r \Longrightarrow a^b \equiv a^{r+(b-r)\bmod \varphi(m)} \pmod m \]又由于 \(m=a^rm'\)，所以 \(\varphi(m) = \varphi(a^r) \varphi(m') \ge \varphi(a^r) =a^{r-1} (a-1) \ge r\)。
所以因为 \(\varphi(m) \ge r\)，故有：
\[ a^r \equiv a^{r+\varphi(m)} \equiv a^{r \bmod \varphi(m) + \varphi(m)} \pmod m \]所以
\[ \begin{aligned} a^b & \equiv a^{r+(b-r) \bmod \varphi(m)} \\ & \equiv a^{r \bmod \varphi(m) + \varphi(m) + (b-r) \bmod \varphi(m)} \\ & \equiv a^{\varphi(m) + b \bmod \varphi(m)} \pmod m \end{aligned} \]即 \(a^b \equiv a^{b \bmod \varphi(m) + \varphi(m)}\)。
-
-
\(a\) 为素数幂的情况，该式成立。
证明：
不妨令 \(a =p^k\)，是否依然有 \(\forall r,a^r \equiv a^{r+\varphi(m)} \pmod m\)？
答案是肯定的，由于第一点可知存在 \(s\) 使得 \(a^s \equiv 1 \pmod m\)，所以 \(p^{\text{lcm} (s,k)} \equiv 1 \pmod m\)，所以令 \(s' = \frac s {\gcd(s,k)}\) 时，我们能有 \(p^{s'k} \equiv 1 \pmod m\)。
此时有关系 \(s' \mid s\) 且 \(s \mid \varphi(m)\)，且 \(r' = \left \lfloor \frac r k \right \rfloor \le r \le \varphi(m)\)，由 \(r',s'\) 与 \(\varphi(m)\) 的关系，依然可以得到 \(a^b \equiv a^{b \bmod \varphi(m) + \varphi(m)} \pmod m\)。
-
\(a\) 为合数的情况，该是成立。
证明：
只证 \(a\) 拆成两个素数的幂的情况，大于两个的用数学归纳法可证。
设 \(a=a_1a_2\)，其中 \(a_i = {p_i}^{k_i}\)，而 \(a_i\) 的循环长度为 \(s_i\)；
则 \(s \mid \text{lcm} (s_1,s_2)\)，由于 \(s_1 \mid \varphi(m),s_2 \mid \varphi(m)\)，那么 \(\text{lcm}(s_1,s_2) \mid \varphi(m)\)，所以 \(s \mid \varphi(m),r=\max \left( \left \lfloor \frac {r_i} {k_i} \right \rfloor \right) \le \max(r_i) \le \varphi(m)\)；
由 \(r,s\) 与 \(\varphi(m)\) 的关系，依然可以得到 \(a^b \equiv a^{b \bmod \varphi(m) + \varphi(m)} \pmod m\)。
证毕。乐。
于是你就可以通过 P5091 【模板】扩展欧拉定理。代码。
10.5 例题
煎蛋题。
P2158 [SDOI2008] 仪仗队
11 莫比乌斯函数
尝试按照《具体数学》的顺序引入莫比乌斯反演。
11.1 引入
莫比乌斯函数 \(\mu(m)\) 是根据 19 世纪的数学家奥古斯特·莫比乌斯命名的，他还发现了著名的莫比乌斯带，\(\mu(m)\) 对所有整数 \(m \ge 1\) 由等式
来定义。
这个式子实际上是一个递归式，左边的 \(\mu(m)\) 和某些满足 \(d \lt m\) 的 \(\mu(d)\) 的值组成的和式。我们有相继带入前面 \(m\)，从而得到前面的值：\(1,-1,-1,0,-1,1,-1,0,0,1,-1,0,\cdots\)。
在 1857 年，Richard Dedekind 和 Joseph Liouville 注意到了如下重要的 “反演原理”。
这就是我们熟知的 莫比乌斯反演。
发现这个东西的证明是容易的：如果 \(g(m) = \sum_{d \mid m} f(d)\)，那么
反之，如果 \(f(m) = \sum_{d \mid m} \mu(d) g(\frac m d)\)，则
于是我们在这里就证明了 莫比乌斯反演。
而莫比乌斯反演还可以被表示成
现在回到最开始的那个式子，发现其实就是 \(\mu * 1 = \epsilon\)，由于我们之前证明的 积性函数的逆元也是积性函数，可以知道 \(\mu\) 必是积性函数。这样一来，只要我们计算出 \(\mu(p^k)\) 就能算出 \(\mu(m)\)。
当 \(m = p^k\) 时，根据最初的式子，我们知道对于所有 \(k \ge 1\)，都有
由此推出
于是，我们就推出了 \(\mu\) 的定义式
如果我们把之前的欧拉反演视作 \(\varphi (m)\) 的递归式，就能应用莫比乌斯反演得到
也就是 \(\mu * \text{id} = \varphi\)，也对狄利克雷卷积逆元有了另一个方面的解释。
接下来我们尝试用容斥的角度来理解莫比乌斯函数。
设 \(g(n) = \sum_{n \mid d} f(d)\)，即 \(g(n)\) 是 \(f\) 在所有 \(n\) 的倍数处的取值和。现在已知 \(g\)，要求 \(f(1)\)。则 \(f(1)\) 等于 \(f\) 在 \(1\) 的倍数处的取值和，减去在质数处的取值和，但是多减去了在两个不同质数乘积处的取值和，一次要加上这些值，但是多加了在三个不同质数乘积处的取值和，以此类推。因此，若 \(n\) 为 \(k\) 个不同质数的乘积，则 \(f(1)\) 会收到 \(g(n)\) 系数为 \((-1)^k\) 的贡献，如下图。
就是对 \(\mathbb N\) 作容斥，得到贡献系数 \(\mu\)。
11.2 线性筛莫比乌斯函数
根据定义是，线性筛莫比乌斯函数是非常容易的。
mu[1]=1;
for(int i=2;i<=maxn;i++){
if(!vis[i]) p[++cnt]=i,mu[i]=-1;
for(int j=1;j<=cnt&&i*p[j]<=maxn;j++){
vis[i*p[j]]=true;
if(i%p[j]==0) break;
mu[i*p[j]]=-mu[i];
}
}
11.3 常见技巧
对于 莫比乌斯反演，我们常常用式子
它把 \(\gcd(i,j)\) 带入了 \(n\)，还将 \(i,j\) 互质的条件转化为枚举 \(\gcd(i,j)\) 的约数 \(d\)，然后对 \(\mu(d)\) 求和。在 \(i,j\) 同样需要枚举的时候，可以先枚举 \(d\) 并计算合法的 \((i,j)\) 对数，这样 \(i,j\) 合法当且仅当 \(d \mid i,d \mid j\)，就把 \(i,j\) 独立开了。
于是我们就可以这样解决下面的式子：
相当于对最大公约数为 \(d\) 的倍数进行了容斥：加上最大公约数为 \(1\) 的倍数的对数，减去最大公约数为 \(p_i\) 的倍数的对数，加上最大公约数为 \(p_ip_j(i \neq j)\) 的倍数的对数，以此类推。得到每个 \(d\) 的贡献系数就是莫比乌斯函数。
而这个式子是非常容易用整除分块完成。
$$ d(ij) = \sum_{x \mid i} \sum_{y \mid j} \left [ x \bot y \right] $$ 考虑单个质因子 $p$，再用中国剩余定理和并。设 $a = v_p(i)$ 即 $i$ 中含有质因子 $p$ 的数量，$b = v_p(j)$，则 $v_p(ij)= a+b$。对于 $ij$ 的约数 $d$，若 $v_p(d) \le a$，则令其对应 $v_p(x) = v_p(d),v_p(y)=0$，即所有都算在 $x$ 上面；若 $v_p(d) \gt a$，则令其对应 $v_p(x)=0,v_p(y) = v_p(d) -a$。容易发现这样会形成双射，因此可证得这一结论。
11.4 例题
接下来，我们来将一些例题，由于比较多，所以单独提出来了。部分例题会用到上面的常见技巧。
前面的一些例题中有些写得比较简略，有些可能是👻乐团中的子问题，如果看不懂可以在那道题题解中找一下，或者直接问我。
P2522 [HAOI2011] Problem b
对于给出的 \(n\) 个询问，每次求有多少个数对 \((x,y)\)，满足 \(a \le x \le b\)，\(c \le y \le d\)，且 \(\gcd(x,y) = k\)，\(\gcd(x,y)\) 函数为 \(x\) 和 \(y\) 的最大公约数。
\(1 \le n,k \le 5 \times 10^4\)，\(1 \le a \le b \le 5 \times 10^4\)，\(1 \le c \le d \le 5 \times 10^4\)。
首先可以对式子做一个二维差分，转换成求下界为 \(1\) 的式子
我们把 \(i,j\) 都 \(\div k\) 就可以得到
这就是一个经典的莫反式子，直接作就可以得到
于是就可以用整除分块在 \(\mathcal O(n + T \sqrt n)\) 时间内完成。代码。
子问题：P3455 [POI2007] ZAP-Queries。
P2257 YY的GCD
给定 \(N, M\)，求 \(1 \leq x \leq N\)，\(1 \leq y \leq M\) 且 \(\gcd(x, y)\) 为质数的 \((x, y)\) 有多少对。
\(T = 10^4\)，\(N, M \leq 10^7\)。
问题就是求
我们用莫反相关知识对其进行推导可以得到
于是我们令 \(T = pd\)，那么可以得到
前面可以直接用整除分块完成，而后面的 \(\sum_{p \mid T \land p \in \mathbb P} \mu(\frac T p)\) 乍一眼是不好完成的。
但细想发现这是可以和莫比乌斯函数一起在线性筛中解决的，所以我们就可以直接预处理出其前缀和，从而做到时间复杂度 \(\mathcal O(n + T \sqrt n)\)。代码。
双倍经验：P2568 GCD。
SP5971 LCMSUM - LCM Sum
\(T\) 次询问，每次询问给出 \(n\)，求
\[\sum_{i=1}^n \text{lcm}(i,n) \]\(1 \le T \le 3 \times 10^5,1 \le n \le 10^6\)。
首先我们考虑把 \(\text{lcm}\) 转化成 \(\gcd\)，那么
我们设
只要我们求出这个，就可以计算上面的式子了。
首先我们可以暴力用莫反计算，也就是
我们发现后面可以被拆成 \(\mu * (\text{id} + 1) = \mu * \text{id} + \mu * 1 = \varphi + \epsilon\)，于是我们就可以得到
还有另外一种方法就是对与 \(n\) 互质的数进行分析，如果 \(x \bot n\) 则一定有 \((n-x) \bot n\)，那么我们把它两两配对，就有 \(\frac {n \varphi(n)} 2\)。
而显然 \(F(1)= 1\)。
现在原式就变成了
而这里的 \(\sum_{d \mid n}d \varphi(d) = 1 * (\text{id} \times \varphi)\) 是可以用线性筛筛出来的，于是总时间复杂度 \(\mathcal O(T+n)\)。代码。
P4318 完全平方数
\(n\) 次询问，每次询问给出一个 \(k\)，查询第 \(k\) 个不是完全平方数的正整数倍数的数。
\(1 \le k \le 10^9,1 \le n \le 50\)。
设 \(f(n)\) 表示 \(\le n\) 的非完全平方数倍数的数，那么我们最后只需要二分答案即可。
怎么求 \(f(n)\) 呢？
发现可以先去掉 \(4,9,\cdots\) 的倍数，而接下来就要加上 \((p_1p_2)^2\) 的倍数，这不难让我们想到莫比乌斯函数。
相当于做容斥，那么我们就有
直接计算，时间复杂度 \(\mathcal O(\sqrt n \log n)\)。代码。
P1829 [国家集训队] Crash的数字表格 / JZPTAB
给定 \(n,m\)，求
\[\sum_{i=1}^n \sum_{j=1}^m \text{lcm}(i,j) \]\(1 \le n,m \le 10^7\)。
后面的式子就和 SP5971 一样了。直接做就是线性的。代码。
P3704 [SDOI2017] 数字表格
Doris 刚刚学习了 fibonacci 数列。用 \(f_i\) 表示数列的第 \(i\) 项，那么
\[f_0=0,f_1=1 \]\[f_n=f_{n-1}+f_{n-2},n\geq 2 \]Doris 用老师的超级计算机生成了一个 \(n\times m\) 的表格，
第 \(i\) 行第 \(j\) 列的格子中的数是 \(f_{\gcd(i,j)}\)，其中 \(\gcd(i,j)\) 表示 \(i,j\) 的最大公约数。
Doris 的表格中共有 \(n\times m\) 个数，她想知道这些数的乘积是多少。
答案对 \(10^9+7\) 取模。
\(1 \leq T\leq 10^3\)，\(1\leq n,m\leq 10^6\)。
稍微推一下式子就是
预处理一下括号内的东西，整除分块之后可以做到单次 \(\mathcal O(\sqrt n)\)。
P5221 Product
给定 \(n\)，求
\[\prod_{i=1}^n \prod_{j=1}^n \frac {\text{lcm}(i,j)} {\gcd(i,j)} \pmod {104857601} \]\(1 \le n \le 10^6\)，时限 \(200ms\)。
有些人喜欢在写完👻之后来写这道题。
话不多说，直接暴力推一下式子就是了。
前面的直接阶乘处理一下，后面的我们枚举 gcd 就可以得到
于是两次整除分块就可以线性完成了，代码相当好些。代码。
P3911 最小公倍数之和
对于 \(a_1,a_2,\cdots,a_n\)，求
\[\sum_{i=1}^n \sum_{j=1}^n \text{lcm}(a_i,a_j) \]\(1 \le n,a_i \le 5 \times 10^4\)。
给定具体的数是不好处理的，所以我们考虑先令 \(c_i = \sum_{j=1}^n \left[a_j =i \right]\)，即 \(i\) 的出现次数，而 \(V\) 为值域。
那么
我们令 \(T = kd\) 就可以得到
于是只要我们令 \(f(T) =\sum_{d \mid T} \mu(d) d,g(T) = \sum_{i=1}^{\left \lfloor \frac V T \right \rfloor} i c_{iT}\)，原式就变成
其中 \(f\) 可以在线性筛中处理出来，而 \(g\) 可以采取暴力枚举倍数和狄利克雷后缀和做到 \(\mathcal O(V \ln V)\) 和 \(\mathcal O(V \ln \ln V)\)。
这样就做完了，查询时可以直接线性完成。代码。（暴力枚举计算 \(g\)）
双倍经验：AGC038C LCMs。
CF1285F Classical?
有 \(n\) 个整数 \(a_1,a_2,\cdots,a_n\)，求 \(\max_{1 \le i \lt j \le n} \left( \text{lcm}(a_i,a_j) \right)\)
\(2 \le n \le 10^5,1\ \le a_i \le 10^5\)。
求完了和，我们现在来求 \(\max\)。
考虑加上一个 \(\gcd(a_i,a_j) = 1\) 的条件。于是我们从大到小枚举每一个数，如果当前加入一个数 \(x\)，存在一个 \(\gt x\) 的数 \(y\) 满足 \(\gcd(x,y)=1\)，那么所有 \(x \le z \le y\) 中间的所有 \(z\) 都是没有用的，因为与他们互质的数比 \(x\) 小了。
容易发现这是可以用单调栈直接维护的，那么我们如何计算出是否有与 \(x\) 互质的数呢？
通过上面有关 \(\mu\) 的容斥原理，我们设 \(cnt_d\) 表示 \(d\) 的倍数的出现次数，那么答案就是
只要这个值不为 \(0\)，我们就可以一直弹栈找到与它互质的那个数。
时间复杂度 \(\mathcal O(\sum_{i=1}^n d(i))\)，也就是 \(\mathcal O(V \log V)\)。
而对于 \(\gcd \neq 1\) 的 \(\text{lcm}\) 如何计算？
我们可以枚举每一个 \(\gcd\) 都这样做一次，这样的时间复杂度是 \(\mathcal O(V \log^2 V)\) 的，但是它并不优秀。
发现如果 \(\gcd(x,y) \neq 1\)，那么对于 \(a \mid x,b \mid y\)，一定存在一组 \((a,b)\) 满足 \(\gcd(a,b)= 1,\text{lcm} (a,b) = \text{lcm}(x,y)\)，于是我们把所有 \(a_i\) 的因数都加到集合中，一起跑上面的算法，时间复杂度还是 \(\mathcal O(V \log V)\)。代码。
P3327 [SDOI2015] 约数个数和
给定 \(n,m\)，求
\[\sum_{i=1}^n \sum_{j=1}^m d(ij) \]\(1 \le T,n,m \le 50000\)。
在常见技巧中，我们列出了整个式子的转化。
于是式子就是
于是直接对这个式子莫反就可以得到
令 \(g(n) = \sum_{i=1}^n \left \lfloor \frac n i \right \rfloor\)，这时可以用整除分块预处理出来的。
于是式子就是
直接整除分块可以做到时间复杂度 \(\mathcal O(n \sqrt n + T \sqrt n)\)。代码。
P4619 [SDOI2018] 旧试题
给定 \(A,B,C\) 求
\[\left( \sum_{i=1}^A \sum_{j=1}^B \sum_{k=1}^C d(ijk) \right) \bmod \left(10^9+7 \right) \]\(1 \le A,B,C \le 10^5,1 \le T \le 10\)。
更为厉害的约数个数和。
多了一个变量，但是按照之前的转化方式是类似的
那么我们就可以对这个东西进行莫反了
我们可以令 \(fa(i) = \sum_{i \mid j} \left \lfloor \frac A j \right \rfloor,fb(i) = \sum_{i \mid j} \left \lfloor \frac B j \right \rfloor,fc(i) = \sum_{i \mid j} \left \lfloor \frac C j \right \rfloor\)，那么这是可以直接枚举倍数在 \(\mathcal O(n \ln n)\) 的时间复杂度内预处理出来的，计算时就变成了
然后你发现这还是 \(O(n^3)\) 的，莫反根本没有用！
真的没有用吗？这个式子在哪些地方是有值的呢？
容易发现我们首先需要保证 \(\mu(u),\mu(v),\mu(w)\) 都不是 \(0\)，并且 \(\text{lcm} (u,v) \le A,\text{lcm}(u,w) \le B,\text{lcm} (v,w) \le C\)。
那么也就是对于任意的两个点，它们可以成为一组就意味着 \(\mu(u) \neq 0,\mu(v) \neq 0,\text{lcm} (u,v) \le \max(A,B,C)\)，然后我们尝试把这样的二元组建立一条边，边权为 \(\text{lcm} (u,v)\)。
这个做法看起来根本不可行，并且直接这样建边是 \(\mathcal O(n^2)\) 的。
但是细想一下，这个条件是相当严的，要求 \(u,v,\gcd(u,v)\) 都不能出现平方因子。于是我们可以先枚举 \(\gcd(u,v)\)，再分别枚举 \(\frac u {\gcd(u,v)}\) 和 \(\frac v {\gcd(u,v)}\) 在枚举的过程中判断 \(\mu(\gcd),\mu(u),\mu(v) \neq 0\) 和 \(\frac u {\gcd(u,v)} \bot \frac v {\gcd(u,v)}\)。
每一次我们枚举了质数的子集，而 \(10^5\) 内只会出现 \(6\) 个不同的质数，所以这样枚举下来是 \(\mathcal O(m)\) 的！
下面给出这一部分的代码实现。
for(int g=1;g<=mx;g++)
for(int i=1;i*g<=mx;i++)
if(mu[i*g])
for(int j=i+1;1ll*i*j*g<=mx;j++)
if(mu[j*g]&&gcd(i,j)==1){
int u=i*g,v=j*g,w=i*j*g;
/*...*/
}
进而，我们发现这样的点对并不多，实测发现只有 760741 个！赢麻了！
那么现在成了怎么找三元环，可以尝试阅读 这篇博客。
也就是我们将无向边定向，每一条边从度数大的点连向度数小的，然后这样去枚举的时间复杂度就是 \(\mathcal O(m \sqrt m)\)，具体可以看上面的链接。
于是我们就做完了！直接枚举就可以了。
for(int x=1;x<=mx;x++){
for(auto i:G[x]) vis[i.fi]=i.se;
for(auto i:G[x]){
int y=i.fi,wxy=i.se;
for(auto j:G[y]){
if(!vis[j.fi]) continue;
int z=j.fi,wyz=j.se,wxz=vis[z],val=mu[x]*mu[y]*mu[z];
/*...*/
}
}
for(auto i:G[x]) vis[i.fi]=0;
}
这样就可以做到 \(O(n \ln n+m \ln m)\)，完全胜利！代码。
存在把上面式子化成 \(\mathcal O(n^2 \log n)\) 的做法，简单来说就是只提出 \(u\)，对后面枚举的 \(i,j\) 都除以 \(u\) 去按照之前我们做过的套路完成，但此做法对后面的推导并无意义，具体可以去看第三篇题解。
三倍经验：CF236B Easy Number Challenge 和 CF235E Number Challenge。
P5518 [MtOI2019] 👻乐团 / 莫比乌斯反演基础练习题
东风谷 早苗(Kochiya Sanae)非常喜欢幽灵乐团的演奏，她想对她们的演奏评分。
因为幽灵乐团有 \(3\) 个人，所以我们可以用 \(3\) 个正整数 \(A,B,C\) 来表示出乐团演奏的分数，她们的演奏分数可以表示为
\[\prod_{i=1}^{A}\prod_{j=1}^{B}\prod_{k=1}^{C}\left(\frac{\text{lcm}(i,j)}{\gcd(i,k)}\right)^{f(type)} \]因为音乐在不同的部分会有不同的听觉感受，所以 \(type\) 会在 \(\{0,1,2\}\) 中发生变化，其中：
\[\begin{aligned} f(0)&=1 \cr f(1)&=i \times j \times k \cr f(2)&=\gcd(i,j,k) \end{aligned}\]因为乐团的歌实在太好听了，导致分数特别高，所以她们的分数要对给定的正整数 \(p\) 取模。
因为有很多歌曲要演奏，所以早苗给出了 \(T\) 组询问。
\[1\leq A,B,C\leq 10^5 \ \ \ \ 10^7 \leq p \leq 1.05\times 10^9\ \ \ \ p\in \{ prime\} \ \ \ \ T =70 \]
再来一个三重循环的题。
众所周知，cdqz 机房有且仅有三个👻，它们是：
-
：有名字 的👻，也是👻乐团的团长！ -
：没名字 的👻。 -
：直接消失的👻，这样的👻常常从你身边飘过！
这三只👻就是 cdqz 机房的真神！/bx
题目名字都告诉你了，还是赶紧来推式子吧。
首先对于 Type = 0/1 是简单的，直接顺着思路往下面推就可以了，也就是这个样子：
对于 Type = 0，我们首先把 \(\text{lcm}\) 给处理掉再把四个量分离出来就是
前面两项是可以直接用阶乘完成的，也就是说我们只需要求
这个时候你可以把上面的指数部分单独提出来用一个函数完成，这样就可以做到线性，是可以过的，也就是后面代码1给出的（把三个 Type 分开了写的）。
但是你推到后面会发现这其实也是 Type = 2 的子问题，而这个式子是可以 \(\mathcal O(\sqrt n)\) 解决的，所以还是再推一步，这一步也是比较显然的，就是你把 \(t\) 提下来再用 \(T\) 代替 \(td\)：
而中间那一块是可以直接用 \(\mathcal O(n \log n)\) 的时间复杂度预处理出来的（直接暴力枚举），所以单次求就是 \(\mathcal O(\sqrt n)\)。
恭喜你，获得了 \(20pts\) 的高分！代码。
现在我们来研究 Type = 1，用同样的思路，写出来的式子也是简单的
把前面两项提出来就是形如
预处理一下 \(\prod_{i=1}^n i^i\) 可以直接用快速幂求出。
而对于后面两项，也就是求
这个东西推导方式就较为类似了，笔者并没有把它推导 \(\mathcal O(\sqrt n)\)，而是把指数单独拿出来处理了。
于是可以把指数部分提出来，就变成了求
然后两次整除分块就可以完成。
很明显和上面 Type = 0 的优化类似，这时可以做到 \(\mathcal O(\sqrt n)\) 的，但是毕竟瓶颈也不在这里，式子也比较复杂，你直接这样作就好了。
预处理一下 \(\prod_{i=1}^n i^{i^2}\) 和 \(\sum_{i=1}^n \mu(i)i^2\) 即可，于是你就有 \(60pts\) 啦！代码。
最后来搞 Type = 2 的情况，按照相同的思路先展开一下
对于前两项，其实就是（后面的 \(\min\) 都代表 \(\min(A,B,C)\)）
像往常一样，我们直接用 \(T\) 代替 \(td\)，而后面的 \(\sum_{d \mid T} d \mu\left( \frac T d \right)\) 也就是 \(\text{id} * \mu = \varphi\)，我们就可以得到
把 \(T\) 提出去就变成了
容易发现这两项都可以用整除分块做出来的，稍微预处理一下 \(\varphi\) 的前缀和就可以了。
而对于后面两项，也就是求
按照之前的经验，我们尝试去枚举 \(\gcd(i,j)\)，令
然后就有
你可以自己尝试推一下，就会发现只能用 \(O(n \log n)\) 的时间复杂度完成单次的查询，根本过不了！
那怎么办呢？
式子中还存在一个 \(\gcd\)，那么我们就尝试去枚举 \(\gcd(i,j,k)\) 于是式子就变成了：
上面的内容和之前推第一部分类似，而现在我们发现其实中间的
其实就是 Type = 0 中我们推的 \(f\) 函数，于是整个式子就变成了
惊奇地发现这两部分前面的那一块其实是一样的，于是它们就抵消啦！
所以我们需要求的就是
和
于是你就可以用两次整除分块线性求出来啦！
这样你就推出了本题的全部式子，写的时候注意细节，在指数上面的东西注意模的是 \(mod-1\) 而不是 \(mod\)，整除分块不要手误了就行。
本题存在对于 Type = 2 的情况求 \(\ln\) 再 \(\exp\) 回去的做法，具体可以去本题题解区。但笔者还是更喜欢莫反做法。
P3700 [CQOI2017] 小 Q 的表格
为了完成任务，小 Q 需要列一个表格，表格有无穷多行，无穷多列，行和列都从 \(1\) 开始标号。为了完成任务，表格里面每个格子都填了一个整数，为了方便描述，小 Q 把第 \(a\) 行第 \(b\) 列的整数记为 \(f(a, b)\)。为了完成任务，这个表格要满足一些条件：
- 对任意的正整数 \(a, b\)，都要满足 \(f(a, b) = f(b, a)\)；
- 对任意的正整数 \(a, b\)，都要满足 \(b \times f(a, a + b) = (a + b) \times f(a, b)\)。
为了完成任务，一开始表格里面的数很有规律，第 \(a\) 行第 \(b\) 列的数恰好等于 \(a \times b\)，显然一开始是满足上述两个条件的。为了完成任务，小 Q 需要不断的修改表格里面的数，每当修改了一个格子的数之后，为了让表格继续满足上述两个条件，小 Q 还需要把这次修改能够波及到的全部格子里都改为恰当的数。由于某种神奇的力量驱使，已经确保了每一轮修改之后所有格子里的数仍然都是整数。为了完成任务，小 Q 还需要随时获取前 \(k\) 行前 \(k\) 列这个有限区域内所有数的和是多少，答案可能比较大，只需要算出答案 \({} \bmod 1,000,000,007\) 之后的结果。
\(1 \le m \le {10}^4\)，\(1 \le a, b, k \le n \le 4 \times {10}^6\)，\(0 \le x \le {10}^{18}\)。
建议盯着式子发呆。
首先，由于 \(f(a,b) = f(b,a)\) 所以我们只需要考虑半张表。
而你把第二个式子换一种写法就是
多写几个就变成了
那么换一种写法就是
我们发现这个东西和辗转相除求 gcd 是一样的，那么我们就可以得到这个很好的性质
那么
于是我们要求的答案其实就是
那么现在就是希望求
我们需要预处理出所有的 \(g(n)\) 才能在查询时以高效的方法完成。
但是按照上一道题我们推式子的方法，发现单次是 \(\sqrt n\) 的，根本过不了！
这时我们考虑在 SP5971 中求的东西
带进去就是
于是这个东西就可以线性递推出来了，只要我们能够动态维护 \(f(d,d)\) 的区间和，那么最后整除分块一下就可以在 \(\mathcal O(\sqrt n)\) 的时间复杂度内求出来了。
但是我们发现上面所说的 \(f(d,d)\) 的区间和需要做到 \(\mathcal O(1)\) 查询，由于 \(m \le 10^4\)，所以我们不难想到用 \(\mathcal O(\sqrt n) - \mathcal O(1)\) 的分块来维护。
具体每个块每个元素维护块内的前缀和，在每一个块末尾维护块间的前缀和，修改时暴力修改即可。
总时间复杂度 \(\mathcal O(m \sqrt n)\)，这道题为我们引入了一类式子新的推导方式。代码。
12 奇妙筛法
已开工。
12.1 杜教筛
12.2 Min-25 筛
Min_25 筛可以在 \(\mathcal O \left( \frac {n^{\frac 3 4}} {\log n} \right)\) 的时间复杂度内解决一类 积性函数 前缀和问题。
说形象点就是 \(1s\) 能跑 \(10^{10}\)，相当优秀。
要求：
- \(f(p)\) 是关于 \(p\) 的低阶多项式，或者可以快速求出。
- \(f(p^c)\) 可以快速求出。
以下 \(p\) 代表质数，\(p_i\) 表示第 \(i\) 小的质数，而 \(p_0 = 1\)。
同时我们钦定 \(lpf(i)\) 表示 \(i\) 的最小质因子。
首先，我们尝试对其进行推导，分成质数和合数两个部分
前面的部分是枚举素数，而后面的部分就是对于每一个合数枚举了其最小的质因子以及其次数。
考虑第一块的质数部分如何求值？
由于我们在前面要求了 \(f(p)\) 是关于 \(p\) 的低阶多项式，也就是
那么
现在问题就变成了求 \(\sum_{1 \le p \le n} p^{k}\) 了。
我们只对一个 \(k\) 拿出来讨论，如果多个就分别做以下的操作就可以了。
我们尝试用一个神秘的 dp，设
（至于它有什么用，你先别急）
我们考虑 \(g(n,j)\) 是如何转移过来的。
不难想到用 \(g(n,j-1)\) 进行转移，也就是减去了 \(lpf(i) = p_j\) 的部分，它可以被表示成
由于 \(i^k\) 是完全积性函数，我们不需要讨论前后是否互质！
而仔细一想，发现我们把质数多算了一次，它一共出现了 \(g \left( p_{j-1},j-1 \right)\) 次，那么实际上减去的就是
其实后面的 \(g(p_{j-1},j-1)\) 就是前 \(j-1\) 个质数的 \(k\) 次方和。
这样一来，我们就有了递推式
发现这是可以依次 dp 下去的，只要我们知道
就可以用质数处理出来。（但是时间复杂度根本过不去啊？你先别急）
回到引入 \(g\) 时的问题，它有什么用呢？
我们发现如果 \(p_x\) 是最大的 \(\le \sqrt n\) 的质数，那么 \(g(n,x)\) 就是我们想要的答案。
所以如果通过 dp 推出 \(g(n,x)\) 的值就完全胜利，而我们把它简记成 \(g(n)\)。（具体细节你先别急）
假设我们现在已经可以快速算出 \(g\) 了，那原式会变成什么样子呢？
发现后面的那一块也是不好处理的，所以我们采用类似的 dp 方式，设 \(S(n,x)\) 表示 \(1 \sim n\) 中 \(lpf(i) \gt p_x\) 的函数值和，容易发现答案就是 \(S(n,0)\)，而我们有
其中 \(sp_x\) 表示前 \(x\) 个素数的答案，而容易发现当 \(p_x \gt \sqrt n\) 是没有意义的，因为只会剩下一些大质数，所以我们的 \(p_x\) 是 \(\sqrt n\) 级别的，完全可以预处理出来。
而后面的 \([e \gt 1]\) 则表示当前的 \(i= p^e\) 的情况。
这样的递归是可以直接求的，时间复杂度为 \(\mathcal O \left( \frac {n^{\frac 3 4}} {\log n} \right)\) 的。具体证明可见 OI-wiki。
现在我们回到对 \(g\) 的求解，具体怎么求？
发现直接处理出 \(1 \sim n\) 的所有值显然不太可能，\(n\) 太大了，但仔细观察两个式子，发现每次需要用到的就是
这不就是类似于整除分块的东西么？
我们在第一章证明过
所以在整个过程中只会用到 \(1,2,\cdots,\sqrt n,\frac n {\sqrt n} ,\cdots ,\sqrt n\) 位置的值，就是整除分块的 \(\left \lfloor \frac n i \right \rfloor\) 那些值。
在第九章我们知道个数是 \(2 \sqrt n\) 的，于是我们只需要对这 \(\mathcal O(\sqrt n)\) 个 \(g\) 值进行 dp 预处理就可以了。
再来讲一讲具体的实现细节，我们用 \(id1[]\) 和 \(id2[]\) 分别记录 \(i \le \sqrt n\) 和 \(i \gt \sqrt n\) 的值在 \(g\) 数组中的下标，而对于 \(i \gt \sqrt n\) 的部分直接维护 \(\frac n i\) 就可以了。
同时，由于 \(1\) 既不是质数也不是合数，所以在上面的推到中它都不应该被算进去，也就是 \(g\) 和 \(S\) 都把 \(1\) 抠掉了，所以最后的答案是 \(S(n,0) + 1\)。
于是就可以通过 P5325 【模板】Min_25 筛 了。代码。
const ll H=1e9+7;
int cnt=0,tot=0;
ll p[N],g1[N],g2[N],s1[N],s2[N],w[N],id1[N],id2[N],B,n,i6,i2;
bool vis[N];
ll qpow(ll a,ll b){
ll res=1;
while(b){if(b&1) res=res*a%H;a=a*a%H,b>>=1;}
return res;
}
ll S1(ll x){x%=H;return x*(x+1)%H*i2%H;}
ll S2(ll x){x%=H;return x*(x+1)%H*(2*x+1)%H*i6%H;}
ll pos(ll x){return x<=B?id1[x]:id2[n/x];}
ll S(ll x,int y){
if(x<=p[y]) return 0;
int k=pos(x),res=((g2[k]-g1[k]-s2[y]+s1[y])%H+H)%H;
for(int i=y+1;i<=cnt&&p[i]*p[i]<=x;i++)
for(ll j=1,nw=p[i];nw<=x;++j,nw*=p[i]){
ll z=nw%H;
res=(res+(z-1)*z%H*(S(x/nw,i)+(j>1))%H)%H;
}
return res;
}
int main(){
ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);
cin>>n,B=sqrt(n);
i2=qpow(2,H-2),i6=qpow(6,H-2);
for(int i=2;i<=B;i++){
if(!vis[i]){
p[++cnt]=i;
s1[cnt]=(s1[cnt-1]+i)%H;
s2[cnt]=(s2[cnt-1]+1ll*i*i%H)%H;
}
for(int j=1;j<=cnt&&i*p[j]<=B;j++){
vis[i*p[j]]=true;
if(i%p[j]==0) break;
}
}
for(ll l=1,r;l<=n;l=r+1){
r=min(n/(n/l),n);
w[++tot]=n/l;
g1[tot]=S1(n/l)-1,g2[tot]=S2(n/l)-1;
if(n/l<=B) id1[n/l]=tot;
else id2[n/(n/l)]=tot;
}
for(int i=1;i<=cnt;i++){//g 数组的 dp，用了滚动数组
for(int j=1;p[i]*p[i]<=w[j]&&j<=tot;j++){
int cur=pos(w[j]/p[i]);
g1[j]=(g1[j]-p[i]*(g1[cur]-s1[i-1])%H+H)%H;
g2[j]=(g2[j]-p[i]*p[i]%H*(g2[cur]-s2[i-1])%H+H)%H;
}
}
cout<<(S(n,0)+1)%H;
return 0;
}
现在我们换几个函数来筛，也就是我们最熟悉的 \(\mu\) 和 \(\varphi\)。
它们在质数处的取值也是相当容易的，\(\mu(p) = -1,\varphi(p)=p-1\)，这样的低阶多项式是可以较为方便计算的。
直接把两个分开写出来就是这样的
const int N=1e6+5;
int cnt=0,tot=0;
ll p[N],w[N],id1[N],id2[N],n,B;
bool vis[N>>2];
int pos(ll x){return x<=B?id1[x]:id2[n/x];}
void Pr(int n){
memset(vis,0,sizeof(vis));
for(int i=2;i<=n;i++){
if(!vis[i]) p[++cnt]=i;
for(int j=1;j<=cnt&&i*p[j]<=n;j++){
vis[i*p[j]]=true;
if(i%p[j]==0) break;
}
}
}
namespace Mu{
ll g[N];
ll S(ll x,int y){
if(x<=p[y]) return 0;
ll res=-g[pos(x)]+y;
for(int i=y+1;i<=cnt&&p[i]*p[i]<=x;i++)
res-=S(x/p[i],i);
return res;
}
void SOL(){
for(int i=1;i<=tot;i++) g[i]=w[i]-1;
for(int i=1;i<=cnt;i++){
for(int j=1;p[i]*p[i]<=w[j]&&j<=tot;j++){
int cur=pos(w[j]/p[i]);
g[j]+=i-g[cur]-1;
}
}
cout<<(S(n,0)+1)<<'\n';
}
}
namespace Phi{
ll g1[N],g2[N],s[N];
ll S(ll x,int y){
if(x<=p[y]) return 0;
int k=pos(x);
ll res=g1[k]-g2[k]-(s[y]-y);
for(int i=y+1;i<=cnt&&p[i]*p[i]<=x;i++){
for(ll j=1,nw=1;nw*p[i]<=x;++j,nw*=p[i])
res+=nw*(p[i]-1)*(S(x/(nw*p[i]),i)+(j>1));
}
return res;
}
void SOL(){
for(int i=1;i<=cnt;i++) s[i]=s[i-1]+p[i];
for(int i=1;i<=tot;i++) g1[i]=1ll*w[i]*(w[i]+1)/2-1,g2[i]=w[i]-1;
for(int i=1;i<=cnt;i++)
for(int j=1;p[i]*p[i]<=w[j]&&j<=tot;j++){
int cur=pos(w[j]/p[i]);
g1[j]-=p[i]*(g1[cur]-s[i-1]);
g2[j]+=i-1-g2[cur];
}
cout<<(S(n,0)+1)<<' ';
}
}
void SOLVE(){
cin>>n,B=sqrt(n),cnt=tot=0;
memset(vis,0,sizeof(vis));
Pr(B);
for(ll l=1,r;l<=n;l=r+1){
r=min(n,n/(n/l));
w[++tot]=n/l;
if(n/l<=B) id1[n/l]=tot;
else id2[n/(n/l)]=tot;
}
Phi::SOL(),Mu::SOL();
}
int main(){
/*2024.4.26 H_W_Y P4213 【模板】杜教筛 Min_25 筛*/
ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);
int _;cin>>_;
while(_--) SOLVE();
return 0;
}
但结果并不理想：（前两个是 Min_25，最下面的是杜教筛）
于是我们开始了常数优化。
优化也是简单的，你把预处理放在一起，就可以跑到 2.95s：
进而，把两个函数一起做（再把不需要 long long 的地方换成 int），于是就可以做到 1.44s！！！
代码。
const int N=1e6+5;
int cnt=0,tot=0,p[N],n,B,w[N],id1[N],id2[N],s[N];
ll g1[N],g2[N];
bool vis[N>>2];
inline int pos(int x){return x<=B?id1[x]:id2[n/x];}
inline void Pr(int n){
for(int i=2;i<=n;i++){
if(!vis[i]) p[++cnt]=i,s[cnt]=s[cnt-1]+i;
for(int j=1;j<=cnt&&i*p[j]<=n;j++){
vis[i*p[j]]=true;
if(i%p[j]==0) break;
}
}
}
inline pll S(int x,int y){
if(x<=p[y]) return {0,0};
int k=pos(x);
pll res={g1[k]-g2[k]-(s[y]-y),-g2[k]+y};
for(int i=y+1;i<=cnt&&p[i]*p[i]<=x;i++)
for(int j=1,nw=1;1ll*nw*p[i]<=x;++j,nw*=p[i]){
pll cur=S(x/(nw*p[i]),i);
if(nw==1) res.se-=cur.se;
res.fi+=1ll*nw*(p[i]-1)*(cur.fi+(j>1));
}
return res;
}
inline void SOLVE(){
cin>>n,B=sqrt(n),tot=0;
for(ll l=1,r;l<=n;l=r+1){
r=n/(n/l);
w[++tot]=n/l;
g1[tot]=1ll*(n/l)*(n/l+1)/2-1,g2[tot]=n/l-1;
if(n/l<=B) id1[n/l]=tot;
else id2[n/(n/l)]=tot;
}
for(int i=1;i<=cnt&&p[i]<=B;i++)
for(int j=1;p[i]*p[i]<=w[j]&&j<=tot;j++){
int cur=pos(w[j]/p[i]);
g1[j]-=1ll*p[i]*(g1[cur]-s[i-1]);
g2[j]+=i-1-g2[cur];
}
pll res=S(n,0);
cout<<(res.fi+1)<<' '<<(res.se+1)<<'\n';
}
int main(){
/*2024.4.26 H_W_Y P4213 【模板】杜教筛 Min_25 筛*/
ios::sync_with_stdio(false);cin.tie(0);cout.tie(0);
int _;cin>>_;Pr(46340);
while(_--) SOLVE();
return 0;
}
12.3 Powerful Number 筛
12.4 扩展埃氏筛
13 阶与原根
模 \(n(n \gt 1)\) 的简化剩余系在模 \(n\) 乘法意义下封闭且构成群。容易证明其满足封闭性和结合律，存在逆元和单位元。将群相关定义应用在其上，得到阶和原根（循环群和生成元）的概念。
——Alex_Wei
13.1 阶
由欧拉定理可知，对于 \(a \in \mathbb Z,m \in \mathbb N^*\)，若 \(a \bot m\)，则有
因此满足同余式 \(a^n \equiv 1 \pmod m\) 的最小正整数 \(n\) 存在，这个 \(n\) 称为 \(a\) 模 \(m\) 的阶，记作 \(\delta_m(a)\) 或 \(\text{ord}_m(a)\)。
阶有一些性质。
\(a,a^2,\cdots,a^{\delta_m(a)}\) 模 \(m\) 两两不同余。
这和我们在 3.5 里面探讨的遗留问题较为相似。
考虑反证，假设存在 \(i \neq j\)，且 \(a^i \equiv a^j \pmod m\)，则有 \(a^{|i-j|} \equiv 1 \pmod p\)。
但是显然有 \(0 \lt |i-j| \lt \delta_m(a)\)，与阶的定义矛盾。得证。
若 \(a^n \equiv 1 \pmod m\)，则 \(\delta_m(a) \mid n\)。
设 \(n = \delta_m(a)q+r,0 \le r \lt \delta_m(a)\)。
若 \(r \gt 0\)，则
与阶的定义矛盾（最小性）。得证。
设 \(m \in \mathbb N^*,a,b \in \mathbb Z,a \bot m,b \bot m\)，则
\[\delta_m(ab)= \delta_m(a)\delta_m(b) \]的充分必要条件是
\[\delta_m(a) \bot \delta_m(b) \]
-
必要性：因为 \(a^{\delta_m(a)} \equiv 1 \pmod m\) 和 \(b^{\delta_m(b)} \equiv 1\) 可知
\[(ab)^{\text{lcm}(\delta_m(a),\delta_m(b))} \equiv 1 \pmod m \]由于上面的性质，我们知道
\[\delta_m(ab) \mid \text{lcm}(\delta_m(a),\delta_m(b)) \]所以我们有
\[\delta_m(a) \delta_m(b) \mid \text{lcm}(\delta_m(a),\delta_m(b)) \]即 \(\gcd(\delta_m(a),\delta_m(b)) = 1\)。
-
充分性：因为
\[(ab) ^{\delta_m(ab)} \equiv 1 \pmod m \]那么
\[1 \equiv (ab)^{\delta_m(ab) \delta_m(b)} \equiv a^{\delta_m(ab) \delta_m(b)} \pmod m \]所以
\[\delta_m(a) \mid \delta_m(ab) \delta_m(b) \]又由于 \(\gcd(\delta_m(a),\delta_m(b)) = 1\)，
\[\delta_m(a) \mid \delta_m(ab) \]同理可得
\[\delta_m(b) \mid \delta_m(ab) \]所以
\[\delta_m(a) \delta_m(b) \mid \delta_m(ab) \]而我们有
\[(ab)^{\delta_m(a) \delta_m(b)} \equiv (a^{\delta_m(a)})^{\delta_m(b)} \times (b^{\delta_m(b)})^{\delta_m(a)} \equiv 1 \pmod m \]故
\[\delta_m(ab) \mid \delta_m(a) \delta_m(b) \]所以
\[\delta_m(ab) = \delta_m(a) \delta_m(b) \]
得证。
设 \(k \in \mathbb N,m \in \mathbb N^*,a \in \mathbb Z,\gcd(a,m) = 1\)，则
\[\delta_m(a^k) = \frac {\delta_m(a)} {\gcd(\delta_m(a),k)} \]
很重要的一个性质。
注意到
而由于 \(a^{\delta_m(a)} \equiv 1 \pmod m\)，
所以
故
得证。
容易发现其实上面的证明方法都挺类似的。/kk
那如何求 \(\delta_m(a)\) 呢？（其中 \(m \bot a\)）
-
直接使用 BSGS 求 \(a^x \equiv 1 \pmod m\) 的最小正整数解。时间复杂度 \(\mathcal O(\sqrt m)\)。
-
根据阶的性质，对于 \(x = k \delta_m(a)\)，必然有 \(a^x \equiv 1 \pmod m\)。
因此，我们考虑先让 \(x =\varphi(m)\)，然后对于 \(\varphi(m)\) 的每个质因子 \(p\)，用 \(x\) 不断试除 \(p\) 知道 \(x\) 无法整除 \(p\) 或者 \(a^{\frac a p} \not \equiv 1 \pmod m\)。
时间复杂度 \(\mathcal O(\sqrt m + \log^2 m)\)，因为需要分解质因数和求解欧拉函数，用 Pollard_rho 可以优化为 \(m^{\frac 1 4}\)。
若 \(m\) 不大，则可以 \(\mathcal O(m)\) 线性处理每个数最小质因子，单次查询即可 \(\log^2m\)。
13.2 原根
原根——数论的神！
设 \(m \in \mathbb N^*,g \in \mathbb Z\)。若 \(\gcd(g,m) = 1\)，且 \(\delta_m(g) = \varphi(m)\)，则称 \(g\) 为模 \(m\) 的原根。
我们需要原根的原因是它可以 生成 模 \(m\) 的缩系，即它的若干次幂取遍了与 \(m\) 互质的数。这使得我们在考虑解与模数互质且模数存在原根的同余方程时可以用原根的若干次幂代替未知数。这一点在求解高次剩余时非常有用。
原根判定定理：设 \(m \ge 3,\gcd(g,m) = 1\)，则 \(g\) 是模 \(m\) 的原根的充要条件是，对于 \(\varphi(m)\) 的每个素因子 \(p\)，都有 \(g^{\frac {\varphi(m)} {p}} \equiv 1 \pmod m\)。
这和我们上面讲的求阶的方法非常像。
-
必要性：显然。
-
充分性：因为 \(\frac {\varphi(m)} {p}\) 的所有因子取遍了 \(\varphi(m)\) 除了其本身的所有因子，所以若存在 \(a^d \equiv 1 \pmod m\)，使得 \(d \lt \varphi(m)\)，则必然存在 \(\varphi(m)\) 的质因子 \(p\) 使得 \(d \mid \frac {\varphi(m)} p\)，这就说明存在 \(a^{\frac {\varphi(m)} p} \equiv 1 \pmod m\)，与假设矛盾。得证
原根个数：若一个数 \(m\) 有原根，则它原根的个数为 \(\varphi(\varphi(m))\)。
若 \(m\) 有原根 \(g\)，则：
所以如果 \(\gcd(\varphi(m),k)=1\)，则有 \(\delta_m(g^k) = \varphi(m)\)，也就是 \(g^k\) 也是模 \(m\) 的原根。
而满足 \(\gcd(\varphi(m),k) = 1\) 且 \(1 \le k \le \varphi(m)\) 的 \(k\) 有 \(\varphi(\varphi(m))\) 个。得证。
原根存在定理：一个数存在原根当且仅当 \(m=2,4,p^a,2p^a\)，其中 \(p\) 是奇素数，\(a \in \mathbb N^*\)。
可以当一个结论背住，具体证明比较复杂，可见 OI-wiki。
最小原根的范围估计：
在 1959 年，中国数学家王元证明了最小原根的级别为 \(\mathcal O(m^{0.25+\epsilon})\)。结合原根存在定理和原根判定定理，我们得以在 \(\mathcal O(m)\) 预处理每个数的最小质因子后单次 \(\mathcal O(m^{0.25} \omega(m) \log m)\) 求出 \(m\) 的最小原根，从而求出 \(m\) 的所有原根。
关于原根，还有一个有趣且重要的性质：
若 \(m\) 存在原根 \(g\)，则 \(n = \varphi(m)\)，则 \(m\) 的简化剩余系与 \(n\) 次单位根同构。
这是因为对于任何一个 \(a \bot m\)，其均可以写成 \(g^k\) 的形式。而任何一个 \(n\) 次单位根 \(\omega_k\) 也可以写成 \(\omega_n^k\) 的形式。
同时，\(g^k \equiv 1 \pmod p\)，且 \(\omega_n^n = 1\)。因此，我们可以将单位圆应用在 \(m\) 的简化剩余系上。
这也是为什么对于特定模数，我们可以使用原根代替 \(n\) 次单位根进行 FFT，所以当需要使用 \(2^k\) 次单位根时，若 \(k\) 不超过 \(\varphi(p)\) 当中质因子 \(2\) 的个数且 \(p\) 存在原根，即可用原根代替 \(2^k\) 次单位根。
于是就得到了 NTT。
13.3 例题
放两道简单题。
P5605 小 A 与两位神仙
这个游戏的规则是这样的：两个神仙先选定一个正整数 \(m\)，保证 \(m\) 是一个奇质数 \(p\) 的正整数次幂。然后进行 \(n\) 轮游戏，每轮中造梦者选定一个正整数 \(x\)，杰瑞米选定一个正整数 \(y\)，保证 \((x, m) = 1, (y, m) = 1\)，即 \(x\) 与 \(m\) 互质，\(y\) 与 \(m\) 互质，接下来询问小 A 是否存在非负整数 \(a\) 使得 \(x^a \equiv y \pmod{m}\)。
神仙们说小 A 只有在每一轮游戏中都回答正确才能回到正常的生活中，不得已之下他只好求助于聪明的你。
\(1\le n\le 2\times 10^4\)，\(3\le m \le 10^{18}\)，\(1 \le x, y < m\) 。
注意到模数 \(p\) 是奇素数的幂，所以存在原根 \(g\)。
我们尝试将 \(x\) 和 \(y\) 表示成 \(g^X\) 和 \(g^Y\) 的形式。
14 高次剩余问题
15 卢卡斯定理
后记
参考文献：
Pollard_Rho 算法 - cjTQX - 博客园 (cnblogs.com)
《具体数学》第四章
