欧几里得&拓展欧几里得

Posted on | In , | |
Words count in article: 2.2k | Reading time ≈ 9

欧几里德算法求最大公约数：

欧几里得算法，又名辗转相除法

1
2
举例：
辗转相除法基于如下原理：两个整数的最大公约数等于其中较小的数和两数的差的最大公约数。例如，252和105的最大公约数是21（252 = 21 × 12；105 = 21 × 5）；因为252 − 105 = 21 × (12 − 5) = 147，所以147和105的最大公约数也是21。在这个过程中，较大的数缩小了，所以继续进行同样的计算可以不断缩小这两个数直至其中一个变成零。这时，所剩下的还没有变成零的数就是两数的最大公约数。

递归公式表示位：

gcd(a,b) = gcd(b,a mod b)

两个数的最大公约数通常写成gcd(a, b)，如果有gcd(a, b)==1，则有a，b互质。

环论定义

在数学中，尤其是高等数学的环论中，最大公约数有一个更加巧妙的定义：a和b的最大公约数g是a和b的线性和中（绝对值）最小的一个，即所有形如ua + vb（其中u和v是整数,可正可负）的数中（绝对值）最小的数。所有ua + vb都是g的整数倍（ua + vb = mg，其中m是整数）。

递归

1
2
3
4
5
6
int gcd(int a, int b)
{
    if(b == 0)
        return a;
    return gcd(b, a%b);
}

迭代

1
2
3
4
5
6
7
8
9
10
int gcd(int a, int b)
{
    while(b != 0)
    {
        int r = b;
        b = a % b;
        a = r;
    }
    return a;
}

拓展欧几里得：

扩展欧几里德算法是用来在已知a, b求解一组p，q使得pa+qb=Gcd(a,b)。

扩展欧几里德常用在求解模线性方程及方程组中。

算法：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
//求解x*a+y*b=gcd(a,b) 中的x, y, gcd(a,b)即exgcd的return 
int exgcd(int a, int b, int &x, int &y)
{
    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的递归实现相比，发现多了下面的x,y赋值过程，这就是扩展欧几里德算法的精髓。 可以这样思考: 对于a’ =b , b’ =a%b 而言，我们求得x, y使得a’ x+b’ y=Gcd(a’, b’) 由于b’ = a % b = a - a / b * b 那么可以得到:

1
2
3
a'x + b'y = gcd(a', b')
bx + (a-a/b*b)y = gcd(a', b') = gcd(a, b)
ay + b(x-a/b*y) = gcd(a, b)

Stein算法

Stein算法由J. Stein 1961年提出，这个方法也是计算两个数的最大公约数。和欧几里德算法 算法不同的是，Stein算法只有整数的移位和加减法，因此对于大素数Stein将更有优势。

描述：

1
2
3
4
5
6
7
8
1. 如果A=0，B是最大公约数，算法结束 
2. 如果B=0，A是最大公约数，算法结束 
3. 设置A1=A、B1=B和C1=1 
4. 如果An和Bn都是偶数，则An+1=An>>1，Bn+1=Bn>>1，Cn+1=Cn<<1 
5. 如果An是偶数，Bn不是偶数，则An+1=An>>1，Bn+1=Bn，Cn+1=Cn(很显然，2不是奇数的约数) 
6. 如果Bn是偶数，An不是偶数，则Bn+1=B>>1，An+1=An，Cn+1=Cn(很显然，2不是奇数的约数) 
7. 如果An和Bn都不是偶数，则An+1=|An-Bn|>>1，Bn+1=min(An,Bn)，Cn+1=Cn //!!!!!
8. n加1，转1

第七条理由：

理由：设在某一次过程中A、B的最大公因子为D,可得A=kD,B=hD;(A>B -> k>h); 所以A-B=|k-h|*D。A,B都为奇数，由奇偶数乘积奇偶性可得D,k,h都为奇数。所以|k-h|一定是偶数，与h互质，所以A-B与B的最大公因数依旧为D。

（偶数乘奇偶数都是偶数，只有奇数相乘才奇数）A为奇数->k和D为奇数，B为奇数->h和D为奇数，两个奇数相减得偶数，所以A-B为偶数，k-h为偶数。

算法：

1
2
3
4
5
6
7
8
9
int gcd(int a, int b)
{
    if(a == 0) return b;
    if(b == 0) return a;
    if(a % 2 == 0 && b % 2 == 0) return 2*gcd(a>>1, b>>1);
    else if(a % 2 == 0) return gcd(a>>1, b);
    else if(b % 2 == 0) return gcd(a, b>>1);
    else return gcd(abs(a-b), min(a, b)); //!!!!!
}

严格意义上来说，负数和正数之间是没有公因子这一说的【因子要大于0】，但是如果非要求他们的正因子【负因子】，可以先把两个数取绝对值，stein之后再加个符号即可。

应用：

求解不定方程

对于不定整数方程pa+qb=c，若 c mod Gcd(p, q)=0,则该方程存在整数解，否则不存在整数解。 上面已经列出找一个整数解的方法，在找到p a+q b = Gcd(p, q)的一组解p0,q0后，p a+q b = Gcd(p, q)的其他整数解满足：

  1. p = p0 + b/Gcd(p, q) * t
  2. q = q0 - a/Gcd(p, q) * t(其中t为任意整数)

至于pa+qb=c的整数解，只需将p a+q b = Gcd(p, q)的每个解乘上 c/Gcd(p, q) 即可。 在找到p a+q b = Gcd(a, b)的一组解p0,q0后，应该是得到p a+q b = c的一组解p1 = p0(c/Gcd(a,b)),q1 = q0(c/Gcd(a,b))，p a+q b = c的其他整数解满足：

  1. p = p1 + b/Gcd(a, b) * t
  2. q = q1 - a/Gcd(a, b) * t(其中t为任意整数)

此处的p 、q就是p a+q b = c的所有整数解。

1
2
3
4
5
6
7
8
9
bool linear_equation(int a,int b,int c,int &x,int &y)
{
    int d=exgcd(a,b,x,y);
    if(c%d)
        return false;
    int k=c/d;
    x*=k; y*=k;    //求得的只是其中一组解
    return true;
}

求解模线性方程（线性同余方程）

同余方程 ax≡b (mod n)对于未知数 x 有解，当且仅当 gcd(a,n) | b。（原因见下一条

且方程有解时，方程有 gcd(a,n) 个解。 求解方程 ax≡b (mod n) 相当于求解方程 ax+ ny= b, (x, y为整数)。

设 d= gcd(a,n)，假如整数 x 和 y，满足 d= ax+ ny(用扩展欧几里德得出)。

如果 d| b，则方程a x0+ n y0= d， 方程两边乘以 b/ d，(因为 d|b，所以能够整除)，得到 a x0 b/ d+ n y0 b/ d= b。

所以 x= x0 b/ d，y= y0 b/ d 为 ax+ ny= b 的一个解，所以 x= x0* b/ d 为 ax≡ b (mod n ) 的解。

ax≡b (mod n)的一个解为 x0= x (b/ d ) mod n，且方程的 d 个解分别为 xi= (x0+ i (n/ d ))mod n {i= 0… d-1}。

设ans=x*(b/d),s=n/d; 方程ax≡b (mod n)的最小整数解为：(ans%s+s)%s;

1
2
3
4
5
6
7
8
9
10
11
12
相关证明：

    证明方程有一解是: x0 = x'(b/d) mod n;
    由 a*x0 = a*x'(b/d) (mod n)
         a*x0 = d (b/d) (mod n)   (由于 ax' = d (mod n))
                 = b (mod n)

    证明方程有d个解: xi = x0 + i*(n/d)  (mod n);
    由 a*xi (mod n) = a * (x0 + i*(n/d)) (mod n)
                             = (a*x0+a*i*(n/d)) (mod n)
                             = a * x0 (mod n)             (由于 d | a)
                             = b
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
首先看一个简单的例子：

5x=4(mod3)

解得x = 2,5,8,11,14.......

由此可以发现一个规律，就是解的间隔是3.

那么这个解的间隔是怎么决定的呢？

如果可以设法找到第一个解，并且求出解之间的间隔，那么就可以求出模的线性方程的解集了.

我们设解之间的间隔为dx.

那么有

a*x = 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是最小的.

设a 和 n的最大公约数为d,那么a 和 n 的最小公倍数为(a*n)/d.

即a*dx = a*n/d;

所以dx = n/d.

因此解之间的间隔就求出来了.
1
2
3
4
5
6
7
8
9
10
11
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=1;i<d;i++)
        printf("%d\n",(x0+i*(n/d))%n);
    return true;
}

用欧几里德算法求模的逆元：

同余方程ax≡b (mod n)，如果 gcd(a,n)== 1，则方程只有唯一解。

在这种情况下，如果 b== 1，同余方程就是 ax=1 (mod n ),gcd(a,n)= 1。

这时称求出的 x 为 a 的对模 n 乘法的逆元。

对于同余方程 ax= 1(mod n )， gcd(a,n)= 1 的求解就是求解方程

ax+ ny= 1，x, y 为整数。这个可用扩展欧几里德算法求出，原同余方程的唯一解就是用扩展欧几里德算法得出的 x 。

Donate? comment?