题目分析：
我记得很久以前有人跟我说NOIP2016的题目出了加强版在清华集训中，但这似乎是一道无关的题目？
由于$k$为素数，那么$lucas$定理就可以搬上台面了。
注意到$inom{i}{j} equiv 0 {mod k}$当且仅当将$i$和$j$用$k$进制表示的时候，有一位上的$i<j$。
位数上的计算用数位DP就没错了。
代码：
1 #include<bits/stdc++.h> 2 using namespace std; 3 4 const int mod = 1000000007; 5 6 int t,k; 7 long long n,m; 8 9 int bn[90],n1,n2,bm[90],nw[90]; 10 11 int f[70][100][2]; //0 0~k-1 1 0~self 12 int sum[70][2],pw[70],dd[70],oo[70],fw[70],yw[70]; 13 14 void init(){ 15 if(n < m) m = n; 16 memset(bn,0,sizeof(bn)); memset(bm,0,sizeof(bm)); n1 = 0,n2 = 0; 17 long long p1 = n,p2 = m; 18 while(p1){bn[++n1] = p1 % k; p1 /= k;} 19 while(p2){bm[++n2] = p2 % k; p2 /= k;} 20 dd[0] = 1;oo[0] = 1; 21 for(int i=1;i<=n1;i++) dd[i] = (dd[i-1] + 1ll*bm[i]*pw[i-1]%mod)%mod; 22 for(int i=1;i<=n1;i++) oo[i] = (oo[i-1] + 1ll*bn[i]*pw[i-1]%mod)%mod; 23 } 24 25 pair<int,int> dfs3(int now){ 26 if(now == 0){return make_pair(0,0);} 27 int ans1 = 0,ans2 = 0; 28 pair<int,int> pt = dfs3(now-1); 29 int cutp = min(bn[now]+1,bm[now]); 30 for(int i=0;i<bn[now];i++){ 31 int cuep = min(i+1,bm[now]); 32 ans1 += (1ll*cuep*sum[now-1][0])%mod; ans1 %= mod; 33 ans1 += (1ll*(bm[now]-cuep)*pw[now-1])%mod*pw[now-1]%mod;ans1%=mod; 34 if(bm[now] > i){ans1 += (1ll*pw[now-1]*dd[now-1])%mod; ans1 %= mod;} 35 else{ans1 += nw[now-1];ans1 %= mod;} 36 ans2 += (1ll*(i+1)*sum[now-1][0])%mod; ans2 %= mod; 37 ans2 += (1ll*(k-i-1)*pw[now-1]%mod*pw[now-1])%mod; ans2 %= mod; 38 } 39 ans1 += (1ll*cutp*pt.second)%mod; ans1 %= mod; 40 ans1 += (1ll*(bm[now]-cutp)*oo[now-1]%mod*pw[now-1])%mod;ans1 %= mod; 41 if(bm[now] > bn[now]){ans1 += (1ll*oo[now-1]*dd[now-1])%mod;ans1%=mod;} 42 else{ans1 += pt.first;ans1 %= mod;} 43 ans2 += (1ll*(bn[now]+1)*pt.second)%mod; ans2 %= mod; 44 ans2 += (1ll*(k-bn[now]-1)*oo[now-1])%mod*pw[now-1]%mod;ans2 %= mod; 45 return make_pair(ans1,ans2); 46 } 47 48 void work(){ 49 memset(f,0,sizeof(f));memset(sum,0,sizeof(sum));memset(nw,0,sizeof(nw)); 50 for(int i=1;i<=n1;i++){ 51 for(int j=0;j<k;j++){ 52 f[i][j][1] = (1ll*j*sum[i-1][0]+sum[i-1][1]+f[i][j][1])%mod; 53 f[i][j][0] += (1ll*(j+1)*sum[i-1][0])%mod; f[i][j][0] %= mod; 54 f[i][j][0] += (1ll*(k-j-1)*((1ll*pw[i-1]*pw[i-1])%mod))%mod; 55 f[i][j][0] %= mod; 56 sum[i][0] += f[i][j][0]; sum[i][0] %= mod; 57 sum[i][1] += f[i][j][1]; sum[i][1] %= mod; 58 } 59 } 60 int ans = 0; 61 for(int now=1;now<=n1;now++){ 62 int ans1 = 0,ans2 = 0; 63 for(int i=0;i<bm[now];i++){ 64 ans1 = (ans1 + 1ll*sum[now-1][0]*(i+1))%mod; 65 ans1 +=(1ll*pw[now-1]*pw[now-1])%mod*(bm[now]-1-i)%mod;ans1%=mod; 66 ans1 += (1ll*pw[now-1]*dd[now-1])%mod; ans1 %= mod; 67 ans2 += (1ll*(i+1)*sum[now-1][0])%mod; ans2 %= mod; 68 ans2 += ((1ll*(k-i-1)*pw[now-1])%mod*pw[now-1])%mod; ans2 %= mod; 69 } 70 ans2 = (ans2+1ll*yw[now-1]*(bm[now]+1))%mod; 71 ans2+=((1ll*pw[now-1]*(k-bm[now]-1))%mod*dd[now-1])%mod;ans2%=mod; 72 ans1 = (1ll*yw[now-1]*bm[now]+fw[now-1]+ans1)%mod; 73 fw[now] = ans1; yw[now] = ans2; 74 } 75 for(int now=1;now<=n1;now++){ 76 for(int i=0;i<bm[now];i++){ 77 nw[now] += (1ll*pw[now-1]*i%mod*pw[now-1])%mod; nw[now] %= mod; 78 nw[now] = (nw[now]+1ll*(k-i)*sum[now-1][0])%mod; 79 } 80 nw[now] += 1ll*bm[now]*pw[now-1]%mod*dd[now-1]%mod; nw[now] %= mod; 81 nw[now] = (nw[now]+1ll*(k-bm[now])*nw[now-1])%mod; 82 } 83 ans += fw[n1];ans-=(m%mod*((m+1)%mod)/2ll)%mod; if(ans < 0) ans += mod; 84 pair<int,int> ans2 = dfs3(n1); 85 ans = ans + (ans2.first-fw[n1]); ans %= mod; ans += mod; ans %= mod; 86 printf("%d ",ans); 87 } 88 89 int main(){ 90 scanf("%d%d",&t,&k); 91 pw[0] = 1; for(int i=1;i<=65;i++) pw[i] = (1ll*pw[i-1]*k)%mod; 92 while(t--){ 93 scanf("%lld%lld",&n,&m); 94 init(); work(); 95 } 96 return 0; 97 }