LOV3 -> OV3

因为一些事情得重新拾起算法了

• Reference(Aka thanks)(可能"Reference"在这里的用法不太对，欢迎提醒我来纠正)
  • Wiki
    • 时间复杂度 - 维基百科，自由的百科全书
    • 快速幂 - OI Wiki
  • 基础
    • [[Binary Exponentiation]]
    • [[Linear Algebra|矩阵乘法]]
  • 视频
    • 509 Luogu P3390【模板】矩阵快速幂 bilibili
  • 题目
    • P1226 【模板】快速幂 - 洛谷 | 计算机科学教育新生态
    • P3390 【模板】矩阵快速幂 - 洛谷 | 计算机科学教育新生态
  • 讨论(想法)
    • 似乎还可以通过剪枝优化，下次()
    • 关于快速幂的名字笔者并无好感，更喜欢称其为“二进制求幂”，非常直观

把基础拿出来讲 -- 矩阵乘法

矩阵乘法:

TL;DR: 前一个矩阵逐行 与 后一个矩阵逐列 相乘求和 -- 结果放入结果矩阵 $Z$ 字排列

设 $A$ 为 $m*p$ 的矩阵，$B$ 为 $p*n$ 的矩阵，称 $m*n$ 的矩阵 $C$ 为 $A$ 与 $B$ 的乘积，则 $C$ 的第 $i$ 行第 $j$ 列的元素为

$c_{ij} = \sum_{k=1}^p a_{ik}b_{kj} = a_{i1}b_{1j} + a_{i2}b_{2j} + \cdots + a_{ip}b_{pj}$

$m$ 为行，$p$ 为列

• Complexity
  • 朴素求幂复杂度为 $\mathcal{O(N)}$，根据主定理(<-朋友告诉的，笔者还不会主定理)，分析得到快速幂的时间复杂度为 $\mathcal{O(logN)}$，
    实际上还可以通过费马小定理加速
    • if isPrime(m) $x^{n \pmod{m - 1}}$
  • 进而来到矩阵快速幂的复杂度: 设矩阵维度为 $m$，矩阵乘法的复杂度为 $m^3$，快速幂的复杂度为 $log_n$，故矩阵快速幂的时间复杂度为 $m^3*log n$，而 $m$ 为常数，而根据 时间复杂度 - 维基百科，自由的百科全书 中所述，在大O表示法中常数不计。矩阵快速幂复杂度即为 $\mathcal{O(log N)}$

初学编程的伙伴可能对于循环中 $i, j$ 不是很敏感，或者容易弄混，这里特意排出这个矩阵，供参考:

一个 $m \times n$ 的矩阵是一个由 $m$ 行 $n$ 列元素排列成的矩形阵列。即形如

$A = \underbrace{ \left. \begin{bmatrix} a_{1 1} & a_{1 2} & \cdots & a_{1 n} \\ a_{2 1} & a_{2 2} & \cdots & a_{2 n} \\ \vdots & \vdots & \ddots & \vdots \\ a_{m 1} & a_{m 2} & \cdots & a_{m n} \\ \end{bmatrix} \right\} }_{n} m \text{.}$

最开始题目所对应的答案(Template Code)

快速幂

Python

def quick_pow(x, n, mod = None):
    res = 1
    while n:
        if n & 1:
            res = res * x % mod
        x = x * x % mod
        # print(x)
        n >>= 1
    return res

print(f'{quick_pow(2, 4) = }') # 16

C++

#include <iostream>
// 快速幂
long long quick_pow(long long a, long long b, long long p){
    long long ans = 1;
    while (b) {
        if (b & 1) ans = ans * a % p;
        a = a * a % p;
        b >>= 1;
    }
    return ans;
}

int main() {
    long long a, b, p;
    std::cin >> a >> b >> p;
    // 2^10 mod 9=7
    std::cout << a << "^" << b << " mod " << p << "=" << quick_pow(a, b, p) << std::endl;
    return 0;
}

矩阵快速幂

C++

#include <cstring>  // memset
#include <iostream> // io

const long long MOD = 1e9 + 7; // 10**9 + 7

struct matrix {
    long long coefficient[105][105];
    matrix() {
        memset(coefficient, 0, sizeof(coefficient)); // 初始化为0
    }
} matrixA, result; // 初始所给的矩阵matrixA 结果矩阵result

int n;
long long k;

// 矩阵乘法 LaTex: 设matrixA为m*p的矩阵 B为p*n的矩阵
// 则称m*n的矩阵C为matrixA与B的乘积 则C的第i行第j列的元素为 $c_{ij} =
// \sum_{k=1}^p a_{ik}b_{kj} = a_{i1}b_{1j} + a_{i2}b_{2j} + \cdots +
// a_{ip}b_{pj}$ m为行 p为列 TL;DR: 前一个矩阵逐行 与 后一个矩阵逐列 相乘求和 --
// 结果放入结果矩阵Z字排列
matrix operator*(matrix &a, matrix &b) { // 重载乘号 实现矩阵乘法
    matrix ans; // 定义矩阵时会自动调用构造函数 也就是 matrix() 使其初始化为0
    for (int i = 1; i <= n; i++) {
        for (int k = 1; k <= n; k++) {
            for (int j = 1; j <= n; j++) {
                // 实现乘积求和的过程 原来的值 + 乘积
                ans.coefficient[i][j] =
                    (ans.coefficient[i][j] +
                     a.coefficient[i][k] * b.coefficient[k][j]) %
                    MOD;
                // 第一维下标尽量放在外层循环 第二维下标放在内层循环
                // 这样可以提高效率 也就是 i j k --> i k j
            }
        }
    }
    return ans;
}

void quickPower(long long exponentiation) { // 快速幂
    for (int i = 1; i <= n; i++) {          // 构造单位矩阵
        result.coefficient[i][i] = 1;       // 主对角线为1
    } // 单位矩阵乘任何矩阵都是原矩阵
    while (exponentiation) {
        if (exponentiation & 1) {
            result = result *
                     matrixA; // c++智能识别变量类型来决定是否调用矩阵乘法函数
                              // 使重载通用
        }
        matrixA = matrixA * matrixA;
        exponentiation >>= 1;
    }
}

int main() {
    // scanf("%d%lld",&n,&k); // k的范围为 10**12 所以用long long读入
    std::cin >> n >> k;
    for (int i = 1; i <= n; i++) { // 选择 index 从 1 开始
        for (int j = 1; j <= n; j++) {
            // 初始矩阵上的系数均为int范围内, 因此%d读入不会爆int,
            // 不过还是用更现代的 cin 吧 scanf("%d",
            // &matrixA.coefficient[i][j]);
            std::cin >> matrixA.coefficient[i][j];
        }
    }
    quickPower(k);
    for (int i = 1; i <= n; i++) {
        for (int j = 1; j <= n; j++) {
            // printf("%d ", result.coefficient[i][j]); //
            // 结果经过mod所以%d不会爆int, 不过避免强制类型转换, 还是改过来
            std::cout << result.coefficient[i][j] << " ";
        }
        fputs("\n", stdout);
    }
    return 0;
}