HDOJ-3037（组合数学）

题意：将不超过m颗的相同的豆子放在n棵不同的树上，每棵树可以为空，求方案数mod p

1 <= n, m <= 1000000000, 1 < p < 100000,p是质数

题解：可以理解为有m颗豆子，在n棵树上放k颗，然后再加一棵树，放m-k颗，于是变成了m颗相同的豆子放在n+1棵不同树上的方案数。

也就是求a[1]+a[2]+a[3]+......+a[n+1]=m，(a[i]>=0)的方案数，但是这种情况并不好计算，我们可以让每一份都加上1，令b[i]=a[i]+1>=1, 则b[1]+b[2]+b[3]+......b[n+1]=m+n+1,

可以用插板法了，m+n+1个元素有m+n个空，分成n+1份就是插n个板子，所以答案就是C(n+m,m) %p

那么问题来了，n,m如此巨大，如何计算这个组合数呢？我们观察到p是质数且p比较小。我们可以用“Lucas定理”来解决！！！

附百度百科：

数论Lucas定理是用来求 c(n,m) mod p的值,p是素数（从n取m组合，模上p）。

描述为:

Lucas(n,m,p)=cm(n%p,m%p)* Lucas(n/p,m/p,p)

Lucas(x,0,p)=1;

而

cm(a,b)=a! * (b!*(a-b)!)^(p-2) mod p

也= (a!/(a-b)!) * (b!)^(p-2)) mod p

这里，其实就是直接求 (a!/(a-b)!) / (b!) mod p

由于 (a/b) mod p = a * b^(p-2) mod p

code：

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <algorithm>
using namespace std;
const int MAXP=100005;
typedef long long ll;

ll f[MAXP];

ll mul(ll x, ll y, ll p)
{
    ll res=0;
    while(y)
    {
        if(y&1) res=(res+x)%p;
        x=(x<<1)%p;
        y>>=1;
    }
    return res;
}

void exgcd(ll a,ll b,ll &d,ll &x,ll &y)
{
    if (!b)
    {
        d=a;
        x=1;
        y=0;
    }
    else
    {
        exgcd(b,a%b,d,y,x);
        y-=a/b*x;
    }
}

void initp(int p)
{
    f[0]=f[1]=1;
    for(int i=2;i<p;i++)f[i]=(f[i-1]*i)%p;
}

ll comb(ll n,ll m,ll p)
{
    if(m>n) return 0;
    ll ans=f[n];
    ll g,x,y;
    exgcd(f[m],p,g,x,y);
    //ans=(ans*(x%p)%p+p)%p;
    ans=mul(ans,(x%p)%p+p,p);
    exgcd(f[n-m],p,g,x,y);
    //ans=(ans*(x%p)%p+p)%p;
    ans=mul(ans,(x%p)%p+p,p);
    return ans;
}

//n取m，如果可能m>n请特判return 0;
ll lucas(ll n,ll m,int p)
{
    ll ans=1;
    initp(p);
    if(m>n) return 0;
    while(n && m && ans)
    {
        //ans=comb(n%p,m%p,p)*ans%p;
        ans=mul(comb(n%p,m%p,p),ans,p);
        n/=p;
        m/=p;
    }
    return ans;
}


void work()
{
    ll n,m,p,ans;
    scanf("%I64d%I64d%I64d",&n,&m,&p);
    ans=lucas(n+m,m,p);
    printf("%I64d\n",ans);
}

int main()
{
    int T;
    scanf("%d",&T);
    while(T--)
        work();
    return 0;
}