zoukankan      html  css  js  c++  java
  • [ACM] POJ 3233 Matrix Power Series (求矩阵A+A^2+A^3...+A^k,二分求和或者矩阵转化)

    Matrix Power Series
    Time Limit: 3000MS   Memory Limit: 131072K
    Total Submissions: 15417   Accepted: 6602

    Description

    Given a n × n matrix A and a positive integer k, find the sum S = A + A2 + A3 + … + Ak.

    Input

    The input contains exactly one test case. The first line of input contains three positive integers n (n ≤ 30), k (k ≤ 109) and m (m < 104). Then follow n lines each containing n nonnegative integers below 32,768, giving A’s elements in row-major order.

    Output

    Output the elements of S modulo m in the same way as A is given.

    Sample Input

    2 2 4
    0 1
    1 1

    Sample Output

    1 2
    2 3

    Source

    POJ Monthly--2007.06.03, Huang, Jinsong


    解题思路:

    第一种方法:

    题意为给定矩阵A(以下代码中用ori表示)。以及k, mod ,求 A+A^2+A^3+......A^k    的和对mod取余。

    一開始用循环k次,递推的做法,超时。

    。。

    看了解题报告,求和的时候要用到二分求和。

    所求的和用s(k)表示。

    当k为偶数时:

    比方 k=6,那么  A+A^2+A^3+A^4+A^5+A^6=    A+A^2+A^3+   A^3*(A+A^2+A^3)   

    s(k)=s(k/2)+A^(n/2) * s(k/2) 即s(k)=(E+A^(n/2))*s(n/2)  (E为单位矩阵)

    当k为奇数时:

    s(k)=s(k-1)+A^k ,  那么k-1为偶数。能够依照上面的二分

    PS:代码要写细致啊,否则一个小错误查半天.....计算两个矩阵相乘时ret.arr[i][j]+=a.arr[i][k]*b.arr[k][j]; 居然写成了ret.arr[i][j]+=a.arr[i][k]*a.arr[k][j]; T T

    代码:

    #include <iostream>
    #include <stdio.h>
    #include <string.h>
    using namespace std;
    const int maxn=31;
    int n,k,mod;
    
    
    struct mat
    {
        int arr[maxn][maxn];
        mat()
        {
            memset(arr,0,sizeof(arr));
        }
    };
    
    
    mat mul(mat a,mat b)
    {
        mat ret;
        for(int i=0;i<n;i++)
            for(int k=0;k<n;k++)
        {
            if(a.arr[i][k])
                for(int j=0;j<n;j++)
            {
                ret.arr[i][j]+=a.arr[i][k]*b.arr[k][j];
                if(ret.arr[i][j]>=mod)
                    ret.arr[i][j]%=mod;
            }
        }
        return ret;
    }
    
    mat add(mat a,mat b)
    {
        mat an;
        for(int i=0;i<n;i++)
            for(int j=0;j<n;j++)
            {
                an.arr[i][j]=a.arr[i][j]+b.arr[i][j];
                if(an.arr[i][j]>=mod)
                    an.arr[i][j]%=mod;
            }
        return an;
    }
    
    mat power(mat p,int k)
    {
        if(k==1) return p;
        mat e;
        for(int i=0;i<n;i++)
            e.arr[i][i]=1;
        if(k==0) return e;
        while(k)
        {
            if(k&1)
                e=mul(p,e);
            p=mul(p,p);
            k>>=1;
        }
        return e;
    }
    
    void output(mat ans)
    {
        for(int i=0;i<n;i++)
            for(int j=0;j<n;j++)
        {
            if(j==n-1)
                cout<<ans.arr[i][j]<<endl;
            else
                cout<<ans.arr[i][j]<<" ";
        }
    }
    
    mat  cal(mat ori,int k)
    {
        if(k==1) return ori;
        if(k&1)
            return add(cal(ori,k-1),power(ori,k));//当k为奇数时,减1变为偶数 S(K)=S(K-1)+ori^K
        else
            return mul(add(power(ori,0),power(ori,k>>1)),cal(ori,k>>1));
            //当K为偶数时,S(K)=(1+ori^(K/2))*S(K/2)
    }
    
    int main()
    {
        while(scanf("%d%d%d",&n,&k,&mod)!=EOF)
        {
            mat ori,ans;
            for(int i=0;i<n;i++)
                for(int j=0;j<n;j++)
            {
                scanf("%d",&ori.arr[i][j]);
                if(ori.arr[i][j]>=mod)
                    ori.arr[i][j]%=mod;
            }
            ans=cal(ori,k);
            output(ans);
        }
        return 0;
    }
    

    另外一种方法:

    解题思路转载自:http://www.cnblogs.com/jiangjing/archive/2013/05/28/3103336.html

    分析:求a^1+..a^n这是矩阵乘法中关于等比矩阵的求法:

    |A  E|

    |0  E|

    当中的A为m阶矩阵。E是单位矩阵,0是零矩阵。而我们要求的是:                                                                             

    A^1+A^2+..A^L。由等比矩阵的性质

    |A  ,  1|                 |A^n , 1+A^1+A^2+....+A^(n-1)|

    |0  ,  1| 的n次方等于|0     ,                   1                   |

    所以我们仅仅须要将A矩阵扩大四倍,变成如上形式的矩阵B。然后开L+1次方就能够得到1+A^1+A^2+....+A^L。

    因为多了一个1。所以最后得到的答案我们还要减去1。

    同理我们把矩阵A变成B:

                                                              |A  E|

                                                              |0  E|

    然后我们就是求B的n+1次幂之后得到的矩阵为|x1   x2|

                                                              |x3   x4|

    右上角的矩阵x2减去单位矩阵E,得到就是要求的矩阵了。


    #include <iostream>
    #include <string.h>
    using namespace std;
    const int maxn=62;
    int n,k,mod;
    
    struct mat
    {
        int arr[maxn][maxn];
        mat()
        {
            memset(arr,0,sizeof(arr));
        }
    };
    
    mat mul(mat a,mat b)
    {
        mat ans;
        for(int i=0;i<n;i++)
            for(int k=0;k<n;k++)
        {
            if(a.arr[i][k])
                for(int j=0;j<n;j++)
            {
                ans.arr[i][j]+=a.arr[i][k]*b.arr[k][j];
                if(ans.arr[i][j]>=mod)
                    ans.arr[i][j]%=mod;
            }
        }
        return ans;
    }
    
    mat power(mat p,int k)
    {
        if(k==1) return p;
        mat e;
        for(int i=0;i<n;i++)
            e.arr[i][i]=1;
        if(k==0) return e;
        while(k)
        {
            if(k&1)
                e=mul(e,p);
            p=mul(p,p);
            k>>=1;
        }
        return e;
    }
    
    void output(mat ans)
    {
        for(int i=0;i<n;i++)
            for(int j=0;j<n;j++)
        {
            if(i==j)
            {
                ans.arr[i][j+n]--;
                while(ans.arr[i][j+n]<0)
                    ans.arr[i][j+n]+=mod;
            }
            if(j==n-1)
                cout<<ans.arr[i][j+n]<<endl;
            else
                cout<<ans.arr[i][j+n]<<" ";
        }
    }
    
    int main()
    {
        while(cin>>n>>k>>mod)
        {
            mat ori,ans;
            for(int i=0;i<n;i++)
                for(int j=0;j<n;j++)
            {
                cin>>ori.arr[i][j];
                if(i==j)
                {
                    ori.arr[i][j+n]=1;
                    ori.arr[i+n][j+n]=1;
                }
            }
            n*=2;
            ans=power(ori,k+1);
            n/=2;
            output(ans);
        }
        return 0;
    }
    


  • 相关阅读:
    BZOJ:4219: 跑得比谁都快 3007: 拯救小云公主
    BZOJ:4816: [Sdoi2017]数字表格
    BZOJ:4333: JSOI2012 智者的考验
    BZOJ:3911: SGU383 Caravans(三角剖分)
    bzoj:2595: [Wc2008]游览计划
    ZOJ3602:Count the Trees
    A Dangerous Maze (II) LightOJ
    Where to Run LightOJ
    Lights inside 3D Grid LightOJ
    Snakes and Ladders LightOJ
  • 原文地址:https://www.cnblogs.com/mthoutai/p/7235218.html
Copyright © 2011-2022 走看看