扩展欧几里得算法

扩展欧几里得算法- -
                                       

因为要用这个算法来计算逆元，所以写一下
扩展欧几里得算法：

欧几里得算法中，计算 x, y 的最大公约数的方法是辗转相除，例如：

gcd (26, 15)

26 % 15 = 1  11
15 % 11 = 1  4
11 % 4 = 2  3
4 % 3 = 1  1
3 % 1 = 3  0 

可知，gcd (26, 15) = 1

如果 gcd(x, y) = r，那么有 ax + by = r，可以看出，上面的步骤实际上是可以直接得出 a, b 的：
null
26 % 15 = 1  11 => 11 = 26 - 15 1 1 -1
15 % 11 = 1  4 => 4 = 15 - 11 = 15 - (26 - 15) = -26 + 2*15 1 -1 2
11 % 4 = 2  3 => 3 = 11 - 4*2 = (26 - 15) - (-26 + 15) * 2 = 3*26 - 5*15 2 3 -5
4 % 3 = 1  1 => 1 = 4 - 3 = (-26 + 2*15) - (3*26 - 5*15) = -4*26 + 7*15 1 -4 7
3 % 1 = 3  0 

在每一轮，我们都可以得到一个模的表达式为：ri = aix + biy
如果不考虑第一轮和第二轮，那么ai 和 bi 可以表示为(qi 为每一轮得到的商)：

ai = ai-2 - qi * ai-1
bi = bi-2 - qi * bi-1

现在来考虑第一轮和第二轮，按照上面的公式，可以认为

a-1 = 1, b-1 = 0
a0 = 0, b0 = 1

有了这两对预设值，上面的两个公式就成立了

求逆元，由上面可以看出，gcd(x, y) = 1 的时候
如果 ax + by = 1，那么 ax = -by + 1 => ax = 1 (mod b)
这时候，a 即是 x 的逆元


1 #include 
2
3  int gcd (int x, int y, int *a1, int *a2, int *b1, int *b2)
4  {
5    int q, r, a, b;
6
7    q = x / y;
8    r = x % y;
9
10    a = *a2 - q*(*a1);
11    b = *b2 - q*(*b1);
12
13    if (0 == r)
14      {
15        return y;
16      }
17
18    *a2 = *a1;
19    *b2 = *b1;
20
21    *a1 = a;
22    *b1 = b;
23
24    return gcd (y, r, a1, a2, b1, b2);
25  }
26
27  int main (int argc, char *argv[])
28  {
29    int x, y, a1, a2, b1, b2, r;
30
31    x = atoi (argv[1]);
32    y = atoi (argv[2]);
33
34    a1 = 0;
35    a2 = 1;
36    b1 = 1;
37    b2 = 0;
38
39    r = gcd (x, y, &a1, &a2, &b1, &b2);
40
41    printf ("%d = (%d) * (%d) + (%d) * (%d)\n",
42            r, a1, x, b1, y);
43
44    return 0;
45  }