HDU 3221 Brute-force Algorithm(指数降幂公式)

题目链接；
HDU 3221 Brute-force Algorithm
题意：
根据递归可以得到$f[1]=a, f[2]=b,f[n]=f[n-1]*f[n-2](n\geq 3)$,给出$a,b,p,n$求$f[n]\% p的值。$
数据范围：$1≤n≤1000000000, 1≤P≤1000000, 0≤a, b<1000000.$
分析：
设菲波那切数列第$n$项为$g[n],g[1]=0,g[2]=1,g[3]=1.$根据递推可以得到:
$f[n]=a^{g(n-1)}*b^{g(n)} \ \% p\ \ (n\geq 3)$
但是考虑到指数$g(n)$会很大，我们需要降幂。
指数降幂公式:

$$A^{x}\ \% \ p=A^{x\%\phi(p)+\phi(p)}\ \% \ p\ （x\geq \phi(p),\phi(p)是p的欧拉函数值）$$
所以用矩阵快速幂求菲波那切数列时的取模数是 $\phi (p)$而不是 $p$。

下面的代码用$G++$提交能 $AC$，用 $C++$就不能。。。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <string>
#include <algorithm>
#include <climits>
#include <cmath>
#include <ctime>
#include <cassert>
#include <bitset>
#define IOS ios_base::sync_with_stdio(0); cin.tie(0);
using namespace std;
typedef long long ll;
const int SIZE = 10;
const int MAX_N = 1000010;bitset<MAX_N> bs;
int prime_cnt, prime[MAX_N], phi[MAX_N];struct Matrix {int row, col;ll data[SIZE][SIZE];Matrix () {}Matrix (int _row, int _col): row(_row), col(_col) {}
};void GetPhi()
{bs.set();prime_cnt = 0;phi[1] = 1;for(int i = 2; i < MAX_N; ++i) {if(bs[i] == 1) {prime[prime_cnt++] = i;phi[i] = i - 1;}for(int j = 0; j < prime_cnt && i * prime[j] < MAX_N; ++j) {bs[i * prime[j]] = 0;if(i % prime[j]) {phi[i * prime[j]] = (prime[j] - 1) * phi[i];}else {phi[i * prime[j]] = prime[j] * phi[i];break;}}}
}Matrix matrix_mul(Matrix a, Matrix b, ll p)
{Matrix res;memset(res.data, 0, sizeof(res.data));res.row = a.row, res.col = b.col;for(int i = 1; i <= res.row; ++i) {for(int j = 1; j <= res.col; ++j) {for(int k = 1; k <= a.col; ++k) {res.data[i][j] += a.data[i][k] * b.data[k][j];if(res.data[i][j] > p) res.data[i][j] = (res.data[i][j] % p + p);}}}return res;
}Matrix matrix_quick_pow(Matrix a, ll b, ll p)
{Matrix res, tmp = a;memset(res.data, 0, sizeof(res));res.row = res.col = a.row;for(int i = 1; i <= res.row; ++i) { res.data[i][i] = 1; }while(b) {if(b & 1) res = matrix_mul(res, tmp, p);tmp = matrix_mul(tmp, tmp, p);b >>= 1;}return res;
}ll quick_pow(ll a, ll b, ll c)
{ll res = 1, tmp = a % c;while(b) {if(b & 1) res = res * tmp % c;tmp = tmp * tmp % c;b >>= 1;}return res;
}int main()
{GetPhi();int T, cases = 0;ll a, b, p, n;scanf("%d", &T); while(T--) {scanf("%lld%lld%lld%lld", &a, &b, &p, &n);if(n == 1) {printf("Case #%d: %lld\n", ++cases, a % p);continue;} else if(n == 2) {printf("Case #%d: %lld\n", ++cases, b % p);continue;} else if(n == 3) {printf("Case #%d: %lld\n", ++cases, a * b % p);continue;} else if(a == 0 || b == 0 || p == 1) {printf("Case #%d: %d\n", ++cases, 0);continue;}Matrix res, tmp;tmp.row = tmp.col = 2;tmp.data[1][1] = 0, tmp.data[1][2] = 1;tmp.data[2][1] = 1, tmp.data[2][2] = 1;res.row = 2, res.col = 1;res.data[1][1] = 0, res.data[2][1] = 1;tmp = matrix_quick_pow(tmp, n - 2, phi[p]);res = matrix_mul(tmp, res, phi[p]);ll ans = quick_pow(a, res.data[1][1], p);ans = ans * quick_pow(b, res.data[2][1], p) % p;printf("Case #%d: %lld\n", ++cases, ans);}return 0;
}