卢卡斯定理

卢卡斯定理:解决一类组合数取模问题

A、B是非负整数，p是质数。AB写成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])  modp同余

即：Lucas(n,m,p)=c(n%p,m%p)*Lucas(n/p,m/p,p) 

这个是单独处理n!的情况，当然C(n,m)就是n!/(m!*(n-m)!)，每一个阶乘都用上面的方法处理的话，就是Lucas定理了，注意这儿的p是素数是有必要的。

Lucas最大的数据处理能力是p在10^5左右，不能再大了，hdu 3037就是10^5级别的！

对于大组合数取模，n,m不大于10^5的话，用逆元的方法，可以解决。对于n,m大于10^5的话，那么要求p<10^5，这样就是Lucas定理了，将n,m转化到10^5以内解。

 

模板:

/*
Lucas定理 C(n,m)%p(p为素数)
C(n,m)与C(a[n],b[n])*C(a[n-1],b[n-1])*C(a[n-2],b[-2])*....*C(a[0],b[0])模p同余
a,b 是n，m在p进制下的数
有的推公式： (C(n%p,m%p,p)*Lcs(n/p,m/p,p))%p;
关键是求 C(n%p,m%p,p) 递归会很慢 for的话会爆掉
这里用一个定理：a/b%p <--> a*x%p  x为b在b%p下的逆元
再来一个定理：x=b^(p-2)   x为b在%p下的逆元  p为素数
然后预处理一下阶乘就ok了 
*/
#include<iostream>
#include<cstdio>
#include<cstring>

#define N 1001

using namespace std;
int f[N];

int Mi(int a,int b,int p)
{
    int res=1;
    while(b)
    {
        if(b&1) res=res*a%p;
        b>>=1;a=a*a%p;
    }return res;
}

int C(int n,int m,int p)
{
    if(m>n)return 0;
    return  f[n]*Mi(f[m]*f[n-m],p-2,p)%p;
}

int Lus(int n,int m,int p)
{
    if(m==0) return 1;
    return (C(n%p,m%p,p)*Lus(n/p,m/p,p))%p;
}

int main()
{
    int n,m,p;
    scanf("%d%d%d",&n,&m,&p);
    f[0]=1;
    for(int i=1;i<=p;i++)
      f[i]=f[i-1]*i%p;
    printf("%d\n",Lus(n,m,p));
}

 

hdu3037 http://acm.hdu.edu.cn/showproblem.php?pid=3037

题目大意：求在n棵树上摘不超过m颗豆子的方案，结果对p取模。

思路：题目可以转换成  x1+x2+……+xn=m 有多少组解，m在题中可以取0~m。

利用插板法可以得出x1+x2+……+xn=m解的个数为C（n+m-1，m）；

则题目解的个数可以转换成求   sum=C（n+m-1，0）+C（n+m-1，1）+C（n+m-1，2）……+C（n+m-1，m）

利用公式C（n，r）=C（n-1，r）+C（n-1，r-1）  == >  sum=C（n+m，m）。

就是要求C（n+m，m）%p。

#include<iostream>
#include<cstdio>
#include<cstring>

#define N 100007

using namespace std;
long long f[N];

long long Mi(long long a,long long b,long long p)
{
    long long res=1;
    while(b)
    {
        if(b&1) res=res*a%p;
        b>>=1;a=a*a%p;
    }return res;
}

long long C(long long n,long long m,long long p)
{
    if(m>n)return 0;
    return  f[n]*Mi(f[m]*f[n-m],p-2,p)%p;
}

long long Lcs(long long n,long long m,long long p)
{
    if(m==0)return 1;
    return (C(n%p,m%p,p)*Lcs(n/p,m/p,p))%p;
}

int main()
{
    long long n,m,p;long long t;        
    cin>>t;
    while(t--)
    {
        cin>>n>>m>>p;
        f[0]=1;
        for(long long i=1;i<=p;i++)
          f[i]=f[i-1]*i%p;
        printf("%lld\n",Lcs(n+m,m,p));
    }
    return 0;
}