本篇為稀疏矩陣求解算法經(jīng)典論著<Direct Methods for Sparse Linear System>的<讀書筆記 4>

Chapter 2 Basic algorithms

上一篇主要介紹了稀疏矩陣的數(shù)據(jù)結(jié)構(gòu)，那么在接下來本章主要包括了其幾個基礎(chǔ)的算法。

Matrix-vector multiplication

對于矩陣-向量乘法操作

,其中 和 是稠密向量， 為稀疏矩陣。

其偽代碼為

這個偽代碼是基于CSC格式的，第一層for循環(huán)實際上是以x向量展開的。以CSR格式進行計算更容易理解。

C代碼為：

int

2.3 Utilities

此部分主要是一些用于內(nèi)存空間處理的應(yīng)用函數(shù)。不展開介紹，比較簡單。這里的函數(shù)僅僅是創(chuàng)建了內(nèi)存空間，但是沒有存儲空間中沒有實際內(nèi)容。

函數(shù)用途

void *cs_malloc (int n, size_t size)	內(nèi)存分配
void *cs_calloc (int n, size_t size)	內(nèi)存分配并清零
void *cs_free (void *p)	清楚分配空間
void *cs_realloc (void *p, int n, size_t size, int *ok)	改變一塊存儲空間p的size
cs *cs_spalloc (int m, int n, int nzmax, int values, int triplet)	創(chuàng)建m*n稀疏矩陣，用于保存nzmax個entries,參數(shù)triplet用于設(shè)定是COO還是CSC格式
int cs_spfree (cs *A, int nzmax)	釋放創(chuàng)建的稀疏矩陣空間
int cs_sprealloc (cs *A, int nzmax)	依據(jù)nzmax,重新對稀疏矩陣A分配空間

2.4 Triplet form

一般情況下，需要從文件中或許原始稀疏矩陣數(shù)據(jù)，常用的數(shù)據(jù)文件，例如mtx文件等，采用的是COO（或者在本書中稱 Triplet form）。程序需要讀入此類文件，獲取數(shù)據(jù)，然后再得到這個數(shù)據(jù)的CSC，或者CSR格式。

如下面的代碼，首先將COO格式的數(shù)據(jù)，存放再Ti，Tj，Tx中，然后再轉(zhuǎn)成CSC格式的數(shù)據(jù)，得到對應(yīng)的p，以及i。如果存在具有相同行和列索引的多個entries，則對應(yīng)的數(shù)值是所有這些重復(fù)entries的和。

cs 函數(shù)用途

int cs_entry (cs *T, int i, int j, double x)	將aij=x的數(shù)據(jù)存入到T中
cs *cs_compress (const cs *T)	將COO轉(zhuǎn)變?yōu)镃SC格式
cs *cs_done (cs *C, void *w, void *x, int ok)	釋放w和x的內(nèi)存空間
double cs_cumsum (int *p, int *c, int n)	計算累加和

2.5 Transpose

稀疏矩陣的轉(zhuǎn)置，與稠密矩陣不同，更像是將矩陣的從CSC格式轉(zhuǎn)變?yōu)镃SR格式的轉(zhuǎn)換。

函數(shù)用途

cs *cs_transpose (const cs *A, int values)	稀疏矩陣A的轉(zhuǎn)置

2.6 Summing up duplicate entries

有限元素方法生成的矩陣是元素（elements）的集合或密集的子矩陣。完備矩陣是元素(elements)的總和。如果兩個元素構(gòu)成同一個條目，則應(yīng)該對它們的值進行求和。

函數(shù)用途

int cs_dupl (cs *A)	對重復(fù)elements求和

2.7 Removing entries from a matrix

CSparse并不要求它的稀疏矩陣不包含數(shù)值上的零項，但它的MATLAB接口要做到這一點。cs_fkeep函數(shù)使用了一個函數(shù)指針fkeep (i ,j ,aij , other)作為參數(shù)，fkeep 用于評估A中的每一個entry。如果fkeep是true，則此entry保持。

函數(shù)用途

int cs_fkeep (cs *A, int (*fkeep) (int, int, double, void *), void *other)	移除零entry
static int cs_nonzero (int i, int j, double aij, void *other)	判斷entry非零
int cs_dropzeros (cs *A)	移除零

2.8 Matrix multiplication

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

，其中 為 ，A為 ，B為 ，取A和B的列，每次計算出C的列。 , , 分別代表矩陣 的列向量，則有 。將A分解為k列，將 分解為k個獨立的entries。

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

定理2.1：

的非零特征是 的非零特征的并集，其中 需要滿足 為非零。如果 表示, 的非零entries的行索引集合，則有

矩陣乘法均需要計算

以及 。函數(shù)用途

int cs_scatter (const cs *A, int j, double beta, int *w, double *x, int mark, cs *C, int nz)	對于給定的j，實現(xiàn)上述兩個公式
cs *cs_multiply (const cs *A, const cs *B)	計算稀疏矩陣乘法C=AB

稀疏矩陣乘法的計算量為

，其中 為浮點運算性能。

2.9 Matrix addition

矩陣加

等價于 。所以函數(shù)cs_add實際上看起來很像cs_multiply。函數(shù)用途

cs *cs_add (const cs *A, const cs *B, double alpha, double beta)	計算C=alpha A+beta B

2.10 Vector permutation

置換矩陣P實際上也是一個稀疏矩陣，它的每一行和列均有只有一個非零元素并且為1。所以置換矩陣P可以用一個稀疏矩陣P來表示。或者用一個長度為n的整型向量p來表示，p即稱為置換向量，其中p[k]=i，意味著

。矩陣P是正交的，所以 。同樣，置換矩陣的逆被定義為pinv[i]=k，如果 。函數(shù)用途

int cs_pvec (const int *p, const double *b, double *x, int n)	計算x=Pb
int cs_ipvec (const int *p, const double *b, double *x, int n)	計算x=PTb
int *cs_pinv (int const *p, int n)	計算p的逆

2.11 Matrix permutation

函數(shù)cs_permute用于計算矩陣的置換，

(在MATLAB中為C=A(p,q))。其中，長度為n的向量q為列轉(zhuǎn)置向量，長度為m的pinv(不是p)為行轉(zhuǎn)置向量，其中A為m*n。如果pinv[i]=k，則A的第i行變?yōu)镃的第k行；同理，A的第j列通過q進行置換。

函數(shù)cs_symperm用于對稱矩陣的置換。如果一個對稱矩陣，僅保存了上三角矩陣或者下三角元素，則利用cs_cs_permute則會使得到C不是一個對稱矩陣，所以需要cs_symperm實現(xiàn)此功能，并且返回值C是與A同樣的三角矩陣形式。

函數(shù)用途

cs *cs_permute (const cs *A, const int *pinv, const int *q, int values)	計算一般矩陣的置換矩陣
cs *cs_symperm (const cs *A, const int *pinv, int values)	計算對稱矩陣的置換矩陣

2.12 Matrix norm

矩陣范數(shù)：

• 1-范數(shù)： ，是最容易計算的，利用cs_norm函數(shù)得到。
• 2-范數(shù)：A的最大奇異值，稀疏矩陣的2-范數(shù)并不重要，此處沒有給計算函數(shù)。
• ∞-范數(shù)：最大的行的和。

函數(shù)用途

double cs_norm (const cs *A)	計算稀疏矩陣的1-范數(shù)

總結(jié)

以上是生活随笔為你收集整理的c++矩阵转置_lt;读书笔记4gt; 稀疏矩阵基础算法的全部內(nèi)容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網(wǎng)站內(nèi)容還不錯，歡迎將生活随笔推薦給好友。