转载

数组和矩阵 4 - 稀疏矩阵

继续是《数据结构算法与应用：C++语言描述》，第四章数组和矩阵的笔记。本小节介绍稀疏矩阵的内容。这也是本章节最后一节内容。

基本概念

如果一个$m/times n$矩阵中有”许多”元素为0，则称该矩阵为 稀疏矩阵(sparse)

对应的非稀疏的矩阵称为 稠密矩阵(dense) 。而实际上，稀疏矩阵和稠密矩阵之间并没有一个精确的界限。

$n/times n$的对角矩阵和三对角矩阵都是稀疏矩阵，二者都有$O(n)$个非0元素和$O(n^2)$个0元素。而对于一个$n/times n$的三角矩阵，它至少有$/frac{n(n-1)}{2}$个0元素，最多有$/frac{n(n+1)}{2}$个非0元素。在本节中，我们规定一个矩阵是稀疏矩阵，则其非0元素的数目应小于$n^2/3$,有些情况下应小于$n^2/5$,因此可以将三角矩阵视为稠密矩阵。

诸如对角矩阵和三对角矩阵这样的稀疏矩阵，其非0区域的结构很有规律，因此可以设计一个很简单的存储结构，该存储结构的大小就等于矩阵非0区域的大小，本小节主要考察具有不规则非0区域的稀疏矩阵。

数组描述

对于下图a中的$4/times 8$矩阵可以按行主次序把非0元素映射到一维数组中，可到到:2,1,6,7,3,9,8,4,5。

为了重建矩阵结构，必须记录每个非0元素所在的行号和列号，所以在把稀疏矩阵的非0元素映射到数组中时必须提供三个域: row(矩阵元素所在的行号)、col(矩阵元素所在列号)和value(矩阵元素的值) 。为此定义了下列所示的模板类 Term :

template<class T> class Term{     int row, col;     T value; };

如果a是一个类型为Term的数组，那么下图a中的稀疏矩阵按行主次序存储到a中所得的结果就如图b所示。
数组和矩阵 4 - 稀疏矩阵

除了存储数组a以外，还必须存储矩阵行数、矩阵列数和非0项的数目。所以存储上图a中的九个非0元素所需要的存储器字节数是 21 sizeof(int)+9 sizeof(T) ,这里每个非0元素都有两个int类型的行号和列号，然后加上总的矩阵行数，列数以及非0数目。

类SparseMatrix

可以定义一个类SparseMatrix，如下所示，用来把稀疏矩阵按行主次序映射到一维数组中。在定义共享成员时，没有定义加法操作符+，因为它会创建一个临时结果，这个临时结果必须复制到所返回的环境才可以使用。由于SparseMatrix的复制构造函数将会复制每一个元素，因此操作符+中的复制代价太大，这里使用Add函数来避免这种情况的发生。

template<class T> class SparseMatrix{ private:     void Append(const Term<T>& t);     int rows, cols;     // 非0元素的数目     int terms;     // 存储非0元素的数组     Term<T> *a;     // 数组a的大小     int MaxTerms; public:     SparseMatrix(int maxTerms = 0);     ~SparseMatrix(){ delete[] a; }     void Transpose(SparseMatrix<T>& b)  const;     void Add(const SparseMatrix<T> &b, SparseMatrix<T>&c)const;     friend std::ostream& operator<<(std::ostream&, const SparseMatrix<T>&);     friend std::istream& operator>>(std::istream&, SparseMatrix<T>&); };

下面给出构造函数，以及输入操作符和输出操作符，两者的时间复杂性都是$/theta(terms)$。

template<class T> SparseMatrix<T>::SparseMatrix(int maxTerms){     if (maxTerms < 1)         throw BadInitializers();     MaxTerms = maxTerms;     a = new Term<T>[maxTerms];     terms = rows = cols = 0; }  // 重载<< template<class T> std::ostream& operator<<(std::ostream& out, const SparseMatrix<T>& x){     // 输出矩阵的特征     out << "rows = " << x.rows << " columns = " << x.rows << std::endl;     out << "nonzeros terms = " << x.terms << std::endl;     // 输出非0元素，每行1个     for (int i = 0; i < x.terms; i++)         out << "a(" << x.a[i].row << ", " << x.a[i].col << ") = " << x.a[i].value << std::endl;     return out; }  // 重载>> template<class T> std::istream& operator>>(std::istream& in, SparseMatrix<T>& x){     // 输入矩阵的特征     std::cout << "Enter number of rows, columns, and terms/n";     in >> x.rows >> x.cols >> x.terms;     if (x.terms > x.MaxTerms)         throw NoMem();     // 输入矩阵元素     for (int i = 0; i < x.terms; i++){         cout << "Enter row, column, and value of term " << (i + 1) << std::endl;         in >> x.a[i].row >> x.a[i].col >> x.a[i].value;     }     return in; }

这里在重载输入操作符>>的时候，如果输入的元素数目大于数组a的大小，则会引发一个异常。一种处理异常的方法是删除数组a，然后使用new重新分配一个更大的数组。

矩阵转置

下面程序给出函数 Tranpose 的代码实现。转置后的矩阵被返回到b中。
首先验证b中是否有足够的空间来存储被转置矩阵的非0元素。如果空间不足，要么重新分配一个更大的数组 b.a ,要么引发一个异常。在下面的程序中是选择引发异常。如果b中有足够的空间来容纳转置矩阵，则创建两个数组 ColSize 和 RowNext 。其中 ColSize[i] 是指矩阵第i列中的非0元素的数目，而 RowNext[i] 则是代表转置矩阵第i行的下一个非0元素在b中的位置。

template<class T> void SparseMatrix<T>::Transpose(SparseMatrix<T>& b)const{     // 把*this的转置结果送入b中     // 验证b有足够空间     if (terms > b.MaxTerms)         throw NoMem();     // 设置转置特征     b.cols = rows;     b.rows = cols;     b.terms = terms;     // 初始化     int *ColSize, *RowNext;     ColSize = new int[cols + 1];     RowNext = new int[rows + 1];     for (int i = 1; i <= cols; i++)         ColSize[i] = 0;     // 计算*this每一列的非0元素数量     for (int i = 0; i < terms; i++)         ColSize[a[i].col]++;     // 给出b中每一行的起始点     RowNext[1] = 0;     for (int i = 2; i <= cols; i++)         RowNext[i] = RowNext[i - 1] + ColSize[i - 1];      // 执行转置操作     for (int i = 0; i < terms; i++){         // 在b中的位置         int j = RowNext[a[i].col]++;         b.a[j].row = a[i].col;         b.a[j].col = a[i].row;         b.a[j].value = a[i].value;     } }

函数Tranpose的时间复杂性是$O(cols+terms)$

矩阵相加

在两个矩阵相加中使用了函数Append，它把一个非0项添加到一个稀疏矩阵的非0项数组的尾部，其时间复杂性是$/theta(1)$。实现代码如下，然后就是两个矩阵相加的实现函数Add，使用两个游标，分别是*this和矩阵b的游标，通过一个while循环来实现相加。

template<class T> void SparseMatrix<T>::Append(const Term<T>& t){     // 把一个非0元素t添加到 *this之中     if (terms >= MaxTerms)         throw NoMem();     a[terms] = t;     terms++; }  template<class T> void SparseMatrix<T>::Add(const SparseMatrix<T>& b, SparseMatrix<T>& c)const{     // 计算 c = (*this) + b     // 验证可行性     if (rows != b.rows || cols != b.cols)         throw SizeMismatch();     // 设置结果矩阵c的特征     c.rows = rows;     c.cols = cols;     c.terms = 0;     // 定义*this 和b 的游标     int ct = 0, cb = 0;     while (ct < terms && cb < b.terms){         // 每一个元素的行主索引         int indt = a[ct].row * cols + a[ct].col;         int indb = b.a[cb].row * cols + b.a[cb].col;         if (indt < indb){             // b的元素在后面             c.Append(a[ct]);             ct++;         }         else{             if (indt == indb){                 // 位置相同                 if (a[ct].value + b.a[cb].value){                     // 仅当和不为0时，才添加到c中                     Term<T> t;                     t.row = a[ct].row;                     t.col = a[ct].col;                     t.value = a[ct].value + b.a[cb].value;                     c.Append(t);                 }                 ct++;                 cb++;             }             else{                 // b的元素在前面                 c.Append(b.a[cb]);                 cb++;             }         }     }     // 复制剩余元素     for (; ct < terms; ct++)         c.Append(a[ct]);     for (; cb < b.terms; cb++)         c.Append(b.a[cb]); }

函数Add的时间复杂性是$O(terms+b.terms)$。而如果用二维数组来描述每个矩阵，则两个矩阵相加耗时$O(rows cols)$，当 terms+b.terms 远小于**rows cols**时，稀疏矩阵的加法执行效率将大大提高。

链表描述

用一维数组来描述稀疏矩阵所存在的缺点是：当我们创建这个一维数组时，必须知道稀疏矩阵中的非0元素总数。
在我们自定义的类SparseMatrix中，当实际非0元素数目多于估计的初始化一维数组时设定的非0元素数目时，会引发一个异常。还有一种做法是可以分配一个更大的、新的数组，然后复制元素，并删除老的数组，但是这种做法会使得算法效率降低，并且也同样需要估计新数组需要多大的问题。

因此，这里就如同线性表一样，除了使用数组描述，还有基于指针的描述，也就是 链表描述 。

描述

链表描述的一种可行方案是把每行的非0元素串接在一起，构成一个链表，如下图所示。

图中每个非阴影节点代表稀疏矩阵中的一个非0元素，它有三个域：col(非0元素所在列号)、value(非0元素的值)和link(指向下一个非阴影节点的指针)。仅当矩阵某行中至少包含一个非0元素才会为该行创建一个链表。在行链表中，每个节点按其col值得升序进行排序。
数组和矩阵 4 - 稀疏矩阵

然后再用一个链表将所有的行链表，即图中阴影链表收集在一起。各个阴影节点按其row值得升序排列，每个阴影节点可以被视为一个行链表的头节点，因此阴影链表可以被视为头节点链表。

链表节点类型

这里分别定义图中非阴影节点 CNode 和阴影节点 HeadNode ，其代码实现如下：

#ifndef LINKMATRIX_H_ #define LINKMATRIX_H_ #include<iostream> #include"ChainList.h"  template<class T> class CNode{ private:     int col;     T value; public:     int operator!=(const CNode<T>& y){         return (value != y.value);     }     void Output(std::ostream& out)const{         out << "column = " << col << ", value= " << value;     } };  template<class T> std::ostream& operator<<(std::ostream& out, const CNode<T>& x){     x.Output(out);     out << std::endl;     return out; }  template<class T> class HeadNode{ private:     int row;     // 行链表     Chain<CNode<T>> a; public:     int operator!=(const HeadNode<T>& y){         return (row != y.row);     }     void Output(std::ostream& out)const{         out << "row = " << row;     } }; template<class T> std::ostream& operator<<(std::ostream& out, const HeadNode<T>& x){     x.Output(out);     out << std::endl;     return out; } #endif

类LinkMatrix

接下来就是定义类LinkMatrix，如下所示。

template<class T> class LinkMatrix{ private:     int rows, cols;     // 头节点链表     Chain<HeadNode<T>> a; public:     LinkMatrix(){}     ~LinkMatrix(){}     void Transpose(LinkMatrix<T>& b)const;     void Add(const LinkMatrix<T>& b, LinkMatrix<T>& c)const;     template<class T>     friend std::ostream& operator<<(std::ostream&, const LinkMatrix<T>&);     template<class T>     friend std::istream& operator>>(std::istream&, const LinkMatrix<T>&); };

重载>>

重载输入操作符>>。首先是要求输入矩阵的维数以及非0元素的个数。然后输入各个非0元素并把它们收集到各行链表中。用变量H代表当前行链表的头节点，如果下一个非0元素不属于当前行链表，则将当前行链表添加到矩阵x的头节点x.a之中；接下来，H被设置为指向一个新的行链表，同时将刚才那个非0元素添加到这个新的行链表之中。如果新的非0元素属于当前行链表，则只需要简单地把它添加到链表H.a中

template<class T> std::istream& operator>>(std::istream& in, const LinkMatrix<T>& x){     // 从输入流中输入矩阵x     // 删除x中所有节点     x.a.Erase();     // 获取矩阵特征     int terms;  // 输入的元素数     cout << "Enter numbers of rows, columns, and terms/n";     in >> x.rows >> x.cols >> terms;     // 虚设第0行     HeadNode<T> H;  // 当前行的头节点     H.row = 0;              // 当前行号     // 输入x的非0元素     for (int i = 1; i <= terms; i++){         // 输入下一个元素         cout << "Enter row, column, and value of term " << i << std::endl;         int row, col;         T value;         in >> row >> col >> value;         // 检查新元素是否属于当前行         if (row > H.row){             // 如果不是第0行，则把当前行的头节点H 添加到头节点链表x.a之中             if (H.row)                 x.a.Append(H);             // 为新的一行准备H             H.row = row;             // 置链表头指针first=0             H.a.Zero();          }         // 添加新元素         CNode<T> *c = new CNode<T>;         c->col = col;         c->value = value;         H.a.Append(*c);     }     // 注意矩阵的最后一行     if (H.row)         x.a.Append(H);     H.a.Zero();     return in; }

重载<<

这里为了输出链表表示的稀疏矩阵，使用了一个链表遍历器依次检查头节点链表中的每个节点。代码的时间复杂性与非0元素的数目呈正比。

template<class T> std::ostream& operator<<(std::ostream& out, const LinkMatrix<T>& x){     // 把矩阵x送至输出流     ChainIterator<HeadNode<T>> p;   // 头节点遍历器     // 输出矩阵的维数     out << "rows = " << x.rows << ",columns = " << x.cols << std::endl;     // 将h指向第一个头节点     HeadNode<T> *h = p.Initialize(x.a);     if (!h){         out << "No non-zero terms/n";         return out;     }     // 每次输出一行     while (h){         out << "row = " << h->row << std::endl;         out << h->a << "/n";        // 输出行链表；         // 下一个头节点         h = p.Next();     }     return out; }

函数Tranpose

对于转置操作，可以采用箱子来从矩阵 this中收集位于同一行的非0元素。 *bin[i] 是结果矩阵b中第i行非0元素所对应的链表。其实现如下所示

template<class T> void LinkMatrix<T>::Transpose(LinkMatrix<T>& b)const{     // 转置 *this，并把结果放入b     b.a.Erase();     // 创建用来收集b中各行元素的箱子     Chain<CNode<T>> *bin;     bin = new Chain<CNode<T>>[cols + 1];     // 头节点遍历器     ChainIterator<HeadNode<T>> p;     // h 指向*this的第一个头节点     HeadNode<T> *h = p.Initialize(a);     // 把*this的元素复制到箱子中     while (h){         int r = h->row;         // 行链表遍历器         ChainIterator<CNode<T>> q;         // 将z指向行链表的第一个节点         CNode<T> *z = q.Initialize(h->a);         // 临时节点         CNode<T> x;         // *this第r行中的元素变成b中第r列的元素         x.col = r;         // 检查*this第r行的所有非0元素         while (z){             x.value = z->value;             bin[z->col].Append(x);             z = q.Next();         }         h = p.Next();     }     // 设置b的维数     b.rows = cols;     b.cols = rows;     // 装配b的头节点链表     HeadNode<T> H;     // 搜索箱子     for (int i = 1; i <= cols; i++){         if (!bin[i].isEmpty()){             // 转置矩阵的第i行             H.row = i;             H.a = bin[i];             b.a.Append(H);             bin[i].Zero();         }     }     H.a.Zero();     delete[] bin; }

其中while循环所需要的时间与非0元素的数目呈线性关系，for循环所需要的时间则与输入矩阵的列数呈线性关系，因此总的时间与这两个量的和呈线性关系。