51nod 1537 分解(矩阵快速幂)

分析:先写出前几项,发现都是有解的.记(1+√2)^n=a+b√2,可以归纳证明,当n为奇数时,m=a^2+1,n为偶数时,m=a^2.写出a的递推式,用矩阵快速幂算一下a即可.

#include<iostream>
#include<cstring>
#include<cstdio>
using namespace std;
const int p=1e9+7;
typedef unsigned long long ull;
ull mat[2][2]={2,1,1,0};
void Copy(ull a[2][2],ull b[2][2]){
    for(int i=0;i<2;i++)
        for(int j=0;j<2;j++)a[i][j]=b[i][j];
}
void multiple(ull result[2][2],ull a0[2][2],ull b0[2][2]){
    ull a[2][2],b[2][2];
    Copy(a,a0);Copy(b,b0);
    for(int i=0;i<2;i++){
        for(int j=0;j<2;j++){
            result[i][j]=0;
            for(int k=0;k<2;k++)
                result[i][j]=(result[i][j]+a[i][k]*b[k][j])%p;
        }
    }
}
void qpow(ull t,ull ans[2][2]){
    ull y[2][2];
    Copy(y,mat);
    ans[0][0]=ans[1][1]=1;
    ans[0][1]=ans[1][0]=0;
    ull k=1;
    while(k){
        if(k&t){
            multiple(ans,ans,y);
        }
        multiple(y,y,y);
        k<<=1;
    }
}
int main(){
    ull n;
    cin>>n;
    if(n==0)cout<<1<<endl;
    else if(n==1)cout<<2<<endl;
    else{
        ull a[2][2];
        qpow(n-1,a);
        ull t=((a[0][0]+a[1][0])*(a[0][0]+a[1][0]))%p;
        if(n%2)cout<<t+1<<endl;
        else cout<<t<<endl;
    }
    return 0;
}