本篇为稀疏矩阵求解算法经典论著<Direct Methods for Sparse Linear System>的<读书笔记 4>

Chapter 2 Basic algorithms

上一篇主要介绍了稀疏矩阵的数据结构，那么在接下来本章主要包括了其几个基础的算法。

1. Matrix-vector multiplication

对于矩阵-向量乘法操作

,其中

和

是稠密向量，

为稀疏矩阵。

其伪代码为

这个伪代码是基于CSC格式的，第一层for循环实际上是以x向量展开的。以CSR格式进行计算更容易理解。

C代码为：

int 

2.3 Utilities

此部分主要是一些用于内存空间处理的应用函数。不展开介绍，比较简单。这里的函数仅仅是创建了内存空间，但是没有存储空间中没有实际内容。

2.4 Triplet form

一般情况下，需要从文件中或许原始稀疏矩阵数据，常用的数据文件，例如mtx文件等，采用的是COO（或者在本书中称 Triplet form）。程序需要读入此类文件，获取数据，然后再得到这个数据的CSC，或者CSR格式。

如下面的代码，首先将COO格式的数据，存放再Ti，Tj，Tx中，然后再转成CSC格式的数据，得到对应的p，以及i。如果存在具有相同行和列索引的多个entries，则对应的数值是所有这些重复entries的和。

cs 

2.5 Transpose

稀疏矩阵的转置，与稠密矩阵不同，更像是将矩阵的从CSC格式转变为CSR格式的转换。

2.6 Summing up duplicate entries

有限元素方法生成的矩阵是元素（elements）的集合或密集的子矩阵。完备矩阵是元素(elements)的总和。如果两个元素构成同一个条目，则应该对它们的值进行求和。

2.7 Removing entries from a matrix

CSparse并不要求它的稀疏矩阵不包含数值上的零项，但它的MATLAB接口要做到这一点。cs_fkeep函数使用了一个函数指针fkeep (i ,j ,aij , other)作为参数，fkeep 用于评估A中的每一个entry。如果fkeep是true，则此entry保持。

2.8 Matrix multiplication

在CSparse中，采用CSC格式，矩阵乘

，其中

为

，A为

，B为

，取A和B的列，每次计算出C的列。

,

,

分别代表矩阵

的列向量，则有

。将A分解为k列，将

分解为k个独立的entries。

矩阵C的非零特征由以下的定理得到：

定理2.1：

的非零特征是

的非零特征的并集，其中

需要满足

为非零。如果

表示,

的非零entries的行索引集合，则有

矩阵乘法均需要计算

以及

。

稀疏矩阵乘法的计算量为

，其中

为浮点运算性能。

2.9 Matrix addition

矩阵加

等价于

。所以函数

cs_add实际上看起来很像cs_multiply。

2.10 Vector permutation

置换矩阵P实际上也是一个稀疏矩阵，它的每一行和列均有只有一个非零元素并且为1。所以置换矩阵P可以用一个稀疏矩阵P来表示。或者用一个长度为n的整型向量p来表示，p即称为置换向量，其中p[k]=i，意味着

。矩阵P是正交的，所以

。同样，置换矩阵的逆被定义为

pinv[i]=k，如果

。

2.11 Matrix permutation

函数cs_permute用于计算矩阵的置换，

(在MATLAB中为C=A(p,q))。其中，长度为n的向量q为列转置向量，长度为m的

pinv(不是p)为行转置向量，其中A为m*n。如果pinv[i]=k，则A的第i行变为C的第k行；同理，A的第j列通过q进行置换。

函数cs_symperm用于对称矩阵的置换。如果一个对称矩阵，仅保存了上三角矩阵或者下三角元素，则利用cs_cs_permute则会使得到C不是一个对称矩阵，所以需要cs_symperm实现此功能，并且返回值C是与A同样的三角矩阵形式。

2.12 Matrix norm

矩阵范数：

• 1-范数：
  ，是最容易计算的，利用

  cs_norm函数得到。

• 2-范数：A的最大奇异值，稀疏矩阵的2-范数并不重要，此处没有给计算函数。
• ∞-范数：最大的行的和。