洛谷P1349 广义斐波那契数列（矩阵快速幂）

最新推荐文章于 2023-05-16 20:51:11 发布

tiany7

最新推荐文章于 2023-05-16 20:51:11 发布

阅读量277

点赞数

分类专栏：数论文章标签：矩阵快速幂洛谷P1349 广义斐波那契数列 acm/oi

本文链接：https://blog.csdn.net/Stagflation/article/details/106804915

版权

数论专栏收录该内容

10 篇文章 1 订阅

订阅专栏

2020.6.17
不得已去学了上一道题，貌似都不是很难的样子。矩阵快速幂真的太好用了，推柿子神器。虽然template并不能用变量赋值，但是我发现一般好像很多时候矩阵的大小在做题之前就是已经确定的，那就可以上template了。

这几天真的不太平。先是咖喱，然后学校突然宣布8.24开学（喵喵喵？一天两万多的确诊坚持开学您是认真的么？甚至还克扣了我们两天假期？）可以线上，线上课程安排也是打乱了的，万一一节课在cdt早上八点，然后一节课在晚上7点那还睡个毛线啊？还让6.26之前回复。甚至线上上课学费也一分不减。好事儿都让学校办了，好话都让学校说了，学校只要不降低教学质量行吧行吧，我一个小老百姓就想赶紧混个大三身份然后下学期去一趟nac退役滚粗，进bat学习学习，赶紧考个研究生然后开始工作。真的是太烦了太烦了，好不容易在一场人类面临的浩劫历经34小时跨过三个地区辗转回国，现在让我再经历一遍谈何容易啊？不想这么多了，矩阵快速幂txdy！

这个矩阵快速幂就是推柿子用的，比如an是从前面某几个状态转移过来的，而递推的操作又会让复杂度接近爆炸，那么这个时候就可以用矩阵快速幂了！！毕竟an只取决于那几个状态，而状态^3跑log次铁定不会超时。就比如这道题：
在这里插入图片描述
这就非常清楚了，按照正常情况我们递推就行了，看到数据范围

这还递推个毛线啊，2^31-1是21亿次运算，常数再小也就是爆炸。这个时候发现转移次数固定，an只取决于前两项，所以就可以构造transition matrix。
这个transition matrix，或者是转置矩阵给出，已知An = An-1 + An-2
那么我们有
在这里插入图片描述

跑出来在乘上系数继续取模就行了，还是蛮简单的。
代码：

#include <bits/stdc++.h>
using namespace std;
#define limit (10 + 5)//防止溢出
#define INF 0x3f3f3f3f
#define inf 0x3f3f3f3f3f
#define lowbit(i) i&(-i)//一步两步
#define EPS 1e-6
#define FASTIO  ios::sync_with_stdio(false);cin.tie(0);
#define ff(a) printf("%d\n",a );
#define pi(a,b) pair<a,b>
#define rep(i, a, b) for(int i = a; i <= b ; ++i)
#define per(i, a, b) for(int i = b ; i >= a ; --i)
#define mint(a,b,c) min(min(a,b), c)
#define MOD 998244353
#define FOPEN freopen("C:\\Users\\tiany\\CLionProjects\\acm_01\\data.txt", "rt", stdin)
typedef long long ll;
typedef unsigned long long ull;
ll read(){
    ll sign = 1, x = 0;char s = getchar();
    while(s > '9' || s < '0' ){if(s == '-')sign = -1;s = getchar();}
    while(s >= '0' && s <= '9'){x = x * 10 + s - '0';s = getchar();}
    return x * sign;
}//快读
void write(ll x){
    if(x < 0) putchar('-'),x = -x;
    if(x / 10) write(x / 10);
    putchar(x % 10 + '0');
}

int n ,k;
int mod;
template <int T>
struct matrix {
    ll mat[T + 1][T + 1];

    matrix() {
        memset(mat, 0 , sizeof(mat));
    }
    ll mul(const matrix<T> &rhs, int x, int y) {
        ll val = 0;
        rep(i, 1, T) {
            (val += ((mat[x][i] * rhs.mat[i][y]) % mod)) %= mod;
        }
        return val % mod;
    }
    matrix<T> operator*(const matrix<T> &rhs) {
        matrix res;
        rep(i, 1, T) {
            rep(j, 1, T) {
                res.mat[i][j] = mul(rhs, i, j);
            }
        }
        return res;
    }
    static matrix<T> norm() {//单位向量
        matrix ans;
        rep(i ,1, T)ans.mat[i][i] = 1;
        return ans;
    }
    void set(int i, int j, ll val){
        mat[i][j] = val;
    }
};
template <int T>
matrix<T> quickPow(matrix<T> base, ll expo){
    matrix ans = matrix<T>::norm();//构造单位矩阵
    while (expo){
        if(expo & 1 )ans = ans * base;
        base = base * base;
        expo >>= 1;
    }
    return ans;
}
int p,q,first,second;
int main() {
#ifdef LOCAL
    FOPEN;
#endif
    matrix<2> ans;
    p = read(), q = read(), first = read(), second = read(),k =read(), mod = read();
    if(k < 3){
        if(k == 1)write(second);
        else write(first);
        return 0;
    }
    ans.set(1,1,p);
    ans.set(1,2,q);
    ans.set(2,1,1);
    matrix res = quickPow(ans, k - 2);//少做两次
    ll ret = res.mat[1][1] * second + res.mat[1][2] * first;
    ret %= mod;
    write(ret);
    return 0;
}