1436 方程的解数 斐波那契数列+快速幂+矩阵快速幂（好题）

1436 方程的解数

基准时间限制：1 秒 空间限制：131072 KB 分值: 80 难度：5级算法题

请计算这个方程组有多少合法的整数解，答案比较大，对m取余后输出。

对于样例，有三组解{1, 1}, {3, 1}, {1,3}。

Input

单组测试数据。

第一行包含四个整数 n, k, l, m (2 ≤ n ≤ 10^18, 0 ≤ k ≤ 10^18, 0 ≤ l ≤ 64, 1 ≤ m ≤ 10^9 + 7)。

Output

对于每一组数据输出答案占一行。

Input示例

2 1 2 10

Output示例

3

题解；公式中只有二进制运算，所以考虑将k拆成二进制数来一位一位的分析，对于k在二进制下的每一位，有且仅有两种状态1,0,；由n个ai组成每一位的方案数是2^n种，因为组成1有很多种方式，但组成0当且仅当n个ai中不存在两个连续的数这一位为1，假设其方案数为x，那么组成1的方案数有2^n-x种，因为题目中还有一个变量L，所以在二进制下实际上是L位（补位0），那么L位中，k有p位为1，q位为0，那么总方案数为q*x*p*(2^n-x)种。分析到这里，就很显然了，问题的关键就是求x，对于n个数，如果第j个数，设方案为xj，如果当前一位为0，那么前一位0和1都可以（xj-1），如果当前这一位为1，那么第j-1个数的这一位只能为0，为第j-2位0和1都可以（xj-2），所以xj=xj-1+xj-2，这就是斐波那契数列，用快速矩阵幂求出xn，最终答案为xn^q*(2^n-xn)^p。（注意此题的特判（巨坑）QAQ）

总结：以后拿到公式题一定不要恐惧，如果和位运算有关，先找性质，比如这道题中的性质就是当k的某一位为0 时，n个数中必定不存在连续的两个数这一位为1；而且，和位运算有关的题目最好是将其中的数字拆成二进制数分析（毕竟位运算就是根据在二进制数下运行的嘛）……

#include <iostream>
#include <cstring>
#include <algorithm>
#include <cstdio>
#define ULL unsigned long long
using namespace std;
struct node
{
	long long a[3][3];
};
int m;
long long n,k,mod;
inline node fast(node re,node opt)
{
	node ans;
	for(int i=1;i<=2;i++)
	    for(int j=1;j<=2;j++)
	    {
	    	ans.a[i][j]=0;
	        for(int k=1;k<=2;k++)
	            ans.a[i][j]=(ans.a[i][j]+re.a[i][k]*opt.a[k][j]%mod)%mod;	    	
		}
	return ans;
}
inline node del(node tmp,long long p)
{
	node ans;
	ans.a[1][1]=ans.a[2][2]=1;ans.a[1][2]=ans.a[2][1]=0;
	while(p)
	{
		if(p&1) ans=fast(ans,tmp);
		tmp=fast(tmp,tmp);
		p>>=1;
	}
	return ans;
}
inline long long prim(long long a,long long b)
{
	long long ans=1;
	while(b)
	{
		if(b&1) ans=(ans*a)%mod;
		a=(a*a)%mod;
		b>>=1;
	}
	return ans;
}
int main()
{
	node tmp;
	scanf("%lld%lld%d%lld",&n,&k,&m,&mod);
	unsigned long long t=1ULL<<m;
	if(mod==1 || ( k>=t && m!=64))//!!!
	{
		printf("0\n");return 0;
	}
	tmp.a[1][1]=tmp.a[1][2]=tmp.a[2][1]=1;tmp.a[2][2]=0;
	tmp=del(tmp,n);
	long long x=(tmp.a[1][1]+tmp.a[2][1])%mod;
	long long y=(prim(2,n)%mod-x%mod+mod)%mod;
	int num1=0,num0=0;
	while(k)
	{
		if(k&1) num1++;
		else num0++;
		k>>=1;
	}
	num0=m-num1;
	long long ans=( prim(x,num0)%mod * prim(y,num1)%mod)%mod;
	cout<<ans<<endl;
	return 0;
}