欧几里得-GCD及扩展

来源:互联网 发布:淘宝网外卖 编辑:程序博客网 时间:2024/06/05 19:07

GCD最大公约数

int gcd(int x,int y){    if(!x || !y) return x>y?x:y;    for(int t;t=x%y;x=y,y=t);    return y;}

快速GCD

int kgcd(int a,int b){    if(a==0) return b;    if(b==0) return a;    if(!(a & 1) && !(b & 1)) return kgcd(a>>1,b>>1)<<1;    else if(!(b & 1)) return kgcd(a,b>>1);    else if(!(a & 1)) return kgcd(a>>1,b);    else return kgcd(abs(a-b),min(a,b));}

扩展GCD
求x,y使得gcd(a,b)=a*x+b*y;

int extgcd(int a,int b,int &x,int &y){    if(b==0) {x=1;y=0;return a;}    int d=extgcd(b,a%b,x,y);    int t=x;    x=y;    y=t-a/b*y;    return d;}

模线性方程 a*x=b(%n)

void modeq(int a,int b,int n) //! n>0{    int e,i,d,x,y;    d=extgcd(a,n,x,y);    if(b%d>0) printf("No answer!\n");    else     {        e=(x*(b/d))%n;        for(i=0;i<d;i++) //!!! here x maybe<0        {            printf("%d-th ans:%d\n",i+1,(e+i*(n/d))%n);        }    }}

模线性方程组
a=b[1](%w[1]);a=b[2](%w[2]);...a=b[k](%w[k]);
其中w,b已知,w[i]>0且w[i]与w[j]互质,求a;(中国剩余定理)

int china(int b[],int w[],int k){    int i,d,x,y,m,a=0,n=1;    for(i=0;i<k;i++) n*=w[i];    for(i=0;i<k;i++)    {        m=n/w[i];        d=extgcd(w[i],m,x,y);        a=(a+y*m*b[i])%n;    }    if(a>0) return a;    else return (a+n);}

不要求w[i]与w[j]互质;
#include <iostream>  #include <cstdio>  #include <cstring>  using namespace std;  typedef __int64 int64;  int64 Mod;  int64 gcd(int64 a, int64 b)  {      if(b==0)          return a;      return gcd(b,a%b);  }    int64 Extend_Euclid(int64 a, int64 b, int64&x, int64& y)  {      if(b==0)      {          x=1,y=0;          return a;      }      int64 d = Extend_Euclid(b,a%b,x,y);      int64 t = x;      x = y;      y = t - a/b*y;      return d;  }    //a在模n乘法下的逆元,没有则返回-1  int64 inv(int64 a, int64 n)  {      int64 x,y;      int64 t = Extend_Euclid(a,n,x,y);      if(t != 1)          return -1;      return (x%n+n)%n;  }    //将两个方程合并为一个  bool merge(int64 a1, int64 n1, int64 a2, int64 n2, int64& a3, int64& n3)  {      int64 d = gcd(n1,n2);      int64 c = a2-a1;      if(c%d)          return false;      c = (c%n2+n2)%n2;      c /= d;     n1 /= d;      n2 /= d;      c *= inv(n1,n2);      c %= n2;      c *= n1*d;      c += a1;      n3 = n1*n2*d;      a3 = (c%n3+n3)%n3;      return true;  }  //求模线性方程组x=ai(mod ni),ni可以不互质  int64 China_Reminder2(int len, int64* a, int64* n)  {      int64 a1=a[0],n1=n[0];      int64 a2,n2;      for(int i = 1; i < len; i++)      {          int64 aa,nn;          a2 = a[i],n2=n[i];          if(!merge(a1,n1,a2,n2,aa,nn))              return -1;          a1 = aa;          n1 = nn;      }      Mod = n1;      return (a1%n1+n1)%n1;  }  int64 a[1000],b[1000];  int main()  {      int i;      int k;      while(scanf("%d",&k)!=EOF)      {          for(i = 0; i < k; i++)              scanf("%I64d %I64d",&a[i],&b[i]);          printf("%I64d\n",China_Reminder2(k,b,a));      }      return 0;  }  


原创粉丝点击