Gauss Fibonacci

Time Limit: 1000/1000 MS (Java/Others)    Memory Limit: 32768/32768 K (Java/Others)
Total Submission(s): 437    Accepted Submission(s): 204

  
  
   
   2 1 4 100
2 0 4 100
  
  

  
  
   
   21
12
  
  
  
  
   
   /*******************************************************************************************
	首先我们看g=k*i+b; 是一个等差数列

	如果能推出f(g)这个函数也是一个等差或者等比数列，就可以得出一个公式

	f(k*i+b)  0<=i<n

	建立一个二分矩阵A=[1 1,1 0]
	
	f(b)=A^b  

	则:

	f(b)=A^b

	f(k+b)=A^k+b

	f(2*k+b)=A^2*k+b
	.
	.
	.
	f((n-1)*k+b)=A^(n-1)*k+b

	我们就可以得出一个等比数列:

	首项：A^b
	
	公比：A^k

	项数：n

	(res是单位矩阵)

	运用等比数列求和公式得出：sum=A^b*(res+A^k+(A^k)^2+(A^k)^3+...+(A^k)^(n-1))

	需要注意的一点是：当b=0的情况
***************************************************************************************/

/***************************************************************************************

第一次这种思路做测试数据都过了，可以TLE了~ 看来主要在
for(i=1;i<n;i++)
	temp=Martix_Add(temp,er_fun(q,i));	
这里超时了，要想办法优化下。

----------想了好久，决定用二分求和------------
例如:A+A^2+A^3+A^4+A^5+A^6=(A+A^2+A^3)+A^3*(A+A^2+A^3)  用递归写

***************************************************************************************/
#include <iostream>
using namespace std;
#define ll __int64
#define N 2
struct Mat
{
	ll martix[N][N];
};

Mat tp,res,A,q,a1,temp;
ll mod,k,b,n;

Mat Martix_Add(Mat a,Mat b)
{
	int i,j;
	Mat c;
	for (i=0;i<N;i++)
	{
		for (j=0;j<N;j++)
		{
			c.martix[i][j]=(a.martix[i][j]+b.martix[i][j])%mod;
		}
	}
	return c;
}
Mat Martix_Mul(Mat a,Mat b)
{
	int i,j,l;
	Mat c;
	for (i=0;i<N;i++)
	{
		for (j=0;j<N;j++)
		{
			c.martix[i][j]=0;
			for (l=0;l<N;l++)
			{
				c.martix[i][j]+=(a.martix[i][l]*b.martix[l][j])%mod;
				c.martix[i][j]%=mod;
			}
		}
	}
	return c;
}

Mat er_fun(Mat e,ll x)  //求矩阵e^x
{
	tp=e;
	e=res;
	while(x)
	{
		if(x&1)
			e=Martix_Mul(e,tp);
		tp=Martix_Mul(tp,tp);
		x>>=1;
	}
	return e;
}

Mat Solve(Mat a,ll p)
{
	if(p==1) return a;
	else if(p&1) return Martix_Add(er_fun(a,p),Solve(a,p-1));
	else return Martix_Mul(Solve(a,p>>1),Martix_Add(er_fun(a,p>>1),res));
}

int main()
{
	int i,j;
	for(i=0;i<N;i++)
	{
		for (j=0;j<N;j++)
		{
			if(i==j)
				res.martix[i][j]=1;
			else
				res.martix[i][j]=0;
		}
	}
	A.martix[0][0]=1;
	A.martix[0][1]=1;
	A.martix[1][0]=1;
	A.martix[1][1]=0;
	while (scanf("%I64d%I64d%I64d%I64d",&k,&b,&n,&mod)!=EOF)
	{
		q=A;
		q=er_fun(q,k);
		temp=res;
		/*for(i=1;i<n;i++)
			temp=Martix_Add(temp,er_fun(q,i));	*/
		temp=Martix_Add(temp,Solve(q,n-1));
		if(b!=0)
		{
			a1=er_fun(A,b);
			temp=Martix_Mul(a1,temp);
		}
		printf("%I64d/n",temp.martix[0][1]%mod);
	}
	return 0;
}