从BZOJ2242看数论基础算法：快速幂，gcd，exgcd，BSGS

LINK

---

其实就是三个板子

1.快速幂

快速幂，通过把指数转化成二进制位来优化幂运算，基础知识

2.gcd和exgcd

gcd就是所谓的辗转相除法，在这里用取模的形式体现出来
(gcd(a,b))，因为b中的a对答案没有贡献，考虑把b变成(b-(b/a)*a)答案是一样的
所以就可以变成了(gcd(b,a%b))，保证大的数在前面，这样当小的数变成0，大的数就是最大公约数

exgcd就是解线性方程(ax+by=c)
有解的条件是(c\%gcd(a,b)=0)
然后考虑gcd的过程对上面的方程进行转化
(ax+by=gcd(a,b)=gcd(b,a\%b)=bx'+(a\%b)y')
这里把(a%b)变成(a-a/b*b)，就变成了
(ax+by=bx'+ay'-(a/b)*b*y')
解得(x=y',y=x'-(a/b)*y')
然后当b变成1的时候，x=1，y=0
递归执行就可以解出来了

那么最小非负解是啥？
首先把x,y乘上(c/gcd(a,b))变成原方程的特解
然后就可以得到通解
(x = x' + b/gcd(a,b)*t)
(y = y' - a/gcd(a,b)*t)
发现这个时候一定满足(ax+by=c),最小非负解也就很好求了

3.BSGS

挺有意思的一个东西，可以求解离散数对问题
(a^x=b(mod c))，给出a，b，c求x

大多数时候会保证c是质数
这个时候我们怎么办？
发现有用的x只有c-1个
考虑把x分解成ir+m的形式
然后式子变成(a^{ir}=b*a^{-m}(mod c))
如果令(r=sqrt{c}),那么枚举i并在hashtable里面查找有没有对应的m就可以了
这样的复杂度是(sqrt{c})的

---

//Author: dream_maker
#include<bits/stdc++.h>
using namespace std;
//----------------------------------------------
//typename
typedef long long ll;
//convenient for
#define fu(a, b, c) for (int a = b; a <= c; ++a)
#define fd(a, b, c) for (int a = b; a >= c; --a)
#define fv(a, b) for (int a = 0; a < (signed)b.size(); ++a)
//inf of different typename
const int INF_of_int = 1e9;
const ll INF_of_ll = 1e18;
//fast read and write
template <typename T>
void Read(T &x) {
  bool w = 1;x = 0;
  char c = getchar();
  while (!isdigit(c) && c != '-') c = getchar();
  if (c == '-') w = 0, c = getchar();
  while (isdigit(c)) {
    x = (x<<1) + (x<<3) + c -'0';
    c = getchar();
  }
  if (!w) x = -x;
}
template <typename T>
void Write(T x) {
  if (x < 0) {
    putchar('-');
    x = -x; 
  }
  if (x > 9) Write(x / 10);
  putchar(x % 10 + '0');
}
//----------------------------------------------
int y, z, P;
int add(int a, int b) {
  return (a += b) >= P ? a - P : a;
}
int sub(int a, int b) {
  return (a -= b) < 0 ? a + P : a;
}
int mul(int a, int b) {
  return 1ll * a * b % P;
}
int fast_pow(int a, int b) {
  int res = 1;
  while (b) {
    if (b & 1) res = mul(res, a);
    b >>= 1;
    a = mul(a, a);
  }
  return res;
}
int gcd(int a, int b) {
  return b ? gcd(b, a % b) : a;
}
void exgcd(int a, int b, ll &x, ll &y) {
  if (!b) {x = 1, y = 0; return;}
  exgcd(b, a % b, y, x);
  y -= a / b * x;
}
void work1() {
  Read(y), Read(z), Read(P);
  Write(fast_pow(y, z));
  putchar('
');
}
void work2() {
  Read(y), Read(z), Read(P);
  int g = gcd(y, P);
  if (z % g) {
    printf("Orz, I cannot find x!
");
  } else {
    ll a, b;
    exgcd(y, P, a, b);
    a *= z / g;
    a = (a % (P / g) + P / g) % (P / g);
    Write(a);
    putchar('
');
  }
}
map<int, int> mp;
void work3() {
  Read(y), Read(z), Read(P);
  y %= P, z %= P;
  if (!y) {
    printf("Orz, I cannot find x!
");  
    return;
  }
  int w = sqrt(P), now = z;
  mp.clear();
  fu(i, 0, w) {
    mp[now] = i;
    now = mul(now, y);
  }
  now = fast_pow(y, w);
  int tmp = now;
  fu(i, 1, w) {
    if (mp.count(tmp)) {
      Write(i * w - mp[tmp]);
      putchar('
');
      return;
    }
    tmp = mul(tmp, now);
  }
  printf("Orz, I cannot find x!
");
}
int main() {
  int T, op;
  Read(T); Read(op);
  if (op == 1) while (T--) work1();
  else if (op == 2) while (T--) work2();
  else while (T--) work3();
  return 0;
}