hdu 1588 矩阵运算的应用

这题弄了两天才做出来，还是去请教了竹教主。

贴个别人的解说吧，自己懒得写了

把斐波那契数列转化为矩阵:
A={1,1}
{1,0};

{f[n+1],f
}
{f
,f[n-1]} = A^n ;最后输出右上角那项
或者用
{f[n+2],f[n+1]}
{f[n+1], f
} = A^(n+1); 最后输出右下角那项

我们用第一个公式
所求即为A^b + A^(k+b) + A^(2*k+b) + ... + A^((n-1)*k+b)
=A^b * ( A^0 + A^k + A^(2*k) + ... + A^((n-1)*k) )
=A^b * ( (A^k)^0 + (A^k)^1 + (A^k)^2 + ...+ (A^k)^(n-1) );
B=A^k;
上式=A^b * ( B^0 + B^1 + B^2 + ... + B^(n-1) );
B^0 + B^1 + B^2 + ... + B^(n-1)用上篇介绍到的递归二分 方法求解
最后输出矩阵的第二项(右上角)即可；

对于求解 B^0 + B^1 + B^2 + ... + B^(n-1)
我们也可以构造矩阵的方法。

我们来设置这样一个矩阵
B I
O I
其中O是零矩阵，I是单位矩阵

将它乘方，得到
B^2 I+B
O I
乘三方，得到
B^3 I+B+B^2
O I
乘四方，得到
B^4 I+B+B^2+B^3
O I
用快速幂求出n方，让第二项再与A^b相乘即可。

/*
* hdu1588/win.cpp
* Created on: 2011-11-30
* Author    : ben
*/
#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <cmath>
#include <ctime>
#include <iostream>
#include <algorithm>
#include <queue>
#include <set>
#include <map>
#include <stack>
#include <string>
#include <vector>
#include <deque>
#include <list>
#include <functional>
#include <numeric>
#include <cctype>
using namespace std;
typedef long long typec;
int MOD;
// 矩阵的最大阶
const int MAX_ORDER = 4;
struct MyMatrix {
int order;
typec num[MAX_ORDER][MAX_ORDER];
MyMatrix(int ord) {
order = ord;
}
void init() {
for (int i = 0; i < order; i++) {
for (int j = 0; j < order; j++) {
num[i][j] = 0;
}
}
}
};
MyMatrix operator*(MyMatrix ma, MyMatrix mb) {
int ord = ma.order;
MyMatrix numc(ord);
numc.init();
int i, j, k;
for (i = 0; i < ord; i++) {
for (j = 0; j < ord; j++) {
for (k = 0; k < ord; k++) {
numc.num[i][j] += ma.num[i][k] * mb.num[k][j];
numc.num[i][j] %= MOD;
}
}
}
return numc;
}
MyMatrix mpow(MyMatrix ma, typec x) {
int ord = ma.order;
MyMatrix numc(ord);
numc.init();
for (int i = 0; i < ord; i++) {
numc.num[i][i] = 1;
}
for (; x; x >>= 1) {
if (x & 1) {
numc = numc * ma;
}
ma = ma * ma;
}
return numc;
}

inline void initA(MyMatrix &m) {
m.num[0][0] = 1;
m.num[0][1] = 1;
m.num[1][0] = 1;
m.num[1][1] = 0;
}

inline void initC(MyMatrix &C, MyMatrix B) {
C.init();
for(int i = 0; i < 2; i++) {
for(int j = 0; j < 2; j++) {
C.num[i][j] = B.num[i][j];
}
}
C.num[0][2] = 1;
C.num[1][3] = 1;
C.num[2][2] = 1;
C.num[3][3] = 1;
}

MyMatrix getBrFromD(MyMatrix &D) {
MyMatrix res(2);
for(int i = 0; i < 2; i++) {
for(int j = 0; j < 2; j++) {
res.num[i][j] = D.num[i][j + 2];
}
}
return res;
}

int main() {
#ifndef ONLINE_JUDGE
freopen("data.in", "r", stdin);
#endif
int k, b, n;
while(scanf("%d%d%d%d", &k, &b, &n, &MOD) == 4 ){
MyMatrix A(2);
initA(A);
MyMatrix B = mpow(A, k);
MyMatrix C(4);
initC(C, B);
MyMatrix D = mpow(C, n);
MyMatrix E = getBrFromD(D);
MyMatrix ans = mpow(A, b) * E;
printf("%d\n", (int)ans.num[1][0]);
}
return 0;
}