拓展欧几里德
例题 同余方程
题目描述
求关于 \(x\) 的同余方程 \(a x \equiv 1 \pmod {b}\) 的最小正整数解。
输入输出格式
输入格式：
一行，包含两个正整数 \(a,b\)，用一个空格隔开。
输出格式：
一个正整数 \(x_0\)，即最小正整数解。输入数据保证一定有解。
输入输出样例
输入样例#1：
3 10
输出样例#1：
7
说明
【数据范围】
对于 100%的数据，\(2 ≤a, b≤ 2,000,000,000\) 。
NOIP 2012 提高组 第二天 第一题
问题转化
-
题目问的是满足 \(ax \mod b = 1\) 的最小正整数 \(x\) 。（ \(a,b\) 是正整数）,但是不能暴力枚举 \(x\) ，会超时。
-
把问题转化一下。观察 \(ax \mod b = 1\) ，它的实质是 \(ax + by = 1\) ：这里 \(y\) 是我们新引入的某个整数，并且似乎是个负数才对。这样表示是为了用扩展欧几里得算法。我们将要努力求出一组 \(x,y\) 来满足这个等式。稍微再等一下——
-
问题还需要转化。扩展欧几里得是用来求 \(ax + by = gcd(a,b)\) 中的未知数的，怎么牵扯到等于 \(1\) 呢？
-
原理是，方程 \(ax + by = m\) 有解的必要条件是 \(m \mod gcd(a,b) = 0\) 。
-
这个简单证一下。
由最大公因数的定义，可知 \(a\) 是 \(gcd(a,b)\) 的倍数，且 \(b\) 是 \(gcd(a,b)\) 的倍数，
若 \(x,y\) 都是整数，就确定了 \(ax + by\) 是 \(gcd(a,b)\) 的倍数，
因为 \(m = ax + by\) ，所以 \(m\) 必须是 \(gcd(a,b)\) 的倍数，
那么 \(m \mod gcd(a,b) = 0\)。
- 可得出在这道题中，方程 \(ax + by = 1\) 的有解的必要条件是 \(1 \mod gcd(a,b) = 0\) 。可怜的 \(gcd(a,b)\) 只能等于 \(1\) 了。这实际上就是 \(a,b\) 互质。
扩展欧几里得
- 前提：知道普通欧几里得算法（辗转相除法）。
(下面字母挺多，希望你耐心地慢慢地读~)
-
我们拿到了一组 \(a,b\) 。设 \(G = gcd(a, b)\) 。那么目标是求出满足 \(ax + by = G\) （①） 的整数 \(x\) 与 \(y\) 。其中 \(x\) 应当是满足条件的最小正整数，即答案，而 \(y\) 是辅助答案。
-
（注意，虽然刚刚已经证明本题的 \(gcd(a,b)\) 必然等于 \(1\) ，但是扩展欧几里得算法本身过程求的是 ax + by = gcd(a,b)$ 的解。它正好非常适合本题，不过我们要按照它求解的过程去做）
-
如果我们先前已经求出了另一组数 \(x_2, y_2\) ，它们满足这么一个式子： \(bx_2 + (a \mod b)y_2 = G\) （②），则此时结合①②一定有：
\(ax + by = bx_2 + (a \mod b)y_2\) （③）
-
可见，在这个“如果”实现的时候，我们的目标变成了“求出满足上式的 \(x\) 和 \(y\) ”。
-
其中 \(a,b,x_2,y_2\) 都已知， \(x,y\) 待求。因为未知数比方程更多，所以没有唯一解。我们先求出一组必然存在的解，最后将在“答案处理”时转为最小解。
-
怎么求呢？取模运算是 \(a \mod b = a - b×(a/b)\) ，所以方程③实际上是：
\(ax + by = bx_2 + (a-b×(a/b))y_2\)
\(\Rightarrow ax + by = bx_2 + ay_2 - b × (a/b)y_2\)
\(\Rightarrow ax + by = ay_2 + b(x_2-(a/b)y_2)\)
- 看上面这个方程，一组必然存在的解出现了：
\(x = y_2, y = x_2 - (a/b)y_2\) （④）
- 可见，我们只要求出 \(x_2,y_2\) ，就能得出正确的 \(x,y\) 。问题是 \(x_2,y_2\) 怎么求。
现在我们手上是 \(b,a \mod b\) 这两个系数，而目标是求出 \(x_2\) 和 \(y_2\) 满足：
\(bx_2 + (a \mod b)y_2 = G\)（②）
把①和②对比一下：
\(ax + by = G\) （①）
原方程中的系数 \(a\) 变成了②中的系数 \(b\) ，原方程中的系数 \(b\) 变成了②中的 \(a \mod b\)而已。
-
所以，把新的方程也看作 \(ax + by = G\) 的形式（只是系数 \(a\) 和 \(b\) 的具体数值改变了）。然后按照上面的一模一样下来（其实都只是推导过程），我们发现，最好有 \(x_3,y_3\) 来支撑 \(x_2, y_2\) 。
-
再一模一样下来，我们又需要 \(x_4,y_4\) 来支撑 \(x_3, y_3\)。
……
这个递归中 \(a, b\) 不断被替代为 \(b, a \mod b\)，这个替换方式与普通欧几里得是一样的，所以最后会出现 \(a_n = G, b_n = 0\) 。
这时要直接返回了，我们需要一组 \(x_n,y_n\) 满足 \(a_nx_n + b_ny_n = G\)（⑤）。
然而该层的 \(a_n = G, b_n = 0\) 。所以只要⑤左边取 \(x_n = 1\) ，这个方程就妥妥的成立了。
（最后一层的 \(y_n\) 建议取 0 。然而由于 \(b = 0\) ，就算返回其它数值，方程也一定成立。但这样的程序容易出错，因为 \(y_n\) 在回溯时滚雪球式增长，容易数值越界。
最后一层结束后，就开始返回，直到最上层。每一层可以轻松地根据下层的 \(x{k+1},y_{k+1}\) 求出当前层的 \(x_k, y_k\) 。
整个过程就是：以辗转相除的方式向下递进，不断缩小系数，保证会出现有确定解的最后一层。
答案处理
一个遗留问题：我们将要求出来的 \(x,y\) 仅仅保证满足 \(ax+by=1\) ，而 \(x\) 不一定是最小正整数解。有可能太大，也有可能是负数。这依然可以解决，那就是， \(x\) 批量地减去或加上 \(b\) ，能保证 \(ax+by=1\)，因为：
\(ax+by=1\)
\(ax+by+k×ba−k×ba=1\)
\(a(x+kb)+(y−ka)b=1\)
1.显然这并不会把方程中任何整数变成非整数。
2.“加上或减去 \(b\) ”不会使 \(x\) 错过任何解。可以这么理解：
已经求出一组整数 \(x,y\) 使得 \(ax+by=1\)，也就是 \((1−ax)/b=y\) 。 \(y\) 是整数，可见目前 \(1−ax\) 是 \(b\) 的倍数。
现在想改变 \(x\) 并使得方程仍然成立。已知 \(a,b\) 互质，假若 \(x\) 的变化量（ \(Δx\)）不是 \(b\) 的倍数，则 \(1−ax\) 的变化量（ \(−a×Δx\) ）也不是 \(b\) 的倍数，这会使得 \(1−ax\) 不再是 \(b\) 的倍数，则 \(y\) 不是整数了。
仅当 \(x\) 的变化量是 \(b\) 的倍数时， \(1−ax\) 能保持自己是 \(b\) 的倍数，此时就出现新的解了。
因此到最后，如果 \(x\) 太小就不断加 \(b\) 直到大于等于 \(0\)$ ，太大则一直减 \(b\) ，直到最小正整数解。也就是这么写：
x = (x % b + b) % b;//括号中取模再加，可以处理负数
代码
推导中的所有 \(x,y\) 共用全局变量 \(long long\) 传递也很方便。
#include <bits/stdc++.h>
#define int long long //全局变量 long long
using namespace std;
int x,y;////目前方程真正的解
void exgcd(int a,int b){////当前目的：求解 ax + by = gcd(a, b) 这么一个方程
if(b==0){ //在 b = 0 时方程还要成立， 使 x = 1, y = 0 ，必然成立
x=1;y=0;//建议返回0。
return;
}
exgcd(b,a%b);//先求下一个方程的解
//现在我们已经拿到了下一个方程的解x, y
int res=x;//暂时存一下x
x=y;
y=res-(a/b)*y;
}
signed main(){
int a,b;
scanf("%lld%lld",&a,&b);
exgcd(a,b);
printf("%lld\n",(x%b+b)%b);//我们求出来的x必然满足方程，但不一定是最小正整数解，所以要进行答案处理
return 0;
}
