【C++算法】稀疏矩阵的一次定位快速转置法

稀疏矩阵，指rows行cols列的矩阵中，非零元素数量非常少的情况。对于这种矩阵，用一个线性空间去表示它会造成不必要的浪费，因此要使用一个三元组表去表示它。

现在假设有如下这个矩阵A：

索引号(i, j)	1	2	3	4	5
1	3	0	0	0	0
2	0	0	0	3	0
3	0	7	0	0	0
4	0	0	0	1	0
5	0	0	0	0	4
6	2	0	0	0	0

此矩阵有5个非零元，那么相应三元组表的结构为：

索引	行号	列号	元素值
0	1	1	3
1	2	4	3
2	3	2	7
3	4	4	1
4	5	5	4
5	6	1	2

有0~4的索引，表明非零元的数目，行号、列号、元素值一一对应了矩阵中非零元的位置。

进一步地，该三元组表用什么数据结构来表示呢？可以使用线性表，也可以使用二维数组，但此处考虑使用结构体更为方便。定义结构体matrixTerm，它有数据成员row、col、value，分别对应三元表的每一行。如果声明一个结构体类型的指针，根据实际的非零元素数目来分配matrixTerm所指空间的大小，就可以通过索引[i]来对应到上述三元表中的每一行位置。

结构体matrixTerm的定义如下：

template <class T>
struct matrixTerm
{
	int row;
	int col;
	T value;
	operator T() const { return value; }
};

本文使用稀疏矩阵类sparseMatrix来标识稀疏矩阵的私有成员以及相应的操作：

template <class T>
class sparseMatrix
{
private:
	int rows;
	int cols; 
	matrixTerm<T>* terms; 
	int size;
public:
	sparseMatrix() : rows(0), cols(0), size(0) {} //默认构造函数
	sparseMatrix(int row, int col) : rows(row), cols(col), size(0) {}     //带参的构造函数
	void transpose();  //转置函数
	void outputAll();  //输出矩阵所有的值
	void updateAll();  //更新矩阵所有的值
	sparseMatrix<T> operator=(const sparseMatrix<T>& m);   //重载=实现类对象的赋值
};

还是以6行5列的矩阵A来举例，需要在转置的实现里声明一个结果矩阵对象result，来储存矩阵A的行数rows、列数cols、非零元数目size、转置后的三元组表terms（注意，新的三元组表必须按照实际转置后矩阵的行主映射来排序）。我们把当前矩阵与结果矩阵分别定义为A、B，分别为：

矩阵A	1	2	3	4	5
1	3	0	0	0	0
2	0	0	0	3	0
3	0	7	0	0	0
4	0	0	0	1	0
5	0	0	0	0	4
6	2	0	0	0	0

矩阵B	1	2	3	4	5	6
1	3	0	0	0	0	2
2	0	0	7	0	0	0
3	0	0	0	0	0	0
4	0	3	0	1	0	0
5	0	0	0	0	4	0

矩阵A的三元组表为termsA：

索引	行号	列号	元素值
0	1	1	3
1	2	4	3
2	3	2	7
3	4	4	1
4	5	5	4
5	6	1	2

矩阵B的三元组表为termsB：

索引	行号	列号	元素值
0	1	1	3
1	1	6	2
2	2	3	7
3	4	2	3
4	4	4	1
5	5	5	4

可以发现，转置操作的主要步骤为：
①把A的行号、列号分别存到B的列号、行号；
②把三元组表termsA的行号和列号分别交换得到termsA的行号和列号；
③对得到的termsB按照行主映射进行排序。

其中最复杂的就是第③步，以下主要分析第③步的实现：

（1）定义两个int类型、cols + 1大小的数组num和position，分别对应矩阵A第[1]列到第[cols]列的位置（第0列没有）。

数组num用于指定矩阵A第i列中非零元的个数，在本例中，有：
num[5] = { φ, 2, 1, 0, 2, 1 };
数组position用于指定第i列的第一个非零元素在三元组表termsB中应在的索引位置。在本例中，有：
position[5] = { φ, 0, 2, φ, 4, 5 };
其中φ表示不存在第0列，或者第i列没有第一个非零元。

首先要获得这num数组的值，具体的实现为：

for (int i = 0; i <= cols; i++)
		num[i] = 0;  //全部初始化为0
for (int i = 0; i < size; i++)  
		num[terms[i].col]++;  //遍历termsA，termsA有size个元素
		//也即矩阵A有size个非零元

num[terms[i].col+]++表明，第terms[i].col列存在元素，把这列对应的num数组中的值加1。
如此之后，便可得到num数组的所有值。

下一步要初始化position数组。
先附上代码：

position[1] = 0;
for (int i = 2; i <= cols; i++)
	position[i] = position[i - 1] + num[i - 1];

这三行代码比较难理解，这也是此算法的重点内容。
第一行position[1] = 0，指出第一列的第一个非零元一定在termsB的第零号元素，这是逻辑上显然成立的。第三行以此为基础进行迭代。
for循环中从第二列开始“遍历”矩阵A，分析一下i = 2时，position[2]为第2列第一个非零元应在的索引位置，它显然等于 第一列第一个非零元所应在的索引位置 + 第一列的非零元个数，这也是逻辑上显然成立的。对于之后的每一个i，也显然都有 此列第一个非零元所应在的索引位置 = 上一列第一个非零元所应在的索引位置 + 上一列的非零元素个数。
因此，公式position[i] = position[i - 1] + num[i - 1]是显然的。可以发现，num数组只是为了求position数组而服务的，后续不会再用到num数组。
如此，这3行代码解析完毕。

（2）之后，仅根据position数组来扫描矩阵A，实现 一边扫描、一边定位，也即“一次定位的快速转置法”。
先给出代码：

//result是结果矩阵
int index;  //用来记录第i列第一个非零元应在的索引位置
for (int i = 0; i < size; i++)
	{                       
		index = position[terms[i].col]; 
		result.terms[index].row = terms[i].col;  
		result.terms[index].col = terms[i].row;
		result.terms[index].value = terms[i].value;
		position[terms[i].col]++;
	}
	*this = result;  //最后，进行赋值运算即可。=重载自己编写即可，且此处必     须使用深复制

我们这里遍历的是矩阵A的三元组表termsA，因此 i 从0到size遍历。
每一次进入for循环都把index赋值为第i列第一个非零元应在的索引位置，在位置找好后，把矩阵B的三元组表进行相应的映射，即可完成该索引的termsB空间初始化。
之后，把index自增1。这点非常重要！很多人不理解这句话是什么意思，这里解释一下。
在矩阵A中，可以看到第一列有两个元素，分别在第一行和第六行。在position[1]中只指出了第一行元素应该对应的索引值，却没有指出第六行元素应该对应的索引值。因此，在进行第一个元素的映射后，此元素就没有用了，该列后面的第一个元素就可以顺理成章地成为“该列的第一个元素” 。并且这个元素所对应的termsB索引值是紧挨着上一元素的，position[terms[i].col]++就这这个道理。

如果以矩阵A第一列来举例，position[1] = 0，赋值后position[1] = 1。之后当i = 5时，col = 1，此时会求得position[1] = 1，并把其映射到索引为1的位置。

同时还要注意：由于termsA是按照行主映射排列的，因此某一列的第一个非零元素一定排在该列其余非零元数索引的前面！

完整的代码为：

template <class T>
void sparseMatrix<T>::transpose()    //一次定位快速转置法
{
	sparseMatrix<T> result;  //结果矩阵
	//使结果矩阵的行数等于本矩阵的列数，列数等于本矩阵的行数，size大小相等
	result.rows = cols;
	result.cols = rows;
	result.size = size;
	result.terms = new matrixTerm<T>[size];

	//快速转置的原理是，先确定当前矩阵每一列的第一个非零元素，
	//然后直接把它放到结果矩阵对应行的第一个非零元素

	//num指当前矩阵各个列中非零元素的个数
	int* num = new int[cols + 1]; 
	//position指当前矩阵的每一列的第一个非零元素在结果矩阵三元表中的恰当位置
	int* position = new int[cols + 1];
	//因为没有第0列，所以均从1开始遍历

	//寻找*this中每一列的非零元的数目
	for (int i = 0; i <= cols; i++)
		num[i] = 0;
	for (int i = 0; i < size; i++)  //遍历termsA
		num[terms[i].col]++;    

	
	//寻找结果矩阵的每一行的起点
	position[1] = 0;    //表示第一列的第一个非零元素在三元组表的第0个位置
	for (int i = 2; i <= cols; i++)
		position[i] = position[i - 1] + num[i - 1];
	//该公式表明，第i列的第一个非零元素，在三元组表中上一列的第一个非零元素所在位置
	//加上上一列的元素数目

	
	int index;  //用来记录第i列第一个非零元应在的索引位置
	//之后扫描*this,将结果矩阵的三元表复制到result
	for (int i = 0; i < size; i++)
	{                         //指向下一个列号为terms[i]的非零元素在结果三元表中的存放位置
		index = position[terms[i].col]++;  //j代表当前三元表第i个元素应该在的位置
		result.terms[index].row = terms[i].col;  
		result.terms[index].col = terms[i].row;
		result.terms[index].value = terms[i].value;
		//position[terms[i].col]++;
	}
	*this = result;
}

个人感觉应该比较详尽了，若有疑问请留言。