稀疏矩阵的三元组顺序表存储及矩阵相乘算法小结
稀疏矩阵的三元组顺序表存储及矩阵相乘算法小结巧若拙（欢迎转载，但请注明出处：http://blog.）
一：稀疏矩阵的三元组顺序表数据结构
typedef int ElemType;
typedef struct
{
int x, y; //该非零元素的行下标和列下标
ElemType e; //该非零元素的值
} Triple;
typedef struct
{
Triple data[MAXSIZE]; //非零元素三元组顺序表
int mu, nu, tu; //矩阵的行数，列数和非零元素的个数
} TSMatrix;
二：三元组顺序表实现矩阵转置
显然，一个稀疏矩阵的转置矩阵仍然是稀疏矩阵。假设a和b是TSMatrix型的变量，分别表示矩阵M和T。那么，如何由a得到b呢？
从分析a和b之间的差异可见只要做到：1。将矩阵的行列值互换；2。将每个三元组中的i和j互换；3。重排三元组之间的次序便可以实现矩阵的转置。
前两条是容易做到的，关键是如何实现第三条，即如何使b.data中的三元组是以T的行序（M的列序）为主序依次排列的。
可以有两种处理方法：
1， 按照b.data中三元组的次序依次在a.data中找到相应的三元组进行转置。即按照M的列序来进行转置。代码如下：
TSMatrix TransposeSMatrix(TSMatrix M)//三元组顺序表 转置矩阵
{
int k, i, col;
TSMatrix T;
T.mu = M.nu;
T.nu = M.mu;
T.tu = M.tu;
if (T.tu > 0)
{
k = 0;
for (col=1; col<=M.nu; col++) //按照T的行序（M的列序）为主序依次排列
for (i=0; i<M.tu; i++)//扫描M的所有元素
if (M.data[i].y == col)
{
T.data[k].x = M.data[i].y;
T.data[k].y = M.data[i].x;
T.data[k++].e = M.data[i].e;
}
}
return T;
}
2，按照a.data中三元组的次序进行转置，并将转置后的三元组置入b中恰当的位置。
即预先确定M中每一列（即T中的每一行）的第一个非零元素在b.data中应有的位置。在此，需要附设num和cpot两个数组，num[col-1]表示矩阵M中第col列中非零元素的个数，cpot[col-1]指示M中第col列中第一个非零元素在b.data中的位置，显然有：
cpot[0] = 0；
cpot[col] = cpot[col-1] + num[col-1] 0<col<a.nu
此算法的时间复杂度为O（mu+tu），在M的非零元素的个数tu和mu*nu等数量级时，其时间复杂度为O（mu.nu），和使用二维数组存储进行转置矩阵运算时时间复杂度相同，故称为快速转置。代码如下：
TSMatrix FastTransposeSMatrix(TSMatrix M) //三元组顺序表 快速转置矩阵
{
TSMatrix T;
int k, i, col;
int num[N] = {0};
int cpot[N] = {0};
T.mu = M.nu;
T.nu = M.mu;
T.tu = M.tu;
if (T.tu > 0)
{
for (i=0; i<M.tu; i++)//扫描M的所有元素
num[M.data[i].y-1]++; //确定矩阵M每一列中非零元素的个数
cpot[0] = 0;
for (col=1; col<M.nu; col++)
cpot[col] = cpot[col-1] + num[col-1];// 确定矩阵M第col列中第一个非零元素在b.data中的位置
for (i=0; i<M.tu; i++)//扫描M的所有元素
{
col = M.data[i].y - 1; //标出该元素所在的列
k = cpot[col]; //k即矩阵M第col列（即T中的第col行）中第一个非零元素在b.data中的位置
T.data[k].x = M.data[i].y;
T.data[k].y = M.data[i].x;
T.data[k].e = M.data[i].e;
cpot[col]++; //矩阵M第col列中第一个非零元素在b.data中的位置向前移动一位
}
}
return T;
}
三：三元组顺序表实现矩阵相乘
压缩存储的矩阵与用二维数组存储的矩阵在进行乘法运算时最大的不同是，在经典（二维数组存储）算法中，不论M（i,k）和N（k,j）的值是否为零，都要进行一次乘法计算，这样造成很大的浪费。而压缩存储则只需在M.data和N.data中找到相应的元素相乘即可。
这里介绍三种算法。
第一种，通过给定的行号i和列号j找出原矩阵的对应的元素，设计一个函数Value()，当在三元组顺序表中找到该元素时，返回其元素值，否则说明该位置元素值为0，返回0。然后利用该函数计算出C的行号i和列号j处元素的值，若该值不为0，则存入C，否则不存入。代码如下：
int Value(TSMatrix M, int i, int j) //计算出矩阵M的i行j列处元素的值
{
int k = 0;
while (k < M.tu && M.data[k].x <= i) //因为是按行序排列，故不必扫描行号比i大的元素
{
if (M.data[k].x == i && M.data[k].y == j)
return M.data[k].e;
k++;
}
return 0;
}
TSMatrix MultSMatrix(TSMatrix A, TSMatrix B) //三元组顺序表低效矩阵乘法
{
TSMatrix C;
int i, j, k, l, p, s;
C.mu = A.mu;
C.nu = B.nu;
C.tu = 0;
if (A.tu * B.tu != 0)
{
p = 0;
for (i=1; i<=A.mu; i++)
{
for (j=1; j<=B.nu; j++)
{
s = 0;
for (k=1; k<=A.nu; k++)
s += Value(A, i, k) * Value(B, k, j);
if (s != 0)
{
C.data[p].x = i;
C.data[p].y = j;
C.data[p++].e = s;
}
}
}
C.tu = p;
}
return C;
}
这种算法的存储空间虽然很少，但时间复杂度为O（A.mu*A.nu*B.nu*（A.tu+B.tu）），比经典的算法时间复杂度还高，因此是一种用时间换空间的方法。
如果我们预先确定矩阵的每一行的第一个非零元素在三元组顺序表中的位置，则可以得到改进的算法。
基本操作是对于A中每个元素A.data[p]，找到N中所有满足条件A.data[p].y==B.data[q].x的元素B.data[q]，求得A.data[p].e和B.data[q].e的乘积，当然这个乘积的值只是乘积矩阵C[i][j]中的一部分。
为了便于操作，我们可以对每个元素设置一个累计乘积值的变量，使其初值为0，然后扫描A，求得相应元素的乘积值并累积到相应的累计乘积值的变量中。
由于C中元素的行号和A中元素的行号一致，又都是以行序为主序排列的，因此可以对C进行逐行处理。代码如下：
void RPos(TSMatrix M, int rpos[])//确定矩阵的每一行的第一个非零元素在三元组顺序表中的位置
{
int num[N] = {0};
int i, row;
for (i=0; i<M.tu; i++)//扫描M的所有元素
num[M.data[i].x-1]++; //确定矩阵M每一行中非零元素的个数
rpos[0] = 0;
for (row=1; row<M.mu; row++)
rpos[row] = rpos[row-1] + num[row-1];// 确定矩阵M第row行中第一个非零元素在M.data中的位置
}
TSMatrix FastMultSMatrix_1(TSMatrix A, TSMatrix B)//三元组顺序表快速矩阵乘法之一
{
TSMatrix C;
int Arpos[N] = {0}, Brpos[N] = {0};//分别存储矩阵A，B的每一行的第一个非零元素在三元组中的位置
int ctemp[N] = {0};//存储正在处理的该行中每一列的乘积值，是一个累积的和值，作为矩阵C在该位置处元素的值
int arow, brow, ccol;
int i, p, q;
int ta, tb;
C.mu = A.mu;
C.nu = B.nu;
C.tu = 0;
if (A.tu * B.tu != 0)//若C是非零矩阵
{
RPos(A, Arpos); //确定矩阵A的每一行的第一个非零元素在三元组顺序表中的位置
RPos(B, Brpos); //确定矩阵B的每一行的第一个非零元素在三元组顺序表中的位置
for (arow=1; arow<=A.mu; arow++)//对A进行逐行处理 ，即对C进行逐行处理
{
for (i=0; i<B.nu; i++)//当前各元素累加器清零
ctemp[i] = 0;
if (arow < A.mu) //处理A中第arow行的每一个非零元素，ta指示下一行第一个非零元素的位置
ta = Arpos[arow];
else //最后一行的 下一行第一个非零元素的位置 当然是 A.tu
ta = A.tu;
for (p=Arpos[arow-1]; p<ta; p++) {
brow = A.data[p].y; //标出该元素的列号，在B中找到与它对应的行号
if (brow < B.mu)
tb = Brpos[brow];
else
tb = B.tu;
for (q=Brpos[brow-1]; q<tb; q++)
{
ccol = B.data[q].y - 1;
ctemp[ccol] += A.data[p].e * B.data[q].e;
}
}
for (ccol=0; ccol<C.nu; ccol++)//得到C的第arow中每一列（个）元素的值
{
if (ctemp[ccol] != 0)//压缩存储该行非零元
{
if (C.tu == MAXSIZE)
{
printf("三元组溢出！\n" );
exit(1);
}
C.data[C.tu].x = arow; //C的行数等于A的行数
C.data[C.tu].y =ccol + 1; //C的列数等于B的列数
C.data[C.tu++].e = ctemp[ccol]; //C的元素值等于A的行数ctemp[ccol]的累积值
}
}
}
}
return C;
}
分析上述算法的时间复杂度有如下结果：确定矩阵的每一行的第一个非零元素在三元组顺序表中的位置的时间复杂度为O（A.tu+A.mu+B.tu+B.mu+），累加器ctemp初始化的时间复杂度为O（A.mu*B.nu），C的所有非零元素的时间复杂度为O（A.tu*B.tu/B.mu），进行压缩存储的时间复杂度为O（A.mu*B.tu），因此总的时间复杂度为O（A.mu*B.nu + A.tu*B.tu/B.mu）。当矩阵为稀疏矩阵时，这种算法的时间复杂度接近O（A.mu*B.nu），比经典算法要快。
还有一种快速矩阵乘法，它是先将矩阵B转置，这样矩阵A与B的列数就相同了，矩阵B的行数即矩阵C的列数。此种算法不需要设置数组ctemp[N]，只需直接计算该行对应列的元素乘积的总和即可得到矩阵C的值。此算法时间复杂度为O（A.mu*B.nu *MAX(ta,tb)），其中MAX(ta,tb)表示A矩阵和B的转置矩阵中每行非零元素个数的最大值。
算法通俗易懂，代码也很简洁。代码如下：
TSMatrix FastMultSMatrix_2(TSMatrix A, TSMatrix B)//三元组顺序表快速矩阵乘法之二
{
TSMatrix C;
int Arpos[N] = {0}, Brpos[N] = {0};//分别存储矩阵A，B的每一行的第一个非零元素在三元组中的位置
int arow, brow, ccol;
int i, p, q, ta, tb;
C.mu = A.mu;
C.nu = B.nu;
C.tu = 0;
if (A.tu * B.tu != 0)//若C是非零矩阵
{
B = FastTransposeSMatrix(B); //求矩阵B的转置矩阵
RPos(A, Arpos); //确定矩阵A的每一行的第一个非零元素在三元组顺序表中的位置
RPos(B, Brpos); //确定矩阵B的每一行的第一个非零元素在三元组顺序表中的位置
for (arow=1; arow<=A.mu; arow++)//对A进行逐行处理 ，即对C进行逐行处理
{
ta = (arow < A.mu) ? Arpos[arow] : A.tu;
for (brow=1; brow<=B.mu; brow++)
{
C.data[C.tu].e = 0;
p = Arpos[arow-1];
tb = (brow < B.mu) ? Brpos[brow] : B.tu;
q = Brpos[brow-1];
while (p < ta && q < tb)
{
if (A.data[p].y < B.data[q].y)
p++;
else if (A.data[p].y > B.data[q].y)
q++;
else
C.data[C.tu].e += A.data[p++].e * B.data[q++].e; //C的元素值等于该行对应列的元素乘积的总和
}
if (C.data[C.tu].e != 0)
{
C.data[C.tu].x = arow;
C.data[C.tu++].y = brow;
}
}
}
}
return C;
}
