质数
试除法判定质数
判定数n是否是质数。
如果d是n的因子，那么n/d也是n的因子，$d\leq\frac{n}{d}$，所以$d\leq\sqrt{n}$，故从1枚举到$\sqrt{n}$。
$O(\sqrt{n})$
bool is_prime(int x) {
if (x < 2)return false;
if (x == 2)return true;
for (int i = 2; i <= x / i; i++) { //等价于i*i<=x,i<=sqrt(x)
if (x % i == 0)return false;
}
return true;
}
试除法分解质因数
算术基本定理：所有的整数都可以分解成若干个质因子的乘积。
$n=P_1^{\alpha_{1}}P_2^{\alpha_{2}}…P_i^{\alpha_{i}}$$，P_i$是质数，$\alpha_i>0,n>1$。
$n$中最多只包含一个大于$\sqrt{n}$的质因子，所以只需枚举$[2,\sqrt{n}]$的质因子，然后特判一下$n$最后是否大于1。
$O(\sqrt{n})$
void divide(int x) {
for (int i = 2; i <= x / i; i++) {
if (x % i == 0) {
int s = 0;
while (x % i == 0)x /= i, s++;
cout << i << ' ' << s << endl; // P_i,alpha_i
}
}
if (x > 1)cout << x << ' ' << 1 << endl;// P_i,alpha_i
}
朴素筛法求素数（埃氏筛）
时间复杂度$O(nloglogn)$，缺点：同一个合数可能被多个质数多次标记。
时间复杂度计算：
如果不加continue优化，也就是把所有数的倍数都筛掉。计算次数为$\frac{n}{2}+\frac{n}{3}+…+\frac{n}{n}=n(\frac{1}{2}+\frac{1}{3}+…+\frac{1}{n})$。
$f(n)=\frac{1}{2}+\frac{1}{3}+…+\frac{1}{n}$是一个调和级数，当$n->\infty$，$f(n)->lnn+C$，也就是$logn$。
加了continue优化，也就是把所有质数的倍数都筛掉。计算次数为$n(\frac{1}{2}+\frac{1}{3}+\frac{1}{5}+…+\frac{1}{n})$，其中的分母都为质数，则$g(n)=\frac{1}{2}+\frac{1}{3}+\frac{1}{5}+…+\frac{1}{n}$趋近于$loglogn$。
int primes[maxn], cnt; // primes[] 存储所有素数
bool st[maxn]; // st[x] 存储x是否被筛掉
void get_primes(int n) {
memset(st, 0, sizeof st);
st[0] = st[1] = true;
for (int i = 2; i <= n; i++) {
if (st[i])continue; //优化
primes[cnt++] = i;
for (int j = 2 * i; j <= n; j += i) {
st[j] = true;
}
}
}
线性筛法求素数（欧拉筛）
每一个数只会被它最小的质因子筛掉，每个数只被筛掉一次，因此是$O(n)$。
注意到线性筛法求素数的同时也得到了每个数的最小质因子。
int primes[maxn], cnt; // primes[] 存储所有素数
bool st[maxn]; // st[x] 存储x是否被筛掉
//int minP[maxn]; // minP[x] 表示x的最小质因子
void get_primes(int n) {
for (int i = 2; i <= n; i++) {
if (!st[i]) {
primes[cnt++] = i;
//minP[i] = i;
}
for (int j = 0; primes[j] <= n / i; j++) { //primes[j]*i<=n //用已经有的质数作为i要乘的数
st[primes[j] * i] = true;
//minP[primes[j] * i] = primes[j];
if (i % primes[j] == 0) {
break;
}
}
}
}
关于代码中的
if (i % primes[j] == 0) {break;}
primes[j]是i的因子，说明i之前被primes[j]筛过了。
由于primes[]里面质数是从小到大的，所以i乘上其他的质数的结果一定会在后面的迭代中被这个primes[j]的倍数筛掉，就不需要在这里先筛一次，所以这里直接break。
约数
试除法求所有约数
$O(\sqrt{n})$
vector<int> get_divisor(int x) {
vector<int> res;//存储x的所有约数
for (int i = 1; i <= x / i; i++) {
if (x % i == 0) {
res.push_back(i);
if (i != x / i)res.push_back(x / i);
}
}
sort(res.begin(), res.end());//把所有约数从小到大排序
return res;
}
约数个数和约数之和
$N=P_1^{\alpha_{1}}P_2^{\alpha_{2}}…P_i^{\alpha_{i}}$
所有约数的个数：$(\alpha_1+1)(\alpha_2+1)…(\alpha_i+1)$
证明：约数$d=P_1^{\beta_{1}}P_2^{\beta_{2}}…P_i^{\beta_{i}}$，其中$\beta_i$取$[0,\alpha_i]$，用乘法原理。
所有约数之和：$(1+P_1+P_1^2+…+P_1^{\alpha_1})(1+P_2+P_2^2+…+P_2^{\alpha_2})…(1+P_1+P_1^2+…+P_i^{\alpha_i})$
证明：约数$d=P_1^{\beta_{1}}P_2^{\beta_{2}}…P_i^{\beta_{i}}$，其中$\beta_i$取$[0,\alpha_i]$。
$$
\begin{aligned}
d=&1+P_2+P_2^2+…+P_2^{\alpha_2}+\\
&P_1+P_1P_2+P_1P_2^2+…+P_1P_2^{\alpha_2}+\\
&P_1^2+P_1^2P_2+P_1^2P_2^2+…+P_1^2P_2^{\alpha_2}+\\
&……\\
&P_1^{\alpha_1}+P_1^{\alpha_1}P_2+P_1^{\alpha_1}P_2^2+…+P_1^{\alpha_1}P_2^{\alpha_2}\\
=&(1+P_1+P_1^2+…+P_1^{\alpha_1})(1+P_2+P_2^2+…+P_2^{\alpha_2})
\end{aligned}
$$
计算：每一个多项式用秦九韶算法。
unordered_map<int, int> ha;//ha[p]=alpha
for (auto i: ha) {
ll p = i.first, t = i.second;
ll tmp = 1;
while (t--)tmp = (tmp * p + 1) % MOD;
ans = ans * tmp % MOD;
}
欧几里得算法（辗转相除法）
最大公约数$(a,b)=(b,a \% b)$
证明：
设$a \% b=a-k*b$，其中$k=\lfloor{a/b}\rfloor$。
若$d$是$a,b$的公约数，则$d|a$且$d|b$，则易知$d|(a-k*b)$，故$d$也是$b,a\%b$的公约数。
若$d$是$b,a \% b$的公约数，则$d|b$且$d|(a-kb)$，则$d|(a-kb+kb)=d|a$，故$d$也是$a,b$的公约数。
因此$a,b$的公约数集合和$b,a \% b$的公约数集合相同，所以它们的最大公约数也相同。一直递归下去直到$b=0$，0和任何数的最大公约数都等于这个数本身。
$O(logn)$
int gcd(int a, int b) {
if (b == 0)return a;
else return gcd(b, a % b);
}
推广到多个数。
$(a,b,c)=((a,b),(b,c))=((b,a \% b),(c,b \% c))=(b,a \% b,c,b \% c)$。
因为$(b,c)=(c,b\%c)$，所以$(b,c,b\%c)=(c,b\%c)$。
所以$(a,b,c)=(a\%b,b\%c,c)$。
更相减损术（辗转相减法）
最大公约数$a\geq{b}, (a,b)=(b,a-b)=(a,a-b)$。
需要高精度时考虑用更相减损术。
欧拉函数
欧拉函数的定义
1~N中与N互质的数的个数被称为欧拉函数，记为$\phi(N)$。
若在算术基本定理中，$N=p_1^{c_{1}}p_2^{c_{2}}…p_i^{c_{i}}$，则：
$\phi(N)=N\frac{p_1-1}{p_1}\frac{p_2-1}{p_2}…\frac{p_i-1}{p_i}$。
证明：
根据欧拉函数的计算式，我们只需要分解质因数，即可顺便求出数$n$的欧拉函数。
$O(\sqrt{n})$。
int phi(int n) { //计算n的欧拉函数
int res = n;
for (int i = 2; i <= n / i; i++) {
if (n % i == 0) {
res = res / i * (i - 1);//先除后乘
while (n % i == 0)n /= i;
}
}
if (n > 1)res = res / n * (n - 1);
return res;
}
筛法求欧拉函数
埃氏筛
根据欧拉函数的计算式，在埃氏筛的过程中，求出1~n的所有数的欧拉函数。
$O(nloglogn)$
void get_eulers(int n) {
phi[1] = 1; //特殊规定
for (int i = 1; i <= n; i++)phi[i] = i;
memset(st, 0, sizeof st);
for (int i = 2; i <= n; i++) {
if (st[i])continue;
phi[i] = i - 1;
for (int j = 2 * i; j <= n; j += i) {
phi[j] = phi[j] / i * (i - 1);
st[j] = true;
}
}
}
线性筛（欧拉筛）
若$p|n$且$p^2|n$，则$\phi(n)=\phi(n/p)*p$。
若$p|n$且$p^2\nmid{n}$，则$\phi(n)=\phi(n/p)*(p-1)$。
证明：
根据欧拉函数以上性质，在线性筛的过程中，求出1~n的所有数的欧拉函数。
$O(n)$
void get_eulers(int n) {
phi[1] = 1; //特殊规定
memset(st, 0, sizeof st);
for (int i = 2; i <= n; i++) {
if (!st[i]) {
primes.push_back(i);
phi[i] = i - 1;
}
for (int j = 0; primes[j] <= n / i; j++) {
st[primes[j] * i] = true;
if (i % primes[j] == 0) {
phi[primes[j] * i] = phi[i] * primes[j];
break;
} else {
phi[primes[j] * i] = phi[i] * (primes[j] - 1);
}
}
}
}
欧拉定理和费马小定理
费马小定理
若$p$为素数，$gcd(a,p)=1$，则$a^{p-1}\equiv1\ (mod\ p)$。
另一个形式：对于任意整数$a$，有$a^p\equiv a\ (mod\ p)$。
欧拉定理
若$gcd(a,m)=1$，则$a^{\phi(m)}\equiv 1\ (mod\ m)$。
证明：
构造一个与$m$互质的数列。设$r_1,r_2,…,r_{\phi(m)}$为模$m$意义下的一个简化剩余系，则$ar_1,ar_2,…,ar_{\phi(m)}$也为模$m$意义下的一个简化剩余系。所以$r_1r_2…r_{\phi(m)}\equiv ar_1\cdot ar_2 \cdot \cdot \cdot ar_{\phi(m)} \equiv a^{\phi(m)}r_1r_2…r_{\phi(m)}\ (mod\ m)$，约去$r_1r_2…r_{\phi(m)}$，即可得$a^{\phi(m)}\equiv 1\ (mod\ m)$。
当m为素数时，由于$\phi(m)=m-1$，代入欧拉定理即可得到费马小定理。
快速幂
求$m^k \% p$，时间复杂度$O(logk)$。
int qpow(int m, int k, int p) {
int res = 1 % p, t = m; //视情况取long long
while (k) {
if (k & 1)res = res * t % p;
t = t * t % p;
k >> 1;
}
return res;
}
快速幂求逆元
乘法逆元的定义：
乘法逆元，是指数学领域群G中任意一个元素a，都在G中有唯一的逆元a‘，具有性质a×a’=a’×a=e，其中e为该群的单位元。
模意义下的乘法逆元：
如果$ax\equiv 1\ (mod\ m)$，则$x$称为$a\ mod\ m$的逆元，记作$a^{-1}$。
$a$存在乘法逆元的充要条件是$a$与模数$m$互质。当模数$m$为质数时，$a^{m-2}$即为$a$的乘法逆元。
证明：
$ax\equiv 1\ (mod\ m)$
$a^{m-1}\equiv 1\ (mod\ m)$（费马小定理）
所以$ax \equiv a^{m-1} \ (mod\ m)$，$x\equiv a^{m-2}\ (mod\ m)$。
int inv = qpow(a, m - 2, m);
if (a % m != 0)cout << t;
线性求逆元
$O(n)$求出$1,2,…,n$中每个数关于$p$的逆元。
逆元的递归形式，使用迭代实现：
$$
i^{-1} \equiv \begin{cases}
1, & \text{if } i = 1, \\
-\lfloor\frac{p}{i}\rfloor (p \bmod i)^{-1}, & \text{otherwises}.
\end{cases} \pmod p
$$
inv[1] = 1;
for (int i = 2; i <= n; i++) {
inv[i] = (long long) (p - p / i) * inv[p % i] % p;
}
使用$p-\lfloor\frac{p}{i}\rfloor$来防止出现负数。
使用以上公式，算2 mod p的逆元的快速方法是，mod - mod/2。
线性求任意n个数的逆元
上面的方法只能求1到n的逆元，如果需要求任意给定n个数$(1\leq a_i<p)$的逆元，需使用下面的方法：
首先计算$n$个数的前缀积，记为$s_i$，然后使用快速幂或扩展欧几里得发计算$s_n$的逆元，记为$sv_n$。
$sv_n$乘上$a_n$，得到$sv_{n-1}$。同理依次计算所有的$sv_i$。
$a_i^{-1}$可以用$sv_i \times s_{i-1}$求得。
在$O(n+logp)$内计算出了n个数的逆元。
s[0] = 1;
for (int i = 1; i <= n; i++)s[i] = s[i - 1] * a[i] % p;
sv[n] = qpow(s[n], p - 2, p);
for (int i = n; i >= 1; i--)sv[i - 1] = sv[i] * a[i] % p;
for (int i = 1; i <= n; i++)inv[i] = sv[i] * s[i - 1] % p;
模运算
模$p$运算下的加、减、乘、除、乘方运算。
$(a+b)\%p=(a\%p+b\%p)\%p$
$(a-b)\%p=(a\%p-b\%p+p)\%p$
$(a\times b)\%p=(a\%p)\times (b\%p)\%p$
$(a/b)\%p=(a*b^{-1}\%p)$
$(a^b)\%p=(a^{b\ mod\ \phi(p) +\phi(p)})\%p$（在$b\geq \phi(p)$成立，通常在指数$b$超级大，结合大数取余使用）（扩展欧拉定理的应用）
扩展欧几里得
贝祖定理
对于任意整数$a$,$b$，存在一对整数$x$，$y$，满足$ax+by=gcd(a,b)$。
证明：
若$b=0$，也就是欧几里得算法的最后一步，显然有一对解$x=1,y=0$，使得$a1+00=gcd(a,0)$。
若$b>0$，则$gcd(a,b)=gcd(b,a\%b)$。假设存在一对整数$x,y$满足$b*x+(a\%b)*y=gcd(b,a\%b)$，因为左式$=bx+(a-b\lfloor a/b \rfloor)y=ay+b(x-\lfloor a/b \rfloor y)$，所以令$x’=y,y’=x-\lfloor a/b \rfloor y$，就得到了$ax’+by’=gcd(a,b)$。
对欧几里得算法的递归过程应用数学归纳，贝祖定理成立。
以上的证明过程就是扩展欧几里得算法。扩展欧几里得算法得到$ax+by=gcd(a,b)$的一组特解$x_0,y_0$，并返回$a,b$的最大公约数$d$。
int exgcd(int a, int b, int &x, int &y) {
if (b == 0) {
x = 1, y = 0;
return a;
}
int d = exgcd(b, a % b, x, y);
int z = x;
x = y, y = z - y * (a / b);
return d;
}
对于一般的方程，$ax+by=c$的通解可以表示为：
$$
x=\frac{c}{d}x_0+k\frac{b}{d} \\
y=\frac{c}{d}y_0-k\frac{a}{d}
$$
通解是所有模$\frac{b}{d}$与$x$同余的整数。因此最小正整数解为$(\frac{c}{d}x_0) \% \frac{b}{d}$（注意要是正余数）。
此方程有解当且仅当$d|c$，其中$d$为$gcd(a,b)$。
线性同余方程
求一个整数$x$，满足$ax\equiv c\ (mod\ b)$，或给出无解。
$ax\equiv c\ (mod\ b)$等价于$ax-c$是$b$的倍数，不妨设为$-y$倍，该方程可以改写为$ax+by=c$。
利用扩展欧几里得算法求解。
求乘法逆元
如果$ax\equiv 1\ (mod\ m)$，则$x$称为$a\ mod\ m$的逆元，记作$a^{-1}$。
$ax+my=1$。
int d = exgcd(a,m,x,y);
//(x%m+m)%m即为所求
中国剩余定理
设$n_1,n_2,…,n_k$是两两互质的整数，对于任意$n$个整数$a_1,a_2,…,a_n$，求方程组的整数解。
$$
\begin{cases}
x &\equiv a_1 \pmod {n_1} \\
x &\equiv a_2 \pmod {n_2} \\
&\vdots \\
x &\equiv a_k \pmod {n_k} \\
\end{cases}
$$
求解步骤：
1.计算所有模数的积$n$。
2.对于第$i$个方程，计算$m_i=\frac{n}{n_i}$，计算$m_i$在模$n_i$意义下的逆元$m_
i^{-1}$，计算$m_im_i^{-1}$（不要对$n_i$取模）。
3.方程组在模$n$意义下的唯一解为：$x=\sum_{i=1}^{k}a_im_im_i^{-1} (mod\ n)$。
证明：
因为$m_i$是除了$n_i$之外所有模数的倍数，所以对于任意$j\neq i$，$a_im_im_i^{-1}\equiv0\ (mod\ n_j)$。又$a_im_im_i^{-1}\equiv a_i\ (mod\ n_i)$。
所以对于每一个方程$i$，代入$x=\sum_{i=1}^{k}a_im_im_i^{-1}$，都能得到$x\equiv a_im_im_i^{-1}\equiv a_i\ (mod\ n_i)$，每一个方程都能成立。
ll a[maxn], r[maxn];//r是模数
ll CRT(int k) {
ll n = 1, ans = 0;
for (int i = 1; i <= k; i++)n = n * r[i];
for (int i = 1; i <= k; i++) {
ll mi = n / r[i], x, y;
//用exgcd求逆元，因为r[i]不一定是素数，所以不要用快速幂求逆元
exgcd(mi, r[i], x, y); //mi*x mod r[i] =1
x = (x % r[i] + r[i]) % r[i];
ans += a[i] * mi * x;
}
ans = (ans % n + n) % n;
return ans;
}
扩展中国剩余定理
扩展：模数不互质的情况
设两个方程分别是 $x\equiv a_1 \pmod {m_1}$、$x\equiv a_2 \pmod {m_2}$；
将它们转化为不定方程：$x=m_1p+a_1=m_2q+a_2$，其中 $p, q$ 是整数，则有 $m_1p-m_2q=a_2-a_1$。
由裴蜀定理，当 $a_2-a_1$ 不能被 $\gcd(m_1,m_2)$ 整除时，无解；
其他情况下，可以通过扩展欧几里得算法解出来一组可行解 $(p, q)$；
则原来的两方程组成的模方程组的解为 $x\equiv b\pmod M$，其中 $b=m_1p+a_1$，$M=\text{lcm}(m_1, m_2)$。
下面证明原来的两方程可以合并为方程$x\equiv b\pmod M$。
对于方程$m_1p-m_2q=a_2-a_1$，求出来一组特解$p_0,q_0$，整个解系为：
$$
p=p_0+k\frac{m_2}{gcd(m_1,m_2)}\\
q=q_0+k\frac{m_1}{gcd(m_1,m_2)}
$$
通过$x=m_1p+a_1$得到$x$的一个特解$x_0=m_1p_0+a_1$。
$x$的通解为：
$$
x=m_1p+a_1\\
=m_1p_0+k\frac{m_1m_2}{gcd(m_1,m_2)}+a_1\\
=x_0+k\frac{m_1m_2}{gcd(m_1,m_2)}\\
=x_0+klcm(m_1,m_2)
$$
所以$x\equiv x_0\ (mod\ lcm(m_1,m_2))$。
多个方程的情况用上面的方法两两合并即可。
$$
\begin{cases}
x &\equiv a_1 \pmod {m_1} \\
x &\equiv a_2 \pmod {m_2} \\
&\vdots \\
x &\equiv a_k \pmod {m_k} \\
\end{cases}
$$
求x的最小非负整数解。
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, x, y);
ll z = x;
x = y, y = z - (a / b) * y;
return d;
}
ll lcm(ll a, ll b, ll d) {
return a / d * b;//先除后乘，防止爆long
}
ll exCRT(int n) {
ll a1, m1, a2, m2, p, q;
cin >> m1 >> a1;
for (int i = 2; i <= n; i++) {
cin >> m2 >> a2;
ll d = exgcd(m1, m2, p, q);
if ((a2 - a1) % d != 0)return -1; //无解
p = p * (a2 - a1) / d;
ll t = m2 / d;
p = (p % t + t) % t;
ll M = lcm(m1, m2, d);
m1 = M;
}
return (a1 % m1);
}
转载请注明来源，欢迎对文章中的引用来源进行考证，欢迎指出任何有错误或不够清晰的表达。可以在下面评论区评论，也可以邮件至 1300773281@qq.com