【洛谷3746】【六省联考2017】组合数问题

题目描述

组合数 CnmC_n^mCnm 表示的是从 n 个互不相同的物品中选出 m 个物品的方案数。举个例子,从 (1;2;3) 三个物品中选择两个物品可以有 (1;2);(1;3);(2;3) 这三种选择方法。根据组合数的定义,我们可以给出计算组合数 CnmC_n^mCnm 的一般公式:

Cnm=n!m!(n−m)!C_n^m = \frac{n!}{m!(n-m)!}Cnm=m!(nm)!n!

其中 n! = 1 × 2 × · · · × n。(特别的,当 n = 0 时, n! = 1 ,当 m > n 时, Cnm=0C_n^m =0Cnm=0

小葱在 NOIP 的时候学习了 CijC_i^jCij 和 k 的倍数关系,现在他想更进一步,研究更多关于组合数的性质。小葱发现, CijC_i^jCij 是否是 k 的倍数,取决于 CijmodkC_i^j mod k Cijmodk是否等于 0,这个神奇的性质引发了小葱对 mod 运算(取余数)的兴趣。现在小葱选择了是四个整数n; p; k; r,小葱现在希望知道

∑i=0infCnkik+rmodp\sum_{i=0}^{\inf} C_{nk}^{ik+r} mod pi=0infCnkik+rmodp的值。

输入输出格式

输入格式:

第一行有四个整数 n; p; k;r,所有整数含义见问题描述。

输出格式:

一行一个整数代表答案。

输入输出样例

输入样例#1:
2 10007 2 0
输出样例#1:
8
输入样例#2:
20 10007 20 0
输出样例#2:
176

说明

• 对于 30% 的测试点, 1 ≤ n; k ≤ 30, p 是质数;

• 对于另外 5% 的测试点, p = 2;

• 对于另外 5% 的测试点, k = 1;

• 对于另外 10% 的测试点, k = 2;

• 对于另外 15% 的测试点, 1 ≤ n ≤ 10^3; 1 ≤ k ≤ 50, p 是质数;

• 对于另外 15% 的测试点, 1 ≤ n × k ≤ 10^6, p 是质数;

• 对于另外 10% 的测试点, 1 ≤ n ≤ 10^9; 1 ≤ k ≤ 50, p 是质数;

• 对于 100% 的测试点, 1 ≤ n ≤ 10^9; 0 ≤ r < k ≤ 50; 2 ≤ p ≤ 2^30 − 1。

题解

一道矩阵加速的题。

对组合数的一个变形。原组合数公式是f[i][j]=f[i-1][j]+f[i-1][j-1]。由于这题是i*k+r,所以设f[i][j]表示i个物品中选x个物品,x满足x%k=j,所以f[i][j]=f[i-1][j]+f[i-1][(j-1+k)%k]。为什么呢??由原组合数公式可知,f[i][a*k+j]=f[i-1][a*k+j]+f[i-1][a*k+j-1]。所以sigmaf[i][a*k+j]=sigmaf[i-1][a*k+j]+sigmaf[i-1][a*k+j]。。因为设g[i][j]=sigma(a从0~一个值)f[a*k+j];(这里g[i][j]含义同上面f[i][j])。由于j-1会出现负的,也就是当j=0时,表示整除,由杨辉三角可知,前一个是k-1,所以f[i][j]=f[i-1][j]+f[i][(j-1+k)%k]。

递推式证完了,下面该说怎么推矩阵了。

像这种二维的递推,第一行矩阵一般把第二维都写出来再根据下一个推矩阵。就是由f[i][0],f[i][1],f[i][2].....f[i][k-1]看怎样推成f[i+1][0],f[i+1][1],f[i+1][2]......f[i+1][k-1],自己手推一下矩阵吧。最后结果就是f[0][r].因为第一行是f[i][0],f[i][1]....

下面就是矩阵加速了。。。不再赘述。

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define ll long long
using namespace std; 
const int maxn=60;
ll n,p,k,r;
struct node{
	ll a[maxn][maxn];
}a,b,c,ans;
node work(node a,node b){
	memset(c.a,0,sizeof(c.a));
	
	for(int i=0;i<k;i++)
	for(int j=0;j<k;j++)
	for(int q=0;q<k;q++)
	c.a[i][j]=(c.a[i][j]+a.a[i][q]*b.a[q][j])%p;
	
	return c; 
}
void qpow(node a,ll b){
	while(b){
		if(b&1) ans=work(ans,a);
		a=work(a,a);
		b>>=1;
	}
}
int main(){
	cin>>n>>p>>k>>r;
	for(int i=0;i<k;i++){
		a.a[i][i]++;
		a.a[(i-1+k)%k][i]++;
	}
	ans=a;
	qpow(a,(ll)(n*k-1));
	printf("%lld\n",ans.a[0][r]);

	return 0;
}


评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值