大整数分解模板

#include<bits/stdc++.h>

using namespace std;

typedef long long LL;
typedef long long ll;

ll ModMul(ll a,ll b,ll n)//快速积取模 a*b%n
{
    ll ans=0;
    while(b)
    {
        if(b&1)
            ans=(ans+a)%n;
        a=(a+a)%n;
        b>>=1;
    }
    return ans;
}
ll ModExp(ll a,ll b,ll n)//快速幂取模 a^b%n
{
    ll ans=1;
    while(b)
    {
        if(b&1)
            ans=ModMul(ans,a,n);
        a=ModMul(a,a,n);
        b>>=1;
    }
    return ans;
}
bool Miller_Rabin(ll n)//Miller-Rabin素数检测算法
{
      ll i,j,a,x,y,t,u,s=10;
      if(n==2)
            return true;
      if(n<2||!(n&1))
            return false;
      for(t=0,u=n-1;!(u&1);t++,u>>=1);//n-1=u*2^t
      for(i=0;i<s;i++)
      {
          a=rand()%(n-1)+1;
          x=ModExp(a,u,n);
          for(j=0;j<t;j++)
          {
              y=ModMul(x,x,n);
              if(y==1&&x!=1&&x!=n-1)
                    return false;
              x=y;
          }
          if(x!=1)
                return false;
      }
      return true;
}


LL Pollard_Rho(LL n,LL c)
{
    LL i=1,j=2,x=rand()%(n-1)+1,y=x;//随机初始化一个基数（p1）
    while(1)
    {
        i++;
        x=(ModMul(x,x,n)+c)%n;//玄学递推
        LL p=__gcd((y-x+n)%n,n);
        if(p!=1&&p!=n)return p;//判断
        if(y==x)return n;//y为x的备份，相等则说明遇到圈，退出
        if(i==j)
        {
            y=x;
            j<<=1;
        }//更新y，判圈算法应用
    }
}
void find(LL n,LL c)//同上，n为待分解数，c为随机常数
{
    if(n==1)return;//1要特判
    if(Miller_Rabin(n))//n为质数
    {
        //保存，根据不同题意有不同写法，在此略去，可以存也可以直接计算
        //printf("%lld\n",n);
        return;
    }
    LL x=n,k=c;
    while(x==n)x=Pollard_Rho(x,c++);//当未分解成功时，换个c带入算法
    find(n/x,k);
    find(x,k);
    //递归操作
}