码迷,mamicode.com
首页 > 其他好文 > 详细

扩展欧几里得算法------扩展欧几里德算法

时间:2014-08-08 18:14:06      阅读:341      评论:0      收藏:0      [点我收藏+]

标签:扩展欧几里得及其应用   扩展欧几里德及其应用   不定方程   同余方程   乘法逆元   

扩展欧几里得算法及其应用

一、扩展欧几里得算法

扩展欧几里得算法:对于不完全为 0 的非负整数 a,b,若gcd(a,b)表示 a,b 的最大公约数,必然存在整数对x,y ,使得 ax+by = gcd(a,b)。

算法过程:

设 a>b,当 b=0时,gcd(a,b)=a。此时满足ax+by = gcd(a,b)的一组整数解为x=1,y=0;当a*b!=0 时,

设 a*x1+b*y1=gcd(a,b);b*x2+(a mod b)*y2=gcd(b,a mod b);

根据欧几里得原理知 gcd(a,b)=gcd(b,a mod b);则:a*x1+b*y1=b*x2+(a mod b)*y2;

即:a*x1+b*y1=b*x2+(a-(a/b)*b)*y2=a*y2+b*x2-(a/b)*b*y2;(a/b:表示a除以b的商)

根据恒等定理得:x1=y2; y1=x2-(a/b)*y2;

按照上面的递归思想,我们不断对gcd进行递归,一定会出现b=0的情况,此时x=1,y=2,递归结束。这样我们就得到了求解一组x1,y1 的方法

将上面的过程写成代码如下:函数exgcd的代码:C/C++代码

可运行代码:

#include<stdio.h>
#include <iostream>
using namespace std;
/*x,y,r定义为全局变量,x,y的结果就是
方程ax+by=gcd(a,b)的一组解*/
int x,y,r;
//递归实现求一组x,y;函数exgcd无返还值
void exgcd(int a,int b)
{
	if(b==0)
	{
		x=1;
		y=0;
		r=a;
	}
	else
	{
		exgcd(b,a%b);
		int temp=x;
		x=y;
		y=temp-(a/b)*y;
	}
}
int main()
{
	int a,b;
	//cin>>a>>b;
	scanf("%d%d",&a,&b);
    exgcd(a,b);
    //cout<<q<<"=("<<x<<")*"<<a<<"+("<<y<<")*"<<b<<endl;
    printf("%d=(%d)*%d+(%d)*%d\n",r,x,a,y,b);
	return 0;
}

注意上面的代码是将x,y,q定义成了全局变量,其中x,y表示方程的一组整数解,q表示a,b的最大公约数。

核心代码:

int x,y,r;
void exgcd(int a,int b)
{
	if(b==0)
	{
		x=1;
		y=0;
		r=a;
	}
	else
	{
		exgcd(b,a%b);
		int temp=x;
		x=y;
		y=temp-(a/b)*y;
	}
}


可运行代码一:递归算法C/C++

#include<stdio.h>
#include <iostream>
using namespace std;
/*递归,a,b代表方程ax+by=gcd(a,b)的a,b
x,y代表方程的一组解,其中x,y并未定义为
全局变量,而是将x,y的地址传过去,exgcd
的返回值时a,b的最大公约数*/
int exgcd(int a,int b,int &x,int &y)
{
    int r,t;
    if(b==0)
    {
        x=1;
        y=0;
        return a;
    }
    r=exgcd(b,a%b,x,y);
    t=x;
    x=y;
    y=t-a/b*y;
    return r;
}
int main()
{
	int m,n,x,y,r;
	scanf("%d%d",&m,&n);
    r=exgcd(m,n,x,y);
    printf("(%d)*%d+(%d)*%d=%d\n",x,m,y,n,r);
	return 0;
}

核心代码:

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

可运行代码二:非递归C/C++

#include<stdio.h>
#include <iostream>
using namespace std;
/*非递归,a,b代表方程ax+by=gcd(a,b)的a,b
x,y代表方程的一组解*/
int exgcd(int a,int b,int &x,int &y)
{
    int x1,y1,x0,y0,r,q;
    x0=1; y0=0;
    x1=0; y1=1;
    x=0; y=1;
    r=a%b;
    q=(a-r)/b;
    while(r)
    {
        x=x0-q*x1; y=y0-q*y1;
        x0=x1; y0=y1;
        x1=x; y1=y;
        a=b; b=r; r=a%b;
        q=(a-r)/b;
    }
    return b;
}
int main()
{
	int m,n,x,y,r;
	scanf("%d%d",&m,&n);
    r=exgcd(m,n,x,y);
    printf("(%d)*%d+(%d)*%d=%d\n",x,m,y,n,r);
	return 0;
}

核心代码:

int exgcd(int a,int b,int &x,int &y)
{
    int x1,y1,x0,y0,r,q;
    x0=1; y0=0;
    x1=0; y1=1;
    x=0; y=1;
    r=a%b;
    q=(a-r)/b;
    while(r)
    {
        x=x0-q*x1; y=y0-q*y1;
        x0=x1; y0=y1;
        x1=x; y1=y;
        a=b; b=r; r=a%b;
        q=(a-r)/b;
    }
    return b;
}

对比①②的两种算法程序,我感觉②较好一些。

方程的多组解:

当求出一组解x,y之后,我们可以求出所有的解,其它的整数解满足

X=x+(b/gcd(a,b))*t  (t为任意整数)      Y=y-(a/gcd(a,b))*t  (t为任意整数)

二、扩展欧几里得应用解不定方程、解线性同余方程、解模的逆元

1.扩展欧几里得解不定方程

不定方程:形如px+qy=c的方程称为不定方程。

使用扩展欧几里得算法解不定方程的办法:

        对于不定整数方程px+qy=c,若 c mod gcd(p,q)=0,则该方程存在整数解,否则不存在整数解。

       上面的欧几里得算法已经详细给出求方程ax+by= gcd(a,b)一个整数解的方法,在找到方程ax+by = gcd(a,b)的一组解x1,y1后,方程ax+by = gcd(a,b)的其他整数解满足:

  x = x1+(b/gcd(a,b))*t  (t为任意整数)      y=y1-(a/gcd(a,b))*t    (t为任意整数)

那么px+qy=c的整数解,只需将ax+by = gcd(a,b)的每个解乘上 c/ gcd(a,b) 即可。

代码如下:C/C++代码

可运行代码:

#include<stdio.h>
#include <iostream>
using namespace std;
/*递归,a,b代表方程ax+by=gcd(a,b)的a,b;
x,y代表方程的一组整数解*/
int exgcd(int a,int b,int &x,int &y)
{
    int r,t;
    if(b==0)
    {
        x=1;
        y=0;
        return a;
    }
    r=exgcd(b,a%b,x,y);
    t=x;
    x=y;
    y=t-a/b*y;
    return r;
}
//扩展欧几里得解不定方程px+qy=c的一组整数解x,y
bool linear_equation(int a,int b,int c,int &x,int &y)
{
    int d=exgcd(a,b,x,y);//a,b的最大公约数
    if(c%d)
        return false;
    else
    {
        int k=c/d;//k为两方程ax+by=gcd(a,b)与px+qy=c相差的系数
        x*=k; y*=k;//求得的只是其中一组解
        return true;
    }
}
int main()
{
	int i,p,q,c,x=0,y=0,r;
	scanf("%d%d%d",&p,&q,&c);
    if(linear_equation(p,q,c,x,y))
        printf("x=%d y=%d\n",x,y);
	return 0;
}

核心代码:

int exgcd(int a,int b,int &x,int &y)
{
    int r,t;
    if(b==0)
    {
        x=1;
        y=0;
        return a;
    }
    r=exgcd(b,a%b,x,y);
    t=x;
    x=y;
    y=t-a/b*y;
    return r;
}
bool linear_equation(int a,int b,int c,int &x,int &y)
{
    int d=exgcd(a,b,x,y);//a,b的最大公约数
    if(c%d)
        return false;
    else
    {
        int k=c/d;//k为两方程ax+by=gcd(a,b)与px+qy=c相差的系数
        x*=k; y*=k;//求得的只是其中一组解
        return true;
    }
}

当求出一组解x,y之后,我们可以求出所有的解,其它的整数解满足

X=x+(b/gcd(a,b))*t  (t为任意整数)    Y=y-(a/gcd(a,b))*t  (t为任意整数)


2.扩展欧几里得算法解线性同余方程

线性同余方程:形如ax≡b (mod n)的方程称为线性同余方程。

性质:当且仅当 b 能够被 a 与 n 的最大公约数d整除(记作 gcd(a,n) | b)时,方程有解,并且方程在[0,n-1]之间恰有gcd(a,n)个解(注意不是方程只有gcd(a,n)个解)。

例如:

①在方程6x ≡ 7 (mod 8)中,d=gcd(6,8)=2,7不能整除2,因此方程无解。

②在方程6x ≡ 10 (mod 8)中,d=gcd(6,8)=2,10能整除2,因此方程有解,并且方程在[0,7]上有2个解,

    分别为x=3,x=7。

③在方程21x ≡ 7 (mod 3)中,d=gcd(21,3)=3,7不能整除3,因此方程无解。

④在方程10x ≡ 5 (mod 5)中,d=gcd(10,5)=5,5能整除5,因此方程有解,并且方程在[0,4]上有5个解,

    分别为x=0,x=1,x=2,x=3,x=4。

⑤在方程30x ≡ 10 (mod 15)中,d=gcd(30,15)=15,10不能整除15,因此方程无解。

⑥在方程15x ≡ 40 (mod 10)中,d=gcd(15,10)=5,40能整除5,因此方程有解,并且方程在[0,14]上有5个解,

    分别为x=0,x=2,x=4,x=6,x=8。

求解:求解线性同余方程ax≡b (mod n)相当于求解方程ax+ny= b, (x, y均为整数)。根据线性同余方程的性质可知,当gcd(a,n) | b时,方程有解,在这里我们讨论有解的情况。

设d= gcd(a,n),则方程ax+ny=d的解可以用扩展欧几里得算法求出一组x0,y0,则有a*x0+n*y0=d,方程两边同时乘以b/d(因为d能整除b)就得到a*x0*b/d+n*y0*b/d=d*b/d=b,所以x=x0*b/d,y=y0*b/d是方程ax+ny=b的一组整数解,所以x=x0*b/d是线性同余方程ax≡b (mod n)的一个解。(设ans=x*(b/d),s=n/d,则方程ax≡b (mod n)的最小整数解为(ans%s+s)%s)

方程ax≡b (mod n)的一个解为 x1= x* (b/ d ) mod n,且方程的 d 个解分别为 Xi= (x1+ i* (n/ d ))mod n {i= 0... d-1}

求解d个解的代码如下:

输入:a b n

输出:No answer!或d个解

#include<stdio.h>
#include<iostream>
using namespace std;
int exgcd(int a,int b,int &x,int &y)
{
    int r,t;
    if(b==0)
    {
        x=1;
        y=0;
        return a;
    }
    r=exgcd(b,a%b,x,y);
    t=x;
    x=y;
    y=t-a/b*y;
    return r;
}
bool modular_linear_equation(int a,int b,int n)
{
    int x,y,x0,i;
    int d=exgcd(a,n,x,y);
    if(b%d)
    {
        printf("No answer!\n");
        return false;
    }
    x0=(x*(b/d)%(n/d)+n/d)%(n/d);//最小整数解
    for(i=0;i<d;i++)
        printf("x=%d\n",(x0+i*(n/d))%n);
    return true;
}
int main()
{
    int a,b,n,x,y;
    scanf("%d%d%d",&a,&b,&n);
    modular_linear_equation(a,b,n);
    return 0;
}

线性同余方程若有解,则有多个解,那么解与解之间有没有关系呢?有什么关系呢?

对上面性质中的举例,我们可以些许发现,当方程有解时,它们的解的间隔是相同的。那么这个间隔与方程的a,b,n有怎样的关系呢?

证明:我们设解之间的间隔为dx,因为ax≡b(mod n)所以a(x+dx)≡b (mod n),两式相减得a*dx (mod n)=0,那么a*dx即是a的倍数也是n的倍数,则a*dx是a,n的公倍数,要求出dx的最小值,我们只需求出a,n的最小公倍数,则此时对应的dx最小。

设d是a和n的最大公约数,则a和n的最小公倍数为(a*n)/d。所以a*dx=a*n/d,则dx=n/d。

求线性同余方程的解的代码如下:

#include<stdio.h>
#include<iostream>
using namespace std;
int exgcd(int a,int b,int &x,int &y)
{
    int r,t;
    if(b==0)
    {
        x=1;
        y=0;
        return a;
    }
    r=exgcd(b,a%b,x,y);
    t=x;
    x=y;
    y=t-a/b*y;
    return r;
}
bool modular_linear_equation(int a,int b,int n)
{
    int x,y,x0,i;
    int d=exgcd(a,n,x,y);
    if(b%d)
        return false;
    x0=x*(b/d)%n;//特解
    for(i=0;i<100;i++)//控制输出多少解
        printf("%d\n",x0+i*(n/d));
    return true;
}
int main()
{
    int a,b,n,x,y;
    scanf("%d%d%d",&a,&b,&n);
    modular_linear_equation(a,b,n);
    return 0;
}

核心代码:

bool modular_linear_equation(int a,int b,int n)
{
    int x,y,x0,i;
    int d=exgcd(a,n,x,y);
    if(b%d)
        return false;
    x0=x*(b/d)%n;//特解
    for(i=0;i<100;i++)//控制输出多少解
        printf("%d\n",x0+i*(n/d));
    return true;
}

3.求模的乘法逆元

乘法逆元:若ax≡1 mod f, 则称x为a关于模f的乘法逆元。也可表示为ax≡1(mod f)。

当a与f互素时,a关于模f的乘法逆元有唯一解(因为方程ax≡1 mod f在  内只有一个解)。如果不互素,则无解。如果f为素数,则从1到f-1的任意数都与f互素,即在1到f-1之间都恰好有一个关于模f的乘法逆元。

代码如下:

#include <stdio.h>
#include<iostream>
using namespace std;
/* 扩展欧几里得算法求乘法逆元*/
int ExtendedEuclid(int f,int d,int *result)
{
    int x1,x2,x3,y1,y2,y3,t1,t2,t3,q;
    x1=y2=1;
    x2=y1=0;
    x3=(f>=d)?f:d;
    y3=(f>=d)?d:f;
    while(1)
    {
	if(y3==0)
		{
		     *result=x3;//两个数不互素则result为两个数的最大公约数,此时返回值为零
		     return 0;
		}
		if(y3==1)
		{
		     *result = y2;//两个数互素则resutl为其乘法逆元,此时返回值为1
		     return 1;
		}
		q=x3/y3;
		t1=x1-q*y1;
		t2=x2-q*y2;
                t3=x3-q*y3;
		x1=y1;
		x2=y2;
		x3=y3;
		y1=t1;
		y2=t2;
		y3=t3;
	}
}
int main()
{
	int x,y,z;
	z = 0;
	printf("输入两个数:\n");
	scanf("%d%d",&x,&y);
	if(ExtendedEuclid(x,y,&z))
        printf("%d和%d互素,乘法的逆元是:%d\n",x,y,z);
	else
        printf("%d和%d不互素,最大公约数为:%d\n",x,y,z);
	return 0;
}

核心代码:

int ExtendedEuclid(int f,int d,int *result)
{
    int x1,x2,x3,y1,y2,y3,t1,t2,t3,q;
    x1=y2=1;
    x2=y1=0;
    x3=(f>=d)?f:d;
    y3=(f>=d)?d:f;
    while(1)
    {
		if(y3==0)
		{
			*result=x3;
			return 0;
		}
		if(y3==1)
		{
			*result = y2;
			return 1;
		}
		q=x3/y3;
		t1=x1-q*y1;
		t2=x2-q*y2;
		t3=x3-q*y3;
		x1=y1;
		x2=y2;
		x3=y3;
		y1=t1;
		y2=t2;
		y3=t3;
	}
}




















扩展欧几里得算法------扩展欧几里德算法,布布扣,bubuko.com

扩展欧几里得算法------扩展欧几里德算法

标签:扩展欧几里得及其应用   扩展欧几里德及其应用   不定方程   同余方程   乘法逆元   

原文地址:http://blog.csdn.net/yanghuaqings/article/details/38440789

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!