本节介绍矩阵快速幂以及矩阵在ACM中的简单应用。本文需要一些线性代数知识作铺垫，这里不再重复。关于快速幂的思想，建议先阅读快速幂。

        首先看这么一道题（洛谷P1307）：

题目描述

        对于Fibonacci数列：1,1,2,3,5,8,13……大家应该很熟悉吧~~~但是现在有一个很“简单”问题：第n项和第m项的最大公约数是多少？

输入输出格式

输入格式：

        两个正整数n和m。（n,m≤109）
        注意：数据很大

输出格式：

        Fn和Fm的最大公约数。
        由于看了大数字就头晕，所以只要输出最后的8位数字就可以了。

输入输出样例

Sample input

4 7

Sample output

1

        要解决本题需要先认识到一个事实（证明略）：

        于是只需要求出n与m的最大公因数，再找它所对应的一位即可。注意在整个过程中对1e8取模，时间复杂度O(gcd(n,m))。
        但是，对于很大的数（比如本题达到1e9的数量级），循环1e9次显然会超时。于是需要借助矩阵这个工具来简化运算过程。
        矩阵是解决线性问题的利器。所谓运用矩阵工具，就是要善于发现题目中的线性关系。注意到斐波那契数列的递推式为：

        这是一个典型的线性关系，立即可以得到（用f(x)表示Fib(x)）：

        发现两侧的非常数矩阵形式不同，不利于递推，于是把它扩展为：

        进行如下定义：

        那么有递推关系：

        这就是斐波那契数列的矩阵递推式，反复利用递推公式可以得到：

        这就是斐波那契数列的矩阵通项公式，问题转化为如何快速求矩阵的幂。结合整数快速幂的思想，仍然可以用倍增方法求解矩阵的幂，将求矩阵乘积次数降至对数级别。
        矩阵快速幂与整数快速幂有以下几点区别：

  1. 矩阵的平方和乘积更复杂，可以采用朴素的行乘列法则暴力求出。
  2. 初始化的值由1改为单位矩阵，即：

        利用矩阵快速幂优化，本题迎刃而解。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
#include<iostream>

#define MOD 100000000
using namespace std;

struct Matrix {//2*2矩阵
    long long Main[2][2]{};

    Matrix(int a, int b, int c, int d) {
        Main[0][0] = a, Main[0][1] = b, Main[1][0] = c, Main[1][1] = d;
    }

    void operator*=(Matrix x) {//重载乘法运算
        long long temp[2][2];
        temp[0][0] = (Main[0][0] * x.Main[0][0] + Main[0][1] * x.Main[1][0]) % MOD;
        temp[0][1] = (Main[0][0] * x.Main[0][1] + Main[0][1] * x.Main[1][1]) % MOD;
        temp[1][0] = (Main[1][0] * x.Main[0][0] + Main[1][1] * x.Main[1][0]) % MOD;
        temp[1][1] = (Main[1][0] * x.Main[0][1] + Main[1][1] * x.Main[1][1]) % MOD;
        Main[0][0] = temp[0][0], Main[0][1] = temp[0][1];
        Main[1][0] = temp[1][0], Main[1][1] = temp[1][1];
    }
};

int gcd(int a, int b) {
    if (a < b)swap(a, b);
    if (b == 0)return a;
    return gcd(b, a % b);
}

int main() {
    int n, m;
    cin >> n >> m;
    Matrix sta(1, 1, 1, 0), ans(1, 0, 0, 1);//ans为单位阵
    int p = gcd(n, m) - 1;//下计算sta的p次方，矩阵快速幂
    while (p > 0) {
        if ((p & 1) == 1)ans *= sta;
        sta *= sta;
        p >>= 1;
    }
    cout << (ans.Main[1][0] + ans.Main[1][1]) % MOD;
    return 0;
}