乘法逆元

在 OI 中，大多数情况下，善良的出题人为了避免高精度等大整数计算，常常会要求输出答案对一个数（大多是质数）取模的情况。
但这衍生了一个问题：若题目中计算需用到除法，而 a≡b(modc) 不能推出⌊a/d⌋≡⌊b/d⌋(modc)，
例如：9≡16(mod7)正确，⌊9/5⌋≡⌊16/5⌋(mod7)错误。

既然除法不能用，我们能否将除法变成乘法呢？
要解决这个问题，就需要用到一个概念：乘法逆元。

乘法逆元，有如下定义：
如果 a * x ≡ 1（mod p） ，那么我们就说在mod p情况下 x 是 a 的逆元。当然啦，反过来， a 也是 x 的逆元。 ~ 有点像倒数哎！~ （逆元小于模数）
因为于此我们可发现 m / a ≡ m /a * a * x ≡ m * x (mod p),问题便迎刃而解！！！

逆元存在性

但是，我们之前一直都没考虑这个问题：乘法逆元一定存在吗？
要知道这个我们需要先知道一个定理

裴蜀定理

裴蜀定理，又称贝祖定理（Bézout's lemma），是一个关于最大公因数的定理。

定义：

设 a,b 是不全为零的整数，则存在整数 x,y，使得 ax+by=gcd(a,b).

证明：

记 a,b 线性组合集为 S={z∣z=ax+by,x∈Z,y∈Z}，那么裴蜀定理说的就是 gcd(a,b)∈S。
实际上，我们可以证明 gcd(a,b) 是 S 中的最小正元素，
| 因为 gcd(a,b) 是 a,b 的最大公约数，那么 a,b 是 gcd(a,b) 的整数倍，
| 即 a=k1 * gcd(a,b), b=k2 * gcd(a,b),
| 那么 z = (k1 * x + k2 * y) * gcd(a,b),一定能被gcd(a,b)整除，
| 设 d 是 S 中最小的正数，则有d≥gcd(a,b)。
| 另设 e 是 S 中任意的一个数，
| e mod d = e - k * d = k1 * x + k2 * y ∈S，
| 又因为 e mod d 的结果是 [0,d−1] 范围内的整数，
| 而 S 中最小的正数是 d，因此 e mod d 只能取 0，即 e 能被 d 整除。
| 也就是说，S 中任意一个整数都能被 d 整除
| 因为 a,b ∈S ,所以a,b都能被d整除，因此d是a,b的公因数，则有d≤gcd(a,b)
| 因此 d = gcd(a,b)

逆元存在性证明

如果 x 是 a 在模 p 意义下的乘法逆元，则 ax≡1(mod p)，那么这个方程可以写成 ax+pk=1 (k∈Z)
还记得上面说的线性组合集的性质吗？对于 a,p 的线性组合集 S={z∣ax+py,x∈Z,y∈Z}，集合中最小的正元素为 gcd(a,p)。若方程 ax+pk=1 有解，说明 gcd(a,p)=1，即 a,p 互质。

存在性结论

我们便可以得出，a 在模 p 意义下的乘法逆元存在的充分必要条件是 a,p 互质。
这也在一定程度上解释了大多数题目模数为质数的原因：（设模数为 p）可以保证在 [1,p−1 中的所有整数都有模 p 意义下的逆元。

逆元性质

我们记a的逆元为inv[a]

1、存在唯一性

存在性在上面已经证明了
下面来看看唯一性
　证明： 假设对于一个数a，有两个不同的逆元 inv1[a] 和 inv2[a]
　　　　不妨先设 inv1[a] < inv2[a] 且 inv2[a] - inv1[a] =k
　　　　由于 a ≠ 0 ，所以 k =0 (mod p).
　　　　所以 inv1[a] 与 inv2[a]的值相等，即只有一个逆元

2、完全积性

即 inv[a] * inv[b]=inv[a * b]
　证明：inv[a] * a ≡ inv[b] * b ≡ 1 (mod p)
　　　　inv[a] * a * inv[b] * b ≡ 1*1 (mod p)
　　　　inv[a] * inv [b] * a * b ≡ 1 (mod p)
　　　　inv[a] * inv [b] = inv[a * b]

3、a*inv[b]=a/b（mod p)

这个性质非常重要：当需要算出 (a/b) mod p 的值时，使用朴素的方法，我们只能在 a 上不断加 p ，直到它能被 b 整除为止。
但是当 a,b,p 都很大的时候，这种方法就GG了，但如果有了逆元，我们就可以 “ 非常方便，快捷 ” 地求解。

逆元求法

1、枚举

太垃圾了，写了估计必炸。。。。。。

2、费马小定理

我们都知道费马小定理：
当 p 为素数时， a^{p-1}=1 (mod p)
那么 a * a^{p-2}=1 (mod p)
所以我们只需要快速幂求出 a^（p-2）就ok了

3、拓欧（扩展欧几里得算法）

原理

前面我们已经证明了方程 ax+by=gcd(a,b) 一定有整数解，但是我们不知道此时 x,y 取什么值。扩展欧几里得就是用来解出 x,y 的值的算法。
定义函数 exgcd(a,b)，它返回 ax+by=gcd(a,b) 的整数解 (x,y)。
如果 b≠0，则有 gcd(a,b)=gcd(b,a mod b)
假设我们已经求出 exgcd(b,a mod b) 的解 (x′,y′)，那么 ax+by=bx′+(a mod b)y′，即
ax+by=bx′+(a−b×⌊a/b⌋)y′
整理得
ax+by=ay′+b(x′−[a/b]y′)
所以 x=y′ ，y=x′−[a/b]y′
问题就成了如何求 exgcd(b,a mod b)
递归下去。当 b=0 的时候可以直接出解。

代码

其实就是辗转相除法的拓展

int gcdEx(int a, int b, int* x, int* y)
{
    if (b == 0)
    {
        *x = 1, * y = 0;
        return a;
    }
    else
    {
        int r = gcdEx(b, a % b, x, y);
        int t = *x;
        *x = *y;
        *y = t - a / b * *y;
        return r;
    }
	}

运算方法

仔细想想，发现其实我们寻找 a * x=1 （ mod p ）的解就相当于寻找方程 a * x+p * y=1 的解。
我们知道 扩展欧几里得 ax+by=c
好的，确定 a 和 p 互质，可以用了。

4、欧拉筛

如果我们要求出 1 ~ p-1 所有数的逆元呢？
还记得逆元是完全积性函数吗？所以对于每个合数 a ，我们把所有它的因子的逆元筛出来再相乘即可。
所以我们可以直接把所有素数筛出来，对它们求逆元，再把它的逆元乘给它的倍数就可以了。