稀疏矩阵
本章描述的函数是用来构造和操作稀疏矩阵,稀疏矩阵有大量的零,只包含少量的非零元素。稀疏矩阵经常出现在偏微分方程的解中。使用专门的数据结构和算法来存储和处理稀疏矩阵是有益的,因为当应用于稀疏矩阵时,密集矩阵算法和结构可能非常慢,并使用大量内存。
头文件gsl_spmatrix.h包含稀疏矩阵函数的原型和相关声明。
43.1 数据类型
所有函数都可用于每种标准数据类型。double的版本有前缀gsl_spmatrix,类似于单精度浮点数组的版本有前缀gsl_spmatrix_float。可用类型的完整列表如下所示,
前缀 |
类型 |
gsl_spmatrix |
double |
gsl_spmatrix_float |
float |
gsl_spmatrix_long_double |
long double |
gsl_spmatrix_int |
int |
gsl_spmatrix_uint |
unsigned int |
gsl_spmatrix_long |
long |
gsl_spmatrix_ulong |
unsigned long |
gsl_spmatrix_short |
short |
gsl_spmatrix_ushort |
unsigned short |
gsl_spmatrix_char |
char |
gsl_spmatrix_uchar |
unsigned char |
gsl_spmatrix_complex |
complex double |
gsl_spmatrix_complex_float |
complex float |
gsl_spmatrix_complex_long_double |
complex long double |
43.2稀疏矩阵存储格式
GSL目前支持三种稀疏矩阵的存储格式:坐标表示(COO)、压缩稀疏列(CSC)和压缩稀疏行(CSR)格式。下面将更详细地讨论这些问题。为了说明不同的存储格式,以下章节将引用这个M乘N的稀疏矩阵,其中M = 4, N = 5:
矩阵中非零元素的个数,在本例中缩写为nnz,等于10。
43.2.1 坐标存储(COO)
坐标存储格式,也称为三元格式,为矩阵的每个非零元素存储三联体(i, j, x)。这个符号表示矩阵A的(i, j)元素是Aij=x.矩阵存储在三个相同长度的数组中,分别表示行索引、列索引和矩阵数据。对于上面的参考矩阵,一种可能的存储格式是:
data |
9 |
7 |
4 |
8 |
-3 |
-1 |
8 |
5 |
6 |
4 |
row |
0 |
1 |
1 |
2 |
0 |
2 |
2 |
3 |
3 |
3 |
col |
0 |
1 |
0 |
1 |
4 |
2 |
3 |
2 |
3 |
0 |
注意,这种表示不是唯一的——坐标三元可以以任何顺序出现,并且仍然表示相同的稀疏矩阵。这三个数组的长度等于矩阵中非零元素的个数nnz,在本例中为10。这种坐标格式对于稀疏矩阵组装,即在稀疏矩阵中添加新元素或改变现有元素的过程是非常方便的。然而,它一般不适合有效地实现矩阵-矩阵乘积,或矩阵分解算法。对于这些应用程序,最好使用下面讨论的一种压缩格式。
为了方便高效的稀疏矩阵组装,除了上面讨论的三个数组外,GSL还将坐标数据存储在一棵平衡的二叉搜索树中,具体来说是一棵AVL树。这允许GSL有效地确定一个条目(i, j)是否已经存在于矩阵中,并将一个现有的矩阵条目替换为一个新的值,而无需搜索未排序的数组。
43.2.2 压缩稀疏列(CSC)
压缩稀疏列存储将稀疏矩阵中每个非零值的列存储在一个连续的内存块中,保持指向该内存块中每个列的开始的指针,并存储每个非零元素的行索引。对于上面的参考矩阵,这些数组如下所示,
data |
9 |
4 |
4 |
7 |
8 |
-1 |
5 |
8 |
6 |
-3 |
row |
0 |
1 |
3 |
1 |
2 |
2 |
3 |
2 |
3 |
0 |
col_ptr |
0 |
3 |
5 |
7 |
9 |
10 |
data和row数组的长度为nnz,与COO存储格式相同。col_ptr数组的长度为N+1,col_ptr[j]给出了data的j列的开始索引。因此,矩阵的第j列被存储在data[col_ptr[j],data[col_ptr[j] + 1],…,data[col_ptr[j+1] - 1]。col_ptr的最后一个元素是nnz。
43.2.3 压缩稀疏行(CSR)
压缩行存储将每一行非零值存储在一个连续的内存块中,保持指针指向块中每一行的开始,并存储每个非零元素的列索引。对于上面的参考矩阵,这些数组如下所示,
data |
9 |
-3 |
4 |
7 |
8 |
-1 |
8 |
4 |
5 |
6 |
row |
0 |
4 |
0 |
1 |
1 |
2 |
3 |
0 |
2 |
3 |
col_ptr |
0 |
2 |
4 |
7 |
10 |
data和col数组的长度为nnz,与COO存储格式相同。row_ptr数组的长度为M+1,而row_ptr[i]给出了data的i行的开始索引。因此,矩阵的第i行存储在data[row_ptr[i]], data[row_ptr[i] + 1],…,数据[row_ptr[i+1] - 1]。row_ptr的最后一个元素是nnz。