HDU 5667 Sequence (矩阵快速幂 + 费马小定理)

Sequence

[align=center]Time Limit: 2000/1000 MS (Java/Others)    Memory Limit: 65536/65536 K (Java/Others)
Total Submission(s): 310    Accepted Submission(s): 99
[/align]

[align=left]Problem Description[/align]

    Holion
August will eat every thing he has found.

    Now
there are many foods,but he does not want to eat all of them at once,so he find a sequence.

fn=⎧⎩⎨⎪⎪1,ab,abfcn−1fn−2,n=1n=2otherwise

    He
gives you 5 numbers n,a,b,c,p,and he will eat fn
foods.But there are only p foods,so you should tell him
fn
mod p.
 

[align=left]Input[/align]

    The
first line has a number,T,means testcase.

    Each
testcase has 5 numbers,including n,a,b,c,p in a line.

    1≤T≤10,1≤n≤1018,1≤a,b,c≤109,p
is a prime number,and p≤109+7.
 

[align=left]Output[/align]

    Output
one number for each case,which is fn
mod p.
 

[align=left]Sample Input[/align]

1
5 3 3 3 233

 

[align=left]Sample Output[/align]

190

 

[align=left]Source[/align]

BestCoder Round #80

题目链接：http://acm.hdu.edu.cn/showproblem.php?pid=5667

题目大意：求公式的结果

题目分析：显然取对数，递推式一眼化成矩阵乘法，然后快速幂，因为矩阵快速幂的结果是a的幂，所以矩阵快速幂取模的时候根据费马小定理要模p-1，不过有可能a模p==0，这个时候答案显然是0，但是有可能指数取模的结果是0，这个时候a^0变成1就不对了

#include <cstdio>
#include <cstring>
#define ll long long
ll n, a, b, c, p;
struct matrix
{
ll m[5][5];
};

matrix multiply(matrix x, matrix y, ll MOD)
{
matrix ans;
memset(ans.m, 0, sizeof(ans.m));
for(int i = 1; i <= 3; i++)
for(int j = 1; j <= 3; j++)
if(x.m[i][j])
for(int k = 1; k <= 3; k++)
ans.m[i][k] = (ans.m[i][k] + x.m[i][j] * y.m[j][k]) % MOD;
return ans;
}

matrix quickmod(matrix a, ll x, ll MOD)
{
matrix ans;
memset(ans.m, 0, sizeof(ans.m));
for(int i = 1; i <= 3; i++)
ans.m[i][i] = 1;
while(x)
{
if(x & 1)
ans = multiply(ans, a, p - 1);
x >>= 1;
a = multiply(a, a, p - 1);
}
return ans;
}

ll qpow(ll x, ll nn)
{
ll res = 1;
while(nn)
{
if(nn & 1)
res = (res * x) % p;
x = (x * x) % p;
nn >>= 1;
}
return res;
}

int main()
{
int T;
scanf("%d", &T);
while(T --)
{
scanf("%I64d %I64d %I64d %I64d %I64d", &n, &a, &b, &c, &p);
if(n == 1)
printf("1\n");
else if(n == 2)
printf("%I64d\n", qpow(a, b));
else
{
if(a % p == 0)
printf("0\n");
else
{
matrix ans;
ans.m[1][1] = (b * c + b) % (p - 1);
ans.m[1][2] = b;
ans.m[1][3] = 1;
matrix mx;
mx.m[1][1] = c; mx.m[1][2] = 1; mx.m[1][3] = 0;
mx.m[2][1] = 1; mx.m[2][2] = 0; mx.m[2][3] = 0;
mx.m[3][1] = b; mx.m[3][2] = 0; mx.m[3][3] = 1;
ans = multiply(ans, quickmod(mx, n - 3, p - 1), p - 1);
printf("%I64d\n", qpow(a, ans.m[1][1]));
}
}
}
}