基础数学算法-Fibonacci 数列(矩阵乘法)

题目-Fibonacci 数列(矩阵乘法)

问题分析

假设求第$n$项, 相当于求矩阵的$A^n$可以用快速幂协助求解
但是还是要求前$n$项和$S_n$, 需要构造矩阵

首先考虑第求$n$项的情况
$$F_n=[f_n,f_{n+1}]$$
$$F_{n+1}=[f_{n+1},f_{n+2}]$$

推导之后发现

$$[f_{n+1},f_{n+2}]=[f_n,f_{n+1}]\times \left[\begin{array}{cc}0& 1\\ 1& 1\end{array}\right]$$

那么矩阵$A=\left[\begin{array}{cc}0& 1\\ 1& 1\end{array}\right]$, 第$n$项 = $F_0\times A^n$, 因为矩阵具有结合律, 可以用快速幂计算

但是题目要求的是$S_n$, 因此需要重新设计矩阵

$$F_n=[f_n,f_{n+1},S_n],F_{n+1}=[f_{n+1},f_{n+2},S_{n+1}]$$

得到等式
[ f n + 1 , f n + 2 , S n + 1 ] = [ f n , f n + 1 , S n ] × [ 0 1 0 1 1 1 0 0 1 ] [f_{n + 1}, f_{n + 2}, S_{n + 1}] = [f_n, f_{n + 1}, S_n] \times \begin{bmatrix} 0& 1 & 0\\ 1& 1 & 1\\ 0& 0 & 1 \end{bmatrix} [fn+1​,fn+2​,Sn+1​]=[fn​,fn+1​,Sn​]× ​010​110​011​ ​
也就是 A = [ 0 1 0 1 1 1 0 0 1 ] A = \begin{bmatrix} 0& 1 & 0\\ 1& 1 & 1\\ 0& 0 & 1 \end{bmatrix} A= ​010​110​011​ ​

实际计算的答案就是$F_n=F_0\times A^n$, 快速幂实现矩阵乘法

算法步骤

主要实现快速幂的乘积使用矩阵乘法计算即可
因为$n\le 2\times 10^9$, 快速幂算法时间复杂度$O(\log n)$, 可以通过

代码实现

#include <bits/stdc++.h>

using namespace std;

typedef long long LL;
const int N = 3;

int n, m;

void mul(int c[], int a[], int b[][N]) {
   
    int tmp[N] = {
   0};
    for (int i = 0; i < N; ++i) {
   
        for (int j = 0; j < N; ++j) {
   
            tmp[i] = (tmp[i] + (LL) a[j] * b[j][i] % m) % m;
        }
    }
    memcpy(c, tmp, sizeof tmp);
}

void mul(int c[][N], int a[][N], int b[][N]) {
   
    int tmp[N][N] = {
   0};
    for (int i = 0; i < N; ++i) {
   
        for (int j = 0; j < N; ++j) {
   
            for (int k = 0; k < N; ++k) {
   
                tmp[i][j] = (tmp[i][j] + (LL) a[i][k] * b[k][j] % m) % m;
            }
        }
    }
    memcpy(c, tmp, sizeof tmp);
}

int main() {
   
    ios::sync_with_stdio(false);
    cin.tie(0);

    cin >> n >> m;
    int f1[N] = {
   1, 1, 1};
    int A[N][N]= {
   {
   0, 1, 0},
                  {
   1, 1, 1},
                  {
   0, 0, 1}};

    n--;
    while (n) {
   
        if (n & 1) mul(f1, f1, A);
        mul(A, A, A);
        n >>= 1;
    }

    cout << f1[2] % m << '\n';

    return 0;
}