欧几里德辗转相除法是最大公约数（greatest common divisor)的求法。
C++代码如下:
  int gcd(int a, int b)
  {
     if(b == 0)
         return a;
     else
         return gcd(b, a%b);   
  }
这个算法就是利用了gcd(a, b) = gcd(b, a mod b)。
证明：
对于a, b的任意公约数r，则r|a, r|b。
那么a mod b = a - [a/b]*b，[a/b]是取整函数
由于r|a， r|b所以r|(a mod b)。
对于b, a mod b的公约数r，则r|b， r|(a mod b)
a = a mod b + [a/b]*b, 所与r|a。这样就证明了
a，b的公约数和b， a mod b完全相同，则最大公约数也相同。

费马小定理
(1)p为素数, 则p|(a^p - a)
(2)p为素数，gcd(a, p) = 1, 则 p|(a^(p-1)-1)

证明：
对于(2)，p是素数，a是整数且gcd(a,p)=1即他们的最大公约数是1。
由于a, 2a, 3a, ……，(p-1)a 模p的余数都不相同。
否则若(i*a) mod p=(j*a) mod p
其中 1 =< i < j <= p-1 则p|(j-i)*a, 而gcd(a,p)=1,
那么p|(j - i),这是不可能的，所以a, 2a, 3a, ……, (p-1)a对p的余数
取遍1，2，……p-1。
所以a*2a*3a*4a*……*(p-1)a对p的余数等于1*2*……(p-1)对p的余数。
这是由于性质：a mod p = x, b mod p = y, 那么(a*b) mod p = (x*y) mod p。
1*2*……(p-1)*a^(p-1) mod p = (p-1)! mod p
((a^(p-1) - 1)*(p-1)! + (p-1)!)mod p = (p-1)! mod p
则(a^(p-1) - 1)*(p-1)! mod p = 0
由于p是素数，gcd(p, (p-1)!)=1
所以a^(p-1)  - 1 mod p = 0
就是 p|(a^(p-1) - 1)。因此(2)得证。
对于(1) 若gcd(a,p) ≠ 1,那么由于p是素数必然p|a 因此p|a*(a^(p - 1) - 1)显然得证。
如果gcd(a,p) = 1, 则由(2)得p|(a^(p-1)-1)因此因此p|(a^(p - 1) - 1)*a显然得证。

欧拉函数是指：对于一个正整数n，小于n且和n互质的正整数的个数，记做：φ(n)
欧拉定理：gcd(a, n) = 1，则a^φ(n) ≡ 1 mod n
对于互质的整数a和n，有a^φ(n) ≡ 1 mod n
证明：
    首先证明下面这个命题：
    对于集合Z = {X1,X2,...,Xφ(n)}，Xi为第i个小于n且和n互质的正整数
    考虑集合
    S = {aX1 mod n,aX2 mod n,...,aXφ(n) mod n}
    则S = Z
    证明:
    a, n互质，Xi与n互质，那么aXi与n互质
    对于Z中两个元素Xi和Xj，如果Xi ≠ Xj
    则aXi mod n ≠ aXj mod n，若aXi mod n = aXj mod n，
    那么a(Xi - Xj) mod n  = 0, a, n互质那么Xi - Xj mod n = 0 由于
     0 < Xi, Xj < n显然不可能。因此S中的元素各不相同，这样S中含有φ(n)个与n互
    质的正整数，且满足每个元素小于n所以，很明显，S = Z        
    既然这样，那么
    (aX1 × aX2×...×aXφ(n))mod n
     = (aX1 mod n × aX2 mod n × ... × aXφ(n) mod n)mod n
     = (x1 × x2 × ... × xφ(n))mod n
     考虑上面等式左边和右边
     最左边等于(a^φ(n) × (x1 × x2 × ... × xφ(n)) mod n
     最右边等于(x1 × x2 × ... × xφ(n))mod n
     得到(a^φ(n) - 1)(x1 × x2 × ... × xφ(n)) mod n = 0
     而(x1 × x2 × ... × xφ(n))mod n 和 n互质,所以 a^φ(n) ≡ 1 mod n
推论：当n为素数时，φ(n) = n - 1,  a^(n - 1) - 1 mod n = 0 就是费马小定理(2)

模p乘法逆元
对于整数a、p，如果存在整数b，满足a*b mod p =1，则说，b是a的模p乘法逆元。
定理：a存在模p的乘法逆元的充要条件是gcd(a,p) = 1
证明：
首先证明充分性
如果gcd(a,p) = 1，根据欧拉定理，a^φ(p) ≡ 1 mod p，因此
显然a^(φ(p)-1) 是a的模p乘法逆元。
再证明必要性
假设存在a模p的乘法逆元为b， a*b ≡ 1 mod p
则a*b = k*p + 1，所以1 = a*b - k*p
因为gcd(a,p) = d
所以d | 1
所以d只能为1

扩展欧几里德定理
如果gcd(a, b) = d, 那么一定存在整数x，y 满足a*x + b*y = d

证明，可以设a = p*d， b = q*d，显然gcd(p, q) = 1
a*x + b*y = p*d*x + q*d*y  若要a*x + b*y = d
必需p*x + q*y = 1，由gcd(p, q) = 1 则p^φ(q) ≡ 1 mod q
令x = p^(φ(q) - 1) 也就是x是p的模q逆元，则p*x - 1 = p^φ(q) - 1 mod q
=  0， 那么y = (p^φ(q) - 1)/q显然是整数， 当然也可以y = q^(φ(p) - 1)，那
么x = (q^φ(p) - 1)/p，注意到x, y有无数组解。

注意到对于gcd(a,b) = d 我们对(a, b)用欧几里德辗转相除会最终得到
(d, 0)此时对于把a =d, b = 0 带入a*x + b*y = d，显然x = 1，y可以为任意值，
这里y可以为任意值就意味着解会有无数个。我们可以用a = d, b = 0的情况逆推出来
任何gcd(a, b) = d 满足a*x + b*y = d的解。如果x0, y0是b*x + (a%b)*y = d 的解，那么对于a*x + b*y = d的解呢？
b*x + (a%b)*y = d => b*x + (a - [a/b]*b)*y = a*y + b*(x - [a/b]*y)，所以a*x + b*y = d的解x1 = y0， y1 = x0 - [a/b]*y0; 这样我们可以程序迭带了。

void extended_gcd(int a, int b, int &x, int &y)
{
  
   if(b == 0)
      {
        x = 1;
        y = 0; //any integer is ok;
       }
       else
       {
         int t;
         extended_gcd(b, a%b, x, y);
         t = x;
         x = y;
         y = t - (a/b)*y;
       }  
}

终于写完了，参考了网上书上的内容。

欧几里德辗转相除法是最大公约数（greatest common divisor)的求法。
C++代码如下:
  int gcd(int a, int b)
  {
     if(b == 0)
         return a;
     else
         return gcd(b, a%b);   
  }
这个算法就是利用了gcd(a, b) = gcd(b, a mod b)。
证明：
对于a, b的任意公约数r，则r|a, r|b。
那么a mod b = a - [a/b]*b，[a/b]是取整函数
由于r|a， r|b所以r|(a mod b)。
对于b, a mod b的公约数r，则r|b， r|(a mod b)
a = a mod b + [a/b]*b, 所与r|a。这样就证明了
a，b的公约数和b， a mod b完全相同，则最大公约数也相同。

费马小定理
(1)p为素数, 则p|(a^p - a)
(2)p为素数，gcd(a, p) = 1, 则 p|(a^(p-1)-1)

证明：
对于(2)，p是素数，a是整数且gcd(a,p)=1即他们的最大公约数是1。
由于a, 2a, 3a, ……，(p-1)a 模p的余数都不相同。
否则若(i*a) mod p=(j*a) mod p
其中 1 =< i < j <= p-1 则p|(j-i)*a, 而gcd(a,p)=1,
那么p|(j - i),这是不可能的，所以a, 2a, 3a, ……, (p-1)a对p的余数
取遍1，2，……p-1。
所以a*2a*3a*4a*……*(p-1)a对p的余数等于1*2*……(p-1)对p的余数。
这是由于性质：a mod p = x, b mod p = y, 那么(a*b) mod p = (x*y) mod p。
1*2*……(p-1)*a^(p-1) mod p = (p-1)! mod p
((a^(p-1) - 1)*(p-1)! + (p-1)!)mod p = (p-1)! mod p
则(a^(p-1) - 1)*(p-1)! mod p = 0
由于p是素数，gcd(p, (p-1)!)=1
所以a^(p-1)  - 1 mod p = 0
就是 p|(a^(p-1) - 1)。因此(2)得证。
对于(1) 若gcd(a,p) ≠ 1,那么由于p是素数必然p|a 因此p|a*(a^(p - 1) - 1)显然得证。
如果gcd(a,p) = 1, 则由(2)得p|(a^(p-1)-1)因此因此p|(a^(p - 1) - 1)*a显然得证。

欧拉函数是指：对于一个正整数n，小于n且和n互质的正整数的个数，记做：φ(n)
欧拉定理：gcd(a, n) = 1，则a^φ(n) ≡ 1 mod n
对于互质的整数a和n，有a^φ(n) ≡ 1 mod n
证明：
    首先证明下面这个命题：
    对于集合Z = {X1,X2,...,Xφ(n)}，Xi为第i个小于n且和n互质的正整数
    考虑集合
    S = {aX1 mod n,aX2 mod n,...,aXφ(n) mod n}
    则S = Z
    证明:
    a, n互质，Xi与n互质，那么aXi与n互质
    对于Z中两个元素Xi和Xj，如果Xi ≠ Xj
    则aXi mod n ≠ aXj mod n，若aXi mod n = aXj mod n，
    那么a(Xi - Xj) mod n  = 0, a, n互质那么Xi - Xj mod n = 0 由于
     0 < Xi, Xj < n显然不可能。因此S中的元素各不相同，这样S中含有φ(n)个与n互
    质的正整数，且满足每个元素小于n所以，很明显，S = Z        
    既然这样，那么
    (aX1 × aX2×...×aXφ(n))mod n
     = (aX1 mod n × aX2 mod n × ... × aXφ(n) mod n)mod n
     = (x1 × x2 × ... × xφ(n))mod n
     考虑上面等式左边和右边
     最左边等于(a^φ(n) × (x1 × x2 × ... × xφ(n)) mod n
     最右边等于(x1 × x2 × ... × xφ(n))mod n
     得到(a^φ(n) - 1)(x1 × x2 × ... × xφ(n)) mod n = 0
     而(x1 × x2 × ... × xφ(n))mod n 和 n互质,所以 a^φ(n) ≡ 1 mod n
推论：当n为素数时，φ(n) = n - 1,  a^(n - 1) - 1 mod n = 0 就是费马小定理(2)

模p乘法逆元
对于整数a、p，如果存在整数b，满足a*b mod p =1，则说，b是a的模p乘法逆元。
定理：a存在模p的乘法逆元的充要条件是gcd(a,p) = 1
证明：
首先证明充分性
如果gcd(a,p) = 1，根据欧拉定理，a^φ(p) ≡ 1 mod p，因此
显然a^(φ(p)-1) 是a的模p乘法逆元。
再证明必要性
假设存在a模p的乘法逆元为b， a*b ≡ 1 mod p
则a*b = k*p + 1，所以1 = a*b - k*p
因为gcd(a,p) = d
所以d | 1
所以d只能为1

扩展欧几里德定理
如果gcd(a, b) = d, 那么一定存在整数x，y 满足a*x + b*y = d

证明，可以设a = p*d， b = q*d，显然gcd(p, q) = 1
a*x + b*y = p*d*x + q*d*y  若要a*x + b*y = d
必需p*x + q*y = 1，由gcd(p, q) = 1 则p^φ(q) ≡ 1 mod q
令x = p^(φ(q) - 1) 也就是x是p的模q逆元，则p*x - 1 = p^φ(q) - 1 mod q
=  0， 那么y = (p^φ(q) - 1)/q显然是整数， 当然也可以y = q^(φ(p) - 1)，那
么x = (q^φ(p) - 1)/p，注意到x, y有无数组解。

注意到对于gcd(a,b) = d 我们对(a, b)用欧几里德辗转相除会最终得到
(d, 0)此时对于把a =d, b = 0 带入a*x + b*y = d，显然x = 1，y可以为任意值，
这里y可以为任意值就意味着解会有无数个。我们可以用a = d, b = 0的情况逆推出来
任何gcd(a, b) = d 满足a*x + b*y = d的解。如果x0, y0是b*x + (a%b)*y = d 的解，那么对于a*x + b*y = d的解呢？
b*x + (a%b)*y = d => b*x + (a - [a/b]*b)*y = a*y + b*(x - [a/b]*y)，所以a*x + b*y = d的解x1 = y0， y1 = x0 - [a/b]*y0; 这样我们可以程序迭带了。

void extended_gcd(int a, int b, int &x, int &y)
{
  
   if(b == 0)
      {
        x = 1;
        y = 0; //any integer is ok;
       }
       else
       {
         int t;
         extended_gcd(b, a%b, x, y);
         t = x;
         x = y;
         y = t - (a/b)*y;
       }  
}

终于写完了，参考了网上书上的内容。