Lucas定理 数位DP

【清华集训 2016】组合数问题

Posted by WHZ0325 on 2019-04-26, Viewed times

题目描述

$t$ 组数据，每组数据给定 $n,m$ 和 $k$，求对于所有的 $0\le i\le n$，$0\le j\le \min(i,m)$ 有多少对 $(i,j)$ 满足 $C_i^j$ 是 $k$ 的倍数。

$1\le n,m\le 10^{18}$，$1\le t,k\le 100$，且 $k$ 是一个质数。

算法分析

问题所求即有多少对 $(i,j)$ 满足 $C_i^j%k=0$，观察到 $k$ 比较小且为质数，可以使用 Lucas 定理求解。

要使结果为 $0$，Lucas 定理的过程中一定有一部分乘积为 $0$，组合数 $C_n^k$ 当且仅当 $k\gt n$ 时值为 $0$，考虑数位 DP，直接求解比较复杂，将至少有一部分乘积为 $0$ 转化为总的方案数减去没有任何一部分为 $0$ 的方案数。

先考虑总的方案数如何计算，分成 $n\le m$ 和 $n\gt m$ 两部分，前半部分显然是 $1+2+\dots +(m+1)=\frac{(m+1)\times(m+2)}{2}$，后半部分显然每个 $n$ 有 $[0,m]$ 共 $m+1$ 种 $i$ 可用，共有 $n-m$ 个大于 $m$ 的 $n$，因此总的方案数为 $\frac{(m+1)(m+2)}{2}+(m+1)(n-m)$。

当 $k=1$ 时答案就是总的方案数，否则 Lucas 定理的过程中最多有 $log_2(10^{18})=60$ 位，则没有任何一部分为 $0$ 的方案数可以设 $f[i][0/1][0/1]$ 表示从左到右前 $i$ 位，$n$ 目前是/否达到上界以及 $m$ 目前是/否达到上界时没有任何一部分为 $0$ 的方案数，设计转移即可。

时间复杂度为 $O(60t)$。

注意由于每一部分的上界是一个余数，因此存在 $k\le n$ 但 $k$ 到达上界而 $n$ 没有的情况，这表明 $f[i][0][1]$ 是有意义的。

代码实现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
#include <cstdio>
#include <algorithm>
typedef long long int ll;
const int mod=(int)1e9+7;
inline int inc(int x,int y) {x+=y;return x>=mod?x-mod:x;}
inline int dec(int x,int y) {x-=y;return x<0?x+mod:x;}
inline int mul(int x,int y) {return (ll)x*y%mod;}
inline int qpow(int n,int k) {
    int ans=1;
    while(k) {
        if(k&1) ans=mul(ans,n);
        n=mul(n,n);k>>=1;
    }
    return ans;
}
inline int calc(int n,int m) {
    m=std::min(n,m);
    return inc(mul(mul(m+1,m+2),qpow(2,mod-2)),(m+1<=n)?mul(m+1,n-m):0);
}
int upn[65],upm[65],cnt=0,f[65][2][2];//当前从高到低处理到第i位，n是否达到上界0/1，m是否达到上界0/1
int main() {
    int t,k;scanf("%d%d",&t,&k);
    while(t--) {
        ll n,m;scanf("%lld%lld",&n,&m);m=std::min(n,m);//保证m不超过n
        int ans=inc(mul(mul((m+1)%mod,(m+2)%mod),qpow(2,mod-2)),(m+1<=n)?mul((m+1)%mod,(n-m)%mod):0);
        if(k==1) {printf("%d\n",ans);continue;}
        cnt=0;while(n||m) {++cnt;upn[cnt]=n%k;n/=k;upm[cnt]=m%k;m/=k;}
        std::reverse(upn+1,upn+1+cnt);std::reverse(upm+1,upm+1+cnt);
        memset(f,0,sizeof(f));f[0][1][1]=1;
        for(register int i=1;i<=cnt;++i) {
            f[i][0][0]=inc(f[i][0][0],mul(f[i-1][0][0],mul(mul(k,k+1),qpow(2,mod-2))));//0~k-1 0~k-1
            f[i][0][0]=inc(f[i][0][0],mul(f[i-1][0][1],calc(k-1,upm[i]-1)));//0~k-1 0~upm-1
            f[i][0][0]=inc(f[i][0][0],mul(f[i-1][1][0],calc(upn[i]-1,k-1)));//0~upn-1 0~k-1
            f[i][0][0]=inc(f[i][0][0],mul(f[i-1][1][1],calc(upn[i]-1,upm[i]-1)));//0~upn-1 0~upm-1
            f[i][0][1]=inc(f[i][0][1],mul(f[i-1][0][1],k-upm[i]));//0~k-1 upm
            f[i][0][1]=inc(f[i][0][1],mul(f[i-1][1][1],upn[i]>=upm[i]?upn[i]-upm[i]:0));//0~upn-1 upm
            f[i][1][0]=inc(f[i][1][0],mul(f[i-1][1][0],std::min(upn[i]+1,k)));//upn 0~k-1
            f[i][1][0]=inc(f[i][1][0],mul(f[i-1][1][1],std::min(upn[i]+1,upm[i])));//upn 0~upm-1
            f[i][1][1]=inc(f[i][1][1],f[i-1][1][1]*(upm[i]<=upn[i]));//upn upm
        }
        ans=dec(ans,inc(inc(inc(f[cnt][0][0],f[cnt][0][1]),f[cnt][1][0]),f[cnt][1][1]));
        printf("%d\n",ans);
    }
    return 0;
}