http://acm.hdu.edu.cn/showproblem.php?pid=4549
f[0] = a^1*b^0%p,f[1] = a^0*b^1%p,f[2] = a^1*b^1%p.....f[n] = a^fib[n-1] * b^fib[n-2]%p。
这里p是质数,且a,p互素,那么我们求a^b%p,当b非常大时要对b降幂。
由于a,p互素,那么由费马小定理知a^(p-1)%p = 1。令b = k*(p-1) + b'。a^b%p = a^(k*(p-1)+b')%p = a^(k*(p-1))%p * a^b'%p
= a^b'%p。
当中p' = b%(p-1)。
那么上题的fib[n-1]和fib[n-2]能够经过矩阵高速幂对(p-1)取模,实现降幂以后再用高速幂求得。
事实上a^b%c较一般的形式是a^(b%phi[c]+phi[c])%c (b>=phi[c]),在c是素数或a与c互素的时候能够简化为a^(b%(c-1))%c。
纳闷了一中午,结果把fib数列弄错了。f[0] = 0, f[1] = 1,f[n] = f[n-1]+f[n-2](n >= 2)。第0项是0,sad....
http://www.narutoacm.com/archives/a-pow-b-mod-m/
对大神表示深深的ym.
#include <stdio.h> #include <iostream> #include <map> #include <set> #include <list> #include <stack> #include <vector> #include <math.h> #include <string.h> #include <queue> #include <string> #include <stdlib.h> #include <algorithm> #define LL long long #define _LL __int64 #define eps 1e-12 #define PI acos(-1.0) #define C 240 #define S 20 using namespace std; const int maxn = 110; const int mod = 1000000007; struct matrix { LL mat[3][3]; void init() { memset(mat,0,sizeof(mat)); for(int i = 0; i < 2; i++) mat[i][i] = 1; } }m; LL a,b,n; matrix mul_matrix(matrix x, matrix y) { matrix ans; memset(ans.mat,0,sizeof(ans.mat)); for(int i = 0; i < 2; i++) { for(int k = 0; k < 2; k++) { if(x.mat[i][k] == 0) continue; for(int j = 0; j < 2; j++) ans.mat[i][j] = (ans.mat[i][j] + x.mat[i][k]*y.mat[k][j])%(mod-1); } } return ans; } matrix pow_matrix(matrix m, LL n) { matrix ans; ans.init(); while(n) { if(n&1) ans = mul_matrix(ans,m); m = mul_matrix(m,m); n >>= 1; } return ans; } LL pow_mod(LL a, LL b) { LL res = 1; a = a%mod; while(b) { if(b&1) res = res*a%mod; a = a*a%mod; b >>= 1; } return res; } int main() { while(~scanf("%lld %lld %lld",&a,&b,&n)) { m.mat[0][0] = m.mat[0][1] = m.mat[1][0] = 1; m.mat[1][1] = 0; if(n == 0) { printf("%lld ",a%mod); continue; } if(n == 1) { printf("%lld ",b%mod); continue; } else if(n == 2) { printf("%lld ",a*b%mod); continue; } matrix ans = pow_matrix(m,n); LL res = pow_mod(a,ans.mat[1][1]) * pow_mod(b,ans.mat[0][1]) %mod; printf("%lld ",res); } return 0; }