hdu4291(循环节+快速幂)

A Short problem

Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 32768/32768 K (Java/Others)
Total Submission(s): 1154    Accepted Submission(s): 449

Problem Description

　　According to a research, VIM users tend to have shorter fingers, compared with Emacs users.
　　Hence they prefer problems short, too. Here is a short one:
　　Given n (1 <= n <= 10 ¹⁸), You should solve for
g(g(g(n))) mod 10 ⁹ + 7
　　where
g(n) = 3g(n - 1) + g(n - 2)
g(1) = 1
g(0) = 0

 

Input

　　There are several test cases. For each test case there is an integer n in a single line.
　　Please process until EOF (End Of File).

 

Output

　　For each test case, please print a single line with a integer, the corresponding answer to this case.

 

Sample Input

  
  

   
   0
1
2
  
  

 

Sample Output

  
  

   
   0
1
42837
  
  

 

Source

2012 ACM/ICPC Asia Regional Chengdu Online

 

Recommend

liuyiding

 

本题题意很简单明确，给定n，要求满足条件的g(g(g(n))) mod 10⁹ + 7

                                                     g(n) = 3g(n - 1) + g(n - 2)  n>=2其中g(1) = 1，g(0) = 0
题目的n很大，常规递推肯定超时。由递推式的形式可以联想到矩阵的快速幂来加速递推公式的求解。

g(n+1) = 3           1       *      g(n)     

g(n)         1            0            g(n-1)
即便采用矩阵的快速幂，在n很大且中间过程的g(n)可能很大，大到无法表示，此题目还是无法解决。

mod 10⁹ + 7，在我看来任何取模的都能找到循环节，本题首先可以找到内层的循环节。此处找循环节是参考大牛的代码http://www.cnblogs.com/kuangbin/archive/2012/09/17/2688852.html

#include<stdio.h>

const long long MOD=222222224;//第一次是MOD=1000000007 找出循环节是222222224
//第二次是MOD=222222224，找出循环节183120
int main()
{
    long long a,b;
    a=1;
    b=3;
    for(int i=1;;i++)
    {
        if(a==0&&b==1)
        {
            printf("%d\n",i);
            break;
        }
        long long c=3*b+a;
        c%=MOD;
        a=b;
        b=c;
    }
    return 0;
}

相当于预处理找到循环节，在解题。

第一次是MOD=1000000007 找出循环节是222222224
第二次是MOD=222222224，找出循环节183120

找循环节暴力找就可以了，矩阵乘法找循环节更加慢。

 

#include<iostream>
#include<cstdio>
#include<cstring>
using namespace std;

const __int64 mod1=183120,mod2=222222224,mod3=1000000007;
//*****************************************************************
//origin存放需计算的矩阵，res存放答案矩阵
struct matrix
{
	__int64 a[2][2];
};

//直接将2个矩阵相乘x*y,返回计算后的矩阵
matrix multiply(matrix &x,matrix &y,__int64 MOD)
{
	matrix temp;
	memset(temp.a,0,sizeof(temp.a));
	for(int i=0;i<2;i++)
	{
		for(int j=0;j<2;j++)
		{
			for(int k=0;k<2;k++)
			{
				temp.a[i][j]+=y.a[i][k]*x.a[k][j];
				temp.a[i][j]%=MOD;
			}
		}
	}
	return temp;
}

//将res初始化为单位矩阵，人为输入origin
void init(matrix &origin,matrix &res)
{
	origin.a[0][0]=3,origin.a[1][0]=origin.a[0][1]=1,origin.a[1][1]=0;
	res.a[0][0]=1,res.a[0][1]=res.a[1][0]=res.a[1][1]=0;

//	res.a[0][0]=3,res.a[1][0]=res.a[0][1]=1,res.a[1][1]=0;
//	origin.a[0][0]=1,origin.a[0][1]=origin.a[1][0]=origin.a[1][1]=0;
}

//矩阵快速幂的计算
void calc(matrix &origin,matrix &res,__int64 n,__int64 MOD)
{
	while(n)
	{
		if(n&1)
			res=multiply(res,origin,MOD);
		n>>=1;
		origin=multiply(origin,origin,MOD);
	}
}
//*****************************************************************


int main()
{
	__int64 n;
	matrix origin,res;
	while(~scanf("%I64d",&n))
	{
		init(origin,res);
		calc(origin,res,n,mod1);
		__int64 tmp=res.a[1][0];

		init(origin,res);
		//origin.a[0][0]=3,origin.a[1][0]=origin.a[0][1]=1,origin.a[1][1]=0;
		calc(origin,res,tmp,mod2);
		tmp=res.a[1][0];

		init(origin,res);
		//origin.a[0][0]=3,origin.a[1][0]=origin.a[0][1]=1,origin.a[1][1]=0;
		calc(origin,res,tmp,mod3);

		printf("%I64d\n",res.a[1][0]);

	}
	return 0;
}