Lucas 定理:n、m是非负整数,p是素数。那么:
1、C(n,m)% p = Lucas(n,m,p)
2、Lucas(n,m,p)= (C(n%p,m%p)% p)*Lucas(n/p,m/p,p)
证明:
我们将n和m表示成p进制数:n = akak-1..a1(p) m = bkbk-1...b1(p)(位数不够高位补0)
通过简单分析,可以得到: C(n,m)和C(a1,b1)*C(a2,b2)*...*C(ak,bk)关于p同余
n = ak*p^(k-1) + ... + a1 m = bk*p^(k-1) + ... + b1
关于为什么同余,可以通过上面的式子观察出来,也可通过精确推导得出,篇幅所限,这里也不在累述。
我们已经知道了,Lucas定理的公式,但是关键是怎么编码来得出结果呢?
也许你会说,那还不简单,直接递归就搞定了,其实没那么简单,还有一个重要的问题没解决,C(ai,bi)(mod p)的求法。
这里还涉及到逆元以及欧拉-费马定理的应用。
关于C(n,m)(mod p)(p为素数)的求法:
1、C(n,m)= n! / (m! * (n - m)!) (组合数的求法)
2、(a / b) (mod p) = (a * x) (mod p) x表示b的逆元 并且 b * x = 1 (mod p)(乘法逆元)
3、b^phi(p) = 1 (mod p)(b和p互质,phi(p)为p的欧拉函数值) (欧拉-费马定理)
通过上述三式,我们知道了C(n,m)(mod p)的求法,用fac[i] = i! % p
C(n,m)(mod p) = fac[n] * pow(fac[m]*fac[n-m],p-2)
code:
1 #include <cstdio> 2 typedef __int64 LL; 3 const int MAXN = 1000005; 4 LL n, m, p; 5 LL fac[MAXN]; 6 7 // 得到阶乘 fac[i] = i! % p 8 void GetFact() 9 { 10 fac[0] = 1; 11 for (LL i = 1; i < MAXN; ++i) 12 fac[i] = fac[i - 1] * i % p; 13 } 14 15 // 快速模幂 a^b % p 16 LL Pow(LL a, LL b) 17 { 18 LL temp = a % p; 19 LL ret = 1; 20 while (b) 21 { 22 if (b & 1) ret = ret * temp % p; 23 temp = temp * temp % p; 24 b >>= 1; 25 } 26 return ret; 27 } 28 29 /* 30 欧拉定理求逆元 31 32 (a / b) (mod p) = (a * x) (mod p) x表示b的逆元 并且 b * x = 1 (mod p) 只有b和p互质才存在逆元 33 34 b * x = 1 (mod p) x是b关于p的逆元 35 36 b^phi(p) = 1 (mod p) 37 38 b * b^(phi(p) - 1) (mod p) = b * x (mod p) 39 40 x = b^(phi(p) - 1) = b^(p - 2) 41 42 (a / b) (mod p) = (a * x) (mod p) = (a * b^(p - 2)) (mod p) 43 44 经过上面的推导,得出: 45 46 (a / b) (mod p) = (a * b^(p - 2)) (mod p) (b 和 p互质) 47 48 */ 49 LL Cal(LL n, LL m) 50 { 51 if (m > n) return 0; 52 return fac[n] * Pow(fac[m] * fac[n - m], p - 2) % p; 53 } 54 55 LL Lucas(LL n, LL m) 56 { 57 if (m == 0) return 1; 58 return Cal(n % p, m % p) * Lucas(n / p, m / p) % p; 59 } 60 61 int main() 62 { 63 int nCase; 64 scanf("%d", &nCase); 65 while (nCase--) 66 { 67 scanf("%I64d %I64d %I64d", &n, &m, &p); 68 printf("%I64d ", Lucas(n, m)); 69 } 70 return 0; 71 }