HDOJ 1005 Number Sequence 矩阵快速幂 代码 + 解析

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

部分题头：

Problem Description

A number sequence is defined as follows:

f(1) = 1, f(2) = 1, f(n) = (A * f(n - 1) + B * f(n - 2)) mod 7.

Given A, B, and n, you are to calculate the value of f(n).

Input

The input consists of multiple test cases. Each test case contains 3 integers A, B and n on a single line (1 <= A, B <= 1000, 1 <= n <= 100,000,000). Three zeros signal the end of input and this test case is not to be processed.

Output

For each test case, print the value of f(n) on a single line.

大意就是给你个递推表达式 f(n) = (A * f(n - 1) + B * f(n - 2)) mod 7，其中 f(1) = 1, f(2) = 1，要你求f(n)。

首先这个n的范围很大，上限是1e8，就别说什么递归了，就是O(n)的一个while递推都会超时。网上不少题解说找规律，这里我觉得用矩阵快速幂来做应该才是正解。思路就是模仿斐波那契数列的矩阵快速幂做法，构造一个包含A,B的矩阵，然后用矩阵快速幂求解，构建出来的矩阵快速幂表达式大概是这样的：

还要注意的是，递推式中还有求余，根据数论里面的知识有：(a * b) % m = (a%m * b%m) % m，修改一下矩阵乘法的代码，在每次两数相乘时先求余，最后得到结果矩阵之后再求一次余即可。

AC代码：

#include <string>
#include <cstring>
#include <algorithm>
#include <cstdio>

using namespace std;

#define M 2

int m = 0;

struct Mat
{

    int data[M][M];
    Mat()
    {
        memset(data, 0, sizeof(data));
    }

    void print()
    {
        for (int i = 0; i < M; ++i)
        {
            for (int j = 0; j < M; ++j)
                printf("%d ", data[i][j]);
            printf("\n");
        }
        printf("\n");
    }
};

Mat multipy(const Mat &ma, const Mat &mb, const int m)
{
    Mat res;
    for (int i = 0; i < M; ++i)
        for (int j = 0; j < M; ++j)
            for (int k = 0; k < M; ++k)
                res.data[i][j] += ma.data[i][k] % m * mb.data[k][j] % m;
    return res;
}

Mat mqp(Mat mat, int n, int m)
{

    Mat res;
    res.data[0][0] = res.data[1][1] = 1;
    while (n)
    {
        if (n & 1)
            res = multipy(res, mat, m);
        mat = multipy(mat, mat, m);
        //mat.print();
        n >>= 1;
    }
    return res;
}

int calc(int n, int a, int b, int m)
{
    if (n == 1 || n == 2)
        return 1;
    --n;

    Mat A;
    A.data[0][0] = 0;
    A.data[0][1] = 1;
    A.data[1][0] = b;
    A.data[1][1] = a;
    //A.print();
    Mat F;
    F.data[0][0] = 1;
    F.data[0][1] = 0;
    F.data[1][0] = 1;
    F.data[1][1] = 0;

    Mat tmp = multipy(mqp(A, n, m), F, m);
    return tmp.data[0][0] % m;
}

int main()
{
    //freopen("1.txt", "r", stdin);
    int a, b, n;
    while (1)
    {
        scanf("%d%d%d", &a, &b, &n);
        if (a == 0 && b == 0 && n == 0)
            break;
        printf("%d\n", calc(n, a, b, 7));
    }

    return 0;
}