矩阵快速幂

16:05:34 2016-05-07

对A^n,我们一般通过连乘(n-1)次,但是我们利用矩阵乘法的结合律做一下简单的改进就能减少连乘的次数，例如,A*A*A*A*A*A  =>  (A*A)*(A*A)*(A*A),可使得连乘次数由5次减少为3次,那么究竟如何利用结合律可以得到最小的连乘次数呢？

答案是:二进制和二分思想的结合

例如,A^21  =>  （A^16）*（A^4）*（A^1）（21=2^4+2^2+2^0）,A^16又可以二分为(A^8)*(A^8),A^8=(A^4)*(A^4)....

从而可将因子数由n降到log(n)

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <string.h>
typedef long long ll;

typedef struct{
    ll m[3][3];
}Mat;

void Ran_Unit_Mat(Mat &a,Mat &e){ 
    int i,j;
    srand((unsigned int)time(NULL));//NULL not null 
    for(i=0;i<3;i++)
        for(j=0;j<3;j++){
            a.m[i][j] = rand()%100;
            e.m[i][j] = (i == j);
        }
}

Mat MatMul(Mat m1,Mat m2){
    Mat rm;
    memset(rm.m,0,sizeof(rm.m)); 
    int i,j,k; 
    for(i=0;i<3;i++)
        for(j=0;j<3;j++)
            for(k=0;k<3;k++)
                rm.m[i][j] += m1.m[i][k]*m2.m[k][j];
    return rm;
}

Mat Matlab(Mat a,Mat e,int n){
    Mat b;
    b = e;
    //利用二进制操作(&1和>>)对n进行分解,同时累乘结果 
    while(n){
        if(n&1)
            b = MatMul(b,a);
        n >>= 1;
        a = MatMul(a,a);
    }
    return b;
}

void MatShow(Mat b){
    int i,j;
    for(i=0;i<3;i++){
        for(j=0;j<3;j++)
            printf("%-12lld",b.m[i][j]);
        printf("\n");
    } 
}

int main(){
    Mat a,b,e;
    int n;
    printf("幂次数:");
    scanf("%d",&n);
    Ran_Unit_Mat(a,e);
    MatShow(a); 
    b = Matlab(a,e,n);
    MatShow(b);
    return 0;
}

参照:http://www.cnblogs.com/yan-boy/archive/2012/11/29/2795294.html