关于数论乘法逆元及相关知识点

在求解a/b%m时，可以转化为(a%(b*m))/b，转化过程如下

令k = (a/b)/m(向下取整), x = (a/b)%m;

a/b = k*m + x (x < m);

a = k*b*m + b*x;

a%(b*m) = b*x;

a%(b*m)/b = x;

得证: a/b%m = a%(b*m)/b;(公式适用于很多情况：m不必是素数，b和m也不必互质)

上面的公式适用于b较小，a需要在线求且较大的时候;

比如：求a/b的值，且另a = x^k, 求a的过程随时会爆精度所以需要对a进行取模，但是求的过程中如果是直接对a进行%mod，之后再进行a/b，会出现错误的。正确的做法是根据上述公式，在求a过程中对(b*m)取模，之后再除以一次b即可。[ a/b%m = (a%(b*m))/b%m ]

所有的除法取模问题都可以用这种方法，但是当b很大的时候，则会出现爆精度问题，所以引出乘法逆元，将除法取模转换为乘法取模。

b存在乘法逆元的充要条件是b与模数m互质。设令c为b的逆元，即b*c ≡ 1(mod m);

解释：满足同余方程，b*c和1对m的模数相同，即b*c对m取余为1(m > 1);

那么 a/b = (a/b)*1 = (a/b)*(b*c)(mod m) = a*c(mod m);

即，除以一个数对m取模等于乘以这个数的逆元对m取模;

三种求逆元的方法&：

1.逆元求解利用扩欧。
2.当m为质数的时候直接使用费马小定理，m非质数时使用欧拉函数。
3.当m为质数的时候，使用神奇的线性方法。

扩展&：

1.利用乘法逆元求解组合数及组合数的其它求法。

2.线性时间求阶乘的逆元。

扩展欧几里得算法.

首先了解一下扩展欧几里得算法,

已知 a 和 b 的最大公约数是 gcd，那么，我们一定能够找到这样的 x 和 y ，使得: a*x + b*y = gcd 这是一个不定方程（其实是一种丢番图方程），有多解是一定的，但是只要我们找到一组特殊的解 x0 和 y0 那么，我们就可以用 x0 和 y0 表示出整个不定方程的通解：
    x = x0 + (b/gcd)*t; ((b/gcd)是b的因子，能表示的数只能多而不会少，下面同理)
    y = y0 – (a/gcd)*t;

现在，我们知道了一定存在 x 和 y 使得 ： a*x + b*y = gcd ， 那么，怎么求出这个特解 x 和 y 呢？只需要在欧几里德算法的基础上加点改动就行了。
我们观察到：欧几里德算法停止的状态是： a = gcd ， b = 0 ，那么，这是否能给我们求解 x, y 提供一种思路呢？这时候，只要 a = gcd 的系数是 1 ，只要 b 的系数是 0 或者其他值（无所谓是多少，反正任何数乘以 0 都等于 0 但是a 的系数一定要是 1），这时，我们就会有： a*1 + b*0 = gcd;
当然这是最终状态，但是我们是否可以从最终状态反推到最初的状态呢？
假设当前我们要处理的是求出 a 和 b的最大公约数，并求出 x 和 y 使得 a*x + b*y= gcd ，而我们已经求出了下一个状态：b 和 a%b 的最大公约数，并且求出了一组x1 和y1 使得： b*x1 + (a%b)*y1 = gcd ， 那么这两个相邻的状态之间是否存在一种关系呢？
我们知道： a%b = a - (a/b)*b（这里的 “/” 指的是整除，例如 5/2=2 , 1/3=0），那么，我们可以进一步得到：
        gcd = b*x1 + (a-(a/b)*b)*y1;
            = b*x1 + a*y1 – (a/b)*b*y1;
            = a*y1 + b*(x1 – a/b*y1);
对比之前我们的状态：求一组 x 和 y 使得：a*x + b*y = gcd，是否发现了什么？
这里：
        x = y1;
        y = x1 – a/b*y1;

上代码：

递归版e_gcd();

还有非递归版的;

前提条件：b与m互质,存在唯一解。

扩展欧几里德算法，我们既可以求出最大公约数，还可以顺带求解出使得： a*x + b*y = gcd 的通解 x 和 y
扩展欧几里德有什么用处呢？
求解形如 a*x +b*y = c 的通解，但是一般没有谁会无聊到让你写出一串通解出来，都是让你在通解中选出一些特殊的解，乘法逆元便是。

a*x ≡ 1(mod m);

这里称 x 是 a 关于 m 的乘法逆元, 这怎么求？可以等价于这样的表达式： a*x + m*y = 1;
看出什么来了吗？没错，当gcd(a , m) != 1 的时候是没有解的, 这也正是 a*x + b*y = c 有解的充要条件： c % gcd(a , b) == 0;
接着乘法逆元讲，一般，我们能够找到无数组解满足条件，但是一般是让你求解出最小的那组解，怎么做？我们求解出来了一个特殊的解 x0 那么，我们用 x0 % m其实就得到了最小的解了。为什么？
可以这样思考：
x 的通解不是 x0 + m*t 吗？
那么，也就是说， a是和 a 关于 m 的逆元是关于 m 同余的，那么根据最小整数原理，一定存在一个最小的正整数，它是 a 关于m 的逆元，而最小的肯定是在（0 , m）之间的，而且只有一个，这就好解释了。
有时候我们得到的特解 x0 是一个负数，还有的时候我们的 m 也是一个负数这怎么办？
当 m 是负数的时候，我们取 m 的绝对值就行了，当 x0 是负数的时候，他模上 m 的结果仍然是负数（在计算机计算的结果上是这样的，虽然定义的时候不是这样的），这时候，我们仍然让 x0 对abs(m) 取模，然后结果再加上 abs(m) 就行了，于是，我们不难写出下面的代码求解一个数 a 对于另一个数 m 的乘法逆元。

当然略微有点长了，比赛时用下面的就足以解决很多。

int ex_gcd(int a, int b, int &d, int &x, int &y){if(b == 0){x = 1; y = 0; d = a;return;}ex_gcd(b, a%b, d, y, x);//不断交错x和yy -= a/b*x;//当递归完成return时，y值return给了x, x值return给了y}int clac(int a, int m){int x, y, d;ex_gcd(a, m, d, x, y);return (m + x%m)%m;}

通过费马小定理或者欧拉函数去求.

1.费马小定理求逆元:

在m是素数的情况下，对任意整数x都有x^m ≡ x(mod m)。 
如果x为整数且无法被m整除，则有x^(m−1) ≡ 1(mod m)。 
所以可以在m为素数的情况下求出一个数的逆元，x * x^(m−2) ≡ 1(mod m)，x^(m−2)即为逆元。

推导过程如下：

然后就可通过快速幂求得逆元。

quickM(int base, int b);

2.欧拉函数求逆元:
令ϕ(m)表示小于等于m且与m互素的正整数的个数(特例(ϕ(1) = 1))。
如果x和m互质，则有x^ϕ(m) ≡ 1(mod m)，即x * x^(ϕ(m)−1) ≡ 1(mod m)，x^(ϕ(m)−1)即为x的逆元。 
在m为质数的情况下，ϕ(m) = m−1，即为费马小定理。

先介绍一下欧拉函数的一些知识(能学一点是一点...)

φ(x) = x*(1-1/p(1))(1-1/p(2))(1-1/p(3))(1-1/p(4))…..(1-1/p(n)) 其中p(1),p(2)…p(n)为x的所有质因数;x是正整数; 注意：每种质因数只有一个。φ(1)=1(唯一和1互质的数，且小于等于1)。

欧拉函数的性质：
(1)p^k型欧拉函数:
若N是质数p(即N=p),φ(N) = φ(p) = p-p^(1-1) = p-1。
若N是质数p的k次幂(即N=p^k)，φ(N)=p^k-p^(k-1) = (p-1)*p^(k-1)。
(2)m*n型欧拉函数:
设n为正整数，以φ(n)表示不超过n且与n互素的正整数的个数，称为n的欧拉函数值。若m,n互质，φ(m*n)=(m-1)*(n-1)=φ(m)*φ(n)。
(3)特殊性质:
若n为奇数时，φ(2n)=φ(n)。
(4)欧拉定理:
对于任何两个互质的正整数x,m(m>2)有:x^φ(m) ≡ 1(mod m)。
(5)扩展(费马小定理)
当m=p且x与素数p互质(即:gcd(x,p)=1)则上式有: x^(p-1) ≡ 1(mod p)。
(6)欧拉函数的延伸:
小于或等于n的数中，与n互质的数的总和为：φ(n) * n/2 (n>1)。

证明就不给出了，下面是模板代码:

求单个数的欧拉函数值:

int Euler(int n){int ans = n;//直接扫一般会超时，所以这里利用任何一个合数都至少//有一个不大于根号n的素因子，所以只需遍历到根号n即可 for(int i = 2; i*i <= n; ++i){if(n%i == 0){ans = ans/i*(i-1);//先除后乘防止中间数据溢出while(n%i == 0) n /= i; }}if(n > 1) ans = ans/n*(n-1);return ans;}

多次用到欧拉函数值，可通过筛选预先打表(时间复杂度O(n·logn)):

void getEuler(){int euler[maxn];for(int i = 1; i <= maxn; ++i)euler[i] = i;for(int i = 2; i <= maxn; i+=2)        euler[i] /= 2; //其实相当于进行了一次euler[i]/2*(2-1);for(int i = 3; i <= maxn; ++i)    {        if(euler[i] == i) //未被筛到, 即是素数, 则用此素数来筛        {            for(int j = i; j <= maxn; j+=i)            {                euler[j] = euler[j]/i*(i-1);            }        }    }} 

Last,

求出模数m的欧拉函数值后，进行一次快速幂 quickM(base,φ(n)-1)求得逆元。

线性时间内求范围内所有整数的逆元.

当m不是很大的时候(数组能够存起来的时候)可以线性时间求逆元。规定m为质数，
首先 1^(−1) ≡ 1(mod m)
然后我们设 m = k*x + r, r < x, 1 < x < p;
再将这个式子放到mod m意义下就会得到
k*x + r ≡ 0 (mod m)
两边同时乘上 x^(−1) * r^(−1)就会得到
k * r^(−1) + x^(−1) ≡ 0 (mod m)
x^(−1) ≡ −k * r^(−1) (mod m)
x^(−1) ≡ -(m/x) * (m%x)^(−1) (mod m)
从头开始扫一遍即可，时间复杂度O(n);

int inv[maxn];inv[1] = 1;for(int i = 2; i < maxn; ++i)    inv[i] = (-p/i + p) * inv[p%i] % p;

扩展：

利用乘法逆元求解组合数及组合数的其它求法.

首先声明C(n, m)是指n个中取m个，网上都没统一写法，真是的...

求组合数常用的公式为C(n, m) = n!/((n-m)! * m!);

当n, m都很小的时候可以利用杨辉三角直接求。
或者通过下面推导公式递归求。
C(n+1, m) = C(n, m) + C(n, m-1);
C(n, m) = C(n-1, m) + C(n-1, m-1);

简单的递归代码:

#include <stdio.h>int calc(int n, int m){if(n < m || n == 0) return 0;if(n == m || m == 0) return 1;return calc(n-1, m) + calc(n-1, m-1);}int main(){int n, m;while(scanf("%d %d", &n, &m)){printf("%d\n", calc(n, m));}return 0;}

再就是最常用的利用乘法逆元进行求解。

//可预先打出指定范围内阶乘的表void init(){fact[0] = 1;for(int i = 1; i <= maxn; ++i)fact[i] = fact[i-1]*i % mod;}

1. (n!/(m!*(n-m)!)) % mod = x %mod ,先对算出n！、m！、(n-m)!对mod取模的余数，就转换为a/b = x%mod;因为m为素数，所以等价于b*x +mod*y = gcd(b,mod); 然后用扩展的欧几里得定理算出 b*x0 +mod*y0 = 1的特解x0，x再进而确定 =(mod + x0%mod) %mod; 则a/b = a*x%mod;

LL e_gcd(LL a, LL b, LL &x, LL &y){if(b == 0){x = 1; y = 0;return a;}LL ans = e_gcd(b, a%b, y, x);y -= x*(a/b);return ans;}LL clac(LL a, LL m){LL x, y;e_gcd(a, m, x, y);return (m + x%m)%m;}LL C(int n, int m){return fact[n]*clac(fact[m]*fact[n-m]%mod, mod)%mod;}

2.如果mod是素数 则b的逆元其实就是b^(mod-2)即 (m!*(n-m)!)的逆元为 (m!*(n-m)!)^(mod-2);

LL quickM(LL a, LL b){LL base = a%mod, ans = 1;while(b){if(b&1) ans = ans*base%mod;base = base*base%mod;b >>= 1; }return ans;}LL C(int n, int m){return fact[n]*quickM(fact[m]*fact[n-m]%mod, mod-2)%mod;}

3.当n和m比较大而mod为素数且比较小(10^5左右)的时候，可以用Lucas定理计算.

Lucas定理:

A、B是非负整数，模数mod = p(p是质数)。A、B写成p进制：
A = a[n]a[n-1]...a[0];
B = b[n]b[n-1]...b[0];
则组合数C(A,B)与C(a[n], b[n])*C(a[n-1], b[n-1])*...*C(a[0], b[0])(mod p)同余
即：Lucas(n, m, p) = c(n%p, m%p)*Lucas(n/p, m/p, p);

模板代码:

LL quickM(LL base, LL b, LL p) {LL ans = 1;while(b) {    if(b&1) ans = (ans * base) % p;    base = (base*base) % p;    b >>= 1;}return ans;}//n,m过大不能打表只能在线求阶乘LL Comb(LL a, LL b, LL p) {if(a < b) return 0;if(a == b) return 1;if(b > a - b) b = a - b;LL ans = 1, ca = 1, cb = 1;for(LL i = 0; i < b; ++i) {    ca = (ca * (a-i)) % p;    cb = (cb * (b-i)) % p;}ans = (ca*quickM(cb, p-2, p)) % p;return ans;}LL Lucas(int n, int m, int p) {LL ans = 1;while(n && m && ans) {ans = (ans*Comb(n%p, m%p, p)) % p;n /= p;m /= p;}return ans;}

4.打表求阶乘数逆元的新方法.

打一个1~n的阶乘的逆元的表，假如n!的逆元叫做f[n]，可以先用费马小定理、扩展欧几里得等 
求出f[n]，再用递推公式求出前面的项。
我们记数字 x 的逆元为 f(x) (mod m)。
　　因为 n! = (n-1)! * n;
　　所以 f(n!) = f((n-1)! * n) = f((n-1)!) * f(n);
　　所以 f((n-1)!) = f(n!) * f(f(n)) = f(n!) * n;  (数的逆元的逆元就是它自身)
这样子我们就可以用后项推出前面的项了。

LL quickM(LL base, LL b){      LL ans = 1;      while(b){          if(b&1) ans = (ans * base) % mod;          base = (base*base) % mod;          b >>= 1;    }      return ans; }void init(){    fact[0] = 1;    for(int i = 1; i <= maxn; ++i)    fact[i] = fact[i-1]*i%mod;    fiv[maxn] = quickM(fact[maxn], mod-2);    for(int i = maxn-1; i >= 0; --i)    {        fiv[i] = fiv[i+1]*(i+1);        fiv[i] %= mod;    }}LL C(int n, int m){if(m > n) return 0;return fact[n]*fiv[n-m]%mod*fiv[m]%mod;}

5.发现一个近乎完美的init(): 阶乘、整数逆元、阶乘逆元，一次O(n)即可。

void init(){fact[0] = fact[1] = 1;fiv[0] = fiv[1] = 1; inv[1] = 1;for(int i = 2; i <= maxn; ++i)  {//递推保存fact阶乘，递推求inv和fiv各个逆元     fact[i] = fact[i-1]*i%mod;    inv[i] = (mod-mod/i)*inv[mod%i]%mod;    fiv[i] = inv[i]*fiv[i-1]%mod;}}

Over~~~

参考博文：

http://blog.csdn.net/yukizzz/article/details/51105009

http://blog.csdn.net/zhjchengfeng5/article/details/7786595

http://blog.csdn.net/acdreamers/article/details/8220787

http://www.cnblogs.com/frog112111/archive/2012/08/19/2646012.html

http://blog.csdn.net/u010582475/article/details/47707739

http://www.cnblogs.com/handsomecui/p/4755455.html

http://www.cnblogs.com/PJQOOO/p/3875545.html

http://blog.csdn.net/pi9nc/article/details/9615359

http://www.cnblogs.com/zxndgv/archive/2011/09/17/2179591.html

http://www.cnblogs.com/vongang/archive/2012/12/02/2798138.html

http://www.lofter.com/tag/%E9%80%86%E5%85%83

http://blog.csdn.net/jeromiewn/article/details/54799318