矩阵乘法

· · 个人记录

矩阵乘法

设A为n∗k矩阵，B为k∗m矩阵,那么它们的乘积C则为一个n∗m矩阵

C_{i,j}=\Sigma^k_{r=1}A_{i,r}\times B_{r,j}

如图

矩阵乘法中满足结合律分配律，但不满足交换律

单位矩阵，都是正方形的，对角线上元素皆为1，其他皆为0，保证乘后矩阵不变

矩阵快速幂

同个矩阵相乘多次，形状不变，因此可以用快速幂计算，重载乘号即可（注意，答案先建成单位矩阵）

#include<algorithm>
#include<iostream>
#include<cstdio>
#include<cstring>
typedef long long ll ;
const int N=105;
const ll P=1e9+7;
using namespace std;
int n;ll k;
struct squ
{
    ll a[N][N];    
    squ()
    {
        memset(a,0,sizeof(a));
    }
    inline void build() //建造单位矩阵，对角线为1，其他为0，乘上单位矩阵，结果不变 
    {    
        for(int i=1;i<=n;i++)a[i][i]=1;
    }
}x;
squ operator *(const squ &x,const squ &y)//重载运算符 
{     
    squ z;
    for(int i=1;i<=n;i++)
        for(int k=1;k<=n;k++)
            for(int j=1;j<=n;j++)
                z.a[i][k]=(z.a[i][k]+x.a[i][j]*y.a[j][k]%P)%P;//（i,k)=第i行*第k列 
    return z;
}
int main()
{
    scanf("%d%lld",&n,&k);
    for(int i=1;i<=n;i++)
     for(int j=1;j<=n;j++)
      scanf("%lld",&x.a[i][j]);
    squ ans;
    ans.build();
    while(k)
    {
        if(k&1==1) ans=ans*x;
        x=x*x;
        k=k>>1;
    }//递推快速幂，与普通的递推快速幂无异，但*不能缩写为*=
    for(int i=1;i<=n;i++)
    {
        for(int j=1;j<=n;j++)
            printf("%d ",ans.a[i][j]);
        printf("\n");
    }
    return 0;
}

矩阵加速

通过构造矩阵，快速计算数列第n项

题目链接

要根据 F[i-1] F[i−2] F[i−3] 推出 F[i]F[i−1]F[i−2]

可以乘上矩阵 \begin{matrix}1&0&1\\1&0&0\\0&1&1\end{matrix}

因此要得到第n项，只需要用单位矩阵*初始矩阵的(n-1)次幂

代码

#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
typedef long long ll;
const ll P=1e9+7;
struct squ
{
    ll a[5][5];   
    squ()
    {
        memset(a,0,sizeof(a));
    }
};
squ operator *(const squ &x,const squ &y)
{     
    squ z;
    for(int i=1;i<=3;i++)
        for(int k=1;k<=3;k++)
            for(int j=1;j<=3;j++)
                z.a[i][k]=(z.a[i][k]+x.a[i][j]*y.a[j][k]%P)%P; 
    return z;
}
int main()
{
    int T,x;
    scanf("%d",&T);
    for(int i=1;i<=T;i++)
    {
        scanf("%d",&x);
        if(x<=3) printf("1\n");
        else 
        {
            squ ans,b;
            for(int i=1;i<=3;i++) ans.a[i][i]=1;//建立单位矩阵 
            b.a[1][1]=b.a[1][3]=b.a[2][1]=b.a[3][2]=1;//建立初始矩阵 
            int k=x-1;
            while(k)
            {
               if(k&1==1) ans=ans*b;
               b=b*b;
               k=k>>1;
            }//x-1次方 
            printf("%lld\n",ans.a[1][1]);//第一列第一行即为答案 
        }
    }
    return 0;
}

斐波那契数列

递推式 f[n]=f[n-1]+f[n-2]

初始矩阵 \begin{matrix}1&1\\1&0\end{matrix}

代码

#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
typedef long long ll;
const ll P=1e9+7;
struct squ
{
    ll a[3][3];   
    squ()
    {
        memset(a,0,sizeof(a));
    }
};
squ operator *(const squ &x,const squ &y)
{     
    squ z;
    for(int i=1;i<=2;i++)
        for(int k=1;k<=2;k++)
            for(int j=1;j<=2;j++)
                z.a[i][k]=(z.a[i][k]+x.a[i][j]*y.a[j][k]%P)%P; 
    return z;
}
int main()
{
    ll n;
    scanf("%lld",&n);
    if(n<=2) printf("1");
    else 
    {
        squ ans,b;
        for(int i=1;i<=2;i++) ans.a[i][i]=1;
        b.a[1][1]=b.a[1][2]=b.a[2][1]=1;
        ll k=n-1;
        while(k)
        {
            if(k&1==1) ans=ans*b;
            b=b*b;
            k=k>>1;
        }
        printf("%lld\n",ans.a[1][1]%P);
    }
    return 0;
}

广义斐波那契

题目链接

递推式 f[n]=p\times f[n-1]+q\times f[n-2]

初始矩阵 \begin{matrix}p&1\\q&0\end{matrix}

代码

#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
typedef long long ll;
ll p,q,a1,a2,n,m;
struct squ
{
    ll a[3][3];   
    squ()
    {
        memset(a,0,sizeof(a));
    }
};
squ operator *(const squ &x,const squ &y)
{     
    squ z;
    for(int i=1;i<=2;i++)
        for(int k=1;k<=2;k++)
            for(int j=1;j<=2;j++)
                z.a[i][k]=(z.a[i][k]+x.a[i][j]*y.a[j][k]%m)%m; 
    return z;
}
int main()
{
    scanf("%lld%lld%lld%lld%lld%lld",&p,&q,&a1,&a2,&n,&m);
    if(n==1) printf("%lld",a1);
    else if(n==2) printf("%lld",a2);
    else 
    {
        squ ans,b;
        ans.a[1][1]=a2;ans.a[1][2]=a1;
        b.a[1][1]=p,b.a[2][1]=q,b.a[1][2]=1;
         ll k=n-2;
        while(k)
        {
            if(k&1==1) ans=ans*b;
            b=b*b;
            k=k>>1;
        }
        printf("%lld\n",ans.a[1][1]);
    }
    return 0;
}