首页 > 编程语言 >算法学习:矩阵快速幂/矩阵加速

算法学习:矩阵快速幂/矩阵加速

时间:2024-08-11 22:39:20浏览次数:6  
标签:Mat int rep 矩阵 st 算法 加速 tem

1.前言

 其实本质上来说,矩阵快速幂或是矩阵加速的题目比较的模板化一些,大体上都是属于我们要先写出来一个递推式子(或者是我们需要递推的式子),然后由于递推的次数过大,1e18之类的,会导致复杂度的飚升,所以我们会用到矩阵来帮我们快速处理。
 另外,从题目的类型上大概是分为两类,一类是线性递推式,再利用矩阵加速处理;另外一类就是图上的问题,一般点数较小,可能100不到,然后对边有限制要求之类的,例如需要走T步,T的范围是[1,1e18],此时我们也可以用矩阵加速处理。
  本篇的主要内容以矩阵加速为主,矩阵快速幂的内容会比较草率。

2.矩阵快速幂

我们来大概了解一下矩阵快速幂,首先要确保我们是会快速幂的。以下贴出来快速幂的代码:

int kuai(int a,int b){
	int tem=1;
	while(b){if(b%2) tem=tem*a; a=a*a,b/=2;}
	return tem;
}

具体了解快速幂请百度,而矩阵的运算和实数的运算在结合律上是一致的,所以我们完全可以套用快速幂的形式写出来矩阵快速幂:Mat是我设的矩阵的结构体

static Mat kuai(Mat A, int b) {
        Mat tem(A.n, A.m, 1, A.st);//参数分别代表,行,列,初始化内容,下标从1或者0开始
        while (b) { if (b % 2) tem = tem * A; b = b / 2, A = A * A; }
        return tem;
    }

另外此处贴上我的矩阵的板子,后面的代码展示中就省略了板子:

struct Mat {
    int a[130][130];
    int n, m, st;
    Mat(int n_, int m_, int fl = 0, int st_ = 1) {
        n = n_, m = m_, st = st_;
        if (!fl) rep(i, st, n) rep(j, st, m) a[i][j] = 0;
        if (fl == INF) rep(i, st, n) rep(j, st, m) a[i][j] = INF;
        if (fl == 1) rep(i, st, n) rep(j, st, m) a[i][j] = (i == j);
    }
    int* operator [] (int x) { return a[x]; }
    Mat operator * (Mat b) {
        Mat res(n, b.m, 0, b.st);//记得赋值
        rep(i, st, n) rep(j, st, b.m) rep(k, st, m) {
            res[i][j] += a[i][k] * b[k][j] % mod;
            res[i][j] %= mod;
        }
        return res;
    }
    Mat operator + (Mat b) {
        Mat res(n, b.m, 0, b.st);//记得赋值
        rep(i, st, n) rep(j, st, b.m) {
            res[i][j] = a[i][j] + b[i][j];
            res[i][j] %= mod;
        }
        return res;
    }
    static Mat kuai(Mat A, int b) {
        Mat tem(A.n, A.m, 1, A.st);
        while (b) { if (b % 2) tem = tem * A; b = b / 2, A = A * A; }
        return tem;
    }
    //指数从0开始的求和
    static Mat kuai_sum(Mat A, int b) {
        Mat tem(A.n, A.m, 1, A.st); Mat sum = tem;
        while (b) {
            if (b % 2) tem = sum + tem * A;
            sum = sum + sum * A; A = A * A, b /= 2;
        }
        return tem;
    }
    void out() {
        rep(i, st, n) { rep(j, st, m) cout << a[i][j] << " "; cout << endl; }
    }
};

3.矩阵加速

第一类:线性递推式加速

学习一个东西是什么,以及怎么用,最好的办法就是去看一道经典的例题。
下面来看一道矩阵加速的板题:洛谷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\) 的值。
对于 \(100\%\) 的数据,\(1\le n < 2^{63}\)。

题解:

这里我们已经有了递推式,\(F_{n-1} + F_{n-2} = F_{n}\),正常我们只需要这样一直递推下去就好了,但是n的数据范围是2的63次方,所以我们不得不考虑别的办法:我们先把这两个式子写成矩阵的形式

\(I_1=\begin{bmatrix}F_{n-1} & F_{n-2} \end{bmatrix}\)

此时我们假设一下,这个矩阵需要乘上一个什么矩阵可以变成

\(I_2=\begin{bmatrix}F_{n} & F_{n-1} \end{bmatrix}\)

在草稿纸上写写画画,我们就发现

\(A=\begin{bmatrix}1 & 1 \\ 1 & 0 \end{bmatrix}\)

当\(I_1*A\)就会变成\(I_2\),而且,这个矩阵都是常数,且任意的一个\(I_{i}\)矩阵乘上它都会满足变成\(I_{i+1}\)
所以我们从\(I_{i}\)一直乘A,可以变成我们想要的\(I_{n}\),此时在利用之前学习的矩阵快速幂,我们就可以在\(log\)的时间里求出来我们最后的\(F_{n}\),另外再提一句,这个A矩阵的构建,是我们要先写出来\(I_1\)(我们本身有的递推式),然后将这个递推式子的下标都+1,这是我们的目标矩阵,根据这两个来倒推出我们的转移矩阵,也就是A矩阵。我们会发现当[1,1]和[2,1]的位置都填1,就是我们斐波那契数列本身的递推式子,然后\(F_{n-1}\)我们直接在[1,2]的位置写1,再次利用我们上一个矩阵本身就有的即可。

此处贴上代码:

int n, m, T;
signed main() {
    read(), kuaidu();
    int st, ed;
    cin >> n;
    if (n == 1 || n == 2) cout << 1 << endl;
    else {
        Mat A(1, 2);
        A[1][1] = 1, A[1][2] = 1;
        Mat base(2, 2);
        base[1][1] = base[1][2] = base[2][1] = 1;
        A = A * Mat::kuai(base, n - 2);
        cout << A[1][1] << endl;
    }
    return 0;
}

第二类:在图上的应用

我们还是接着来看一道例题:洛谷P2233

题目描述

在长沙城新建的环城公路上一共有 \(8\) 个公交站,分别为 A、B、C、D、E、F、G、H。公共汽车只能够在相邻的两个公交站之间运行,因此你从某一个公交站到另外一个公交站往往要换几次车,例如从公交站 A 到公交站 D,你就至少需要换 \(3\) 次车。

图片来自洛谷

Tiger 的方向感极其糟糕,我们知道从公交站 A 到公交 E 只需要换 \(4\) 次车就可以到达,可是 tiger 却总共换了 \(n\) 次车,注意 tiger 一旦到达公交站 E,他不会愚蠢到再去换车。现在希望你计算一下 tiger 有多少种可能的乘车方案。

输入格式
仅有一个正整数 \(n\),表示 tiger 从公交车站 A 到公交车站 E 共换了 \(n\) 次车。
\(4\le n\le10^7\)。

题解:

首先这个题我们很显然的能够知道,转化成图上就是有8个点,8个点只有互相相邻的可以抵达,我们建立一个邻接矩阵,其中给相邻的能到达的赋值为1,(除了E这个点,因为到达了就不会再去别的车站了),E车站向其他车站的赋值是0。因为这里n的数值很大,所以我们还是可以进行矩阵加速。
这里稍微浅扯一下关于这里为什么图上的邻接矩阵我们依旧可以进行矩阵加速,先看矩阵乘法的代码实现:

Mat operator * (Mat b) {
        Mat res(n, b.m, 0, b.st);//记得赋值
        rep(i, st, n) rep(j, st, b.m) rep(k, st, m) {
            res[i][j] += a[i][k] * b[k][j] % mod;
            res[i][j] %= mod;
        }
        return res;
    }

这里我们中间的式子,可以抽象的理解成是\(i\)到\(k\)的路径方案数乘以\(k\)到\(j\)的路径方案数,加到\(res[i][j]\)里面,此刻\(res[i][j]\)代表的就是我们将两个矩阵相乘后(意义上是步数相加)的方案数,原本这个邻接矩阵代表着我们走1步的方案数,所以再进行矩阵快速幂即可。

signed main() {
    read(), kuaidu();
    int t;
    cin >> n;
    Mat A(7, 7);
    rep(i, 0, 7) {
        int x = i, y = (i + 1) % 8, z = (i - 1 + 8) % 8;
        A[x][z] = 1;
        A[x][y] = 1;
    }
    rep(i, 0, 7) A[4][i] = 0;
    Mat dp = Mat::kuai(A, n);
    // dp.out();
    // A.out();
    cout << dp[0][4] << endl;
    return 0;
}

到此为止,我们算是将矩阵加速的基础题过了一遍,本人剩下的训练内容来自洛谷的矩阵题单,可以自行搜索。
我在这里面的题当中抽出来几道讲讲:

4.习题

4.1 P4159 [SCOI2009] 迷路

题目描述

该有向图有 \(n\) 个节点,节点从 \(1\) 至 \(n\) 编号,windy 从节点 \(1\) 出发,他必须恰好在 \(t\) 时刻到达节点 \(n\)。
现在给出该有向图,你能告诉 windy 总共有多少种不同的路径吗?
答案对 \(2009\) 取模。
注意:windy 不能在某个节点逗留,且通过某有向边的时间为该边长度。

输入格式

第一行包含两个整数,分别代表 \(n\) 和 \(t\)。
第 \(2\) 到第 \((n + 1)\) 行,每行一个长度为 \(n\) 的字符串,第 \((i + 1)\) 行的第 \(j\) 个字符 \(c_{i, j}\) 是一个数字字符,若为 \(0\),则代表节点 \(i\) 到节点 \(j\) 无边,否则代表节点 \(i\) 到节点 \(j\) 的边的长度为 \(c_{i, j}\)。

- 对于 \(100\%\) 的数据,保证 \(2 \leq n \leq 10\),\(1 \leq t \leq 10^9\)。

题解:

这里我们知道,如果有向边的长度是1的话,就直接变成模板题了,但是长度是权值,我们不能这么做。细心一点,会发现因为输入的是字符串,所以长度最大就是9,边权的范围是在1-9之间。我们可以像网络流那样来进行拆点,每一个点\(i\)拆成9个点,\(i_1,i_2,...,i_9\),而我们每一次从一个点走,都是从\(i_1\)走,举个例子,我们从\(i\)到\(j\),权值是3,我们可以建立\(i_1 -> j_3\),权值是1的边,然后再从\(j_3->j_2->j_1\),每一条边权值都是1,这样用拆点的方法解决了权值的问题,剩下的就和正常的板子题一样了。

int fan(int x, int y) {
    return (x - 1) * 10 + y;
}
signed main() {
    read(), kuaidu();
    cin >> n >> m;
    Mat A(100, 100);
    rep(i, 1, n) {
        string S; cin >> S;
        int len = S.size();
        rep(j, 0, len - 1) {
            if (S[j] == '0') continue;
            int x = fan(i, 1), y = fan(j + 1, (S[j] - '0'));
            A[x][y] = 1;
        }
    }
    rep(i, 1, n) {
        rep(j, 2, 10) A[fan(i, j)][fan(i, j - 1)] = 1;
    }
    A = Mat::kuai(A, m);
    cout << A[fan(1, 1)][fan(n, 1)] << endl;
    return 0;
}

4.2 P3216 [HNOI2011] 数学作业

题目描述

给定正整数 \(n,m\),要求计算 \(\text{Concatenate}(n) \bmod \ m\) 的值,其中 \(\text{Concatenate}(n)\) 是将 \(1 \sim n\) 所有正整数 顺序连接起来得到的数。

例如,\(n = 13\),\(\text{Concatenate}(n) = 12345678910111213\)。小 C 想了大半天终于意识到这是一道不可能手算出来的题目,于是他只好向你求助,希望你能编写一个程序帮他解决这个问题。

对于 \(100\%\) 的数据,\(1\le n \le 10^{18}\),\(1\le m \le 10^9\)。

题解:

这里我们先用递推的形式,大概能写出来,当i是个位数的时候,我们只需要\(F_{i-1}\)乘10,再加上\(i\)就好了,但是我们发现当i是两位数的时候,我们就是乘以100而不是10,以此类推,当我们是k位数时,我们需要乘\(10^k\),所以这里我们没有办法像之前那样,找到转移矩阵一口气直接乘到最后就好了,我们可以先写出来当i是个位数时的转移矩阵,然后需要每当位数+1的时候我们的转移矩阵也更新,将10变成100,再变成1000...,就可以了,时间复杂度会比之前的log多一个lg出来。
(顺便一提,这个题需要开unsigned long long,并且不要使用std的pow函数,返回类型是int,还有判断位数的时候不要使用log10函数,要手写!!不要问我怎么知道的qwq)

int kuai(int a, int b) {
    int tem = 1;
    while (b) {
        if (b % 2) tem = tem * a;
        a = a * a, b /= 2;
    }
    return tem;
}
int fan(int x) {
    int cnt = 0;
    while (x) {
        x /= 10;
        cnt++;
    }
    return cnt;
}
signed main() {
    read(), kuaidu();
    // cout << log10(124551) << endl;
    cin >> n >> mod;
    Mat base(3, 3);
    base[1][1] = 10, base[2][1] = 1, base[2][2] = 1, base[3][2] = 1, base[3][3] = 1;
    Mat dp(1, 3);
    dp[1][1] = 0, dp[1][2] = 1, dp[1][3] = 1;

    int k = fan(n) - 1;
    int las = n - kuai(10, k);
    int q = 1;
    while (1) {
        base[1][1] = q * 10 % mod;
        q *= 10;
        if (q <= n) { dp = dp * Mat::kuai(base, q - q / 10); }
        else break;
    }
    dp = dp * Mat::kuai(base, las + 1);
    cout << dp[1][1] << endl;

    return 0;
}

关于矩阵加速的内容先写这么多,之后还会再更一篇矩阵的题,是2024杭电7的循环图。就是因为这个题才让我开始补矩阵的-_-

标签:Mat,int,rep,矩阵,st,算法,加速,tem
From: https://www.cnblogs.com/advisedy/p/18354016

相关文章

  • 人工智能算法工程师(高级)课程11-自然语言处理之NLP的语言模型-seq2seq模型,seq+注意
    大家好,我是微学AI,今天给大家介绍一下人工智能算法工程师(高级)课程11-自然语言处理之NLP的语言模型-seq2seq模型,seq+注意力,word2vec与代码详解。本课程面向高级人工智能算法工程师,深入讲解自然语言处理(NLP)中的关键语言模型技术,包括seq2seq模型及其增强版加入注意力机制......
  • haproxy负载均衡之-调度算法详解
    HAProxy的调度算法分为静态调度算法、动态调度算法和其他调度算法静态算法:按照事先定义好的规则轮询公平调度,不关⼼后端服务器的当前负载、链接数和响应速度等,且⽆法实时修改权重,只能靠重启HAProxy⽣效。动态算法:基于后端服务器状态进⾏调度适当调整,⽐如优先调度⾄当前负载较......
  • 2-梯度下降算法
    梯度下降算法只能保证找到的是局部最优,不是全局最优平常我们经过大量实验,发现局部最优点不是很多,所以可以使用梯度下降算法。但是还要提防鞍点下面进行实现梯度下降算法点击查看代码importnumpyasnpimportmatplotlib.pyplotaspltx_data=[1.0,2.0,3.0]y_......
  • 05-组件生命周期及diff算法
    组件生命周期及diff算法生命周期回调函数理解旧新生命周期总结重要的勾子即将废弃的勾子DOM的diffing算法key的作用经典面试题:生命周期回调函数理解1.组件从创建到死亡它会经历一些特定的阶段。2.React组件中包含一系列勾子函数(生命周期回调函数),会在特......
  • 【C++算法】双指针
    移动零题目链接:移动零https://leetcode.cn/problems/move-zeroes/description/算法原理这类题是属于数组划分、数组分开题型代码步骤:使用cur遍历数组当cur所指的元素等于0时,cur向后面移动当cur所指的元素不等于0时,dest向后面移动,cur所指元素与dest移动后所指的元素交换当......
  • 第九天:K-Means算法
    K-Means算法简介K-Means算法是一种广泛使用的聚类算法,旨在将数据集分成K个预定义的簇。每个簇的中心是簇中所有点的均值,称为质心。K-Means算法的目标是最小化每个数据点到其所属簇的质心的距离的平方和。算法原理K-Means算法的工作原理可以分为以下几个步骤:初始化:随机......
  • 算法题系列5
    题目描述模拟一套简化的序列化传输方式,请实现下面的数据编码与解码过程1.编码前数据格式为[位置,类型,值],多个数据的时候用逗号分隔,位置仅支持数字,不考虑重复等场景;类型仅支持:Integer/String/Compose(Compose的数据类型表示该存储的数据也需要编码)2.编码后数据参考图示,数据......
  • 用电量预测 | 基于BiLSTM双向长短期记忆神经网络算法的用电量预测附matlab完整代码
    用电量预测|基于BiLSTM双向长短期记忆神经网络算法的用电量预测附matlab完整代码数据收集:收集历史用电量数据,包括时间戳和相应的用电量值。选择模型:选择合适的模型进行预测,可以根据数据特点和需求选择合适的模型。训练模型:使用历史数据训练模型,并根据评估指标来调整......
  • 【WSN覆盖优化】基于鱼鹰优化算法OOA求解无线传感器节点2D覆盖优化问题附Matlab代码
    鱼鹰优化算法(OspreyOptimizationAlgorithm,OOA)是一种基于鱼鹰捕鱼行为的启发式优化算法,可用于解决优化问题。在无线传感器网络(WSN)中,覆盖优化是一个关键问题,涉及到最大化网络覆盖范围并减少节点数量。以下是一个简单的示例框架,展示如何基于OOA算法求解无线传感器节点的二......
  • 什么是算法
    1.概述算法指的是为了完成某一事情(或者解决某一问题),而经过特定步骤的处理之后得到结果的一种手段具有明确的步骤/顺序/可行性计算机擅长做固定的运算,例如求和等计算型的处理,通过合适的算法(合适的处理策略),可以大大降低运算所需要的时间2.举例为了对数字进行排序......