题解 P1082 【同余方程】
学委
2018-08-06 21:01:29
*2020-02-27 更新*
一步一步来。
## 问题转化
题目问的是满足 $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$ 在回溯时滚雪球式增长,容易数值越界。下面的代码在最后一层令 $y = 7$,开了long long,通过了此题。)
最后一层结束后,就开始返回,直到最上层。每一层可以轻松地根据下层的 $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$,也就是 $\frac{1-ax}{b} = y$。$y$ 是整数,可见目前 $1-ax$ 是 $b$ 的倍数。
现在想改变 $x$ 并使得方程仍然成立。已知 $a,b$ 互质,假若 $x$ 的**变化量**($\Delta x$)不是 $b$ 的倍数,则 $1-ax$ 的变化量($-a×\Delta x$)也不是 $b$ 的倍数,这会使得 $1-ax$ 不再是 $b$ 的倍数,则 $y$ 不是整数了。
仅当 $x$ 的变化量是 $b$ 的倍数时,$1-ax$ 能保持自己是 $b$ 的倍数,此时就出现新的解了。
___
因此到最后,如果 $x$ 太小就不断加 $b$ 直到大于等于 $0$,太大则一直减 $b$,直到最小正整数解。也就是这么写:
```cpp
x = (x % b + b) % b;//括号中取模再加,可以处理负数
```
## 代码
推导中的所有 $x,y$ 共用全局变量 long long x, y,传递也很方便。
```cpp
#include<bits/stdc++.h>
using namespace std;
long long x, y;//目前方程真正的解
void exgcd(long long a, long long b)
{
//当前目的:求解 ax + by = gcd(a, b) 这么一个方程
if(b == 0) //a, b不断改变的过程中,b最终必然会成为0
{
//在 b = 0 时方程还要成立? 使 x = 1, y = 0 ,必然成立
x = 1;
y = 7; //建议返回0。不过y = 7能AC,证明了最后一个等式不受最后一个y影响
return;
}
exgcd(b, a % b);//把下一层系数传进去(先求下一个方程的解 )
//现在我们已经拿到了下一个方程的解x, y
long long tx = x;//暂时存一下x,别丢了
x = y;
y = tx - a / b * y;
}
int main()
{
long long a, b;
cin >> a >> b;
exgcd(a, b);
x = (x % b + b) % b;//我们求出来的x必然满足方程,但不一定是最小正整数解,所以要进行答案处理
printf("%lld\n", x);
return 0;
}
```
求这样一个满足 $ax \mod b = 1$ 的 $x$ 有什么用呢?一个重要用途如下。
这个 $x$ 就是:$a$ 在模 $b$ 意义下的乘法逆元。
在模 $b$ 意义下,如果想要**除以 $a$** 就非常麻烦。这时候**乘以 $a$ 的逆元**等效于**除以 $a$。**
***
假设我们已经算出了:
```cpp
ans = (a * b * c) % p;
```
现在要从 $ans$ 中把 $b$ 给除掉,如何处理 $ans$?
如果像下面这样直接除会出错(举个实例会很直观):
```cpp
ans = ans / b % p;
```
所以我们就求出 $b$ 的逆元 $x$(满足 $bx \mod p = 1$),然后直接算这个:
```cpp
ans = ans * x % p;
```
原理可以这么理解:
$abcx \equiv ac(bx) \equiv ac \pmod p$
***
更新记录:
2020-02-01 证明略有改动;加入了逆元一种作用的说明。
2020-02-03 加点补充说明减少误解。
2020-02-27 语句修改。