扩展欧几里得算法&同余方程&模m乘法逆元详解

复习：求最大公约数算法

int gcd(int a, int b)
{
	return b ? gcd(b, a % b) : a;
}

首先介绍扩展欧几里得定理：

对于两个不全为0的整数a,b，必存在一组解x,y，使得ax+by=gcd(a,b)。

换句话说，形如ax+by的最小正整数等于gcd(a,b)。

实现代码如下：（一般题目都要用64位）

（复杂度：O(log max(a,b))）

typedef long long LL;

///求整数x和y，使得ax+by=gcd(a,b)，且|x|+|y|最小(|x|<=b,|y|<=a)
///注意：即使a,b在int范围内，中间计算过程可能超出int范围

void exgcd(LL a, LL b, LL& d, LL& x, LL& y)
{
	if (!b) d = a, x = 1, y = 0;///此时ax+by=ax,gcd(a,b)=gcd(a,0)=a,故ax=a,所以x=1,y可以取任意整数,这里为计算方便取y=0,你也可以取y=1,同样可以得到一组不同的解！
	else exgcd(b, a % b, d, y, x), y -= x * (a / b);///注意后面这一步这么写是由于前面的递归使得x和y交换了位置
	
	/*
	上面公式的推导过程:
	设x,y表示第一次递归时的值,x',y'表示第二次递归时的值
	由于
	        gcd(a,b)=gcd(b,a%b)
	所以
	        ax+by=gcd(a,b)等价于bx'+(a%b)y'=gcd(b,a%b)=gcd(a,b)
	即
	        bx'+(a-(a/b)*b)y'=gcd(a,b)
	变形整理得
	        ay'+b*(x'-(a/b)*y')=gcd(a,b)
	所以
	        x=y',y=x'-(a/b)*y'
	*/
}

至于一般情况：ax+by=c，gcd(a,b)|c 的解，请参考《数论概论》P24-25

这里补充下为何通解是x=x0+kb，y=y0-ka（当然你也可以写成x=x0-kb，y=y0+ka）：

首先把a化为a/gcd(a,b)，b化为b/gcd(a,b)，这样ax+by=1，gcd(a,b)=1（注意这里a，b实际上写作a'和b'，为简单起见省掉撇号）

由于

ax+by = ax0+by0

所以

a(x-x0) = -b(y-y0) = P

所以

a|P，b|P

因为a，b互素，所以

ab|P

即

P=kab

所以

x-x0=kb，y-y0=-ka

整理得到

x=x0+kb，y=y0-ka

变形1：同余方程

扩展欧几里得算法除了解决线性不定方程ax+by=c之外，

还可以求解同余方程a≡c(mod m)，b≡c(mod n)

将其化为a=c+mx，b=c+ny，作差得

mx-ny=a-b，从而化为线性不定方程解决

变形2：模m乘法逆元

定义：对于整数a，m，如果存在整数b，满足ab ≡ 1(mod m)，则说，b是a的模m乘法逆元。

定理：a存在模m的乘法逆元的充要条件是gcd(a,m) = 1

充分性：
因为

gcd(a,m) = 1

根据欧拉定理，有

a^φ(m) ≡ 1(mod m)

因此

a * a^(φ(m)-1) mod m = 1

所以存在a的模m乘法逆元，即a^(φ(m)-1)

必要性：
假设存在a模m的乘法逆元为b，则

ab ≡ 1 (mod m)

所以

ab = km +1

所以

1 = ab - km

由欧几里得定理，有

gcd(a,m) = 1

由定理知：

对于ax + by = 1，可以看出x是a模b的乘法逆元，y是b模a的乘法逆元。

反过来，要计算a模b的乘法逆元，就相当于求ax + by = 1的x的最小正整数解，从而化为线性不定方程解决。

（当然，你也可以求a^(φ(b)-1) mod b，但即使我们用φ函数公式和逐次平方法（复杂度也是O(log b)），由于取模的次数远大于扩展欧几里得算法和代码量略大，所以此法不实用）

****转载请注明出处：http://blog.csdn.net/synapse7/article/details/9901195****