Trudbot's Blog

Appearance

来聊聊如何将快速幂的思想应用到矩阵乘法上， 以及矩阵快速幂的应用。

矩阵快速幂 ​

矩阵乘法 ​

在线性代数中学过， n行x列的矩阵A与x行m列的矩阵B是可以相乘的， 结果为一个n行m列的矩阵R， 且Rij=∑k=1xAik∗Bkj.

而对于方阵Mnn, 又有幂的概念， Mn=M∗M∗M∗...∗M, 即n个M矩阵相乘.

由此我们可以定义矩阵类型， 并实现矩阵乘法。

c++
#define N 2
struct matrix {
    int m[N][N];

    matrix() {
        memset(m, 0, sizeof(m));
    }

    matrix operator * (const matrix& b) const {
        matrix res;
        for (int i = 0; i < N; i ++) 
            for (int j = 0; j < N; j ++) 
                for (int k = 0; k < N; k ++)
                    res.m[i][j] += m[i][k] * b.m[k][j];
        return res;
    }
};

快速幂思想 ​

快速幂的细节不多说， 简单说快速幂的思想， 就是： 计算an时， 可以先计算an2, 再用an2∗an2得到an。

可以发现， 计算an只比an2多了一次（若n为奇数还要多一次), 相比于暴力做法(n / 2次)优秀很多， 而时间复杂度也从O(n)优化为O(logn)。

矩阵快速幂 ​

我们可以套用快速幂的模板， 写出矩阵快速幂的代码。

cpp
matrix qpow(matrix a, int p) {
    matrix res;
    for (int i = 0; i < N; i ++) res.m[i][i] = 1; //单位矩阵
    while (p) {
        if (p & 1) res = res * a;
        a = a * a;
        p >>= 1;
    }
    return res;
}

到此， 矩阵快速幂的实现就完成了。

矩阵快速幂的应用 ​

典中典——斐波那契数列 ​

假设有这样一个问题：对于整数n， n∈(1,1e9)， 在1s内求出斐波那契数列的第n项模p的结果。

显然常规O(n)做法是行不通的。

斐波那契数列我们都知道存在关系式fn=fn−1+fn−2, 那我们能否通过构造多项式， 将fn和f1联系起来， 得到一个公式呢？倘若这样， 才有可能在O(n)内完成。

将fn和fn−1设为一个行向量， 则

(fn fn−1)=(fn−1+fn−2  fn−1)

是否能通过(fn−1 fn−2)来表示(fn fn−1)呢？

可以发现

(fn fn−1)=(fn−1+fn−2 fn−1)=(fn−1 fn−2)∗(1 01 1)

将fn−1,fn−2看作自变量， 则1 01 1即为(fn,fn−1)的系数矩阵.

由此递推可得

(fn fn−1)=(f1 f0)∗(1 01 1)n−1

所以我们的主要任务就是算出此系数矩阵的n-1次幂，也就可以用上述提到的矩阵快速幂， 这样即使n是109也不会超时。

整数的连续幂次和 ​

对于这样一个问题：求a1+a2+...+ak模p的结果， k∈(1,1e9), a为正整数.

我们不妨定义fn=a1+a2+...+an，则观察可知fn=a∗fn−1+a

以此我们可以构造得到关系式

(fn 1)=(fn−1 1)∗(a 0a a)

问题解决。