zoukankan      html  css  js  c++  java
  • 【题解】【模板】矩阵级数

    【题解】矩阵乘法

    给定(n)阶矩阵,求(Sigma_{i=1}^{k} A^i),对(mod)取模

    此题可以直接分治,但我用纯数学办法做出来了QAQ

    在我的另一篇博客里,有这道题的分治做法。

    我们直接对题目要求什么设未知数吧,(sum_i=Sigma_{j=0}^{i} A^j),联系无穷级数的套路,从(0)开始,可能有比较好的性质(实际上是为了初始化方便)。构造一个矩阵把它们联系到一起。

    [left[egin{matrix} A^i\ S_i\ end{matrix} ight] ]

    考虑增广

    [left[egin{matrix} ? end{matrix} ight] imesleft[egin{matrix} A^i\ S_i\ end{matrix} ight]=left[egin{matrix} A^{i+1}\ S_{i+1}\ end{matrix} ight] ]

    根据矩阵乘法法则,直接待定系数法/观察法解出

    [left[egin{matrix} ? end{matrix} ight]=left[egin{matrix} A&0\ I&I\ end{matrix} ight] ]

    其中(left[egin{matrix} I end{matrix} ight])是纯量阵。

    那么

    [{left[egin{matrix} A&0\ I&I\ end{matrix} ight]} imes {left[egin{matrix} A^i\ S_i\ end{matrix} ight]} = {left[egin{matrix} A^{i+1}\ S_{i+1}\ end{matrix} ight]} ]

    所以

    [{left[egin{matrix} A&0\ I&I\ end{matrix} ight]}^{k+1} imes {left[egin{matrix} I\ 0\ end{matrix} ight]} = {left[egin{matrix} A^{k+1}\ S_{k+1}\ end{matrix} ight]} ]

    由于

    [S_n=Sigma_{i=0}^{n}A^i ]

    所以

    [Sigma_{i=1}^{k} A^i=S_{k+1}-I ]

    对于

    [{left[egin{matrix} A&0\ I&I\ end{matrix} ight]}^{k+1} ]

    考虑矩阵快速幂,而最后的答案是

    [{left[egin{matrix} S_{k+1} end{matrix} ight]} ]

    直接把

    [{left[egin{matrix} A^{k+1}\ S_{k+1}\ end{matrix} ight]} ]

    下面的(S_{k+1})输出出来就好了。记得(A)(S)都是矩阵,所以可以直接输出

    考场代码

    #include<bits/stdc++.h>
    
    using namespace std;typedef long long ll;
    #define DRP(t,a,b) for(register int t=(a),edd=(b);t>=edd;--t)
    #define RP(t,a,b)  for(register int t=(a),edd=(b);t<=edd;++t)
    #define ERP(t,a)   for(register int t=head[a];t;t=e[t].nx)
    #define midd register int mid=(l+r)>>1
    #define TMP template < class ccf >
    TMP inline ccf qr(ccf b){
        register char c=getchar();register int q=1;register ccf x=0;
        while(c<48||c>57)q=c==45?-1:q,c=getchar();
        while(c>=48&&c<=57)x=x*10+c-48,c=getchar();
        return q==-1?-x:x;}
    TMP inline ccf Max(ccf a,ccf b){return a<b?b:a;}
    TMP inline ccf Min(ccf a,ccf b){return a<b?a:b;}
    TMP inline ccf Max(ccf a,ccf b,ccf c){return Max(a,Max(b,c));}
    TMP inline ccf Min(ccf a,ccf b,ccf c){return Min(a,Min(b,c));}
    TMP inline ccf READ(ccf* _arr,int _n){RP(t,1,_n)_arr[t]=qr((ccf)1);}
    //----------------------template&IO---------------------------
    const int maxn=31<<1|1;
    ll n,yyb,k;
    struct mtx{
        ll data[maxn][maxn];
        mtx(){memset(data,0,sizeof data);}
        inline ll *operator [](int x){return data[x];}
        inline mtx operator *=(mtx x){return *this=(*this)*x;}
        inline mtx operator ^=(int x){return *this=(*this)^x;}
        inline void unis(){
    	RP(t,1,n)
    	    data[t][t]=1;
        }
        inline mtx operator * (mtx x){mtx ret;
    	RP(t,1,n)RP(i,1,n)RP(k,1,n)
    	    (ret[t][i]+=data[t][k]*x[k][i])%=yyb;
    	return ret;
        }
        inline mtx operator ^ (int x){
    	mtx base=*this,ret;ret.unis();
    	for(register ll t=x;t;t>>=1,base*=base)
    	    if(t&1)ret*=base;
    	return ret;
        }
    }zsy,psj,mona;
    
    
    int main(){
    #ifndef ONLINE_JUDGE
        freopen("t1.in","r",stdin);
        freopen("t1.out","w",stdout);
    #endif
        n=qr(1ll);
        k=qr(1ll)+1ll;
        yyb=qr(1ll);
        RP(t,1,n){
    	RP(i,1,n){
    	    zsy[t][i]=qr(1)%yyb;
    	}
        }
        RP(t,1,n)
    	zsy[t+n][t]=zsy[t+n][t+n]=1;
        n<<=1;
        zsy^=k;
        psj.unis();
        n>>=1;
        RP(t,1,n) zsy[t+n][t]=(zsy[t+n][t]-1+yyb)%yyb;
        RP(t,1,n){
    	RP(i,1,n){
    	    cout<<zsy[t+n][i]<<' ';
    	}
    	cout<<endl;
        }
        return 0;
    }
    /*
      矩阵快速幂板子题
      
      群论板子题
      
      看一下那个<线性代数>的矩阵分割和线性代数的那个方法就会做了
      
         A 0
      I=
         I I
        
    */
    
    
  • 相关阅读:
    什么样的代码称得上是好代码?
    九年程序人生 总结分享
    Docker入门 第一课 --.Net Core 使用Docker全程记录
    阿里云 Windows Server 2012 r2 部署asp.net mvc网站 平坑之旅
    Visual studio 2015 Community 安装过程中遇到问题的终极解决
    Activiti6.0 spring5 工作流引擎 java SSM流程审批 项目框架
    java 进销存 库存管理 销售报表 商户管理 springmvc SSM crm 项目
    Leetcode名企之路
    24. 两两交换链表中的节点
    21. 合并两个有序链表
  • 原文地址:https://www.cnblogs.com/winlere/p/10384725.html
Copyright © 2011-2022 走看看