说在前面
默认了解一些基本定义，如整除、取模、质数等，仅有算法的思想和实现，没有且不做证明
如果需要更详细的说明、了解，也许你需要：基础数论，OI-Wiki
一些表示方法
整数：\(\mathbf{Z}\)
属于：\(a \in \mathbf{Z}\)（\(a\) 属于整数）
存在：$ \exists$
任意，任取：\(\forall i\)
整除：\(a \mid b\)
不整除：\(a \nmid b\)
绝对值：\(|a|\)
向上取整：$\left \lceil \frac{a}{b} \right \rceil $
向下取整：$\left \lfloor \frac{a}{b} \right \rfloor $
组合数：$\binom{m}{n} $
内积：$\left \langle a,b \right \rangle $
范数：$\left | a,b \right | $
最大公约数：\(\gcd(a, b)\)
互质：\(\gcd(a, b) = 1\)
同余：\(a \equiv b \pmod m\)
不同余：\(a \not \equiv b \pmod m\)
求和：\(\sum_{i = 1}^{n}\)
连乘：\(\prod_{i = 1}^{n}\)
单位函数：\(\varepsilon(n) = [n = 1]\)（当 \(n\) 为真时，函数值为 \(1\)，否则为 \(0\)）
最大公约数（\(\gcd\)）
欧几里得算法 / 辗转相除法
原理
\(\gcd(a, b) = \gcd(b, a \bmod b)\)
实现
int gcd(int a, int b) {
return b == 0 ? a : gcd(b, a % b);
}
标准库
C++14：__gcd(a, b)
C++17：std::gcd(a, b)，需要头文件 <byneric>
时间复杂度
\(O(n)\)（例如求斐波那契数列相邻两项的最大公约数）
更替减损术
对于大整数更优的时间复杂度
原理
\(\gcd(a, b) = \gcd(a - b, b)\)
实现
int gcd(int a, int b) {
if (a == b) return a;
else if (a > b) a -= b;
else b -= a;
return gcd(a, b);
}
优化：stein 算法
OI-Wiki：
若 \(2\mid a,2\mid b\)，\(\gcd(a,b) = 2\gcd\left(\dfrac a2, \dfrac b2\right)\)。
否则，若 \(2\mid a\)（\(2\mid b\) 同理），因为 \(2\mid b\) 的情况已经讨论过了，所以 \(2 \nmid b\)。因此 \(\gcd(a,b)=\gcd\left(\dfrac a2,b\right)\)。
实现
int stein(int a, int b) {
if (!a) return b;
if (!b) return a;
if ((a & 1) && !(b & 1)) return stein(a >> 1, b >> 1) << 1;
else if (!(a & 1)) return stein(a >> 1, b);
else if (!(b & 1)) return stein(a, b >> 1);
return stein(abs(a - b), min(a, b));
}
大整数（高精）实现
你需要：高精减法，高精大小比较，高精左、右移（或乘除），高精判断奇偶
ll gcd(ll a, ll b) {
int atimes = 0, btimes = 0;
while (a % 2 == 0) {
a >>= 1;
atimes++;
}
while (b % 2 == 0) {
b >>= 1;
btimes++;
}
while (a != b) {
while (a % 2 == 0) a >>= 1;
while (b % 2 == 0) b >>= 1;
if (a < b) swap(a, b);
a -= b;
}
return a << min(atimes, btimes);
}
最小公倍数（LCM）
原理
\(\gcd(a, b) \times \operatorname{lcm}(a, b) = a \times b\)，所以只需先求出最大公约数即可。
扩展欧几里得算法（EXGCD）
用于求解 \(ax + bt = \gcd(a, b)\) 的一组解
实现
int Exgcd(int a, int b, int &x, int &y) {
if (!b) {
x = 1, y = 0;
return a;
}
int d = Exgcd(b, a % b, x, y);
int t = x;
x = y, y = t - (a / b) * y;
return d;
}
欧拉函数
定义
小于等于 \(n\) 的书中与 \(n\) 互质的数的数量，表示为 \(\varphi(n)\)
性质
- 欧拉函数是极性函数
即当 \(\gcd(a, b) = 1\) 时，\(\varphi(a \times b) = \varphi(a) \times \varphi(b)\)
1.5 当 \(n\) 为奇数时，\(\varphi(2n) = \varphi(n)\)
2. \(n = \sum_{d \mid n} \varphi(d)\)
3. 若 \(n = p^k\)，则有 \(\varphi(n) = p^k - p^{k - 1}\)
4. 设 \(n = \prod_{i = 1}^s p_i^{k_i}\)，其中 \(p_i\) 为质数，那么有 \(\varphi(n) = n \times \prod_{i = 1}^s \frac{p_i - 1}{p_u}\)
实现
int euler_phi(int n) {
int ans = n;
for (int i = 2; i * i <= n; i++)
if (n % i == 0) {
ans = ans / i * (i - 1);
while (n % i == 0) n /= i;
}
if (n > 1) ans = ans / n * (n - 1);
return ans;
}
筛法求欧拉函数
vector<int> pri;
bool not_prime[N];
int phi[N];
void pre(int n) {
phi[1] = 1;
for (int i = 2; i <= n; i++) {
if (!not_prime[i]) {
pri.push_back(i);
phi[i] = i - 1;
}
for (int pri_j : pri) {
if (i * pri_j > n) break;
not_prime[i * pri_j] = 1;
if (i % pri_j == 0) {
phi[i * pri_j] = phi[i] * pri_j;
break;
}
phi[i * pri_j] = phi[i] * phi[pri_j];
}
}
}
欧拉定理
若 \(\gcd(a, m) = 1\)，则有 \(a^{\varphi(m)} \equiv 1 \pmod m\)
扩展欧拉函数
\(a^b\equiv \begin{cases} a^{b\bmod\varphi(p)},\,&\gcd(a,\,p)=1\\ a^b,&\gcd(a,\,p)\ne1,\,b<\varphi(p)\\ a^{b\bmod\varphi(p)+\varphi(p)},&\gcd(a,\,p)\ne1,\,b\ge\varphi(p) \end{cases} \pmod p\)
筛法
素数筛：埃氏筛法
定义
对于一个质数 \(p\)，\(x \times p\) 一定为合数，可以进行标记
实现
vector<int> prime;
bool is_prime[N];
void Eratosthenes(int n) {
is_prime[0] = is_prime[1] = 0;
for (int i = 2; i <= n; ++i) is_prime[i] = 1;
for (int i = 2; i * i <= n; ++i) {
if (is_prime[i]) {
prime.push_back(i);
if ((long long)i * i > n) continue;
for (int j = i * i; j <= n; j += i)
is_prime[j] = 0;
}
}
}
时间复杂度：\(O(n \log \log n)\)（\(n \ln \ln \sqrt{n} + o(n)\)）
优化
- 可以只筛奇数
- 可以使用
bitset
实现：
bitset<N> vis;
void Prime(int n) {
vis.set();
vis[0] = vis[1] = 0;
for (int i = 2; i * i <= n; i++) {
if (vis[i]) {
for (int j = i << 1; j <= n; j += i) vis[j] = 0;
}
}
}
- 分块筛选
实现：
int count_primes(int n) {
const int S = 10000;
vector<int> primes;
int nsqrt = sqrt(n);
vector<char> is_prime(nsqrt + 1, 1);
for (int i = 2; i <= nsqrt; i++) {
if (is_prime[i]) {
primes.push_back(i);
for (int j = i * i; j <= nsqrt; j += i) is_prime[j] = 0;
}
}
int result = 0;
vector<char> block(S);
for (int k = 0; k * S <= n; k++) {
fill(block.begin(), block.end(), 1);
int start = k * S;
for (int p : primes) {
int start_idx = (start + p - 1) / p;
int j = max(start_idx, p) * p - start;
for (; j < S; j += p) block[j] = 0;
}
if (k == 0) block[0] = block[1] = 0;
for (int i = 0; i < S && start + i <= n; i++) {
if (block[i]) result++;
}
}
return result;
}
素数筛：线性筛法 / 欧拉筛法
原理
每个合数只标记一次
实现
vector<int> pri;
bool not_prime[N];
void pre(int n) {
for (int i = 2; i <= n; ++i) {
if (!not_prime[i]) pri.push_back(i);
for (int pri_j : pri) {
if (i * pri_j > n) break;
not_prime[i * pri_j] = 1;
if (i % pri_j == 0) break;
}
}
}
时间复杂度
\(O(n)\)
筛法求约数个数
定理
若 \(n=\prod_{i=1}^m p_i^{c_i}\)，则 \(d_i=\prod_{i=1}^m (c_i+1)\)。
实现
vector<int> pri;
bool not_prime[N];
int d[N], num[N];
void pre(int n) {
d[1] = 1;
for (int i = 2; i <= n; ++i) {
if (!not_prime[i]) {
pri.push_back(i);
d[i] = 2, num[i] = 1;
}
for (int pri_j : pri) {
if (i * pri_j > n) break;
not_prime[i * pri_j] = 1;
if (i % pri_j == 0) {
num[i * pri_j] = num[i] + 1;
d[i * pri_j] = d[i] / num[i * pri_j] * (num[i * pri_j] + 1);
break;
}
num[i * pri_j] = 1;
d[i * pri_j] = d[i] * 2;
}
}
}
筛法求约数和
实现
vector<int> pri;
bool not_prime[N];
int g[N], f[N];
void pre(int n) {
g[1] = f[1] = 1;
for (int i = 2; i <= n; ++i) {
if (!not_prime[i]) {
pri.push_back(i);
g[i] = i + 1, f[i] = i + 1;
}
for (int pri_j : pri) {
if (i * pri_j > n) break;
not_prime[i * pri_j] = 1;
if (i % pri_j == 0) {
g[i * pri_j] = g[i] * pri_j + 1;
f[i * pri_j] = f[i] / g[i] * g[i * pri_j];
break;
}
f[i * pri_j] = f[i] * f[pri_j];
g[i * pri_j] = 1 + pri_j;
}
}
}
分解质因数
朴素算法
实现
vector<int> breakdown(int N) {
vector<int> result;
for (int i = 2; i * i <= N; i++) {
if (N % i == 0) {
while (N % i == 0) N /= i;
result.push_back(i);
}
}
if (N != 1) result.push_back(N);
return result;
}
Pollard Rho 算法
略。不会。没学过。
裴蜀定理
定义
设 \(a,b\) 是不全为零的整数，对任意整数 \(x,y\)，满足 \(\gcd(a,b)\mid ax+by\)，且存在整数 \(x,y\), 使得 \(ax+by=\gcd(a,b)\)。
费马小定理
定义
若 \(p\) 为素数，\(\gcd(a, p) = 1\)，则 \(a^{p - 1} \equiv 1 \pmod{p}\)。
对于任意整数 \(a\)，有 \(a^p \equiv a \pmod p\)
乘法逆元
定义
如果一个线性同余方程 \(ax \equiv 1 \pmod b\)，则 \(x\) 称为 \(a \bmod b\) 的逆元，记作 \(a^{-1}\)。
实现
扩展欧几里得
#define ll long long
ll y, p;
ll exgcd(ll a, ll b, ll &x, ll &y) {
if (b == 0) {
x = 1, y = 0;
return a;
}
ll ans = exgcd(b, a % b, x, y), t = x;
x = y, y = t - (a / b) * y;
return ans;
}
ll mmi(ll y, ll p) {
ll x = 0, k = 0, sum = exgcd(y, p, x, k);
if (sum == 1) return (x % p + p) % p;
else return -1;
}
快速幂（费马小定理）
#define ll long long
ll qpow(ll a, ll n, ll p) {
ll ans = 1;
while (n) {
if (n & 1) ans = ans * a % p;
a = a * a % p;
n >>= 1;
}
return ans;
}
ll mmi(ll a, ll p) {
return qpow(a, p - 2, p);
}
中国剩余定理（CRT）
定义
求解如下形式的一元线性同余方程组（其中 \(n_1, n_2, \cdots, n_k\) 两两互质）：
\(\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}\)
过程
- 计算所有模数的积 n；
- 对于第 i 个方程：
2.1 计算 \(m_i=\frac{n}{n_i}\)；
2.2 计算 \(m_i\) 在模 \(n_i\) 意义下的逆元 \(m_i^{-1}\)；
2.3 计算 \(c_i=m_im_i^{-1}\)（不要对 \(n_i\) 取模）。 - 方程组在模 n 意义下的唯一解为：\(x=\sum_{i=1}^k a_ic_i \pmod n\)。
实现
void exgcd(int a, int b, int &x, int &y) {
if (b == 0) {
x = 1, y = 0;
return;
}
exgcd(b, a % b, y, x);
y -= a / b * x;
}
int CRT(int k, int *a, int *r) {
int n = 1, ans = 0;
for (int i = 1; i <= k; i++) n = n * r[i];
for (int i = 1; i <= k; i++) {
int m = n / r[i], b, y;
exgcd(m, r[i], b, y);
ans = (ans + a[i] * m * b % n) % n;
}
return (ans % n + n) % n;
}
int solve() {
return CRT(n, b, a);
}
exCRT
用于求解模数不互质的情况
实现
#include<bits/stdc++.h>
using namespace std;
const int N = 100005;
__int128 n, m[N], r[N];
__int128 exgcd(__int128 a,__int128 b,__int128 &x,__int128 &y) {
if(b==0) {
x=1, y=0;
return a;
}
__int128 d, x1, y1;
d = exgcd(b, a%b, x1, y1);
x = y1, y = x1-a/b*y1;
return d;
}
__int128 EXCRT(__int128 m[], __int128 r[]) {
__int128 m1, m2, r1, r2, p, q;
m1 = m[1], r1 = r[1];
for(int i=2; i<=n; i++) {
m2 = m[i], r2 = r[i];
__int128 d = exgcd(m1,m2,p,q);
if((r2-r1)%d) {
return -1;
}
p=p*(r2-r1)/d;
p=(p%(m2/d)+m2/d)%(m2/d);
r1 = m1*p+r1;
m1 = m1*m2/d;
}
return (r1%m1+m1)%m1;
}
int main() {
scanf("%lld",&n);
for(int i = 1; i <= n; ++i)
scanf("%lld%lld",m+i,r+i);
printf("%lld\n", EXCRT(m,r));
return 0;
}
Lucas 定理
定义
对于质数 p，有
实现
// 需要先预处理出fact[]，即阶乘
inline ll C(ll m, ll n, ll p)
{
return m < n ? 0 : fact[m] * inv(fact[n], p) % p * inv(fact[m - n], p) % p;
}
inline ll lucas(ll m, ll n, ll p)
{
return n == 0 ? 1 % p : lucas(m / p, n / p, p) * C(m % p, n % p, p) % p;
}
exLucas 定理
用于求解 \(p\) 不为质数的情况
实现
LL calc(LL n, LL x, LL P) {
if (!n) return 1;
LL s = 1;
for (LL i = 1; i <= P; i++)
if (i % x) s = s * i % P;
s = Pow(s, n / P, P);
for (LL i = n / P * P + 1; i <= n; i++)
if (i % x) s = i % P * s % P;
return s * calc(n / x, x, P) % P;
}
LL multilucas(LL m, LL n, LL x, LL P) {
int cnt = 0;
for (LL i = m; i; i /= x) cnt += i / x;
for (LL i = n; i; i /= x) cnt -= i / x;
for (LL i = m - n; i; i /= x) cnt -= i / x;
return Pow(x, cnt, P) % P * calc(m, x, P) % P * inverse(calc(n, x, P), P) %
P * inverse(calc(m - n, x, P), P) % P;
}
LL exlucas(LL m, LL n, LL P) {
int cnt = 0;
LL p[20], a[20];
for (LL i = 2; i * i <= P; i++) {
if (P % i == 0) {
p[++cnt] = 1;
while (P % i == 0) p[cnt] = p[cnt] * i, P /= i;
a[cnt] = multilucas(m, n, i, p[cnt]);
}
}
if (P > 1) p[++cnt] = P, a[cnt] = multilucas(m, n, P, P);
return CRT(cnt, a, p);
}