大组合数:Lucas定理

最近碰到一题,问你求mod (p1*p2*p3*……*pl) ,其中n和m数据范围是1~1e18 , l ≤10 , pi ≤ 1e5为不同的质数,并保证M=p1*p2*p3*……*pl ≤ 1e18 。

要解决这个问题首先需要Lucas定理 或者 C!解法。

Lucas定理

我们令n=sp+q , m=tp+rq , r ≤ p

那么,然后你只要继续对调用Lucas定理即可。

代码可以递归的去完成这个过程,其中递归终点为t = 0

伪代码,时间O(logp(n)*p)

int Lucas (ll n , ll m , int p) {
        return m == 0 ? 1 : 1ll*comb (n%p , m%p , p) * Lucas (n/p , m/p , p) %  p ; 
}
//comb()函数中,因为q , r ≤ p , 所以这部分暴力完成即可。

 Lucas定理证明:

证明资料:http://www.cut-the-knot.org/arithmetic/algebra/LucasTheorem.shtml

首先你需要这个算式:,然后

(1 + x) n Ξ (1 + x) sp+q  Ξ ( (1 + x)p)s • (1 + x) q Ξ (1 + xp) s • (1 + x)   (mod p) ;

.

所以,通过左右系数比较,你就会发现当i=t , j=r  (及xtp+r的系数等式)成立。

 

 题目:http://acm.hdu.edu.cn/showproblem.php?pid=5446

#include<bits/stdc++.h>
using namespace std;
typedef long long ll ;
const int M = 1e5 + 10 ;
ll n , m ;
int k ;
int prime[15] ;
int rem[15] ;

ll MM ;

int F[15] ;
//int Finv[15][M] , F[15][M] , inv[15][M] ;

ll mul (ll a , ll b , ll mod) {
        ll tmp = 0 ;
        while (b) {
                if (b & 1) tmp = (1ll*tmp+a)%mod ;
                b >>= 1 ;
                a = (a+a)%mod ;
        }
        return tmp ;
}

int Fermat (ll a , int b) {
        int ret = 1 ;
        int mod = b ;
        b -= 2 ;
        while (b) {
                if (b & 1) ret = mul(ret , a , mod) ;
                b >>= 1 ;
                a = mul(a , a , mod) ;
        }
        return ret ;
}

int fact (int n , int p) {
        int ret = 1 ;
        for (int i = 1 ; i <= n ; i ++) ret = 1ll*ret*i%p ;
        return ret ;
}

int comb (int n , int m , int p) {
        if (m < 0 || m > n) return 0 ;
        return 1ll* fact (n , p) * Fermat(fact (m , p) , p) * Fermat (fact(n-m , p) , p) % p ;
}

int Lucas (ll n , ll m , int p) {
        return m == 0 ? 1 : 1ll*comb (n%p , m%p , p) * Lucas (n/p , m/p , p) %  p ; 
}

void solve () {
        MM = 1 ;
        for (int i = 1 ; i <= k ; i ++) {
                rem[i] = Lucas (n , m , prime[i]) ;
                MM *= 1ll*prime[i] ;
                //cout << "rem[i] = " << rem[i] << endl;
        }
        ll sum = 0 ;
        for (int i = 1 ; i <= k ; i ++) {
                ll tmp = MM/prime[i] ;
                ll ans = mul (rem[i] , Fermat(tmp , prime[i]) , MM) ;
                sum = (sum + mul (ans , tmp , MM) ) % MM ;
        }
        printf ("%I64d\n" , sum);
}

int main () {
        int T ;
        scanf ("%d" , &T) ;
        while (T --) {
                scanf ("%I64d%I64d%d\n" , &n , &m , &k) ;
                for (int i = 1 ; i <= k ; i ++) {
                        scanf ("%d" , &prime[i]) ;
                }
                solve () ;
        }
        return 0 ;
}

 

posted @ 2015-09-16 16:29  92度的苍蓝  阅读(699)  评论(0编辑  收藏  举报
http://images.cnblogs.com/cnblogs_com/Running-Time/724426/o_b74124f376fc157f352acc88.jpg