HDU 1588 Gauss Fibonacci ★(矩阵 && 求和)

题目大意：f(i)为fibonacci函数. g(i)=k*i+b，给定k,b,n,M求sum(f(g(i)) mod M for 0<=i<n 和这道题很像，不过这道题有个地方不好想点儿. 我们把sum展开 = f(b) + f(k+b) + f(2k+b) + …… + f((n-1)k+b). 到这儿我们发现这个和函数不像上题那么好看出递推规律……那么这题的难点就在这儿…… 自己推的时候就隐约觉得他们存在某种关系，但就是想不出来，后来看了解题报告才明白了…… 通常时候我做fibonacci递推时是用 ( f(n), f(n-1) ) = (1, 1, 1, 0)^(n-1) * (f(1), f(o)) ，这样的话就不好把f(n)表示成一个形式比较好的具体表达式，现在我们用另一种方法表示f(n)，把上面那个矩阵递推关系变形一下得：f(n) = (1, 1, 1, 0) ^ n 的右上角元素。设A = (1, 1, 1, 0)，我们就大胆的把f(n)表示成A^n (事实上这个正确性也好说明，因为求和用f(n)时都是加法，而加法是对应位置相加的，并不影响其他位置，也不会受其他位置影响，所以可以这样表示) 那么就有 sum = A^b + A^(k+b) + A^(2k+b) + …… + A^((n-1)k + b). 再变一下形式就很显然了： sum = A^b * (A^0 + A^k + A^2k + A^3k + …… + A^(n-1)k ) 如果我们定义SUM(A, n) = A^1 + A^2 + …… + A^n，那么sum = A^b * ( E + SUM(A^k, n) ). (别忘了sum实际上是矩阵的右上角元素)

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#define MID(x,y) ((x+y)>>1)
using namespace std;
typedef long long LL;

const int MAX = 4;
struct Mat{
    int row, col;
    LL mat[MAX][MAX];
};
//initialize square matrix to unit matrix
Mat unit(int n){
    Mat A;
    A.row = A.col = n;
    memset(A.mat, 0, sizeof(A.mat));
    for (int i = 0; i < n; i ++)
        A.mat[i][i] = 1;
    return A;
}
//return T(A)
Mat transpose(Mat A){
    Mat T;
    T.row = A.col;
    T.col = A.row;
    for (int i = 0; i < A.col; i ++)
        for (int j = 0; j < A.row; j ++)
            T.mat[i][j] = A.mat[j][i];
    return T;
}
//return (A+B)%mod
Mat add(Mat A, Mat B, int mod){
    Mat C = A;
    for (int i = 0; i < C.row; i ++)
        for (int j = 0; j < C.col; j ++)
            C.mat[i][j] = (A.mat[i][j] + B.mat[i][j]) % mod;
    return C;
}
//return A*B%mod
Mat mul(Mat A, Mat B, int mod){
    Mat C;
    C.row = A.row;
    C.col = B.col;
    for (int i = 0; i < A.row; i ++){
        for (int j = 0; j < B.col; j ++){
            C.mat[i][j] = 0;
            for (int k = 0; k < A.col; k ++)
                //注意这里要保证乘法不溢出，否则还需要设计特殊的乘法模
                C.mat[i][j] += A.mat[i][k] * B.mat[k][j];
            C.mat[i][j] %= mod;
        }
    }
    return C;
};
//return A^n%mod
Mat exp_mod(Mat A, int n, int mod){
    Mat res = unit(A.row);
    while(n){
        if (n & 1){
            res = mul(res, A, mod);
        }
        A = mul(A, A, mod);
        n >>= 1;
    }
    return res;
}
//return S = A^1 + A^2 + …… + A^n   (binary search)
Mat Sum(Mat A, int n, int mod){
    if (n == 1)
        return A;
    if (n & 1){
        return add(Sum(A, n-1, mod), exp_mod(A, n, mod), mod);
    }
    else{
        return mul(Sum(A, n/2, mod), add(exp_mod(A, n/2, mod), unit(A.col), mod), mod);
    }
}
int main(){
    int k, b, n, M;
    while(cin >> k >> b >> n >> M){
        Mat A;
        A.col = A.row = 2;
        A.mat[0][0] = 1;    A.mat[0][1] = 1;    A.mat[1][0] = 1;    A.mat[1][1] = 0;
        Mat Ab = exp_mod(A, b, M);
        Mat Ak = exp_mod(A, k, M);
        Mat E = unit(2);
        Mat SUM = Sum(Ak, n-1, M);
        Mat res = mul(Ab, add(E, SUM, M), M);
        cout << res.mat[0][1] << endl;
    }
	return 0;
}

举杯独醉，饮罢飞雪，茫然又一年岁。 ------AbandonZHANG