当前位置: 代码迷 >> 综合 >> 【组合数学】【Lucas】HDU6372 sacul
  详细解决方案

【组合数学】【Lucas】HDU6372 sacul

热度:19   发布时间:2023-09-27 06:48:01.0

分析:

原来倒着读题目就是题解啊。。。。出题人很皮。。。

恩。。lucas定理在组合数学中还算是比较常见的了,之前在CQOI2018就遇到一道可以用Lucas骗分的(运气好还能卡过)的题。。。

其实lucas定理本身还是比较简单的(看着隔壁的EXLucas瑟瑟发抖。。。)
C(n,m)%p=(C(np,mp)?C(n%p,m%p))%pC(n,m)%p=(C(np,mp)?C(n%p,m%p))%p

具体的证明请自行百度。。。(或者跟我一样不会证明直接用)

观察这个式子,发现其实就是把n和m按p进制分解,每一位分别求组合数,再乘起来。
而这个式子为0的条件也很容易发现

由于组合数满足:如果n<mC(n,m)=0n<m,则C(n,m)=0
所以n和m在p进制下,只要任意一位m比n大,那么值就为0.

现在看回这道题,海绵宝宝矩阵中,只有当C(n,m)>0C(n,m)>0时,这个位置才是1,否则就是0.

可以把这个矩阵看做一个图,若C(i,j)=1C(i,j)=1则说明从到j有一条单向边。换言之,可以将这个矩阵的每一位(i,j)(i,j)看作从ii走一步到 j 的方案数。那么这个矩阵的k次幂就表示:从ii走k步到 j 的方案数。。。。

(卧槽这不就是之前那场的HDU6331 Walking Plan么。。。。)

然后题目非常恶心,它求的其实就是这个图中,从i出发,到达j的步数不超过k的方案数。

现在考虑怎么求答案
(i,j)(i,j)存在的条件是,C(i,j)%p>0C(i,j)%p>0,根据lucas定理,可以按p进制下每一位独立来看。因为若从i走到j,必须每一位都满足i>ji>j,它只与自己这一位有关,无法限制其他位的情况。因此每一位是互相独立的。所以只要知道某一位的方案数,然后利用乘法原理再组合起来,形成一个幂的形式,就是整个的方案数。

现在分析某一位的方案怎么求:
假设要求走x步的方案数:
其实就是[0,p-1]这些数组成长度为x+1的不上升序列的方案数。。。(因为要求每一步都不能变大)

可以利用隔板法,把[0,p-1]中的每个数看作一个空隙,在这些空隙中放x+1个隔板,隔板可以重合。每种隔板的摆放方案,就对应了一个不上升序列。(每个隔板到0的距离就表示这一位的值)

所以方案数就是C(p+x+1?1,x+1)=C(p+x,x+1)=C(p+x,p?1)C(p+x+1?1,x+1)=C(p+x,x+1)=C(p+x,p?1)

总的答案就是iki=1jnj=1Cj(p+i,p?1)∑i=1i≤k∑j=1j≤nCj(p+i,p?1)
后面那部分可以用等比数列求和公式算出来。前面那部分枚举就好了。

#include<cstdio>
#include<cstring>
#include<algorithm>
#include<vector>
#include<cmath>
#define SF scanf
#define PF printf
#define EPS 1e-6
#define MAXN 200010
#define MAXM 1300000
#define MOD 1000000007
using namespace std;
typedef long long ll;
int t,n;
bool isprime[MAXM];
int primes[MAXN],cnt;
void linearsieve(){for(int i=2;i<=MAXM&&cnt<100000;i++){if(isprime[i]==0)primes[++cnt]=i;    for(int j=1;i*primes[j]<=MAXM;j++){int x=i*primes[j];isprime[x]=1;if(i%primes[j]==0)break;  }}
}
ll fac[MAXM+MAXN],ifac[MAXM+MAXN];
ll fsp(ll x,int y){ll res=1;while(y){if(y&1)res=res*x%MOD;x=x*x%MOD;y>>=1;  }return res;
}
void prepare(){fac[0]=1;for(int i=1;i<=1500000;i++)fac[i]=fac[i-1]*i%MOD;ifac[1500000]=fsp(fac[1500000],MOD-2);for(int i=1500000;i>0;i--)ifac[i-1]=ifac[i]*i%MOD;
}
ll C(int n,int m){return fac[n]*ifac[m]%MOD*ifac[n-m]%MOD;    
}
ll solve(ll x,int n){return (fsp(x,n)-1+MOD)%MOD*(fsp(x-1,MOD-2))%MOD*x%MOD;
}
int main(){linearsieve();prepare();SF("%d",&t);int c,n,k;while(t--){SF("%d%d%d",&c,&n,&k);ll ans=0;ll p=primes[c];for(int i=1;i<=k;i++){ll x=C(p+i,p-1);if(x==1)ans+=n;elseans+=solve(x,n);ans%=MOD;}PF("%lld\n",ans);}
}