acm数论之欧几里得gcd

叁歲伎倆 2021-12-15 04:07 282阅读 0赞

## 1.欧几里得定理 ##

1.  *同余定理的公式：(a+b)%mod=(a%mod+b%mod)%mod* 
2.   (a\*b)%mod=(a%mod\*b%mod)%mod
    
     
3.   扩展欧几里得也有自己的一个公式：a\*x+b\*y=gcd(a,b)

##### **1.1定义：** #####

##### gcd(a,b)=gcd(b,a%b) #####

##### 1.2用法： #####

##### 获得最大公约数与最小公倍数\----gcd(a,b)\*lcm(a,b)=a\*b; #####

int gcd(int a,int b)   //最大公约数
    {
        if(b==0)  return a;
        else return gcd(b,a%b);
    }     
    int lcm(int a,int b)  //最小公倍数
    {
        return a/gcd(a,b)*b;    //防止溢出
    }

**2.扩展欧几里得**  
2.1定义：  
ax+by=gcd(a,b)=d，在已知a,b的情况下求解出一组x,y

2.2推导过程：  
因为a%b=a-(a/b)b,  则有 d=bx1+\[a-(a/b)b\]y1=bx1+ay1-(a/b)by1=ay1+b(x1-a/by1),  故 x=y1,y=x1-a/by1

2.3性质：  
**1.若通过扩展欧几里得求出一组特解(x0,y0)，那么有ax0+by0=d.**

则方程的通解为，其中k为任意整数   x0=d/exgcd()\*x，y0=d/exgcd()\*y；  
      x=x0+k\*(b/d)--------x = x0 +( b / gcd( a, b ) ) \* k  
      y=y0-k\*(a/d)--------y = y0 - ( a / gcd( a, b ) ) \* k  
**2.已知ax+by=d的解，对于ax+by=c的解**，c为任意正整数，只有当d|c时才有解  
    其通解为  
    x=(c/d)x0+k(b/d)  
    y=(c/d)y0-k(a/d)

//扩展欧几里得处理ax+by=gcd(a,b)=d解的问题 
    int exgcd(int a,int b,int &x,int &y)  //得到gcd(a,b)的值 
    {
        if(b==0)
        {
            x=1,y=0;
            return a;
         } 
         int r=exgcd(b,a%b,x,y);
         int t=x; x=y; y=t-a/b*y;
         return r;//gcd(a,b)
    }

**2.4常见用法：**

（1）求形如ax+by=c的通解，或从中选取某些特解  
（2）求乘法逆元  
（3）求解线性同余方程

**写的超级好的博客：[http://www.cnblogs.com/wkfvawl/p/9350867.html][http_www.cnblogs.com_wkfvawl_p_9350867.html]**

**解题思路**：

1**.先化成不定方程ax+by=c**

**2.a,b参数取正数，判断a,b是否为负数，为负数就a,b,c都乘上-1，另外注意若a,b中只有一个是负数，那另一个不需要乘-1，因为a,b本身是含未知数参数的。**

**3.r=exgcd(a,b,x,y),然后判断是否有解采用if(c%r==0)**

**求特解：x0=x\*c/r,   y0=y\*c/r;**

**求通解：x=x0+b/r\*c , y=y0-a/r\*c （其中 t 为整数）**

**求最小解：int s=b/r;    minx= (x0%s+s)%s;//最小解**

**2.4.1例题实例：**

**问题描述：对于 ax+by=c 的不定方程求通解或特解？**

**设 r=gcd(a,b)，**

**若 c%r!=0 这此方程无整数解**

**若 c%r==0，特解为x1=x0\*c/r , y1=y0\*c/r，通解 x=x1+b/r\*t , y=y1-a/r\*t （其中 t 为整数）**

**不定方程的最小解：**

**r=exgcd(a,b,x,y);**

**x=x\*c/r;**  
 **int s=b/r;**  
**minx= (x%s+s)%s;//最小解**

**例题1：**

题目链接：[https://cn.vjudge.net/problem/HDU-2669][https_cn.vjudge.net_problem_HDU-2669]

题目大意：0<a, b<=2^31，X>=0，给定a,b求满足X\*a + Y\*b = 1的最小正数X，以及与之对应的Y。

解题思路：求形如ax+by=c的最小解，x=x0+k\*(b/d)， y=y0-k\*(a/d)，从通解可看出，其中b / d 都正数，x 和 t成正比，当x0为正数的时候 t == 0，有最小值；当x0为负数的时候，x最小就是0，与横轴相交，所以可以求出 x0 + b / t \* t = 0时候的 t，因为此时 t是整数，不一定是准确相交的那一个点，所以当前的x可能是负数，如果为负数，t++，取下一个整数点。

ac代码：

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

1 #include<iostream>
     2 using namespace std;
     3 //扩展欧几里得处理ax+by=gcd(a,b)=d解的问题 
     4 int exgcd(int a,int b,int &x,int &y)  //得到gcd(a,b)的值 
     5 {
     6     if(b==0)
     7     {
     8         x=1,y=0;
     9         return a;
    10      } 
    11      int r=exgcd(b,a%b,x,y);
    12      int t=x; x=y; y=t-a/b*y;
    13      return r;//gcd(a,b)
    14 }
    15 int main()
    16 {
    17     int a,b;
    18     while(~scanf("%d%d",&a,&b))
    19     {
    20         int x,y;
    21         int d=exgcd(a,b,x,y);
    22         if(1%d)//有余数 
    23         {
    24             cout<<"sorry"<<endl; 
    25             continue;
    26         }
    27         else
    28         {
    29             //通过exgcd函数得到特解 
    30             x=x/d,y=y/d;
    31             if(x>0)//由于a,b为正数，所以k=0取最小 
    32             cout<<x<<" "<<y<<endl;
    33             else
    34             {
         //当x0为负数的时候，x最小就是0，与横轴相交, 但x可能是负数，如果为负数，取下一个整数点
    35                 int k=(-1)*x/b*d;
    36                 int kx=x+b/d*k;
    37                 int ky=y-k*(a/d);
    38                 if(kx<0)
    39                 {
    40                     kx=x+b/d*(k+1);
    41                     ky=y-(k+1)*(a/d);
    42                 }
    43                 cout<<kx<<" "<<ky<<endl;                
    44             }
    45         }    
    46     }
    47     return 0; 
    48  }

**例题2：**

题目链接：[https://cn.vjudge.net/problem/ZOJ-3593][https_cn.vjudge.net_problem_ZOJ-3593]

题意：一个人要从A走到B  只能走a布、b步、(a+b)步，可以往左或右走， 问 最少要走几步，不能走到的输出-1

题目思路：可以设走x步a米的，走y部b米的，得到式子：ax + by = B - A;

通过扩展欧几里德可得x和y

x、y的解集：X = x + b/gcd(a,b) \* t;   Y = y - a /gcd(a,b) \* t;

x与y最接近时为答案，这样的话c可以取的更多，那么步数会最小

它们

（1）X、Y同号时，即方向相同，可以走a+b来减少步数，取max(X,Y)（这里可以画X、Y的坐标图来帮助理解）

（2）X、Y异号时，即方向相反，步数绝对值相加，abs(X) + abs(Y)

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

1 #include<iostream>
     2 #include<cstdio>
     3 #include<cmath>
     4 using namespace std;
     5 const long long maxn = 0x3f3f3f3f;
     6 //#define maxn 0x3f3f3f3f
     7 long long exgcd(long long a,long long b,long long &x,long long &y)
     8 {
     9     if(b==0)
    10     {
    11         x=1,y=0;
    12         return a;
    13     }
    14     long long r=exgcd(b,a%b,x,y);
    15     long long t=x;x=y;y=t-a/b*y;
    16     return r;
    17 }
    18 long long Abs(long long a)
    19 {
    20     if(a<0)
    21     a=-1*a;
    22     return a;
    23 }
    24 int main()
    25 {
    26     long long a,b,A,B;
    27     long long ans;
    28     int T;
    29     cin>>T;
    30     while(T--)
    31     {
    32         cin>>A>>B>>a>>b;
    33         long long x,y;
    34         long long h=exgcd(a,b,x,y);
    35         if((B-A)%h!=0)
    36             printf("-1\n");
    37         else{
    38             x=(B-A)/h*x;//注释修改
    39             y=(B-A)/h*y;//注释修改
    40             long long minx,miny;
    41             long long k;
    42             b=b/h;
    43             a=a/h;
    44             k=(y-x)/(b+a);//注释修改
    45             ans=maxn*maxn;   //注释修改  这里定义成maxn数据定义小了
    46             for(int i=k-1;i<=k+1;i++)
    47             {
    48                 minx=x+b*i,miny=y-a*i;
    49                 long long m;
    50                 if(Abs(minx+miny)==Abs(minx)+Abs(miny))//说明x与y同号
    51                     m=max(Abs(minx),Abs(miny));
    52                 else
    53                     m=Abs(minx)+Abs(miny);
    54                 if(m<ans)
    55                  ans=m;
    56 //                if((minx>0&&miny>0)||(minx<0&&miny<0))//x与y同号
    57 //                {
    58 //                m=Abs(min(minx,miny))+Abs(minx-miny);//这句有问题的  测试这个存在数据Abs(min(minx,miny))+Abs(minx-miny)！=max(Abs(minx),Abs(miny));
    59 //                }
    60 //                 else
    61 //                 m=Abs(minx)+Abs(miny);
    62 //                 if(m<ans)
    63 //                 ans=m;
    64 
    65             }
    66             cout<<ans<<endl;
    67         }
    68     }
    69 }

ac代码

wa的心得或者可以说是经验：

1.maxn定义不同而wa:

const long long maxn = 0x3f3f3f3f;   //正确  
//\#define maxn 0x3f3f3f3f   错误

2.除法运算要考虑是否倍数除，可以倍数除最好倍数除

x=(B-A)/h\*x;  y=(B-A)/h\*y;//正确

x=x/h\*(B-A);  y=y/h\*(B-A);  //错误

3.ans=maxn\*maxn;   //注释修改  这里定义成maxn数据定义小了

3.前一个为正确写法，后一个为错误写法

分析：两者的不同在于在for循环里面minx,miny，后一种不能保证i=k-1~k+1这里面取得答案，因为多了一个参数。

1 b=b/h;
    2 a=a/h;
    3 k=(y-x)/(b+a);//注释修改
    4 for(int i=k-1;i<=k+1;i++)
    5 minx=x+b*i,miny=y-a*i;

1 b=b;
    2 a=a;
    3 k=(y-x)/(b+a)*h;
    4 for(int i=k-1;i<=k+1;i++)
    5 minx=x+b*i/h,miny=y-a*i/h;

**2.4.2例题实例：求乘法逆元**

**问题描述：求解一个数 a 对于另一个数 m 的乘法逆元：a\*x=1(mod m) 也即：a\*x%m=1求x   也即a\*x+m\*y=1 **

**当gcd(a , m) != 1 ,a\*x + m\*y = 1没有解**

**当1 % gcd(a , m) == 0，a\*x + m\*y = 1有解**

**x0=x\*1/exgcd();**

**最小的解:x0 % m**

**x 的通解: x0 + m\*t **

#include<iostream>
    #include<cmath>
    using namespace std;
    //扩展欧几里得处理ax+by=gcd(a,b)=d解的问题 
    int exgcd(int a,int b,int &x,int &y)  //得到gcd(a,b)的值 
    {
        if(b==0)
        {
            x=1,y=0;
            return a;
         } 
         int r=exgcd(b,a%b,x,y);
         int t=x; x=y; y=t-a/b*y;
         return r;//gcd(a,b)
    }
    //求解一个数 a 对于另一个数 m 的乘法逆元：a*x=1(mod m) 也即：a*x%m=1求x   也即a*x+m*y=1 
    int inversenum(int a,int m)
    {
        int x,y,ans,gcd;
        gcd=exgcd(a,m,x,y);
        if(1%gcd!=0)//无解 
        return -1;
        x=x*1/gcd;
        m=abs(m);
        ans=x%m;
        if(ans<=0)
        ans+=m;
        return ans;
    }
    int main()
    {
        int a,m;
        cin>>a>>m;
        cout<<inversenum(a,m)<<endl;
     }

##### 例题：ZOJ3609 Modular Inverse #####

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

1 #include<iostream>
     2 #include<cmath>
     3 using namespace std;
     4 //扩展欧几里得处理ax+by=gcd(a,b)=d解的问题 
     5 int exgcd(int a,int b,int &x,int &y)  //得到gcd(a,b)的值 
     6 {
     7     if(b==0)
     8     {
     9         x=1,y=0;
    10         return a;
    11      } 
    12      int r=exgcd(b,a%b,x,y);
    13      int t=x; x=y; y=t-a/b*y;
    14      return r;//gcd(a,b)
    15 }
    16 //求解一个数 a 对于另一个数 m 的乘法逆元：a*x=1(mod m) 也即：a*x%m=1求x   也即a*x+m*y=1 
    17 int inversenum(int a,int m)
    18 {
    19     int x,y,ans,gcd;
    20     gcd=exgcd(a,m,x,y);
    21     if(1%gcd!=0)//无解 
    22     return -1;
    23     x=x*1/gcd;
    24     m=abs(m);
    25     ans=x%m;
    26     if(ans<=0)
    27     ans+=m;
    28     return ans;
    29 }
    30 int main()
    31 {
    32     int T;
    33     cin>>T;
    34     while(T--){
    35     int a,m;
    36     cin>>a>>m;
    37     if(inversenum(a,m)==-1)
    38     cout<<"Not Exist"<<endl;
    39     else
    40     cout<<inversenum(a,m)<<endl;
    41     }    
    42  }

**2.4.2例题实例：求线性同余方程**

**问题描述：**关于 x 的模方程 ax%b=c (或者是(x+mt)%L=(y+nt)%L)的解，方程转换为 ax+by=c 其中 y 一般为非正整数

解得 x1=x0\*c/r，通解为 x=x1+b/r\*t

设 s=b/r （已证明 b/r 为通解的最小间隔），则 x 的最小正整数解为 (x1%s+s)%s

**例题：poj 1061**

根据题意：可写出(x+mt)%L=(y+nt)%L，求可取的最小t

(m-n)\*t=(y-x)%L;    (m-n)\*t+L\*b=y-x;   (m-n)\*a+L\*b=y-x;  看成a\*m1+b\*n1=t;(m1=m-n,n1=L,t=y-x)

注意：数据类型应该为long long

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

1 #include<iostream>
     2 #include<cmath>
     3 using namespace std;
     4 //扩展欧几里得处理ax+by=gcd(a,b)=d解的问题 
     5 long long t;
     6 long long exgcd(long long a,long long b,long long &x,long long &y)  //得到gcd(a,b)的值 
     7 {
     8     if(b==0)
     9     {
    10         x=1,y=0;
    11         return a;
    12      } 
    13      long long r=exgcd(b,a%b,x,y);
    14      long long t=x; x=y; y=t-a/b*y;
    15      return r;//gcd(a,b)
    16 }
    17 long long solve(long long m,long long n)
    18 {
    19     long long x,y;
    20     long long r=exgcd(m,n,x,y);
    21     if(t%r)
    22     {
    23         return -1;
    24     }
    25     x=x*t/r;
    26     long long s=n/r;
    27     return (x%s+s)%s;//不定方程的最小解x 
    28 }
    29 int main()
    30 {
    31     long long x,y,m,n,L;
    32     cin>>x>>y>>m>>n>>L;
    33     //不定式a*(m-n)+b*L=y-x;   看成a*m+b*n=t; 
    34     t=y-x;
    35     m=m-n;
    36     n=L;
    37     if(m<0)//整个不定式a,b取正数，由于b>0，b不需要变(b前面有系数),所以t需要改变 
    38     {
    39         m=-1*m;
    40         t=-1*t;
    41     }
    42     if(solve(m,n)==-1)
    43     cout<<"Impossible"<<endl;
    44     else
    45     {
    46         cout<<solve(m,n)<<endl;
    47     }
    48     return 0;
    49  }

**例题：p2421荒岛野人**

题目链接：[https://www.luogu.org/problemnew/show/P2421][https_www.luogu.org_problemnew_show_P2421]

**题目思路：**暴力+exgcd

首先判断两个野人是否会在同一个洞穴：(Ci+Pi\*t)%M==(Cj+Pj\*t)%M    (p\[i\]-p\[j\])\*t=(C\[j\]-C\[i\])%M     (P\[i\]-P\[j\])\*a+M\*b=C\[j\]-C\[i\]

暴力求M,从1递增，满足条件就停下来，check(M)

check(M)的作用是判断n个野人之间没有住在同一个洞穴，即余数不同或者是相同余数数字大于两野人最小寿命。

题目代码：

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

1 #include<iostream>
     2 #include<algorithm> 
     3 using namespace std;
     4 int c[200],p[200],l[200],n;
     5 int gcd(int a,int b){
         return a%b==0?b:gcd(b,a%b);}
     6 //扩展欧几里得处理ax+by=gcd(a,b)=d解的问题 
     7 int exgcd(int a,int b,int &x,int &y)  //得到gcd(a,b)的值 
     8 {
     9     if(b==0)
    10     {
    11         x=1,y=0;
    12         return a;
    13      } 
    14      int r=exgcd(b,a%b,x,y);
    15      int t=x; x=y; y=t-a/b*y;
    16      return r;//gcd(a,b)
    17 }
    18 bool go(int i,int j,int b){
    19       int a=p[i]-p[j],f=c[j]-c[i],g,x,y;
    20       if(a<0)
    21       {
    22           a=a*-1,f=-1*f;
    23       }
    24       g=exgcd(a,b,x,y);
    25       if(f%g==0){
    26           x=x*f/g;
    27           b=b/g;
    28           x=(x%b+b)%b;//最小解
    29           if(x<=min(l[i],l[j]))
    30           return 0;
    31       }
    32       return 1;
    33 }
    34 bool check(int m){
    35       int i,j,k;
    36       for(i=1;i<n;i++)
    37         for(j=i+1;j<=n;j++)
    38           if(!go(i,j,m))return 0;
    39       return 1;
    40 }
    41 int main()
    42 {     int m=0,i,j,k;
    43       scanf("%d",&n);
    44       for(i=1;i<=n;i++)scanf("%d%d%d",&c[i],&p[i],&l[i]),m=max(m,c[i]);
    45       while(m){
    46           if(check(m)){
    47             printf("%d\n",m);
    48             return 0;
    49         }
    50           m++;
    51       }
    52       return 0;
    53 }

**题目集应用：**[http://www.cnblogs.com/yzxverygood/category/1232095.html][http_www.cnblogs.com_yzxverygood_category_1232095.html]

**第二题：hdu 1576** [https://cn.vjudge.net/problem/HDU-1576][https_cn.vjudge.net_problem_HDU-1576]

题目思路：求（A/B)%9973，已知条件有：n=A%9973,gcd(B,9973)=1,A%B=0;输入n,B，求（A/B)%9973 的值？

设（A/B)%9973=x,即（A/B)=9973h+x;  A=9973\*h\*B+x\*B ;  n=A%9973即n=(9973\*h\*B+x\*B)%9973，即n=(x\*B)%9973，即x\*B-9973y=n;求最小的正数x

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

**第三题：UVA12169 **  [https://cn.vjudge.net/problem/UVA-12169][https_cn.vjudge.net_problem_UVA-12169]

题目大意：有3个整数 x\[1\], a, b 满足递推式x\[i\]=(a\*x\[i-1\]+b)mod 10001。由这个递推式计算出了长度为2T的数列，现在要求输入x\[1\],x\[3\],......x\[2T-1\], 输出x\[2\],x\[4\]......x\[2T\]. T<=100,0<=x<=10000. 如果有多种可能的输出，任意输出一个结果即可。

题目思路：第一种直接暴力方法，枚举a,b，找到合适的a,b（即满足式子，又因为我们只知道x\[1\],x\[3\]..，所以只要它满足所有的x\[i\]=(a\*((a\*x\[i-2\]+b)mod 10001)+b)mod 10001）然后输出相对于应的x\[2\],x\[4\]...

ac代码：

![ContractedBlock.gif][] ![ExpandedBlockStart.gif][]

1 #include<iostream>
     2 using namespace std;
     3 #define maxn 100005
     4 int a[maxn];
     5 int main()
     6 {
     7     int n;
     8     scanf("%d",&n);
     9     for(int i=0;i<n;i++)
    10     scanf("%d",&a[i]);
    11     int flag=0;
    12     for(int i=1;i<=10000;i++)
    13     {
    14     for(int j=1;j<=10000;j++)
    15     {
    16         flag=0;//i,j为a,b   筛选满足式子的a,b 
    17         for(int h=1;h<n;h++)
    18         {
    19             if(a[h]!=((i*((i*a[h-1]+j)%10001)+j)%10001))
    20             {
    21                 flag=1;
    22                 break;
    23             }
    24         }
    25         if(flag==0)
    26         {
    27             for(int h=0;h<n;h++)
    28             {
    29                 printf("%d\n",(i*a[h]+j)%10001);
    30             }
    31             break;
    32         }
    33     }
    34     if(flag==0)
    35     break;
    36     }
    37     return 0;
    38 }

暴力法

第二种方法，比第一种方法快，我们枚举a,利用简化式子(a+1)\*b mod 10001 = (x\[3\]-a\*a\*x\[1\]) mod 10001。这里就变成了同模方程，扩展欧几里得即可解答,求b。然后由a,b求相对于应的x\[2\],x\[4\]...

相应博客：[https://www.cnblogs.com/pach/p/6057160.html][https_www.cnblogs.com_pach_p_6057160.html]

*  ZOJ3609
 *  ZOJ3593
 *  POJ1061
 *  HDU1576
 *  HDU2669
 *  **UVA12169**

1.  *同余定理的公式：(a+b)%mod=(a%mod+b%mod)%mod* 
2.   
    
     *(a\*b)%mod=(a%mod\*b%mod)%mod*
3.   
    
    *扩展欧几里得也有自己的一个公式：a\*x+b\*y=gcd(a,b)*

转载于:https://www.cnblogs.com/Aiahtwo/p/10894280.html

[http_www.cnblogs.com_wkfvawl_p_9350867.html]: http://www.cnblogs.com/wkfvawl/p/9350867.html
[https_cn.vjudge.net_problem_HDU-2669]: https://cn.vjudge.net/problem/HDU-2669
[ContractedBlock.gif]: https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif
[ExpandedBlockStart.gif]: /images/20211214/e9d575f0c298469b8c99abe212a21b99.png
[https_cn.vjudge.net_problem_ZOJ-3593]: https://cn.vjudge.net/problem/ZOJ-3593
[https_www.luogu.org_problemnew_show_P2421]: https://www.luogu.org/problemnew/show/P2421
[http_www.cnblogs.com_yzxverygood_category_1232095.html]: http://www.cnblogs.com/yzxverygood/category/1232095.html
[https_cn.vjudge.net_problem_HDU-1576]: https://cn.vjudge.net/problem/HDU-1576
[https_cn.vjudge.net_problem_UVA-12169]: https://cn.vjudge.net/problem/UVA-12169
[https_www.cnblogs.com_pach_p_6057160.html]: https://www.cnblogs.com/pach/p/6057160.html