基本数论

说一下一些基本数论算法的实现和简洁证明。

持续更新中。

  1. \(exgcd\)

    目的:求解方程 \(ax+by=(a,b)\) 的一组特解。

    做法:由于他叫exgcd,我们就按求 \(\gcd\) 方法试一下。 \[ \begin{aligned} 不妨设 (a,b)=1,a>b\\ a\cdot x+b\cdot y=&1\\ b\cdot y+(a\bmod b+\lfloor\frac{a}{b}\rfloor b)\cdot x=&1\\ 提取b\\ b\cdot(y+\lfloor \frac{a}{b}\rfloor x)+(a\bmod b)x=&1\\ 问题化归,求解\ bx_1+(a\bmod b)y_1=&1\\ \end{aligned} \] \(\mathcal{Code}\)

    1
    2
    3
    4
    5
    inline 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;
    }

    为啥能求逆元? \[ \begin{aligned} ax+by&=1\\ ax+by&=1\pmod{b}\\ ax&=1\pmod{b} \end{aligned} \] 显然 \(x\)\(a\)\(\bmod b\) 意义下的逆元

  2. \(excrt\)

    目的:合并若干线性同余方程组。 \[ \begin{cases} x\equiv a_1\pmod{p_1}\\ x\equiv a_2\pmod{p_2}\\ \ \ \ \ \vdots\\ x\equiv a_n\pmod{p_n} \end{cases} \] 问题本质等价,合并两个方程。

    \(x\equiv a_1\pmod{p_1},x\equiv a_2\pmod{p_2}\)

    若存在解则解一定能表达成 \(x=x_0+k\cdot lcm(p_1,p_2)\),\(x_0\) 是一个特解。

    显然 \(p_1\mid lcm(p_1,p_2),p_2\mid lcm(p_1,p_2)\)

    现在目的找到这样的一个特解 \(k\)\[ k=k_1\cdot p_1+a_1,k=k_2\cdot p_2+a_2\Rightarrow k_1\cdot p_1+a_1=k_2\cdot p_2+b_2\\ \Rightarrow k_1\cdot p_1-k_2\cdot p_2=b_2-b_1 \] 可以用 \(exgcd\) 求一组 \(k_1,k_2\) 特解。

    ·

  3. \(Miller\ rabin\)

    1. 首先费马小定理判断一波,\(p\in \text{prime} \Leftrightarrow \forall (a,p)=1,a^{p-1}\equiv1\pmod{p}\)

    2. 如果 \(a^{2}\equiv 1\pmod{p}\),如果 \(p\) 是质数,那么必定有 \(p\mid(a-1)\),或者 \(p\mid(a+1)\)

      否则 \(p\) 是合数。

      \(a^{p-1}\equiv 1\pmod{p}\)

      此时如果 \(p\)是质数,\(a^{\frac{p-1}{2}}\equiv 1\ or\ -1\pmod{p}\)

      如果是 \(-1\) 的时候,继续不下去,只得暂时假装他是素数。

      否则继续检查 \(a^{\frac{p-1}{2}}\equiv 1\pmod{p}\)

      本质上我们就是通过看 \(a^b\equiv 1\pmod{p}\) 这种形式来看 \(p\) 是否为质数。

      1
      2
      3
      4
      5
      6
      7
      8
      9
      10
      11
      12
      13
      14
      15
      16
      17
      18
      19
      20
      21
      22
      23
      24
      25
      inline bool mr(int x,int p)
      {
      int b=x-1;
      while(!(b&1)) b>>=1;
      int k=qpow(p,b,x);
      while(b<=x-1)
      {
      int muls=k*k%x;
      if(muls==1&&(k!=1&&k!=x-1)) return 0;
      k=muls;b<<=1;
      }
      if(k!=1) return 0;
      return 1;
      }
      int p[]={2,3,5,7,11,13,15,17,19,23};
      inline bool judge(int x)
      {
      if(!(x&1)) return 0;
      for(int i=0;i<=8;i++)
      {
      if(x==p[i]) return 1;
      if(mr(x,p[i])==0) return 0;
      }
      return 1;
      }

      \(p\) 的值域是long long 范围时,乘法可能会需要快速乘。

  4. 找阶&原根。

    这个东西说实话只要真正了解定义就能如同行云流水般解决了。

    1. 最小的使得 \(a^x\equiv 1\pmod{p}\) 的正整数 \(x\),被称为 \(a\bmod p\) 的阶。

      由于 \(a^{\phi(p)}\equiv 1\pmod{p}\),所以 \(ord_p(a)\mid\phi(p)\)

      证明很显然,如果不整除,可以用类似辗转相除的方法推出,\(\exists \ x<ord_p(a),a^x\equiv 1\pmod{p}\)

      求法可以把 \(\phi(p)\) 的质因子一个一个消去,来看是否有 \(a^{\frac{\phi(p)}{k}}\equiv 1\pmod{p}\)

      复杂度 \(O(pollard-\rho)+\log(n)\)

      其实,如果我们知道阶的大小,也能很好知道阶为这个数的数有多少个。

      若阶为 \(d\)

      那么 \(ord_{n}(x)=d\)\(x\) 一定有 \(\phi(d)\) 个。

      考虑我们随便找一个 \(ord_n(x)=d\),那么 \(x,x^2,cdots,x^{d-1},x^{d}\) 这样的一个集合。

      显然我们应该在指数上与 \(d\) 互质的数。

      所以 \(ord_x(x^k)\) 当且仅当 \((k,d)=1\),此时也就是 \([1,d]\)\(d\) 互质的数,也就是 \(\phi(d)\) 个。

    2. 原根

      满足 \(ord_p(a)=\phi(p)\)\(a\) 被叫做 \(p\) 的一个原根。

      原根性质:若 \(g\)\(p\) 其中一个原根,那么 \(1,g,g^2\cdots g^{\phi(p)-1}\) 互不相同。

      证明可以用反证法。

      由于 如果 \((k,p)=1\),有 \(0,k,2k,\cdots,(p-1)k\bmod p\) 一定互不相同。

      也是反证,如果 \(xk\equiv yk\pmod{p}\),那么 \((x-y)k\equiv 0\pmod{p}\),与 \((k,p)=1\) 矛盾。

      所有 \(p\)\(\phi(\phi(p))\) 个原根。

      原根有很多,随机找可以在 \(O(\frac{p}{\phi(\phi(p))})\) 的期望次找到。

      \(\phi(p)=\prod p_i^{\alpha_i}\)

      那么判断一个数 \(x\) 是不是原根就看 是否 \(\forall \ k,x^\frac{\phi(p)}{p_k}\not \equiv 1\pmod{p}\)

      复杂度 \(O(pollard-\rho)+\frac{p}{\phi(\phi(p))}\log(n)\)

  5. \(pollard-\rho\)

    \(\mathcal(O(\sqrt p)\)\(p\)\(n\) 最小质因子的复杂度,找到 \(n\) 的一个因数。

    考虑对于最小质因子 \(p\),我们随机枚举,枚举到 \(p\) 的倍数的概率为 \(\frac{1}{p}\)

    根据生日悖论,我们随机枚举,发生碰撞的概率则为 \(\frac{1}{\sqrt p}\)

    似乎问题马上就解决了。

    但是考虑我们并不知道 \(p\) 的具体值,判断碰撞也只能使用 \(\gcd(\mid a-b\mid,n)\) 这么来判断。

    那该怎么办?

    \(pollard\) 发现了一种伪随机数列 \(a_{i+1}=a_i^2+c\bmod m\)

    可以证明该数列的形状为一条链加上一个环。

    如果我们认为该数列是纯随机的那,环和链的长度均为期望 \(\sqrt m\) 级别的,根据生日悖论。


    现在我们的任务就是给定 \(c,a_0\) 的条件下找到这个环。

    但是,注意这样一点,我们并不能得到 \(a\) 的值,因为我们甚至都得不到 \(p\) 的值。

    我们接下来做的都是用 \(\bmod n\) 的数列,反应 \(\bmod p\) 的数列。

    找环的过程可以使用 \(a_i,a_{2i}\) 看什么呢时候冲突。

    设链长为 \(l\),环长为 \(m\)

    显然当 \(i=km>l\) 时,\(a_i=a_{l+(i-l)}=a_{l+(i-l+km)}=a_{i+km}=a_{2i}\)

    也就是说,碰撞一定会产生。

    在考虑 \(\bmod{p}\) 的数列,环和链长均为期望 \(\sqrt p\) 级别的,并且一定会发生碰撞,那么找到 \(p\) 的次数,一定是期望 \(\sqrt p\) 次。

    那么一个复杂度大概为 \(\mathcal{O(\log_2 n\log n\sqrt p)}\)\(pollard-\rho\) 算法出现了。


    考虑在 \(a_i,a_{2i}\) 找环过程,我们是对每一组 \(\mid a_i-a_{2i}\mid\) 都和 \(n\)\(\gcd\) 了。

    这让我们有将所有 \(\mid a_i-a_{2i}\mid\) 乘起来一起做 \(\gcd\)

    考虑 \(\gcd(a,b)=\gcd(a,b\bmod a)\),我们可以让 \(\mid a_{i}-a_{2i}\mid\)\(\bmod n\) 乘积。

    那就可以分块,每 \(128\) 分成一块,看 \(\gcd\),这样我们可以认为一次找因数的复杂度被均摊成 \(\mathcal{O(\sqrt p)}\)

    如果 \(\sqrt p\) 比较小时,\(128\) 会太大,所以可以考虑倍增。

    复杂度可以认为是 \(\mathcal{O(\log_2 n\sqrt p)}\)。由于每次找到因数后 \(n\) 的最小因数会改变,所以复杂度可能会更优,但是具体怎么证明我不是很会,在具体应用时你可以认为前面那个 \(\log_2\) 因子会十分的小,甚至忽略。

    \(\mathcal{O(\log_2 n \sqrt p)}\) 是期望复杂度上界。


    实现,边界。

    1. 如果发现 \(a_i=a_{2i}\) 立马停止,重新换 \(a_0,c\) 重新 \(pollard-\rho\)
    2. 如果 \(n=4\) 那么 \(a_i=a_{2i}\) 据说不可避免,可以特判掉。
    3. 如果 \(n\) 为质数,立马停止。

    大概代码:

    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    14
    15
    16
    17
    18
    19
    20
    21
    22
    23
    24
    25
    26
    27
    28
    29
    30
    31
    32
    33
    34
    35
    36
    37
    38
    39
    40
    41
    inline ll nxt(ll a,ll c,ll mod){
    ll tmp=mul(a,a,mod)+c;
    return tmp>=mod?tmp-mod:tmp;
    }
    inline ll gcd(ll a,ll b){
    if(!b) return a;
    return gcd(b,a%b);
    }
    inline ll pollard(ll n){
    if(n==4) return 2;
    ll a1=1ll*rand()*rand()%n,c=1ll*rand()*rand()%n;
    ll a2=a1;
    a1=nxt(a1,c,n),a2=nxt(nxt(a2,c,n),c,n);
    int lim=1;
    while(a1!=a2){
    ll v=1;
    for(int i=1;i<=lim&&a1!=a2;i++,a1=nxt(a1,c,n),a2=nxt(nxt(a2,c,n),c,n)){
    ll backv=v;
    v=mul(v,abs(a1-a2),n);
    if(!v) return gcd(abs(a1-a2),n);
    }
    ll g=gcd(v,n);
    if(g>1) return g;
    if(lim<128) lim<<=1;
    }
    return pollard(n);
    }
    vector<ll> f;
    inline void factor(ll n){
    if(prime(n)){
    f.push_back(n);
    return ;
    }
    ll fac=pollard(n);
    if(prime(fac)){
    f.push_back(fac);
    return factor(n/fac);
    }
    factor(fac);
    factor(n/fac);
    }


基本数论
https://proton-z.github.io/2022/08/30/basic-number-theory/
作者
Imitators
发布于
2022年8月30日
许可协议