hdu 5667 Sequence 矩阵快速幂

                    1          n=1

f（n）={     a^b     n=2

                    a^b*f（n-1）^c *f（n-2）  其他

给了你几个数:n,a,b,c,你需要输出f（n）模p后的数值 ,p是质数

思路：列几项后发现f（n） =（a^b）^指数，（a^b）是常数，对f（n）取

以（a^b）为底 的对数后，令F（n）=logf（n），F（n）=c*F（n-1）+F（n-2）n-2+1。

类比斐波那契数列，可以快速求出F（n），             F（n+2）         c    1      1                F（n+1）

                                                                             F（n+1） =      1     0      0      *        F（  n  ）

                                                                             1                       0     0      1                1

因为F（n）可能很大所以用费马小定理x^(p-1)%p=1  （x与p互质），既算F（n）的过程中对p-1

取模。最后用快速幂算a^(  b*F( n )  ) %p  。然后就愉快的WA了 QAQ。

理由很简单，用费马小定理条件为a与p互质，但是a与p不一定互质。因为p为质数所以a如果与p不互质

a一定为p的倍数，此时a的任意次幂对p取模都为0，直接输出答案。a如果与p互质则可以用费马小定理。

所以加个特判就过了。感觉把这道题做复杂了 23333333333 。

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cmath>
#include <cstring>
#include <cstdlib>
#include <string>
#include <queue>
#include <map>
#include <vector>
using namespace std;
const int N=100005;
//const int mod=1e9+7;
typedef long long ll;



typedef vector<ll > vec;
typedef vector<vec > mat;

 ll p;
ll n,a,b,c;



mat mul(mat &A,mat &B){
    mat C(A.size(), vec( B[0].size() ) );
    for(int i=0;i<A.size();i++){
        for(int k=0;k<B.size();k++){
            for(int j=0;j<B[0].size();j++){
                C[i][j]=( C[i][j]+A[i][k]*B[k][j] )%(p-1);
            }
        }
    }
    return C;
}

mat pow(mat A,ll n){
    mat B( A.size(),vec(A.size() ,0 ) );
    for(int i=0;i<A.size();i++){
        B[i][i]=1;
    }
    while(n>0){
        if(n&1)B=mul(B,A);
        A=mul(A,A);
        n>>=1;
    }
    return B;
}


ll  mod_pow(ll x, ll n1, ll mod){
    ll res=1;
    while(n1>0){
        if(n1&1)
        res=res*x%mod;
        x=x*x%mod;
        n1>>=1;
    }
    return res;
}




int main()
{
    int  t;
    scanf("%d",&t);
    while(t--){
        scanf("%I64d%I64d%I64d%I64d%I64d",&n,&a,&b,&c,&p);
        if( a%p==0 ){printf("0\n");continue;}
        mat A( 3,vec( 3,0 ) );
        A[0][0]=c; A[0][1]=1;A[0][2]=1;
        A[1][0]=1; A[1][1]=0;A[1][2]=0;
        A[2][0]=0; A[2][1]=0;A[2][2]=1;
        A=pow(A,n-1);
        ll logfn=A[1][0]+A[1][2];
        ll tmp=( b*logfn )%(p-1);
        printf("%I64d\n",  mod_pow(a,tmp,p)  );
    }
	return 0;
}