acm数论之欧几里得gcd

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

解题思路：

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

题目大意：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代码：

#include<iostream>
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)
}
int main()
{
    int a,b;
    while(~scanf("%d%d",&a,&b))
    {
        int x,y;
        int d=exgcd(a,b,x,y);
        if(1%d)//有余数 
        {
            cout<<"sorry"<<endl; 
            continue;
        }
        else
        {
            //通过exgcd函数得到特解 
            x=x/d,y=y/d;
            if(x>0)//由于a,b为正数，所以k=0取最小 
            cout<<x<<" "<<y<<endl;
            else
            {//当x0为负数的时候，x最小就是0，与横轴相交, 但x可能是负数，如果为负数，取下一个整数点
                int k=(-1)*x/b*d;
                int kx=x+b/d*k;
                int ky=y-k*(a/d);
                if(kx<0)
                {
                    kx=x+b/d*(k+1);
                    ky=y-(k+1)*(a/d);
                }
                cout<<kx<<" "<<ky<<endl;                
            }
        }    
    }
    return 0; 
 } 

 

例题2：

题目链接：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)

#include<iostream>
#include<cstdio>
#include<cmath>
using namespace std;
const long long maxn = 0x3f3f3f3f;
//#define maxn 0x3f3f3f3f
long long exgcd(long long a,long long b,long long &x,long long &y)
{
    if(b==0)
    {
        x=1,y=0;
        return a;
    }
    long long r=exgcd(b,a%b,x,y);
    long long t=x;x=y;y=t-a/b*y;
    return r;
}
long long Abs(long long a)
{
    if(a<0)
    a=-1*a;
    return a;
}
int main()
{
    long long a,b,A,B;
    long long ans;
    int T;
    cin>>T;
    while(T--)
    {
        cin>>A>>B>>a>>b;
        long long x,y;
        long long h=exgcd(a,b,x,y);
        if((B-A)%h!=0)
            printf("-1\n");
        else{
            x=(B-A)/h*x;//注释修改
            y=(B-A)/h*y;//注释修改
            long long minx,miny;
            long long k;
            b=b/h;
            a=a/h;
            k=(y-x)/(b+a);//注释修改
            ans=maxn*maxn;   //注释修改  这里定义成maxn数据定义小了
            for(int i=k-1;i<=k+1;i++)
            {
                minx=x+b*i,miny=y-a*i;
                long long m;
                if(Abs(minx+miny)==Abs(minx)+Abs(miny))//说明x与y同号
                    m=max(Abs(minx),Abs(miny));
                else
                    m=Abs(minx)+Abs(miny);
                if(m<ans)
                 ans=m;
//                if((minx>0&&miny>0)||(minx<0&&miny<0))//x与y同号
//                {
//                m=Abs(min(minx,miny))+Abs(minx-miny);//这句有问题的  测试这个存在数据Abs(min(minx,miny))+Abs(minx-miny)！=max(Abs(minx),Abs(miny));
//                }
//                 else
//                 m=Abs(minx)+Abs(miny);
//                 if(m<ans)
//                 ans=m;

            }
            cout<<ans<<endl;
        }
    }
}

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这里面取得答案，因为多了一个参数。

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

b=b;
a=a;
k=(y-x)/(b+a)*h;
for(int i=k-1;i<=k+1;i++)
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

#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 T;
    cin>>T;
    while(T--){
    int a,m;
    cin>>a>>m;
    if(inversenum(a,m)==-1)
    cout<<"Not Exist"<<endl;
    else
    cout<<inversenum(a,m)<<endl;
    }    
 } 

 

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

#include<iostream>
#include<cmath>
using namespace std;
//扩展欧几里得处理ax+by=gcd(a,b)=d解的问题 
long long t;
long long exgcd(long long a,long long b,long long &x,long long &y)  //得到gcd(a,b)的值 
{
    if(b==0)
    {
        x=1,y=0;
        return a;
     } 
     long long r=exgcd(b,a%b,x,y);
     long long t=x; x=y; y=t-a/b*y;
     return r;//gcd(a,b)
}
long long solve(long long m,long long n)
{
    long long x,y;
    long long r=exgcd(m,n,x,y);
    if(t%r)
    {
        return -1;
    }
    x=x*t/r;
    long long s=n/r;
    return (x%s+s)%s;//不定方程的最小解x 
}
int main()
{
    long long x,y,m,n,L;
    cin>>x>>y>>m>>n>>L;
    //不定式a*(m-n)+b*L=y-x;   看成a*m+b*n=t; 
    t=y-x;
    m=m-n;
    n=L;
    if(m<0)//整个不定式a,b取正数，由于b>0，b不需要变(b前面有系数),所以t需要改变 
    {
        m=-1*m;
        t=-1*t;
    }
    if(solve(m,n)==-1)
    cout<<"Impossible"<<endl;
    else
    {
        cout<<solve(m,n)<<endl;
    }
    return 0;
 } 

 

例题：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个野人之间没有住在同一个洞穴，即余数不同或者是相同余数数字大于两野人最小寿命。

题目代码：

#include<iostream>
#include<algorithm> 
using namespace std;
int c[200],p[200],l[200],n;
int gcd(int a,int b){return a%b==0?b:gcd(b,a%b);}
//扩展欧几里得处理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)
}
bool go(int i,int j,int b){
      int a=p[i]-p[j],f=c[j]-c[i],g,x,y;
      if(a<0)
      {
          a=a*-1,f=-1*f;
      }
      g=exgcd(a,b,x,y);
      if(f%g==0){
          x=x*f/g;
          b=b/g;
          x=(x%b+b)%b;//最小解
          if(x<=min(l[i],l[j]))
          return 0;
      }
      return 1;
}
bool check(int m){
      int i,j,k;
      for(i=1;i<n;i++)
        for(j=i+1;j<=n;j++)
          if(!go(i,j,m))return 0;
      return 1;
}
int main()
{     int m=0,i,j,k;
      scanf("%d",&n);
      for(i=1;i<=n;i++)scanf("%d%d%d",&c[i],&p[i],&l[i]),m=max(m,c[i]);
      while(m){
          if(check(m)){
            printf("%d\n",m);
            return 0;
        }
          m++;
      }
      return 0;
}

 

题目集应用：http://www.cnblogs.com/yzxverygood/category/1232095.html

第二题：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

#include<iostream>
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)
}
int main()
{
    int T,n,B;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d",&n,&B);
        int x,y;
        int d=exgcd(B,9973,x,y);
        x=x*n/d;//最小解 
        int s=9973/d;
        x=(x%s+s)%s;
        cout<<x<<endl;
    }
    return 0; 
 }

 

第三题：UVA12169   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代码：

#include<iostream>
using namespace std;
#define maxn 100005
int a[maxn];
int main()
{
    int n;
    scanf("%d",&n);
    for(int i=0;i<n;i++)
    scanf("%d",&a[i]);
    int flag=0;
    for(int i=1;i<=10000;i++)
    {
    for(int j=1;j<=10000;j++)
    {
        flag=0;//i,j为a,b   筛选满足式子的a,b 
        for(int h=1;h<n;h++)
        {
            if(a[h]!=((i*((i*a[h-1]+j)%10001)+j)%10001))
            {
                flag=1;
                break;
            }
        }
        if(flag==0)
        {
            for(int h=0;h<n;h++)
            {
                printf("%d\n",(i*a[h]+j)%10001);
            }
            break;
        }
    }
    if(flag==0)
    break;
    }
    return 0;
}

第二种方法，比第一种方法快，我们枚举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

 

• 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)