1 int gcd(int x, int y) {
 2     return y == 0 ? x : gcd(y, x % y);
 3 }
 4 
 5 int extgcd(int a, int b, int &x, int &y) {
 6     int d = a;
 7     if (b) {
 8         d = extgcd(b, a % b, y , x);
 9         y -= (a / b) * x;
10     } else {
11         x = 1;
12         y = 0;
13     }
14     return d;
15 }

//求解逆元
int mod_inverse(int a, int m) {
    int x, y;
    extgcd(a, m, x, y);
    return (m + x % m) % m;
}

//求欧拉函数值 O(√n)
int euler_phi(int n) {
    int res = n;
    for (int i = 2; i * i <= n; i++) {
        if (n % i == 0) {
            res -= res / i;
            while (n % i == 0) n /= i;
        }
    }
    if (n != 1) res -= res / n;
    return res;
}

 1 // O(maxn)时间内筛出欧拉函数值的表
 2 int euler[maxn];
 3 void euler_pji2() {
 4     for (int i = 0; i < maxn; i++) euler[i] = i;
 5     for (int i = 2; i < maxn; i++) {
 6         if (euler[i] == i) {
 7             for (int j = i; j < maxn; j += i) {
 8                 euler[j] -= euler[j] / i;
 9             }
10         }
11     }
12 }

 1 // 线性同余方程组
 2 // a?x ￡ b?(mod m?)
 3 // 返回一个(b, m)的数对
 4 pair<int, int> linear_congruence(const vector<int>& A, const vector<int>& B, const vector<int>& M) {
 5     //由于最开始没有任何限制，所以先把解设为表示所有整数的x￡0(mod 1)
 6     int x = 0, m = 1;
 7     
 8     for (int i = 0; i < A.size(); i++) {
 9         int a = A[i] * m, b = B[i] - A[i] * x, d = gcd(M[i], a);
10         if (b % d != 0) return make_pair(0, -1); // 无解
11         int t = b / d * mod_inverse(a / d, (M[i] / d) % (M[i] / d));
12         x += m * t;
13         m *= M[i] / d;
14     }
15     return make_pair(x % m, m);
16 }

// 中国剩余定理
// f(x) ￡ 0(mod n) 等价于 f(x) ￡ 0 (mod pa) (pa | n)

 1 //n! O(log n)
 2 int fact[maxp];
 3 //分解n! ￡ b pa, 返回b mod p
 4 int mod_fact(int n, int p, int &a) {
 5     a = 0;
 6     if (n == 0) return 1;
 7     int res = mod_fact(n / p, p, a);
 8     a += n / p;
 9     if (n / p % 2 != 0) return res * (p - fact[n % p]) % p;
10     return res * fact[n % p] % p;
11 }

1 //求组合数
2 int mod_comb(int n, int k, int p) {
3     if (n < 0 || k < 0 || n < k) return 0;
4     int e1, e2, e3;
5     int a1 = mod_fact(n, p , e1), a2 = mod_fact(k, p, e2), a3 = mod_fact(n - k, p, e3);
6     if (e1 > e2 + e3) return 0;
7     return a1 * mod_inverse(a2 * a3 % p, p) % p;
8 }