数据结构-稀疏矩阵的乘法（C语言）

引言

稀疏矩阵是什么？人们无法给出确切的定义，它只是一个凭人们的直觉来了解的概念。

假设在m×n的矩阵中，有t个元素不为零。令σ=，称σ为矩阵的稀疏因子。通常认为σ<=0.05时称为稀疏矩阵。

三元组顺序表

// 稀疏矩阵的三元组顺序表储存表示
typedef struct
{
    int i, j; // 非零元的行下标和列下标
    int e;
} Triple;

typedef struct
{
    Triple data[MAXSIZE + 1]; // 非零元三元组表，data[0]未用
    int ru, cu, nu;           // 矩阵的行数，列数，非零元个数
} TSMatrix;

在讲矩阵乘法之前我们先讨论矩阵的转置，因为其算法很相近。

矩阵转置

转置运算时一种最简单的矩阵运算，对于一个m×n的矩阵M，它的转置矩阵T是一个n×m的矩阵，且T(i,j)=M(j,i),1<=i<=n,1<=j<=m。

对于一个三元组顺序表储存的矩阵M，

1	2	12
1	3	9
3	1	-3
3	6	14
4	3	24
5	2	18
6	1	15
6	4	-7

它的转置矩阵T应该是

1	3	-3
1	6	15
2	1	12
2	5	18
3	1	9
3	4	24
4	6	-7
6	3	14

由图可知，M和转置后的T都是按顺序排列的，我们应该从M中找出列为1(也就是在T中行为1)的元素，由于是顺序表，我们不用担心T中同一行列数的顺序，每找到一个，它就是T中的一个元素，并且M.data[i].i=T.data[cnt].j ; M.data[i].j=T.data[cnt],i ; M.data[i].e=T.data[cnt].e ; cnt++;代码如下：

void TransposeSMatrix(TSMatrix M)
{
    T.ru=M.ru,T.cu=M.cu,T.nu=M.nu;
    int cnt=1;
    for(int j=1;j<=T.cu;j++)  // 以列为单位开始遍历
    {
        for(int i=1;i<=M.nu;i++)  //遍历M的每一个元素
        {
            if(M.data[i].j==j)
            {
                T.data[cnt].i=j;
                T.data[cnt].j=M.data[i].i;
                T.data[cnt].e=M.data[i].e;
                cnt++;
            }
        }
    }
    for(int i=1;i<=T.nu;i++)
    {
        printf("%d %d %d\n",T.data[i].i,T.data[i].j,T.data[i].e);
    }
}

矩阵转置的改进

上面算法的时间复杂度为O(cu*nu),我们在每次开始新一列的时候都会遍历一遍M矩阵，T矩阵是从1开始的，但是如果我们知道M矩阵中每一列元素的第一个非零元在T中的位置以及每一列非零元的个数，问题就解决了！

此时我们只需遍历M矩阵的每一个元素，就知道对应列在T矩阵的位置，定义T，此时T矩阵不是从1开始的

我们有cpot[M.cu]表示M矩阵某一列第一个非零元在T中的位置，num[MAXSIZE]表示某一列非零元的个数。

int num[MAXSIZE];
void TransposeSMatrix_(TSMatrix M)
{
    T.ru = M.ru, T.cu = M.cu, T.nu = M.nu;
    int cnt = 1;
    int cpot[T.cu];
    for (int i = 1; i <= M.nu; i++)
    {
        num[M.data[i].j]++;
    }
    cpot[1] = 1;
    for (int i = 2; i <= M.cu; i++)
    {
        cpot[i] = cpot[i - 1] + num[i - 1];
    }
    for (int i = 1; i <= M.nu; i++)
    {
        T.data[cpot[M.data[i].j]].i=M.data[i].j;//M矩阵中第j列的元素在T中的位置
        T.data[cpot[M.data[i].j]].j=M.data[i].i;
        T.data[cpot[M.data[i].j]].e=M.data[i].e;
        cpot[M.data[i].j]++;
    }
    for (int i = 1; i <= T.nu; i++)
    {
        printf("%d %d %d\n", T.data[i].i, T.data[i].j, T.data[i].e);
    }
}

矩阵乘法

根据线性代数的知识，我们知道两个矩阵相乘：

Q=M×T

其中，M是m1×n1矩阵，T是m2×n2矩阵，当n1=m2时有：

for (i = 1; i <= m1; i++)
{
    for (j = 1; j <= n2; j++)
    {
        Q[i][j] = 0;
        for (k = 1; k <= n1; k++)
            Q[i][j] += M[i][k] * T[k][j];
    }
}

此时算法的时间复杂度时O(m1*n2*n1)。

在三元组顺序表中该乘法变现为：

void MultSMatrix(TSMatrix M, TSMatrix T)
{
    if (M.cu == T.ru)
    {
        Q.ru = M.ru, Q.cu = T.cu;
        int cnt = 1;
        for (int i = 1; i <= M.nu; i++)
        {
            for (int j = 1; j <= T.nu; j++)
            {
                if (M.data[i].j == T.data[j].i)
                {
                    Q.data[cnt].i = M.data[i].i;
                    Q.data[cnt].j = T.data[i].j;
                    Q.data[cnt].e = (M.data[i].e) * (T.data[j].e);
                    cnt++;
                }
            }
        }
        Q.nu = cnt - 1;
        for (int i = 1; i <= Q.nu; i++)
            printf("%d %d %d\n", Q.data[i].i, Q.data[i].j, Q.data[i].e);
    }
}

举个例子，比如M=

1	1	3
1	4	5
2	2	-1
3	1	2

T=

1	2	2
2	1	1
3	1	-2
3	2	4

在做矩阵乘法时，我们遍历到了M和T中的每一个元素，当M中元素的列数和T中元素的行数相同时，他们就可以相乘，此时Q矩阵的行数等于M矩阵的行数，Q矩阵的列数等于T矩阵的列数，Q矩阵的值等于M和T的值的乘积。

此时算法的时间复杂度是O(nu*nu)。

矩阵乘法的改进

我们想一想刚才的转置，我们在相乘时是不是也可以保存一些值来减少我们的计算量呢，其实思路和转置一摸一样！

我们如果知道M矩阵每一列元素在T中的位置，我们就不需要费那么大劲去找了！

我们有cpos[M.cu]为M矩阵中某一列在T矩阵中对应行的位置，num[MAXSIZE]为M矩阵某一列在T矩阵中对应行的个数

int num[MAXSIZE];
void MultSMatrix_(TSMatrix M, TSMatrix T)
{
    if (M.cu == T.nu)
    {
        int cpos[T.ru + 1];
        int cnt = 1;
        for (int i = 1; i <= T.nu; i++)
        {
            num[T.data[i].i]++;
        }
        cpos[1] = num[1];
        for (int i = 2; i <= T.ru; i++)
        {
            cpos[i] = cpos[i - 1] + num[i - 1];
        }
        for (int i = 1; i <= M.nu; i++)
        {
            int tmp = cpos[M.data[i].j];
            if (tmp && tmp < T.nu)
            {
                for (int j = 1; j <= num[M.data[i].j]; j++)
                {
                    Q.data[cnt].i = M.data[i].i;
                    Q.data[cnt].j = T.data[tmp].j;
                    Q.data[cnt].e += (M.data[i].e) * (T.data[tmp].e);
                    tmp++;
                }
                cnt++;
            }
        }
        Q.nu = cnt - 1;
        for (int i = 1; i <= Q.nu; i++)
            printf("%d %d %d\n", Q.data[i].i, Q.data[i].j, Q.data[i].e);
    }
}

对于M矩阵每一列在T矩阵对应行的乘积运算应该是累加的！

此时的时间复杂度是O(nu+ru+nu*num),其中num为矩阵T每一行非零元的个数。