矩阵快速幂与应用

矩阵快速幂

矩阵快速幂很好理解，和普通的快速幂差不多，只不过将普通的乘法换成了矩阵的乘法。

---

模板

const int MAX = 3;
const int mod = 10007;
typedef struct
{
    int m[MAX][MAX];
} Matrix;
Matrix I = {1, 0, 0,
            0, 1, 0,
            0, 0, 1};//单位矩阵
Matrix P;//基数矩阵
Matrix mul(Matrix a, Matrix b)//矩阵乘法
{
    Matrix c;
    for (int i = 0; i < MAX; i++)
        for (int j = 0; j < MAX; j++)
        {
            c.m[i][j] = 0;
            for (int k = 0; k < MAX; k++)
                c.m[i][j] += (a.m[i][k] * b.m[k][j]) % mod;
            c.m[i][j] %= mod;
        }
    return c;
}
Matrix quickpow(int n)//矩阵快速幂
{
    Matrix m = P, b = I;
    while (n >= 1)
    {
        if (n & 1)
            b = mul(b, m);
        n = n >> 1;
        m = mul(m, m);
    }
    return b;
}

---

应用

矩阵快速幂常用于类斐波那契数列的相关问题。

HDU 1757 A Simple Math Problem

x < 10 f(x) = x.
x >= 10 f(x) = a0 * f(x-1) + a1 * f(x-2) + a2 * f(x-3) + …… + a9 * f(x-10);
求f(k)%mod

构造递推矩阵

P：										S：
a0 a1 a2 a3 a4 a5 a6 a7 a8 a9           f(n-1)          f(n)
1  0  0  0  0  0  0  0  0  0            f(n-2)          f(n-1)
0  1  0  0  0  0  0  0  0  0            f(n-3)          f(n-2)
0  0  1  0  0  0  0  0  0  0            f(n-4)          f(n-3)
0  0  0  1  0  0  0  0  0  0            f(n-5)          f(n-4)
0  0  0  0  1  0  0  0  0  0      x     f(n-6)     =    f(n-5)
0  0  0  0  0  1  0  0  0  0            f(n-7)          f(n-6)
0  0  0  0  0  0  1  0  0  0            f(n-8)          f(n-7)
0  0  0  0  0  0  0  1  0  0            f(n-9)          f(n-8)
0  0  0  0  0  0  0  0  1  0            f(n-10)         f(n-9)

#include <bits/stdc++.h>
using namespace std;
const int MAX = 10;
typedef struct
{
    int m[MAX][MAX];
} Matrix;
Matrix I, P;
int k, mod;
int f[10];
Matrix mul(Matrix a, Matrix b)
{
    Matrix c;
    for (int i = 0; i < MAX; i++)
        for (int j = 0; j < MAX; j++)
        {
            c.m[i][j] = 0;
            for (int k = 0; k < MAX; k++)
                c.m[i][j] += (a.m[i][k] * b.m[k][j]) % mod;
            c.m[i][j] %= mod;
        }
    return c;
}
Matrix quickpow(int n)
{
    Matrix m = P, b = I;
    while (n >= 1)
    {
        if (n & 1)
            b = mul(b, m);
        n = n >> 1;
        m = mul(m, m);
    }
    return b;
}
int main()
{
    memset(I.m, 0, sizeof(I));
    memset(P.m, 0, sizeof(P));
    for (int i = 0; i < 10; i++)
    {
        I.m[i][i] = 1;
        f[i] = i;
    }
    for (int i = 1; i < 10; i++)
        P.m[i][i - 1] = 1;
    while (~scanf("%d%d", &k, &mod))
    {
        for (int i = 0; i < 10; i++)
            scanf("%d", &P.m[0][i]);
        if (k < 10)
        {
            printf("%d\n", f[k] % mod);
            continue;
        }
        Matrix S = quickpow(k - 9);
        int sum = 0;
        for (int i = 0; i < 10; i++)
            sum += (S.m[0][i] * f[9 - i]) % mod;
        sum %= mod;
        printf("%d\n", sum);
    }
    return 0;
}

fibonacci 数列的前 N 项和

大意求第 a 个到第 b 个之间所有 fibonacci 的数和。 这个有点小技巧在里面，可以自己推导一下：
F(3) = F(1) + F(2)
F(4) = F(2) + F(3) = 1 * F(1) + 2 * F(2)
F(5) = F(3) + F(4) = 2 * F(1) + 3 * F(2)
F(6) = F(4) + F(5) = 3 * F(1) + 5 * F(2)
F(7) = F(5) + F(6) = 5 * F(1) + 8 * F(2)
F(8) = F(6) + F(7) = 8 * F(1) + 13 * F(2)
S(3) = 2 * F(1) + 2 * F(2)
S(4) = 3 * F(1) + 4 * F(2)
S(5) = 5 * F(1) + 7 * F(2)
S(6) = 8 * F(1) + 12 *F(2)
S(7) = 13 *F(1) + 20 *F(2)
不难发现， S(n) = F(n + 2) - F(2) 若求 a 到 b 个之间的和：如下 因此题目就转换为了求 F(b + 2) - F(a + 2 - 1)
也可用如下递推矩阵

1  1  1         s(n-1)         s(n) = s(n-1) + f(n-1) + f(n-2)
0  1  1    x    f(n-1)    =    f(n) = f(n-1) + f(n-2)
0  1  0         f(n-2)         f(n-1) = f(n-1)