矩阵乘法

Posted CHADLZX

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了矩阵乘法相关的知识,希望对你有一定的参考价值。

p1353:

一个由自然数组成的数列按下式定义:
对于i <= k:ai = bi 
对于i > k: ai = c[1]*a[i-1] + c[2]*a[i-2] + ... + c[k]*a[i-k]
其中bj 和 cj (1<=j<=k)是给定的自然数。写一个程序,给定自然数m <= n,
计算a[m] + a[m+1] + a[m+2] + ... + a[n], 并输出它除以给定自然数p的余数的值

矩阵乘法应用入门题,设矩阵A={a1,a2,a3...an-1,an,Sn},然后再构造一个B矩阵,使A*B=A’,A’={a2,a3,a4...an,an+1,Sn+1};用快速幂优化,求Sn,Sm-1,答案就出来了;

下面代码,因为是入门,所以也就没写什么struct Matrix;

技术分享
 1 #include<iostream>
 2 #include<cstdio>
 3 #include<cstring>
 4 #include<string>
 5 #include<cstdlib>
 6 #include<ctime>
 7 #include<vector>
 8 #include<algorithm>
 9 #include<queue>
10 using namespace std;
11 #define LL long long
12 const int maxn=20;
13 LL k,mod,m,n;
14 int a[maxn],b[maxn][maxn];
15 int y[maxn][maxn],c[maxn];
16 LL find(LL q){
17     if(q<=0)return 0;
18     memcpy(y,b,sizeof(y));
19     int ans[maxn][maxn],t[maxn][maxn];
20     memset(ans,0,sizeof(ans));
21     memset(t,0,sizeof(t));
22     for(int i=1;i<=k;i++)ans[i][i]=1;
23     while(q!=0){
24         if(q%2==1){
25             memset(t,0,sizeof(t));
26             for(int i=1;i<=k;i++)
27                 for(int j=1;j<=k;j++)
28                     for(int l=1;l<=k;l++)
29                         t[i][j]=((LL)t[i][j]+((LL)ans[i][l]*y[l][j])%mod)%mod;
30             memcpy(ans,t,sizeof(t));
31         }
32         q/=2;
33         memset(t,0,sizeof(t));
34         for(int i=1;i<=k;i++)
35             for(int j=1;j<=k;j++)
36                 for(int l=1;l<=k;l++)
37                     t[i][j]=(t[i][j]+((LL)y[i][l]*y[l][j])%mod)%mod;
38         memcpy(y,t,sizeof(t));
39     }
40     LL sum=0;
41     for(int i=1;i<=k;i++)sum=(sum+(LL)a[i]*ans[k][i]%mod)%mod;
42     return sum;
43 }
44 int main(){
45     cin>>k;
46     for(int i=1;i<=k;i++)cin>>a[i];
47     for(int i=1;i<k;i++)b[i][i+1]=1;
48     for(int i=1;i<=k;i++)cin>>b[k][k-i+1];
49     b[k+1][1]=1;
50     b[k+1][k+1]=1;
51     k++;
52     cin>>m>>n>>mod;
53     printf("%d\n",(find(n)-find(m-1)+mod)%mod);
54     return 0;
55 }
View Code

p1354:

HH有个一成不变的习惯,喜欢饭后百步走。所谓百步走,就是散步,就是在一定的时间
内,走过一定的距离。
但是同时HH又是个喜欢变化的人,所以他不会立刻沿着刚刚走来的路走回。
又因为HH是个喜欢变化的人,所以他每天走过的路径都不完全一样,他想知道他究竟有多
少种散步的方法。
现在给你学校的地图(假设每条路的长度都是一样的都是1),问长度为t,从给定地
点A走到给定地点B共有多少条符合条件的路径;

这道题的原型是二维动归入门,然后这道题经历了以下变迁: 二维dp入门->二维dp(图论形式)->二维dp(矩阵乘法优化)->二维dp(矩阵乘法优化)+点边转化;

二维dp入门:f[t][i][j]=f[t-1][i-1][j]+f[t-1][i][j-1]+f[t-1][i+1][j]+f[t-1][i][j-1];很水;

二维dp(图论形式):建成图,然后f[t][i]=sum(f[t-1][j]+d[j][i]);

二维dp(矩阵乘法优化):t很大的时候,上矩阵乘法,建一个邻接矩阵,利用矩阵乘法转移;

这道题:由于不能走上一步走的路,点边转化一下,拆成两条有向边,其他和上面一样;

写得很难受,从来没在图上搞这么难受的设计,尤其是在不知道有没有效果的情况下;

技术分享
 1 #include<iostream>
 2 #include<cstdio>
 3 #include<cstring>
 4 #include<string>
 5 #include<cstdlib>
 6 #include<ctime>
 7 #include<vector>
 8 #include<algorithm>
 9 #include<queue>
10 using namespace std;
11 #define LL long long
12 const int maxn=125;
13 const int mod=45989;
14 struct node{
15     int x,y,next;
16 }e[maxn];
17 int linkk[maxn],len=-1;
18 int n,m,t,A,B;
19 struct Matrix{
20     int v[maxn][maxn];
21     Matrix(){memset(v,0,sizeof(v));}
22     friend Matrix operator*(const Matrix& a,const Matrix& b){
23         Matrix c;
24         for(int i=0;i<=len;i++)
25             for(int j=0;j<=len;j++)
26                 for(int k=0;k<=len;k++)
27                     c.v[i][j]=(c.v[i][j]+a.v[i][k]*b.v[k][j]%mod)%mod;
28         return c;
29     }
30     friend Matrix operator^(Matrix a,int b){
31         Matrix c;
32         for(int i=0;i<=len;i++)c.v[i][i]=1;
33         while(b){
34             if(b%2)c=c*a;
35             b/=2;
36             a=a*a;
37         }
38         return c;
39     }
40 }a,b;
41 void insert(int x,int y){
42     e[++len].y=y;e[len].next=linkk[x];linkk[x]=len;e[len].x=x;
43 }
44 void init(){
45     memset(linkk,-1,sizeof(linkk));
46     cin>>n>>m>>t>>A>>B;
47     int x,y;A++,B++;
48     for(int i=1;i<=m;i++){
49         cin>>x>>y;
50         x++,y++;
51         insert(x,y);insert(y,x);
52     }
53     for(int i=linkk[A];i!=-1;i=e[i].next)a.v[1][i]++;
54     for(int i=0;i<=len;i++)
55         for(int j=0;j<=len;j++)
56             if(e[i].y==e[j].x)if(i!=(j^1))b.v[i][j]++;
57     a=a*(b^(t-1));
58     LL sum=0;
59     for(int i=linkk[B];i!=-1;i=e[i].next)sum=(sum+a.v[1][i^1])%mod;
60     cout<<sum%mod<<endl;
61 }
62 int main(){
63     init();
64 }
View Code

p1355:

FJ的N(2 <= N <= 1,000,000)头奶牛选择了接力跑作为她们的日常锻炼项目。至于进行接力跑的地点,自然是在牧场中现有的T(2 <= T <= 100)条跑道上。 

  农场上的跑道有一些交汇点,每条跑道都连结了两个不同的交汇点I1_i和I2_i(1 <= I1_i <= 1,000; 1 <= I2_i <= 1,000)。每个交汇点都是至少两条跑道的端点。奶牛们知道每条跑道的长度length_i(1 <= length_i <= 1,000),以及每条跑道连结的交汇点的编号。并且,没有哪两个交汇点由两条不同的跑道直
接相连。你可以认为这些交汇点和跑道构成了一张图。

  为了完成一场接力跑,所有N头奶牛在跑步开始之前都要站在某个交汇点上(有些交汇点上可能站着不只1头奶牛)。当然,她们的站位要保证她们能够将接力棒顺次传递,并且最后持棒的奶牛要停在预设的终点。

  你的任务是,写一个程序,计算在接力跑的起点(S)和终点(E)确定的情况下,奶牛们跑步路径可能的最小总长度。显然,这条路径必须恰好经过N条跑道.

 

这道题和上道题的转移方程差不多:f[t][i]=min(f[t-1][j]+d[j][i]);

但不同的是上一题是sum,这一题是min。

和上题一样,建立矩阵,但矩阵的乘法规则换成c[i][j]=min(a[i][k]+b[k][j]);

然后就可以转移了,初始矩阵^n表示走n步后的最短长度;

在看题解之前还真不知道矩阵乘法的规则可以这么改变;

技术分享
#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
#include<cstdlib>
#include<ctime>
#include<vector>
#include<algorithm>
#include<queue>
using namespace std;
#define LL long long
const int maxn=1010;
const int inf=1000000000;
struct node1{
    int x,y,v,next;
}e[maxn];
int linkk[maxn],len=0,n,N,m,S,E,top=0;
struct node{int x,y,v;}b[maxn];
bool flag[maxn];
int q[maxn],f[maxn];
void insert(int x,int y,int v){
    e[++len].y=y;e[len].x=x;
    e[len].next=linkk[x];
    linkk[x]=len;
    e[len].v=v;
}
void init(){
    scanf("%d%d%d%d",&n,&m,&S,&E);
    int x,y,v;
    for(int i=1;i<=m;i++){
        scanf("%d%d%d",&v,&x,&y);b[i].x=x,b[i].y=y,b[i].v=v;
        if(!flag[x])q[++top]=x,flag[x]=1;
        if(!flag[y])q[++top]=y,flag[y]=1;
    }
    sort(q+1,q+top+1);
    for(int i=1;i<=maxn;i++)f[q[i]]=i;
    for(int i=1;i<=m;i++){
        insert(f[b[i].x],f[b[i].y],b[i].v);
        insert(f[b[i].y],f[b[i].x],b[i].v);
    }
}
struct Mat{
    LL v[202][202];
    Mat(){memset(v,0,sizeof(v));}
    Mat operator*(Mat& x){
        Mat c;
        memset(c.v,10,sizeof(c.v));
        for(int i=1;i<=top;i++)
            for(int j=1;j<=top;j++)
                for(int k=1;k<=top;k++)
                {
                    c.v[i][j]=min(c.v[i][j],v[i][k]+x.v[k][j]);
                }
        return c;
    }
    Mat operator^(int b){
        Mat c;
        bool flag1=1;
        while(b){
            if(b%2)
            {
                if(flag1)c=*this,flag1=0;
                else c=c*(*this);
            }
            b/=2;
            *this=(*this)*(*this);
        }
        return c;
    }
}a;
void work(){
    memset(a.v,10,sizeof(a.v));
    for(int x=1;x<=top;x++){
        for(int i=linkk[x];i;i=e[i].next){
            a.v[x][e[i].y]=e[i].v;
        }
    }
    a=a^(n);
    cout<<a.v[f[S]][f[E]]<<endl;
}
int main(){
    init();
    work();
}
View Code

矩阵乘法有这么一段话:如果一个向量是其他向量的线性组合,便可以考虑用矩阵乘法;

这其实是在考验建模能力;

以上是关于矩阵乘法的主要内容,如果未能解决你的问题,请参考以下文章

C语言实现矩阵乘法

大型矩阵的 CUDA 矩阵乘法中断

将 SSE 矩阵向量乘法代码转换为 AVX

C++ 乘法大矩阵

疯子的算法总结 矩阵乘法 (矩阵快速幂)

需要帮助使用 MPI 调试并行矩阵乘法