数论进阶:Lucas定理

Lucas定理

Lucas定理详述

内容

对于一个组合数$C_m^n$，取模一个素数$p$的值，假设$n=a\times p+b$,$m=c\times p+d$有：

$C_n^m≡C_a^c\times C_b^d(\mod p)$

证明

先证明一个一会儿要用到的定理：

$(1+x)^p≡1+x^p(\mod p)$

这里可以用费马小定理加以证明:

由于$(1+x)^p≡(1+x)(\mod p)$

且$x^p≡x(\mod p)$即$1+x^p≡1+x(\mod p)$

于是得证：$(1+x)^p≡1+x^p(\mod p)$

然后设$n=a\times p+b$,$m=c\times p+d$

于是：$(1+x)^n≡(1+x)^{a\times p+b}(\mod p)$

即：

$(1+x)^n≡(1+x)^{a\times p}\times (1+x)^b(\mod p)$

由之前证明的定理得：

$(1+x)^{a\times p}≡(1+x^p)^a(\mod p)$

即：

$(1+x)^{a\times p+b}≡(1+x^p)^a\times (1+x)^ b(\mod p)$

由二项式定理得：

$(1+x)^{a\times p+b}≡\displaystyle \sum^{i\to a}{i=0}(C_a^i\times x^{p\times i})\times \displaystyle \sum^{j\to b}{j=0}(C_b^j\times x^j)(\mod p)$

考虑等号两边$x^m$的系数：

由二项式定理得左边$x^m$的系数为：$C_n^m$

将右边的式子展开，显然当i,j满足$p\times i+j==m$时为$x^m$的系数,也就是$i==c,j==d$时，也就是：

$C_a^c\times C_b^d$

则能得到等式$C_n^m≡C_a^b\times C_b^d(\mod p)$

于是定理得到了证明。

代码

#include <iostream>
#include <cmath>
#include <cstring>
#include <cstdio>
#include <algorithm>
using namespace std;
typedef long long ll;


//计算n！%mod
ll fact(ll n,ll mod)
{
    ll res=1;
    for (int i = 1; i <= n; ++i) {
        res=(res*i)%mod;
    }
    return res;
}


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


//用exgcd计算t模mod的逆元
ll inv(ll t,ll mod)
{
    ll d,x,y;
    exgcd(t,mod,x,y,d);
    return d==1?(x%mod+mod)%mod:-1;
}


//计算C(n,m)%mod:C(n,m)%mod=(n!/[m!*(n-m)!])%mod=n!%mod*inv(m!,mod)%mod*inv((n-m)!,mod)%mod
ll comb(ll n,ll m,ll mod)
{
    if(m<0||m>n)return 0;
    return fact(n,mod)*inv(fact(m,mod),mod)%mod*inv(fact(n-m,mod),mod)%mod;
}


//lucas定理:C(n,m)%mod=C(n/mod,m/mod)*C(n%mod,m%mod)%mod
ll lucas(ll n,ll m,ll mod)
{
    return m?lucas(n/mod,m/mod,mod)*comb(n%mod,m%mod,mod)%mod:1;
}

拓展Lucas定理

考虑如果p不是素数的情况，我们考虑利用唯一分解定理将p分解为若干素数的乘积，于是可以将式子分解为一组方程组，即：

$C_n^m≡r_1(\mod p_1)$

$C_n^m≡r_2(\mod p_2)$

$C_n^m≡r_3(\mod p_3)$

$C_n^m≡r_k(\mod p_k)$

于是这个方程组又让我们想到了中国剩余定理（CRT），于是我们就可以使用Lucas定理把$r_1…r_k$求出来，再用中国剩余定理求解$C_n^m$。

代码

#include <iostream>
#include <cmath>
#include <cstring>
#include <cstdio>
#include <algorithm>
#include <vector>
using namespace std;
typedef long long ll;
#define MAXN 15
ll t;
ll n,m;
ll k;
ll a[MAXN];
ll r[MAXN];

//唯一分解定理将mod分解为素数
ll fj(ll n)
{
    int c=0;
    for(int i=2;i<=n;i++)
    {
        while(n%i==0)
        {
            a[c++]=i;
            n/=i;
        }
    }
    return c;
}

//快速乘取模
ll mul(ll a,ll b,ll mod)
{
    ll res=0;
    while(b)
    {
        if(b&1)res=(res+a)%mod;
        a=(a+a)%mod;
        b>>=1;
    }
    return res;
}
//计算n！%mod
ll fact(ll n,ll mod)
{
    ll res=1;
    for (int i = 1; i <= n; ++i) {
        res=(res*i)%mod;
    }
    return res;
}
void exgcd(ll a,ll b,ll &x,ll &y,ll &d)
{
    if(!b)
    {
        d=a;
        x=1;
        y=0;
    } else{
        exgcd(b,a%b,y,x,d);
        y-=x*(a/b);
    }
}
//用exgcd计算t模mod的逆元
ll inv(ll t,ll mod)
{
    ll d,x,y;
    exgcd(t,mod,x,y,d);
    return d==1?(x%mod+mod)%mod:-1;
}
//计算C(n,m)%mod:C(n,m)%mod=(n!/[m!*(n-m)!])%mod=n!%mod*inv(m!,mod)%mod*inv((n-m)!,mod)%mod
ll comb(ll n,ll m,ll mod)
{
    if(m<0||m>n)return 0;
    return fact(n,mod)*inv(fact(m,mod),mod)%mod*inv(fact(n-m,mod),mod)%mod;
}
//lucas定理:C(n,m)%mod=C(n/mod,m/mod)*C(n%mod,m%mod)%mod
ll lucas(ll n,ll m,ll mod)
{
    return m?lucas(n/mod,m/mod,mod)*comb(n%mod,m%mod,mod)%mod:1;
}

//中国剩余定理求解C(n,m)%mod=C(n,m)%(a1*a2*a3...)(ai为素数)
ll china(ll n,ll *d,ll *m)
{
    ll M = 1, ret = 0;
    for(int i = 0; i < n; i ++) M *= m[i];
    for(int i = 0; i < n; i ++)
    {
        ll w = M / m[i];
        ret = (ret + mul(w * inv(w, m[i]), d[i], M)) % M;
    }
    return (ret + M) % M;
}
int main()
{
    cin>>t;
    while(t--)
    {
        cin>>n>>m>>k;
        ll num=fj(k);//分解模k
        for (int i = 0; i <num; ++i) {
            r[i]=lucas(n,m,a[i]);
        }
        cout<<china(num,r,a)<<endl;
    }
    return 0;
}