UVA 11149-Power of Matrix （等比矩阵求和）

又一个矩阵的知识点，这个没什么好说的，按着公式的来。

这里放出两个模板，后面是白皮书的。

书上更加简洁些~

#include<cstdio>
#include<iostream>
#include<cstring>
using namespace std;

int n,k;
struct matrix {
    int m[55][55];
    matrix () {
        memset(m,0,sizeof(m));
    }
};

matrix mul(matrix a,matrix b) {
    matrix tmp;
    for(int i = 1;i <= n;i++)
        for(int j = 1;j <= n;j++)
            if(a.m[i][j]) {
                for(int k = 1;k <= n;k++)
                    tmp.m[i][k] = (tmp.m[i][k] + a.m[i][j] * b.m[j][k]) % 10;
            }
    return tmp;
}

void pi(matrix base) {
    for(int i = 1;i <= n;i++) {
        for(int j = 1;j < n;j++)
            printf("%d ",base.m[i][j]);
        printf("%d\n",base.m[i][n]);
    }
    printf("\n");
}

matrix powmul(matrix a,int k) {
    matrix tmp;
    for(int i = 1;i <= n;i++)
        tmp.m[i][i] = 1;
    while(k) {
        if(k & 1)
            tmp = mul(tmp,a);
        a = mul(a,a);
        k >>= 1;
    }
    return tmp;
} 

matrix add(matrix a,matrix b) {
    matrix tmp;
    for(int i = 1;i <= n;i++)
        for(int j = 1;j <= n;j++)
            tmp.m[i][j] = (a.m[i][j] + b.m[i][j]) % 10;
    return tmp;
}

matrix getsum(matrix e,int k) {
    matrix s,a,b;
    if(k == 1)
        return e;
    a = getsum(e,k/2);
    s = add(a,mul(a,powmul(e,k/2))); //等比矩阵求和模板
    if(k & 1)
        s = add(s,powmul(e,k));
    return s;
}

int main(){
    while(scanf("%d%d",&n,&k) != EOF && n) {
        matrix base,ans;
        for(int i = 1;i <= n;i++)
            for(int j = 1;j <= n;j++) {
                int x;
                scanf("%d",&x);
                base.m[i][j] = x%10; //一开始就要取余，否则会溢出！
            }
        ans = getsum(base,k);
        pi(ans);
    }
    return 0;
}

再研究一下白皮书的方法~，感觉这个更加简洁。

#include<cstdio>
#include<iostream>
#include<cstring>
using namespace std;

int n,k;
struct matrix {
    int m[85][85];
    matrix () {
        memset(m,0,sizeof(m));
    }
};

matrix mul(matrix a,matrix b) {
    matrix tmp;
    for(int i = 1;i <= 2*n;i++)
        for(int j = 1;j <= 2*n;j++)
            if(a.m[i][j]) {
                for(int k = 1;k <= 2*n;k++)
                    tmp.m[i][k] = (tmp.m[i][k] + a.m[i][j] * b.m[j][k]) % 10;
            }
    return tmp;
}

void pi(matrix base) {
    for(int i = 1;i <= 2*n;i++) {
        for(int j = 1;j < 2*n;j++)
            printf("%d ",base.m[i][j]);
        printf("%d\n",base.m[i][2*n]);
    }
    printf("\n");
}

matrix powmul(matrix a,int k) {
    matrix tmp;
    for(int i = 1;i <= 2*n;i++)
        tmp.m[i][i] = 1;
    while(k) {
        if(k & 1)
            tmp = mul(tmp,a);
        a = mul(a,a);
        k >>= 1;
    }
    return tmp;
} 

void getsum(matrix a,int k) {
    matrix b;
    for(int i = 1;i <= n;i++) {
        for(int j = 1;j <= n;j++) {
            b.m[i][j] = a.m[i][j];
        }
        b.m[n+i][i] = b.m[n+i][n+i] = 1;
    }
    b = powmul(b,k+1);
    //pi(b);
    for(int i = 1;i <= n;i++) {
        for(int j = 1;j <= n;j++) {
            int t = b.m[n+i][j] % 10;
            if(i == j) t = (t + 10 - 1) % 10;
            printf("%d%c",t,j == n ? '\n' : ' ');
        }
    }
}

int main(){
    while(scanf("%d%d",&n,&k) != EOF && n) {
        matrix base,ans;
        for(int i = 1;i <= n;i++)
            for(int j = 1;j <= n;j++) {
                int x;
                scanf("%d",&x);
                base.m[i][j] = x%10;
            }
        getsum(base,k);
        printf("\n");
        //pi(ans);
    }
    return 0;
}

数组要开到2*n，相应的快速幂也是

void getsum(matrix a,int k) {
    matrix b;
    for(int i = 1;i <= n;i++) {
        for(int j = 1;j <= n;j++) {
            b.m[i][j] = a.m[i][j];
        }
        b.m[n+i][i] = b.m[n+i][n+i] = 1;
    }
    b = powmul(b,k+1);
    for(int i = 1;i <= n;i++) {
        for(int j = 1;j <= n;j++) {
            int t = b.m[n+i][j] % 10;
            if(i == j) t = (t + 10 - 1) % 10;
            printf("%d%c",t,j == n ? '\n' : ' ');
        }
    }
}