Luogu P12847 [蓝桥杯 2025 国 A] 斐波那契数列 题解 [ 蓝 ] [ 矩阵加速 ] [ 扩展欧拉定理 ]

斐波那契数列：比较巧妙的矩阵加速。

普通的矩阵是 \((+，\times)\) 的，显然无法处理此类乘积的值。而如果修改成 \((\times,\times)\) 的矩阵，则 \(\times\) 对 \(\times\) 不存在分配律，无法使用矩阵加速。因此可以从本题的特殊性质入手。

另记一下 \((\oplus,\otimes)\) 能使用矩阵加速的条件：

• \(\oplus\) 满足交换律。
• \(\otimes\) 满足交换律、结合律，且 \(\otimes\) 对 \(\oplus\) 具有分配律。

注意到该斐波那契的初值为 \(2,3\)，且递推式是前两个元素的乘积，因此该数列中的所有数的质因子只能是 \(2,3\)。这启发我们从这两个质因子的角度来考虑问题。下面以质因子 \(2\) 为例：

进一步发现，\(\log_2 f_i + \log_2 f_{i+1}=\log_2 f_{i+2}\)。所以每个质因子的 \(\log\) 就是普通的斐波那契数列，可以使用矩阵加速递推。注意矩阵加速的模数不能直接取 \(998244353\)，根据扩展欧拉定理，降幂的时候幂的模数应该为 \(\phi(998244353) = 998244352\)。在算出每个质因子的个数后快速幂求得答案即可。

因为扩展欧拉定理需要分两种情况：指数小于 \(\phi\)、指数大于等于 \(\phi\)，而判断指数的大小不太好做。所以当 \(n\) 较小的时候，直接暴力递推；当 \(n\) 较大的时候，采用扩展欧拉定理即可。时间复杂度 \(O(B^3\log n)\)，其中 \(B = 3\)。

#include <bits/stdc++.h>
#define fi first
#define se second
#define eb(x) emplace_back(x)
#define pb(x) push_back(x)
#define lc(x) (tr[x].ls)
#define rc(x) (tr[x].rs)
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
typedef long double ldb;
using pi = pair<int, int>;
const ll mod = 998244353, phi = 998244352, N = 1e6 + 5;
ll n, ans1, ans2, fx[N];
ll qpow(ll a, ll b)
{
    ll res = 1;
    while(b)
    {
        if(b & 1) res = (res * a) % mod;
        b >>= 1;
        a = (a * a) % mod;
    }
    return res;
}
struct Matrix{
    ll a[3][3];
    Matrix(){ memset(a, 0, sizeof(a)); }
    Matrix operator * (const Matrix & t) const{
        Matrix res;
        for(int i = 0; i < 3; i++)
            for(int k = 0; k < 3; k++)
                for(int j = 0; j < 3; j++)
                    res.a[i][j] = (res.a[i][j] + a[i][k] * t.a[k][j]) % phi;
        return res;
    }
}s1, s2, dp;
Matrix matqpow(Matrix a, ll b)
{
    Matrix res;
    res.a[0][0] = res.a[1][1] = res.a[2][2] = 1;
    while(b)
    {
        if(b & 1) res = res * a;
        b >>= 1;
        a = a * a;
    }
    return res;
}
void construct_mat()
{
    s1.a[0][0] = s1.a[0][2] = s2.a[0][1] = s2.a[0][2] = 1;
    dp.a[0][1] = dp.a[1][0] = dp.a[1][1] = dp.a[0][2] = dp.a[1][2] = dp.a[2][2] = 1;
}
int main()
{
    //freopen("sample.in", "r", stdin);
    //freopen("sample.out", "w", stdout);
    ios::sync_with_stdio(0);
    cin.tie(0);
    cout.tie(0);
    cin >> n;
    if(n <= 1e6)
    {
        fx[1] = 2, fx[2] = 3; ans1 = 1;
        for(int i = 3; i <= n; i++) fx[i] = (fx[i - 1] * fx[i - 2]) % mod;
        for(int i = 1; i <= n; i++) ans1 = (fx[i] * ans1) % mod;
        cout << ans1;
        return 0;
    }
    construct_mat();
    s1 = s1 * matqpow(dp, n - 2);
    s2 = s2 * matqpow(dp, n - 2);
    ll x = qpow(2, s1.a[0][2] + phi), y = qpow(3, s2.a[0][2] + phi);
    cout << x * y % mod;
    return 0;
}