acm模板-数论

发布于 2021-10-18  574 次阅读

数论

位运算

  1. $ a + b = a \& b + a | b $
  2. 二进制逆序
a = ((a & 0xAAAA) >> 1) | ((a & 0x5555) << 1);
a = ((a & 0xCCCC) >> 2) | ((a & 0x3333) << 2);
a = ((a & 0xF0F0) >> 4) | ((a & 0x0F0F) << 4);
a = ((a & 0xFF00) >> 8) | ((a & 0x00FF) << 8);
  1. 要求集合中不能有两个相邻的元素
if ((mask >> 1) & mask) continue;
  1. 在限定必须不取某些元素的前提下枚举子集
// mask的第x位为0表示x必须不在子集中(原集合中不含这个元素)：
for (int mask1 = mask; mask1 >= 0; mask1 = (mask1 - 1) & mask)
  1. int __builtin_ffs(int x) ：返回x的二进制末尾最后一个1的位置，位置的编号从1开始（最低位编号为1）。当x为0时返回0
  2. int __builtin_clz(unsigned int x) ：返回x的二进制的前导0的个数。当x为0时，结果未定义。
  3. int __builtin_ctz(unsigned int x) ：返回x的二进制末尾连续0的个数。当x为0时，结果未定义。
  4. int __builtin_clrsb(int x) ：当x的符号位为0时返回x的二进制的前导0的个数减一，否则返回x的二进制的前导1的个数减一。
  5. int __builtin_popcount(unsigned int x) ：返回x的二进制中1的个数。
    1. int __builtin_parity(unsigned int x) ：判断x的二进制中1的个数的奇偶性。

这些函数都可以在函数名末尾添加 l 或 ll （如 builtin_popcountll ）来使参数类型变为 ( unsigned ) long 或 ( unsigned ) long long （返回值仍然是 int 类型）。 例如，我们有时候希望求出一个数以二为底的对数，如果不考虑 0 的特殊情况，就相当于这个数二进制的位数 -1 ，而一个数 n 的二进制表示的位数可以使用 32-__builtin_clz(n) 表示，因此 31-builtin_clz(n) 就可以求出 n 以二为底的对数。

拓展欧几里得算法

求解$ ax + by = gcd(a, b) $

LL exgcd(LL a, LL b, LL &x, LL &y) {
  LL d = a;
  if (b == 0)
    x = 1, y = 0;
  else {
    d = exgcd(b, a % b, y, x), y -= a / b * x;
  }
  return d;
}

基础快速幂

LL ksm(LL n, LL m, LL p = MOD) {
  LL res = 1;
  while (m > 0) {
    if (m & 1) res = n * res % p;
    m >>= 1;
    n = n * n % p;
  }
  return res % p;
}

多幂数快速幂

对于给定x和n个幂数 $ a_i $ 和模数p，求 $ x^{a_i} \mod p $

解法：设 $ S = \sqrt p + 1 $ ，预处理S以内快速幂，那么 $ x^y = x^{y \mod S} * x^{\lfloor \frac{y}{S} \rfloor *S} $

快速乘

LL mul(LL a, LL b, LL mod) {
  LL res = a * b - mod * (LL)((long double)a / mod * b);
  if (res < 0) return res + mod;
  if (res < mod) return res;
  return res - mod;
}

素数判定

Pollard-Rho算法

int pri[] = {0, 19260817, 1145141923};
bool isprime(LL n) {
  if (n < 2) return 0;
  for (int i = 1; i <= 2; ++i)
    if (pri[i] == n) return 1;
  int cnt = 0;
  LL res = n - 1;
  while (!(res & 1)) cnt++, res >>= 1;
  for (int i = 1; i <= 2; ++i) {
    LL x = ksm(pri[i], res, n), t;
    for (int j = 1; j <= cnt; ++j) {
      t = mul(x, x, n);
      if (t == 1 && x != 1 && x != n - 1) return 0;
      x = t;
    }
    if (x != 1) return 0;
  }
  return 1;
}
LL random(LL n) { return 1LL * rand() * rand() * rand() * rand() % n; }
LL PoLLardRho(LL n) {
  if (n == 4) return 2;
  while (1) {
    LL c = random(n - 1) + 1;
#define f(x) (mul(x, x, n) + c)
    LL s = f(0), t = f(s), v = 1, cnt = 0, o = 1;
    while (s != t) {
      v = mul(v, abs(s - t), n);
      if ((++cnt) % 114 == 0) {
        LL d = __gcd(n, v);
        if (d > 1 && d < n) return d;
      }
      if (cnt == o) {
        LL d = __gcd(n, v);
        if (d > 1 && d < n) return d;
        v = 1, cnt = 0, o *= 2;
      }
      s = f(s), t = f(f(t));
    }
    LL d = __gcd(n, v);
    if (d > 1 && d < n) return d;
  }
}
void dfs(LL x, LL& ans) {
  if (x < ans) return;
  if (isprime(x)) return ans = max(ans, x), void();
  LL d = PoLLardRho(x);
  dfs(d, ans);
  while (x % d == 0) x /= d;
  dfs(x, ans);
}
int main() {
  srand(time(0));
  LL n, ans;
  read(n);
  ans = 2;
  dfs(n, ans);
  if (ans == n) puts("Prime");
}

筛法

线性筛

const LL p_max = 1E6 + 100;
LL pr[p_max], p_sz;
void get_prime() {
  static bool vis[p_max];
  for (int i = 2; i < p_max; ++i) {
    if (!vis[i]) pr[p_sz++] = i;
    for (int j = 0; j < p_sz; ++i) {
      if (pr[j] * i >= p_max) break;
      vis[pr[j] * i] = 1;
      if (i % pr[j] == 0) break;
    }
  }
}

欧拉函数

求欧拉函数值

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;
}

筛法求欧拉函数

void pre() {
  memset(is_prime, 1, sizeof(is_prime));
  int cnt = 0;
  is_prime[1] = 0;
  phi[1] = 1;
  for (int i = 2; i <= 5000000; i++) {
    if (is_prime[i]) {
      prime[++cnt] = i;
      phi[i] = i - 1;
    }
    for (int j = 1; j <= cnt && i * prime[j] <= 5000000; j++) {
      is_prime[i * prime[j]] = 0;
      if (i % prime[j])
        phi[i * prime[j]] = phi[i] * phi[prime[j]];
      else {
        phi[i * prime[j]] = phi[i] * prime[j];
        break;
      }
    }
  }
}

裴蜀定理

设a, b是不全为零的整数，则存在整数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 $

欧拉定理

若$ gcd(a, m) = 1 $，则$ a ^ {\phi (m)} \equiv 1 \pmod m $

拓展欧拉定理

威尔逊定理

对于素数$ p $， 有$ (p - 1)! \equiv -1 (\mod p) $

线性求逆元

inv[1] = 1;
for (int i = 2; i <= n; ++i) {
  inv[i] = (LL)(p - p / i) * inv[p % i] % p;
}

线性求任意n个数的逆元（复杂度$ O(n + \log p) $）

s[0] = 1;
for (int i = 1; i <= n; ++i) s[i] = s[i - 1] * a[i] % p;
sv[n] = ksm(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;

中国剩余定理(余数为a，模数为r)

LL CRT(int k, LL* a, LL* r) {
  LL n = 1, ans = 0;
  for (int i = 1; i <= k; i++) n = n * r[i];
  for (int i = 1; i <= k; i++) {
    LL m = n / r[i], b, y;
    exgcd(m, r[i], b, y);
    ans = (ans + a[i] * m * b % mod) % mod;
  }
  return (ans % mod + mod) % mod;
}

二次剩余欧拉判别准则

$ \left(\frac{n}{p}\right)\equiv n^{\frac{p-1}{2}}\pmod p $

若n是二次剩余，当且仅当$ n^{\frac{p-1}{2}}\equiv 1\pmod p $

若n是非二次剩余，当且仅当$ n^{\frac{p-1}{2}}\equiv -1\pmod p $

BSGS

求$ a^x \equiv b \pmod p $

LL BSGS(LL a, LL b, LL p) {
  a %= p;
  if (!a && !b) return 1;
  if (!a) return -1;
  static map<LL, LL> mp;
  mp.clear();
  LL m = sqrt(p + 1.5);
  LL v = 1;
  for (int i = 1; i <= m; ++i) {
    v = v * a % p;
    mp[v * b % p] = i;
  }
  LL vv = v;
  for (int i = 1; i <= m; ++i) {
    auto it = mp.find(vv);
    if (it != mp.end()) return i * m - it->second;
    vv = vv * v % p;
  }
  return -1;
}
LL exBSGS(LL a, LL b, LL p) {
  a %= p, b %= p;
  if (a == 0) return b > 1 ? -1 : b == 0 && p != 1;
  LL c = 0, q = 1;
  while (1) {
    LL g = __gcd(a, p);
    if (g == 1) break;
    if (b == 1) return c;
    if (b % g) return -1;
    ++c;
    b /= g, p /= g;
    q = a / g * q % p;
  }
  static map<LL, LL> mp;
  mp.clear();
  LL m = sqrt(p + 1.5);
  LL v = 1;
  for(int i = 1; i <= m; ++i) {
    v = v * a % p;
    mp[v * b % p] = i;
  }
  for(int i = 1; i <= m; ++i) {
    q = q * v % p;
    auto it = mp.find(q);
    if (it != mp.end()) return i * m - it->second + c;
  }
  return -1;
}

杜教筛求欧拉函数和及莫比乌斯函数前缀和

#include <algorithm>
#include <cstdio>
#include <cstring>
#include <map>
using namespace std;
const int maxn = 2000010;
long long T, n, pri[maxn], cur, mu[maxn], sum_mu[maxn];
bool vis[maxn];
map<long long, long long> mp_mu;
long long S_mu(long long x) {  //求mu的前缀和
  if (x < maxn) return sum_mu[x];
  if (mp_mu[x]) return mp_mu[x];  //如果map中已有该大小的mu值，则可直接返回
  long long ret = (long long)1;
  for (long long i = 2, j; i <= x; i = j + 1) {
    j = x / (x / i);
    ret -= S_mu(x / i) * (j - i + 1);
  }
  return mp_mu[x] = ret;  //路径压缩，方便下次计算
}
long long S_phi(long long x) {  //求phi的前缀和
  long long ret = (long long)0;
  long long j;
  for (long long i = 1; i <= x; i = j + 1) {
    j = x / (x / i);
    ret += (S_mu(j) - S_mu(i - 1)) * (x / i) * (x / i);
  }
  return (ret - 1) / 2 + 1;
}
int main() {
  scanf("%lld", &T);
  mu[1] = 1;
  for (int i = 2; i < maxn; i++) {  //线性筛预处理mu数组
    if (!vis[i]) {
      pri[++cur] = i;
      mu[i] = -1;
    }
    for (int j = 1; j <= cur && i * pri[j] < maxn; j++) {
      vis[i * pri[j]] = true;
      if (i % pri[j])
        mu[i * pri[j]] = -mu[i];
      else {
        mu[i * pri[j]] = 0;
        break;
      }
    }
  }
  for (int i = 1; i < maxn; i++)
    sum_mu[i] = sum_mu[i - 1] + mu[i];  //求mu数组前缀和
  while (T--) {
    scanf("%lld", &n);
    printf("%lld %lld\n", S_phi(n), S_mu(n));
  }
  return 0;
}

Min_25筛求素数和

#include <bits/stdc++.h>

using namespace std;
using LL = long long;
const int N = 1000000 + 10;
int prime[N], id1[N], id2[N], flag[N], ncnt, m;
LL g[N], sum[N], a[N], T, n;

int ID(LL x) { return x <= T ? id1[x] : id2[n / x]; }

LL calc(LL x) { return x * (x + 1) / 2 - 1; }

LL init(LL x) {
  T = sqrt(x + 0.5);
  for (int i = 2; i <= T; i++) {
    if (!flag[i]) prime[++ncnt] = i, sum[ncnt] = sum[ncnt - 1] + i;
    for (int j = 1; j <= ncnt && i * prime[j] <= T; j++) {
      flag[i * prime[j]] = 1;
      if (i % prime[j] == 0) break;
    }
  }
  for (LL l = 1; l <= x; l = x / (x / l) + 1) {
    a[++m] = x / l;
    if (a[m] <= T)
      id1[a[m]] = m;
    else
      id2[x / a[m]] = m;
    g[m] = calc(a[m]);
  }
  for (int i = 1; i <= ncnt; i++)
    for (int j = 1; j <= m && (LL)prime[i] * prime[i] <= a[j]; j++)
      g[j] = g[j] - (LL)prime[i] * (g[ID(a[j] / prime[i])] - sum[i - 1]);
}

LL solve(LL x) {
  if (x <= 1) return x;
  return n = x, init(n), g[ID(n)];
}

int main() {
  memset(g, 0, sizeof(g));
  memset(a, 0, sizeof(a));
  memset(sum, 0, sizeof(sum));
  memset(prime, 0, sizeof(prime));
  memset(id1, 0, sizeof(id1));
  memset(id2, 0, sizeof(id2));
  memset(flag, 0, sizeof(flag));
  ncnt = m = 0;
  scanf("%LLd", &n);
  printf("%LLd\n", solve(n));
}
Hello, world!
最后更新于 2021-10-18