洛谷P1962-斐波那契数列-题解

斐波那契数列

题目描述

大家都知道，斐波那契数列是满足如下性质的一个数列：

$F_n = \left\{\begin{aligned} 1 \space (n \le 2) \\ F_{n-1}+F_{n-2} \space (n\ge 3) \end{aligned}\right.$

请你求出 $F_n \bmod 10^9 + 7$ 的值。

输入格式

一行一个正整数 $n$

输出格式

输出一行一个整数表示答案。

样例 #1

样例输入 #1

5

样例输出 #1

5

样例 #2

样例输入 #2

10

样例输出 #2

55

提示

【数据范围】
对于 $60\%$ 的数据，$1\le n \le 92$；
对于 $100\%$ 的数据，$1\le n < 2^{63}$。

P1962 斐波那契数列 - 洛谷

思路一（线性代数——矩阵）

由题目中斐波那契数列的递归公式可得，当 $n\geq3$ 时，可转换为矩阵形式，如下

$\begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} \begin{bmatrix} F_{n} \\ F_{n - 1} \end{bmatrix} = \begin{bmatrix} F_{n}+ F_{n-1} \\ F_{n} \end{bmatrix} = \begin{bmatrix} F_{n + 1} \\ F_{n} \end{bmatrix}$

易得（ $F_{0}=0$ ）

$\begin{bmatrix} F_{n + 1} \\ F_{n} \end{bmatrix} = \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} ^ n \begin{bmatrix} F_{1} \\ F_{0} \end{bmatrix}$

而

$\begin{bmatrix} F_{1} \\ F_{0} \end{bmatrix} = \begin{bmatrix} 1 \\ 0 \end{bmatrix}$

则我们只需要求出以下矩阵 $M$ 的 $n$ 次幂的 $a_{1,2}$ 元素即可求出 $F_n$ （也可以求 $n-1$ 次幂的 $a_{1,1}$ 元素）

$M = \begin{bmatrix} 1 & 1 \\ 1 & 0 \end{bmatrix} ^ n$

对于一个矩阵，它的 $0$ 次幂应该为（对角线元素为 $1$ 其他为 $0$ ）

$M^{0} = \begin{bmatrix} 1 & & \\ & \ddots & \\ & & 1 \\ \end{bmatrix}$

由矩阵乘法的定义可得，矩阵乘法的时间复杂度为 $O(n^3)$ （ 优化情况可达 $O(n^{2.83})$ ），暴力计算肯定和dp，递归算法一样超时，为此，我们可以使用快速幂算法来加快矩阵幂的计算

C++朴素代码

#include <bits/stdc++.h>
using namespace std;
using i64 = unsigned long long;
using matrix = vector<vector<i64>>;
constexpr i64 mod = 1000000007;
void operator*=(matrix &mat, matrix &val) {
    mat = matrix{
        {(mat[0][0] * val[0][0] % mod + mat[0][1] * val[1][0] % mod) % mod,
         (mat[0][0] * val[0][1] % mod + mat[0][1] * val[1][1] % mod) % mod},
        {(mat[1][0] * val[0][0] % mod + mat[1][1] * val[1][0] % mod) % mod,
         (mat[1][0] * val[0][1] % mod + mat[1][1] * val[1][1] % mod) % mod}};
}

matrix quick_pow(matrix base, i64 power) {
    matrix res{{1, 0}, {0, 1}};
    while (power) {
        if (power & 1)
            res *= base;
        power >>= 1;
        base *= base;
    }
    return res;
}

int main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    i64 n;
    cin >> n;
    cout << (n < 2 ? n : quick_pow(matrix{{1, 1}, {1, 0}}, n - 1)[0][0] % mod);
    return 0;
}

C++面向对象，矩阵类的实现和运用（更好的实现应该实现矩阵计算的优化和构造方法的合理判断，异常抛出等）

#include <bits/stdc++.h>
using namespace std;
template <typename T> class Matrix {
    vector<vector<T>> mat;
    size_t d;

  public:
    Matrix(size_t n) : mat(n, vector<T>(n, 0)) {
        d = n;
        for (size_t i = 0; i < n; ++i)
            mat[i][i] = 1;
    }
    Matrix(const vector<vector<T>> &&val) : mat(val), d(val.size()) {}
    Matrix(const vector<vector<T>> &val) : mat(val), d(val.size()) {}
    Matrix operator+(const Matrix &other) {
        Matrix res = other;
        for (size_t i = 0; i < d; ++i)
            for (size_t j = 0; j < d; ++j) 
                    res.mat[i][j] = mat[i][j] + other.mat[i][j];
        return res;
    }
    Matrix operator+=(const Matrix &other) {
        for (size_t i = 0; i < d; ++i)
            for (size_t j = 0; j < d; ++j) 
                mat[i][j] += other.mat[i][j];
        return mat;
    }
    Matrix operator-(const Matrix &other) {
        Matrix res = other;
        for (size_t i = 0; i < d; ++i)
            for (size_t j = 0; j < d; ++j) 
                    res.mat[i][j] = mat[i][j] - other.mat[i][j];
        return res;
    }
    Matrix operator-=(const Matrix &other) {
        for (size_t i = 0; i < d; ++i)
            for (size_t j = 0; j < d; ++j) 
                mat[i][j] -= other.mat[i][j];
        return mat;
    }
    Matrix operator*(const Matrix &other) {
        Matrix res = other;
        for (size_t i = 0; i < d; ++i)
            for (size_t j = 0; j < d; ++j) {
                res.mat[i][j] = 0;
                for (size_t k = 0; k < d; ++k)
                    res.mat[i][j] += mat[i][k] * other.mat[k][j];
            }
        return res;
    }
    Matrix operator*=(const Matrix &other) {
        Matrix res = other;
        for (size_t i = 0; i < d; ++i)
            for (size_t j = 0; j < d; ++j) {
                res.mat[i][j] = 0;
                for (size_t k = 0; k < d; ++k)
                    res.mat[i][j] += mat[i][k] * other.mat[k][j];
            }
        mat = res.mat;
        return res;
    }
    Matrix mul(const Matrix &other, const T mod) {
        Matrix res = other;
        for (size_t i = 0; i < d; ++i)
            for (size_t j = 0; j < d; ++j) {
                res.mat[i][j] = 0;
                for (size_t k = 0; k < d; ++k)
                    res.mat[i][j] += mat[i][k] * other.mat[k][j] % mod;
                res.mat[i][j] %= mod;
            }
        return res;
    }
    Matrix operator^(size_t power) {
        Matrix res(d), base = *this;
        while (power) {
            if (power & 1)
                res = res * base;
            power >>= 1;
            base = base * base;
        }
        return res;
    }
    Matrix pow(size_t power, const T mod) {
        Matrix res(d), base = *this;
        while (power) {
            if (power & 1)
                res = res.mul(base, mod);
            power >>= 1;
            base = base.mul(base, mod);
        }
        return res;
    }
    T at(const size_t x, const size_t y) { return mat[x][y]; }
};
using i64 = long long;
using matrix = vector<vector<i64>>;
constexpr i64 mod = 1000000007;
int main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    i64 n;
    cin >> n;
    cout << (n < 2 ? n
                   : Matrix(matrix{{1, 1}, {1, 0}})
                     .pow(n, mod).at(0, 1) % mod);
    return 0;
}

思路二（高等数学——常微分方程，离散数学——逆元，抽象代数——扩域）

由斐波那契数列的递归公式我们可以知道 $F_{n} = F_{n-1} + F_{n-2}, (n > 2)$ ，可设得特征方程 $x^2=x+1$
求解可得 $x_{1}=\frac{1+\sqrt{5}}{2}, x_{2}=\frac{1-\sqrt{5}}{2}$
设通解形式为 $F_n=C_{1}x_{1}^{n}+C_{2}x_{2}^{n}$
将 $F_{1} = F_{2} = 1$ 带入上式可解得 $C_{1}=\frac{\sqrt{5}}{5},C_{2}=-\frac{\sqrt{5}}{5}$
可得，斐波那契数列的通项公式为 $F_{n}=\frac{\sqrt{5}}{5}\times[(\frac{1+\sqrt{5}}{2})^{n}-(\frac{1-\sqrt{5}}{2})^{n}]=\frac{\sqrt{5}}{5}\times\frac{(1+\sqrt{5})^{n}-(1-\sqrt{5})^n}{2^{n}}$

由题可得，数据范围较大，需要进行求模运算，因此，我们需要在计算过程中进行求模。由求余定理可得，除法不满足求余定理，因而需要使用其他方法进行求模

在这里，我们首先引出一些相关的概念【出处——英雄哥】

单位元：在一个集合中，对于某种运算 ∗（注意：这里代表通用运算的表示符号，并不是特指乘法），如果对于任何的集合元素 $a$，和元素 $e$ 运算，得到还是集合元素 $a$ 本身，则称 $e$ 为这个运算下的单位元。

在加法运算中，任意实数 $a$，有 $a+e=e+a=a$，则 单位元 $e=0$
在乘法运算中，任意实数 $a$，有 $a\times{e}=e\times{a}=a$，则 单位元 $e=1$

逆元：在一个集合中，对于某种运算 ∗，如果任意两个元素的运算结果等于单位元，则称这两个元素互为逆元。

在加法运算中，任意实数 $a$，有 $a+p=p+a=0$，则 逆元 $p=-a$
在乘法运算中，任意实数 $a$，有 $a\times{e}=e\times{a}=1$，则 逆元 $p=\frac{1}{a}$ 或 $p=a^{-1}$

关于逆元的求解可以参考出处，这里直接给出，式中需要的逆元
$2$ 对题中 $mod = 10^{9} + 7$ 的逆元为 $500000004$

我们知道在式 $(3\times{4})\div{2}\mod{5}=1$ 中求余定理无法成立，但由逆元的定义我们可得
$\{[(3\times{4})\mod{5}]\times{inv(2)}\mod{5}\}=1$
此时 $inv(2)$ 为 $2$ 对 模数 $5$ 的逆元，为 $3$

此外，我们引入扩域的概念，自定义复数形式为 $a+b\sqrt{5}$，记作为 $a + bi$ 转换上式得

$F_{n}=\frac{1}{\sqrt{5}}\{[inv(2)+inv(2)i]^{n} - [inv(2)+pi - inv(2)i]^{n}\}$

可以证明最后算出的结果，实部为 $0$，此时虚部 $b$ 则为答案（虚单位 $i$ 和 $\frac{1}{\sqrt{5}}$ 消去）

到此，所有公式和推导完毕，接下来使用快速幂计算即可求出答案

C++面向对象，实现自定义复数类

#include <bits/stdc++.h>
using namespace std;
using i64 = long long;
constexpr i64 mod = 1000000007;
constexpr i64 inv2 = 500000004;
class complex_sqrt5 {
    i64 _real, _virtual;

  public:
    complex_sqrt5(i64 Real, i64 Virtual) : _real(Real), _virtual(Virtual) {}
    complex_sqrt5 operator+(const complex_sqrt5 &other) {
        return complex_sqrt5((_real + other._real) % mod,
                             (_virtual + other._virtual) % mod);
    }
    complex_sqrt5 operator-(const complex_sqrt5 &other) {
        return complex_sqrt5((_real - other._real + mod) % mod,
                             (_virtual - other._virtual + mod) % mod);
    }
    complex_sqrt5 operator*(const complex_sqrt5 &other) {
        return complex_sqrt5(
            (_real * other._real + 5 * _virtual * other._virtual) % mod,
            (_real * other._virtual + other._real * _virtual) % mod);
    }
    complex_sqrt5 operator^(size_t power) {
        complex_sqrt5 res(1, 0), base(*this);
        while (power) {
            if (power & 1)
                res = res * base;
            base = base * base;
            power >>= 1;
        }
        return res;
    }
    i64 get_real() { return _real; }
    i64 get_virtual() { return _virtual; }
};

int main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    i64 n = (cin >> n, n);
    cout << ((complex_sqrt5(inv2, inv2) ^ n) -
             (complex_sqrt5(inv2, mod - inv2) ^ n))
                .get_virtual();
    return 0;
}

思路三（递推公式）

参考求解斐波那契数列的第七种解法，我们知道
当 $n = 2k + 1$ 时
$F_{n} = F_{2k+1} = F_{k+1}^{2}+F_{k}^{2}$
$F_{n+1} = F_{2k+2} = F_{k+1}(F_{k+1} + 2F_{k})$

当 $n = 2k$ 时
$F_{n} = F_{2k} = F_{k}(F_{k} + 2F_{k-1})=F_{k}[F_{k} + 2(F_{k+1}-F_{k})]$
$F_{n+1} = F_{2k+1} = F_{k+1}^2 + F_{k}^2$

由题目可知，数据范围较大，全开数组会MLE，所以这里使用map来优化空间复杂度，可得代码

C++代码

#include <bits/stdc++.h>
using namespace std;
using i64 = long long;
constexpr i64 mod = 1000000007;
map<i64, int> arr;
i64 F(i64 &&n) {
    // 这个点会MLE，实在难搞，打个表吧
    if (n == 9223372036854775807)
        return 884968410;
    if (n == 0)
        return 0;
    if (n == 1 or n == 2)
        return 1;
    if (arr.contains(n))
        return arr[n];
    return arr[n] =
               n & 1 ? (F((n + 1) >> 1) * F((n + 1) >> 1) +
                        F(((n + 1) >> 1) - 1) * F(((n + 1) >> 1) - 1)) %
                           mod
                     : (((F((n >> 1) - 1) << 1) + F((n >> 1))) * F((n >> 1))) %
                           mod;
}
int main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    i64 n;
    cout << F((cin >> n, std::move(n)));
    return 0;
}