HDU4549M-斐波那契数列(矩阵快速幂，二分幂)

M斐波那契数列

Time Limit: 3000/1000 MS (Java/Others)    Memory Limit: 65535/32768 K (Java/Others)
Total Submission(s): 519    Accepted Submission(s): 145

Problem Description

M斐波那契数列F[n]是一种整数数列，它的定义如下：

F[0] = a
F[1] = b
F[n] = F[n-1] * F[n-2] ( n > 1 )

现在给出a, b, n，你能求出F[n]的值吗？

 

 

Input

输入包含多组测试数据；
每组数据占一行，包含3个整数a, b, n（ 0 <= a, b, n <= 10^9 ）

 

 

Output

对每组测试数据请输出一个整数F[n]，由于F[n]可能很大，你只需输出F[n]对1000000007取模后的值即可，每组数据输出一行。

 

 

Sample Input

0 1 06 10 2

 

 

Sample Output

060

分析：我们写出该数列的前几项：a, b, ab, ab^2, a^2*b^3, a^3*b^5。。。。我们可以发现，F(n)=a^f(n)*b^f(n-1)，这里f(n)是斐波那契数列（f(0)=0,f(1)=1,f(2)=1,f(n)=f(n-1)+f(n-2)）,所以，我们可以先用矩阵快速幂算出f(n)和f(n-1)，然后再用二分快速求幂计算出F(n)。但是，我们会发现一个问题，如果n很大的话，则f(n)太大，无法表示。所以我们想到了费马小定理，描述如下：如果a和p互素，那么a^(p-1)mod p=1。所以，我们要可以用费马小定理对a^f(n)中指数f(n)缩减为F(n)=a^(f(n)%M) * b^(f(n-1)% M)，其中M=p-1，在这里，M=1000000007-1。

#include<stdio.h>#include<string.h>#define M 1000000007void fun(__int64 a[][2],__int64 b[][2]){    __int64 c[2][2]={0},i,j,k;    for(i=0;i<2;i++){        for(j=0;j<2;j++){            for(k=0;k<2;k++){                c[i][j]+=a[i][k]*b[k][j];                c[i][j]%=(M-1);            }        }    }    for(i=0;i<2;i++){        for(j=0;j<2;j++){            a[i][j]=c[i][j];        }    }}__int64 getx(__int64 x, __int64 y){    __int64 sum=1;    while(y>0){        if(y&1)            sum=sum*x%M;        x=x*x%M;        y=y/2;    }    return sum;}int main(){    __int64 x,y,n,a[2][2],b[2][2],t,sum;    while(scanf("%I64d %I64d %I64d",&x,&y,&n)!=EOF){if(n==0){printf("%I64d\n",x%M);continue;}else if(n==1){printf("%I64d\n",y%M);continue;}        a[0][0]=1,a[0][1]=1,a[1][0]=1,a[1][1]=0;        b[0][0]=1,b[0][1]=0,b[1][0]=0,b[1][1]=1;        t=n;        while(t>0){            if(t&1)                fun(b,a);            fun(a,a);            t=t/2;        }        sum=getx(x,b[1][1])*getx(y,b[1][0])%M;        printf("%I64d\n",sum);    }    return 0;}

对于以上计算斐波那契的过程，我们可以如下结论：我们定义f(0)=a,f(1)=b,f(2)=a+b...f(n)=f(n-1)+f(n-2),这样，我么也可以通过矩阵快速幂迅速求得 f(n)矩阵连乘为 A*B^(n-1)，其中 A 为初始矩阵 {{b+a,b},{b,a}}, B={{1,1},{1,0}},通过快速幂即可求得。这很容易推导的出来，假设A={{f(n+1),f(n)},{f(n),f(n-1)}},则A*B={{f(n+1)+f(n),f(n+1)},{f(n+1),f(n)}}=F(n+1)，所以，当A初始值为{{b+a,b},{b,a}}时，我们可以用F(n)=A*B^(n-1)快速计算斐波那契的第n项。代码如下：

#include<stdio.h>void fun(__int64 b[][2],__int64 a[][2]){//矩阵连乘快速幂     __int64 c[2][2]={0},i,j,k;    for(i=0;i<2;i++){        for(j=0;j<2;j++){            for(k=0;k<2;k++){                c[i][j]=c[i][j]+b[i][k]*a[k][j];            }        }    }    for(i=0;i<2;i++){        for(j=0;j<2;j++){            b[i][j]=c[i][j];        }    }}int main(){    __int64 a[2][2],b[2][2]={{1,1},{1,0}},c[2][2]={{1,0},{0,1}},i,j,x,y,n,t;    while(scanf("%I64d %I64d %I64d",&x,&y,&n)!=EOF){        a[0][0]=x+y,a[0][1]=y,a[1][0]=y,a[1][1]=x;        b[0][0]=1,b[0][1]=1,b[1][0]=1,b[1][1]=0;        c[0][0]=1,c[0][1]=0,c[1][0]=0,c[1][1]=1;        t=n-1;        while(t>0){            if(t%2==1){                fun(c,b);            }            fun(b,b);            t=t/2;        }        fun(a,c);        printf("%I64d\n",a[1][0]);    }}