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

012

 

Sample Output

0142837

 

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，找出循环节183120int 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初始化为单位矩阵，人为输入originvoid 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;}