模板之数论大全1

来源:互联网 发布:天刀明月心捏脸数据女 编辑:程序博客网 时间:2024/05/09 21:26

1.扩展的欧几里德定理

//拓展欧几里得定理,求ax+by=gcd(a,b)的一组解(x,y),d=gcd(a,b)void gcd(int a,int b,int &d,int &x,int &y){    if(!b){d=a;x=1;y=0;}    else{gcd(b,a%b,d,y,x);y-=x*(a/b);}}

用法1:

求ax+by=c的整数解。ax+by=gcd(a,b)=g的一组解为(x0,y0),则ax+by=c的一组解为(x0*c/g,y0*c/g)。当c不是g的倍数时无整数解

若(x1,y1)是ax+by=c的一组解,则其任意整数解为(x1+k*bb,y1-k*aa),其中aa=a/gcd(a,b),bb=b/gcd(bb),k为任意整数

2.求模乘法的逆

//求得a在模n条件下的逆int inv(int a,int n){   int d,x,y;   gcd(a,n,d,x,y);   return d==1?(x+n)%n:-1;}

或则:

v=pow_mod(a,n-m-1,n);//n为素数,pow_mod(a,n-1,n)=1,费马小定理。所以a^m*a^(n-m-1)=a^(n-1)=1(mod n).

3.快速幂求解a^b

//快速幂求a^bint pow_mod(int a,int b){    int s=1;    while(b)    {        if(b&1)            s=(s*a)%mod;        a=(a*a)%mod;        b=b>>1;    }    return s;}

4.求解模方程a^x=b(mod n)。

1.n为素数。无解返回-1

//求解模方程a^x=b(mod n)。n为素数,无解返回-1int log_mod (int a,int b,int n){    int m,v,e=1,i;    m=ceil(sqrt(n+0.5));    v=inv(pow_mod(a,m),n);    map<int,int>x;    x[1]=0;    for(i=1;i<m;i++)    {        e=(e*a)%n;        if(!x.count(e))x[e]=i;    }    for(i=0;i<m;i++)    {        if(x.count(b))return i*m+x[b];        b=(b*v)%n;    }    return -1;}

2.n不是素数。

//hdu 2815 Mod Tree#include <cstdio>#include <cstring>#include <cmath>#include <map>#include <iostream>#include <algorithm>using namespace std;#define LL __int64LL gcd(LL a,LL b){    return b==0?a:gcd(b,a%b);}//拓展欧几里得定理,求ax+by=gcd(a,b)的一组解(x,y),d=gcd(a,b)void gcd_mod(LL a,LL b,LL &d,LL &x,LL &y){    if(!b){d=a;x=1;y=0;}    else{gcd_mod(b,a%b,d,y,x);y-=x*(a/b);}}//求解模方程d*a^(x-c)=b(mod n)。d,a和n互质,无解返回-1LL log_mod (LL a,LL b,LL n,LL c,LL d){    LL m,v,e=1,i,x,y,dd;    m=ceil(sqrt(n+0.5));     //x=i*m+j    map<LL,LL>f;    f[1]=m;    for(i=1;i<m;i++)  //建哈希表,保存a^0,a^1,...,a^m-1    {        e=(e*a)%n;        if(!f[e])f[e]=i;    }    e=(e*a)%n;//e=a^m    for(i=0;i<m;i++)//每次增加m次方,遍历所有1<=f<=n    {        gcd_mod(d,n,dd,x,y);//d*x+n*y=1-->(d*x)%n=1-->d*(x*b)%n==b        x=(x*b%n+n)%n;        if(f[x])        {            LL num=f[x];            f.clear();//需要清空,不然会爆内存            return c+i*m+(num==m?0:num);        }        d=(d*e)%n;    }    return -1;}int main(){    LL a,b,n;    while(scanf("%I64d%I64d%I64d",&a,&n,&b)!=EOF)    {        if(b>=n)        {            printf("Orz,I can’t find D!\n");            continue;        }        if(b==0)        {            printf("0\n");            continue;        }        LL ans=0,c=0,d=1,t;        while((t=gcd(a,n))!=1)        {            if(b%t){ans=-1;break;}            c++;            n=n/t;            b=b/t;            d=d*a/t%n;            if(d==b){ans=c;break;}//特判下是否成立。        }        if(ans!=0)        {            if(ans==-1){printf("Orz,I can’t find D!\n");}            else printf("%I64d\n",ans);        }        else        {            ans=log_mod(a,b,n,c,d);            if(ans==-1)printf("Orz,I can’t find D!\n");            else printf("%I64d\n",ans);        }    }    return 0;}/*    求解模方程a^x=b(mod n),n不为素数。拓展Baby Step Giant Step    模板题。        方法:    初始d=1,c=0,i=0;    1.令g=gcd(a,n),若g==1则执行下一步。否则由于a^x=k*n+b;(k为某一整数),则(a/g)*a^k=k*(n/g)+b/g,(b/g为整除,若不成立则无解)令n=n/g,d=d*a/g,b=b/g,c++则d*a^(x-c)=b(mod n),接着重复1步骤。    2.通过1步骤后,保证了a和d都与n互质,方程转换为d*a^(x-c)=b(mod n)。由于a和n互质,所以由欧拉定理a^phi(n)==1(mod n),(a,n互质)可知,phi(n)<=n,a^0==1(mod n),所以构成循环,且循环节不大于n。从而推出如果存在解,则必定1<=x<n。(在此基础上我们就可以用Baby Step Giant Step方法了)    3.令m=ceil(sqrt(n)),则m*m>=n。用哈希表存储a^0,a^1,..,a^(m-1),接着判断1~m*m-1中是否存在解。    4.为了减少时间,所以用哈希表缩减复杂度。分成m次遍历,每次遍历a^m长度。由于a和d都与n互质,所以gcd(d,n)=1,所以用拓展的欧几里德定理求得d*x+n*y=gcd(d,n)=1,的一组整数解(x,y)。则d*x+n*y=1-->d*x%n=(d*x+n*y)%n=1-->d*(x*b)%n=((d*x)%n*b%n)%n=b。所以若x*b在哈希表中存在,值为k,则a^k*d=b(mod n),答案就是ans=k+c+i*m。如果不存在,则令d=d*a^m,i++后遍历下一个a^m,直到遍历a^0到a^(m-1)还未找到,则说明不解并退出。*/



5.中国剩余定理

互质情况

//中国剩余定理,求得M%A=a,M%B=b,...中的M,其中A,B,C...互质  int china(int a[]){      int i,j,k,d,ans=0,x,y,M=1;    for(i=0;i<n;i++)          M=M*x[i];    for(i=0;i<n;i++)      {          m=M/x[i];          gcd(x[i],m,d,x,y);          ans=(ans+y*m*a[i])%M;      }      return (ans+M)%M;  }

不互质的情况

int China(int n)  {      int m1,r1,m2,r2,flag=0,i,d,x,y,c,t;      scanf("%d%d",&m1,&r1);      flag=0;      for(i=1;i<n;i++)      {          scanf("%d%d",&m2,&r2);          if(flag)continue;          gcd(m1,m2,d,x,y);//d=gcd(m1,m2);x*m1+y*m2=d;          c=r2-r1;          if(c%d)//对于方程m1*x+m2*y=c,如果c不是d的倍数就无整数解          {              flag=1;              continue;          }          t=m2/d;//对于方程m1x+m2y=c=r2-r1,若(x0,y0)是一组整数解,那么(x0+k*m2/d,y0-k*m1/d)也是一组整数解(k为任意整数)                  //其中x0=x*c/d,y0=x*c/d;          x=(c/d*x%t+t)%t;//保证x0是正数,因为x+k*t是解,(x%t+t)%t也必定是正数解(必定存在某个k使得(x%t+t)%t=x+k*t)          r1=m1*x+r1;//新求的r1就是前i组的解,Mi=m1*x+M(i-1)=r2-m2*y(m1为前i个m的最小公倍数);对m2取余时,余数为r2;                      //对以前的m取余时,Mi%m=m1*x%m+M(i-1)%m=M(i-1)%m=r          m1=m1*m2/d;      }      if(flag)return -1;     else return r1;  }  


6.素数

普通筛选法

const int maxn=10000001;const int maxc=700000;int vis[maxn];//vis[i]=0时,i为素数或者1,否则为合数int prime[maxc];//筛素数void sieve(int n){    int m=(int)sqrt(n+0.5);//避免浮点误差    memset(vis,0,sizeof(vis));    for(int i=2;i<=m;i++)if(!vis[i])        for(int j=i*i;j<=n;j+=i)vis[j]=1;}//生成素数表,存在prime数组中,返回素数个数int primes(int n){    sieve(n);    int t=0;    for(int i=2;i<=n;i++)if(!vis[i])        prime[t++]=i;    return t;}

线性筛法

:在O(n)时间复杂度内找出所有的素数

void init()//预处理,找出所有1e7以内的素数,以减少查找1e14范围数的因子的时间  {           //现行筛素数的方法,时间复杂度为O(n)      memset(check,false,sizeof(check));      int i,j;      tot=0;      for(i=2;i<=1e7;i++)      {          if(!check[i])prime[tot++]=i;          for(j=0;j<tot;j++)          {              if(i*prime[j]>1e7)break;              check[i*prime[j]]=true;              if(i%prime[j]==0)break;          }      }      //printf("%d\n",tot);      //for(i=0;i<20;i++)      //    printf("prime[%d]:%d\n",i,prime[i]);  } 




7.欧拉phi函数

,求不超过n且与n互质的正整数个数

int euler_phi(int n)//求单个欧拉函数{    int m=(int)sqrt(n+0.5);    int i,ans=n;    for(i=2;i<=m;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;}

int phi[maxn];  void euler_phi()  {      int i,j,k;      //欧拉函数,phi[i]表示不超过i的与i互质的整数个数      for(i=2;i<maxn;i++)phi[i]=0;      phi[1]=1;      for(i=2;i<maxn;i++)          if(!phi[i])          for(j=i;j<=maxn;j+=i){              if(!phi[j])phi[j]=j;              phi[j]=phi[j]/i*(i-1);          }  } 

8.求三角形两边向量积

,面积S=1/2*|a×b|

int cross(node a,node b,node c)//向量积{    return (a.x-c.x)*(b.y-c.y)-(b.x-c.x)*(a.y-c.y);//向量ac×bc}

9.求凸包

,res数组存凸包上的点。

int convex(int n){    sort(e,e+n,cmp);//根据x从小到大排序,x相等时根据y从小到大排序    int m=0,i,j,k;    //求得下凸包,逆时针    for(i=0;i<n;i++)    {        while(m>1&&cross(res[m-1],e[i],res[m-2])<=0)m--;        res[m++]=e[i];    }    k=m;    //求得上凸包    for(i=n-2;i>=0;i--)    {        while(m>k&&cross(res[m-1],e[i],res[m-2])<=0)m--;        res[m++]=e[i];    }    if(n>1)m--;//起始点重复。    return m;//返回凸包边界结点数。}


10.三分法

for(i=0;i<100;i++){    m1=l+(r-l)/3;    m2=r-(r-l)/3;    if(find(m1)<find(m2)) r=m2;    else l=m1;}



原创粉丝点击