稀疏矩陣演算法
1. 稀疏矩陣的運算
import java.util.*;
//稀疏矩陣演算法。
//稀疏態昌矩陣演算法是為了在大型矩陣中非零元素少時，減少存貯空間，並提高矩陣運算速度的。
//但本例中的矩陣只是為了演示演算法，都比較小，時間告備和空間效率提升可以忽襪閉毀略。
public class SparseMatrix{
public static void main(String[] args){
TripleSMatrix tsm=new TripleSMatrix(7,4);
//tsm.printTriple();
tsm.printMatrix();
TripleSMatrix tsm2=new TripleSMatrix(7,4);
System.out.println("矩陣a:");
tsm.printMatrix();
System.out.println("矩陣b:");
tsm2.printMatrix();
int[][] matrixSum=addSMatrix(tsm,tsm2);
System.out.println("矩陣a+矩陣b:");
for(int i=0;i <matrixSum.length;i++){
for(int j=0;j <matrixSum[i].length;j++){
System.out.print(" "+matrixSum[i][j]);
}
System.out.println("");
}
}
public static int[][] addSMatrix(TripleSMatrix t1,TripleSMatrix t2){ //計算兩個三元組表示的矩陣之和，返回結束數組
if(t1.rows!=t2.rows||t1.columns!=t2.columns){
System.out.println("這兩個矩陣不能相加");
return null;
}
int[][] c=new int[t1.rows][t2.columns];
int i=1,j=1;
while(i <=t1.nonzeroElements||j <=t2.nonzeroElements){
if(t1.triple[i][0] <t2.triple[j][0]&&i <=t1.nonzeroElements){
c[t1.triple[i][0]-1][t1.triple[i][1]-1]=t1.triple[i][2];
i++;
}else if(t2.triple[j][0] <t1.triple[i][0]&&j <=t2.nonzeroElements){
c[t2.triple[j][0]-1][t2.triple[j][1]-1]=t2.triple[j][2];
j++;
}else{
if(t1.triple[i][1] <t2.triple[j][1]&&i <=t1.nonzeroElements){
c[t1.triple[i][0]-1][t1.triple[i][1]-1]=t1.triple[i][2];
i++;
}else if(t1.triple[i][1]>t2.triple[j][1]&&j <=t2.nonzeroElements){
c[t2.triple[j][0]-1][t2.triple[j][1]-1]=t2.triple[j][2];
j++;
}else{
c[t1.triple[i][0]-1][t1.triple[i][1]-1]=t1.triple[i][2]+t2.triple[j][2];
i++;j++;
}
}
}
return c;
}
}
//下面的類定義不一定是最好的，比如其中的屬性大多是包訪問許可權，可以改進。
class TripleSMatrix{ //定義了一個三元組的類。
int[][] triple=new int[2001][3]; //三元組數組，假設稀疏矩陣的值都是整數。最多可以有2000個非零元素。第零行沒有用。
int rows,columns,nonzeroElements; //稀疏矩陣的行列數和非零元素個數。
TripleSMatrix(int rows,int columns){ //構造方法，rows是稀疏矩陣的行數，columns是稀疏矩陣的列數。
Scanner input=new Scanner(System.in);
System.out.println("請輸入稀疏矩陣三元組");
System.out.println("以行 列 值的形式輸入，如：1 2 4表示第1行第2列元素的值為4，當輸入的行為999時結束：");
int count=1;
int i=0,j,v; //i行j列,值v
while(i!=999&&input.hasNext()){
i=input.nextInt();
j=input.nextInt();
v=input.nextInt();
if(i>rows||i <1||j>columns||j <1){
System.out.println("剛才的行，列值錯，將被忽略");
continue;
}
triple[count][0]=i;
triple[count][1]=j;
triple[count][2]=v;
count++;
}
this.rows=rows;
this.columns=columns;
this.nonzeroElements=count-1;
sortTriple(triple,1,count); //對輸入的三元組排序。
}
static void sortTriple(int[][] triple,int first,int end){ //對三元組排序方法，按行排，行一樣按列排。
Arrays.sort(triple,first,end,new Comparator <int[]>(){
public int compare(int[] t1,int[] t2){
if(t1[0]>t2[0]) return 1;
if(t1[0] <t2[0]) return -1;
if(t1[0]==t2[0]) return t1[1]-t2[1];
return 0; //沒有用的一個語句，但沒有它編譯通不過。
}
});
}
public void printMatrix(){ //列印出當前三元組表示的稀疏矩陣。
int row=1,column=1; //row當前要列印的行，column當前要列印的列。
for(int t=1;t <=nonzeroElements;t++){
while(triple[t][0]>row){ //三元組中的行比當前行大
if(column!=1){ //前面列印的行沒有列印完，繼續列印完
for(;column <=columns;column++) System.out.print(" "+0);
column=1; //新的一行列從1開始。
}else{ //當前行全為0
for(int i=1;i <=columns;i++){
System.out.print(" "+0);
}
}
System.out.println(""); //換行
row++; //下一行
}
for(;column <triple[t][1];column++){ //當前列印的列小於三元組中的列，前面要補零。
System.out.print(" "+0);
}
System.out.print(" ".substring(0,6-(String.valueOf(triple[t][2])).length())+triple[t][2]); //列印三元組對應的元素。
column++;
}
if(column!=1){ //前面列印的行沒有列印完，繼續列印完
for(;column <=columns;column++) System.out.print(" "+0);
System.out.println("");
column=1;
row++ ;
}
for(;row <=rows;row++){ //三元組中沒有對應的值了，矩陣後面的元素全為0
for(column=1;column <=columns;column++){
System.out.print(" "+0);
}
System.out.println("");
}
}
public void printTriple(){ //列印三元組
for(int i=1;i <=nonzeroElements;i++){
for(int j=0;j <3;j++){
System.out.print(triple[i][j]+" ");
}
System.out.println("");
}
}
}
為了節省存儲空間並且加快處理速度，需要對這類矩陣進行壓縮存儲，壓縮存儲的原則是：不重復存儲相同元素；不存儲零值元素。稀疏矩陣，有三元組表示法、帶輔助行向量的二元組表示法（也即行邏輯鏈表的順序表），十字鏈表表示法等。演算法基本思想：num[col]:第col列的非零元素個數；cpot[col]:第col列第一個非零元在b.data中的恰當位置；在轉置過程中，指示該列下一個非零元在b.data中的位置。
3. 稀疏矩陣的運算
M AT L A B中對滿矩陣的運算和函數同樣可用在稀疏矩陣中.結果是稀疏矩陣還是滿矩陣,
這取決於運算符或者函數及下列的操作數:當函數用一個矩陣作為輸入參數,輸出參數為一個標量或者一個給定大小的向量時,輸出參數的格式總是返回一個滿陣形式,如命令s i z e.
當函數用一個標量或者一個向量作為輸入參數,輸出參數為一個矩陣時,輸出參數的格式也總是返回一個滿矩陣,如命令e y e.還有一些特殊的命令可以得到稀疏矩陣,如命令s p e y e.
對於單參數的其他函數來說,通常返回的結果和參數的形式是一樣的,如d i a g.
對於雙參數的運算或者函數來說,如果兩個參數的形式一樣,那麼也返回同樣形式的結果.在兩個參數形式不一樣的情況下,除非運算的需要,均以滿矩陣的形式給出結果.
兩個矩陣的組和[A B],如果A或B中至少有一個是滿矩陣,則得到的結果就是滿矩陣.
表達式右邊的冒號是要求一個參數的運算符,遵守這些運算規則.
表達式左邊的冒號不改變矩陣的形式. 假設有:
這是一個5×5的單位滿矩陣和相應的稀疏矩陣.
(a) C = 5*B,結果為:
這是一個稀疏矩陣.
(b) D = A + B,給出的結果為:
這是一個滿矩陣.
(c) x = B h,結果為:
這是一個滿向量.
有許多命令可以對非零元素進行操作.
命令集8 9矩陣的非零元素
n n z ( A )求矩陣A中非零元素的個數.它既可求滿矩陣也可求稀疏矩陣.
s p y ( A )畫出稀疏矩陣A中非零元素的分布.也可用在滿矩陣中,在
這種情況下,只給出非零元素的分布.
s p y ( A , c s t r , s i z e )用指定的顏色c s t r(見表1 3 - 1 )和在s i z e規定的范圍內畫出稀疏
矩陣A中非零元素的分布.
n o n z e r o s ( A )按照列的順序找出矩陣A中非零的元素.
s p o n e s ( A )把矩陣A中的非零元素全換為1.
s p a l l o c ( m , n ,產生一個m×n階只有n z m a x個非零元素的稀疏矩陣.這樣可以
n z m a x )有效地減少存儲空間和提高運算速度.
n z m a x ( A )給出為矩陣A中非零元素分配的內存數.不一定和n n z ( A )得
到的數相同;參見s p a r s e或者s p a l l o c.
i s s p a r s e ( A )如果矩陣A是稀疏矩陣,則返回1;否則返回0.
s p f u n ( f c n , A )用A中所有非零元素對函數f c n求值,如果函數不是對稀疏矩
陣定義的,同樣也可以求值.
s p r a n k( A )求稀疏矩陣A的結構秩.對於所有的矩陣來說,都有
s p r a n k ( A）≥r a n k ( A ). 用下面的命令定義稀疏矩陣:
創建一個大矩陣:
Big=kron(A, A)
這個矩陣B i g是什麼樣子呢?
K r o n e c k e r張量積給出一個大矩陣,它的元素是矩陣A的元素之間可能的乘積.因為參量都是稀疏矩陣,所以得到的矩陣也是一個稀疏矩陣.可以用命令 w h o s和i s s p a r s e來確認一下.
查看矩陣B i g的結構圖,可輸入s p y ( B i g ),結構圖如右圖所示. 從圖中可以看出B i g是一個塊雙對角矩陣. MATLAB中有四個基本稀疏矩陣,它們是單位矩陣,隨機矩陣,對稱隨機矩陣和對角矩陣.
命令集9 0單位稀疏矩陣
s p e y e ( n )生成n×n的單位稀疏矩陣.
s p e y e ( m , n )生成m×n的單位稀疏矩陣.
命令speye(A) 得到的結果和s p a r s e ( e y e ( A ) )是一樣的,但是沒有涉及到滿陣的存儲.
命令集9 1隨機稀疏矩陣
s p r a n d ( A )生成與A有相同結構的隨機稀疏矩陣,且元素服從均勻分布.
s p r a n d ( m , n , d e n s )生成一個m×n的服從均勻分布的隨機稀疏矩陣,有d e n s×m×
n個非零元素,0≤d e n s≤1.參數d e n s是非零元素的分布密度.
s p r a n d ( m , n , d e n s ,生成一個近似的條件數為1 /rc,大小為m×n的隨機稀疏矩
r c )陣.如果rc=rc是一個長度為l≤l ( m i n (m, n) )的向量,那麼
矩陣將rci作為它l個奇異值的第一個,其他的奇異值為0.
s p r a n d n ( A )生成與A有相同結構的隨機稀疏矩陣,且元素服從正態分布.
s p r a n d n ( m , n , d e n s ,生成一個m×n的服從正態分布的隨機稀疏矩陣,和sprand
r c )一樣.
s p r a n d s y m ( S )生成一個隨機對稱稀疏矩陣.它的下三角及主對角線部分與S的結構相同,矩陣元素服從正態分布.
s p r a n d s y m ( n , d e n s )生成一個m×n的隨機對稱稀疏矩陣.矩陣元素服從正態分布,分布密度為d e n s.
s p r a n d s y m ( n , d e n s ,生成一個近似條件數為1 /rc的隨機對稱稀疏矩陣.元素以0r c )對稱分布,但不是正態分布.如果rc=rc是一個向量,則矩陣有特徵值rci.也就是說,如果rc是一個正向量,則矩陣是正定矩陣.
s p r a n d s y m ( n , d e n s ,生成一個正定矩陣.如果k= 1,則矩陣是由一正定對稱矩r c , k )陣經隨機J a c o b i旋轉得到的,其條件數正好等於1 /rc;如果k= 2,則矩陣為外積的換位和,其條件數近似等於1 /rc.
s p r a n d s y m ( S , d e n s ,生成一個與矩陣S結構相同的稀疏矩陣,近似條件數為1 /rc.
r c , 3 )參數d e n s被忽略,但是這個參數在這個位置以便函數能確認最後兩個參數的正確與否. (a) 假設有矩陣A:
輸入R a n d o m = s p r a n d n ( A ),可得到隨機稀疏矩陣:
矩陣中隨機數的位置和矩陣A中非零元素的位置相同.
(b) 對於( a )中的矩陣A,輸入:
B = s p r a n d s y m ( A )
結果為:
這是一個用矩陣A的下三角及主對角線部分創建的對稱矩陣,在非零元素的位置用隨機數作為元素值.
用命令s p d i a g s可以取出對角線元素,並創建帶狀對角矩陣.假設矩陣A的大小為m×n,
在p個對角線上有非零元素.B的大小為m i n (m×n)×p,它的列是矩陣A的對角線.向量d的長度為p,其整型分量給定了A的對角元:
di0 主對角線上的對角線
命令集9 2對角稀疏矩陣
[ B , d ] = s p d i a g s ( A )求出A中所有的對角元,對角元保存在矩陣B中,它們的下標保存在向量d中.
s p d i a g s ( A , d )生成一個矩陣,這個矩陣包含有矩陣A中向量d規定的對角元.
s p d i a g s ( B , d , A )生成矩陣A,用矩陣B中的列替換d定義的對角元.
A = s p d i a g s ( B , d , m , n )用保存在由d定義的B中的對角元創建稀疏矩陣A.
例11 . 4給出了如何使用s p d i a g s命令來解普通微分方程組. 在許多實際應用中要保留稀疏矩陣的結構,但是在計算過程中的中間結果會減弱它的稀疏性,如L U分解.這就會導致增加浮點運算次數和存儲空間.為了避免這種情況發生,在第稀疏矩陣1 2 9
M AT L A B中用命令對矩陣進行重新安排.這些命令都列在下面的命令集9 3中.通過h e l p命令
可以得到每個命令更多的幫助信息,也可見h e l p d e s k.
命令集9 3矩陣變換
c o l m m d ( A )返回一個變換向量,使得矩陣A列的秩為最小.
s y m m m d ( A )返回使對稱矩陣秩為最小的變換.
s y m r c m ( A )矩陣A的C u t h i l l - M c K e e逆變換.矩陣A的非零元素在主對角線附近.
c o l p e r m ( A )返回一個矩陣A的列變換的向量.列按非零元素升序排列.有時這是L U因式分解前有用的變換:lu(A(:, j)).如果A是一個對稱矩陣,對行和列進行排序,這有利於C h o l e s k y分解:chol(A(j, j)).
r a n d p e r m ( n )給出正數1,2,. . .,n的隨機排列,可以用來創建隨機變換矩陣.
d m p e r m ( A )對矩陣A進行D u l m a g e - M e n d e l s o h n分解,輸入help dmperm可得更多信息. 創建一個秩為4的變換矩陣,可輸入:
一旦運行p e r m = r a n d p e r m ( 4 ),就會得到:
給出的變換矩陣為:
如果矩陣A為:
輸入命令:
運行結果為:
有兩個不完全因式分解命令,它們是用來在解大線性方程組前進行預處理的.用h e l p d e s k命令可得更多信息.命令集9 4不完全因式分解c h o l i n c ( A , o p t )進行不完全C h o l e s k y分解,變數o p t取下列值之一:
d r o p t o l指定不完全分解的舍入誤差,0給出完全分解.
m i c h o l如果m i c h o l = 1,則從對角線上抽取出被去掉的元素.
r d i a g用s q r t ( d r o p t o l*n o r m ( X ( : , j ) ) )代替上三角分
解因子中的零元素,j為零元素所在的列.
[ L , U , P ]=返回矩陣X的不完全分解得到的三個矩陣L,U和P,變數o p t取
l u i n c ( X , o p t )下列值之一:
d r o p t o l指定分解的舍入誤差.
m i l u改變分解以便從上三角角分解因子中抽取被去掉的列元素.
u d i a g用d r o p t o l值代替上三角角分解因子中的對角線上的零元素.
t h r e s h中心極限.
解稀疏線性方程組既可用左除運算符解,也可用一些特殊命令來解.
命令集9 5稀疏矩陣和線性方程組
s p p a r m s ( k e y s t r , o p )設置稀疏矩陣演算法的參數,用help spparms可得詳細信息.
s p a u g m e n t ( A , c )根據[ c*l A; A' 0 ]創建稀疏矩陣,這是二次線性方程組的最
小二乘問題.參見7 . 7節.
s y m b f a c t ( A )給出稀疏矩陣的C h o l e s k y和L U因式分解的符號分解因子.
用help symbfact可得詳細信息.
稀疏矩陣的范數計算和普通滿矩陣的范數計算有一個重要的區別.稀疏矩陣的歐幾里德范數不能直接求得.如果稀疏矩陣是一個小矩陣,則用n o r m ( f u l l ( A ) )來計算它的范數;但是對於大矩陣來說,這樣計算是不可能的.然而M AT L A B可以計算出歐幾里德范數的近似值,在計算條件數時也是一樣.
命令集9 6稀疏矩陣的近似歐幾里德范數和條件數
n o r m e s t ( A )計算A的近似歐幾里德范數,相對誤差為1 0-6.
n o r m e s t ( A , t o l )計算A的近似歐幾里德范數,設置相對誤差t o l,而不用預設時的1 0-6.
[ n r m , n i t ] =計算近似n r m范數,還給出計算范數迭代的次數n i t.
n o r m e s t ( A )
c o n d e s t ( A )求矩陣A條件數的1 -范數中的下界估計值.
[ c , v ]=求矩陣A的1 -范數中條件數的下界估計值c和向量v,使得
c o n d e s t ( A , t r )| |Av| | = ( | |A| | | |v| | ) / c.如果給定t r,則給出計算的過程.t r= 1,
給出每步過程;t r=-1,給出商c / r c o n d ( A ). 假設給出:
用n o r m A p p r o x = n o r m e s t ( S p r s )計算出:
用t h e N o r m = n o r m ( f u l l ( S p r s ) )得:
為了找到它們之間的差別,計算d i f f e r e n c e = t h e N o r m - n o r m A p p r o x,得:
在許多應用中,n o r m e s t計算得到的近似值是一個很好的近似歐幾里德范數,它的計算步數要比n o r m要少得多;可參見7 . 6節.
用e t r e e命令來找到稀疏對稱矩陣的消元樹,用向量f來描述消元樹,還可用e t r e e p l o t命令畫出來.元素fi是矩陣的上三角C h o l e s k y分解因子中i行上第1非零元素的列下標.如果有非零元素,則fi= 0.消元樹可以這樣來建立:
節點i是fi的孩子,或者如果fi= 0,則節點i是樹的根節點.
命令集9 7矩陣的消元樹
e t r e e ( A )求A的消元樹向量f,這個命令有可選參數;輸入h e l p
e t r e e獲取幫助.
e t r e e p l o t ( A )畫出向量f定義的消元樹圖形.
t r e e p l o t ( p , c , d )畫出指針向量p的樹圖形,參數c和d分別指定節點的顏色和分支數.e t r e e p l o t可以調用這個命令.
t r e e l a y o u t顯示樹的結構,t r e e p l o t可以調用這個命令. 假設有對稱稀疏矩陣B:
運行命令b t r e e = e t r e e ( B ),結果為:
開始的數字2不難理解,它是矩陣的第1列上第1個非零元素的行數,它決定了在C h o l e s k y分解因子的第1行第2列處有一個非零元素.當縮減第1列的元素時就得到第2列的數字5.B在縮減後,在( 5 , 2 )位置的元素是非零的,這樣消元樹向量中第2個元素的值為5.
s p y ( c h o l ( B ) )給出了C h o l e s k y分解因子的結構圖,如圖9 - 2所示:
圖9-2 Cholesky分解結構圖
圖9-3 矩陣B的消元樹
這個向量消元樹可以這樣來建立:上三角中只有一行有非零元素,節點8,因此這就是樹
唯一的根.節點1是節點2的孩子,節點2和3是節點5的孩子,而節點5是節點6的孩子.節點4和6是節點7的孩子,而節點7又是節點8的孩子,即根的孩子.
命令e t r e e p l o t ( B )給出了樹的結構圖,如圖9 - 3所示.
消元樹的形狀取決於列和行序,它可以用來分析消元過程.
用g p l o t命令可以畫出坐標和矩陣元素間的聯系圖形.必須在n×2的矩陣中給出n個坐標,
矩陣的每一行作為一個點.這樣就創建出點點之間連接的n×n矩陣,如果點4連接到點8,則(4, 8)的值為1.由於是一個大矩陣,而且非零元素較少,所以它應該被建成稀疏矩陣.
這個圖可以說明網路問題,如傳遞問題.它還包含有線性方程組中未知量之間的相關信息.
命令集9 8網路圖形
g p l o t ( A , K )如果矩陣A的a(i, j)不為0,則將點ki連接到點kj.K是一個n×
2的坐標矩陣,A是一個n×n的關聯矩陣.
g p l o t ( A , K , s t r )用字元串s t r給定的顏色和線型畫出的同上圖形.字元串s t r的
取值參見表1 3 - 1.
[ X , A ] = u n m e s h ( E )求網格邊界矩陣E的L a p l a c e矩陣A和網格點的坐標矩陣X.
例七
假設有下面的坐標矩陣K和關聯矩陣A:
矩陣A在稀疏化後,用命令g p l o t ( A , K )畫出圖9 - 4,給出了點(0, 1)和點(4, 1)之間所有可能的路徑.
4. c++ 十字鏈表實現稀疏矩陣的轉置
mplate <class Type>SparseMatrix<Type>
SparseMatrix<Type>::SparseMatrix::FastTranspose()
{
int *RowSize = new int[Cols]; // 統計各列非零元素個數
int *RowStart = new int[Cols]; // 預計轉置後各行存放位置
SparseMatrix b; // 存放轉置結果
b.Rows = cols; b.cols = Rows; b.Terms = Terms;
if (Terms>0) //nonzero matrix
{
for (int i=0;i<Cols;i++) RowSize[i] = 0; //initialize
for (i=0;i<terms;i++) RowSize[smArray[i].col]++;
//RowStart[i]=starting position of row i in b
RowStart[0]=0;
for (i=1;i<Cols;i++) RowStart[i] = RowStart[i-1]+
RowSize[i-1];
for (i=0;i<Terms;i++) //move from a to b
{
int j=Rowstart[smArray[i].col];
b.smArray[j].row=smArray[i].col;
b.smArray[j].col= smArray[i].row;
b.smArray[j].value= smArray[i].value;
RowStart[smArray[i].col]++;
} //end of for
}
delete [ ]Rowsize;
delete [ ]Rowstart;
return b;
}
另外,站長團上有產品團購,便宜有保證
5. 稀疏矩陣的乘法的演算法思想
/*Multiplicate part*/
//C = A * B
/*演算法分析：
首先，由於樓主沒有給出輸入函數，也沒有對三元組的稀疏矩陣的數據結構做完整的說明，
所以我只能猜測這個稀疏矩陣是以行為主序存儲的。（後面的乘法函數佐證了我的猜測，
但是如果我不幸猜錯了，還請樓主告知）
另一個猜測是樓主的程序中設定矩陣和數組的下標都是從1開始，而非從0開始。
接下來，我們說一下普通矩陣的乘法，這個在線性代數裡面有定義，無需贅言。
我要說的是稀疏矩陣的乘法也是用這個公式來計算，但卻有一個問題——效率。
我們舉一個例子來說明：
假設我們已知A:m*n和B:n*l,要計算C:m*l，那麼C(i,j)的計算公式就是：
C(i,j) = A(i,1)*B(1,j) + A(i,2)*B(2,j) + …… + A(i,n)*B(n,j) 公式1
如果A、B都是普通矩陣（並且以二維數組存儲），那麼直接用循環語句就可以完成。
如果A、B都是稀疏矩陣（並且亂序存儲，即是說沒有以行序或列序存儲），那麼計算就會很麻煩。
我們需要首先遍歷整個A矩陣去查找是否存在A(i,1)（若有則取其值，若無則其值為0），然後再去
遍歷B矩陣查找B(1,j)，並將兩者相乘；接著又是A(i,2)和B(2,j)，以此類推。
這樣做當然是可行的，但是顯然效率太低了，一種改進的方法（就是樓主的程序中所用的方法）如下：
首先我們要優化稀疏矩陣的存儲，不能亂序存儲，而是以行序或列序為主序來存儲，
比如這里我們以行序為主序，以列序為次序。
具體來說，就是將稀疏矩陣的第一行中的非零元素排列在前面，然後才是第二行、第三行……
在各行的非零元素中又以列序來排列，這樣存儲的稀疏矩陣就是有序的。
接下來，在真正用公式來計算之前，我們要做一個預處理。
這個預處理是對矩陣B所做的，其實就是要計算得到矩陣B所對應的pos數組。
那麼這個pos數組代表什麼意思呢？
我們從代碼中不難看出pos數組的長度是矩陣B的行數（也就是B->m）加1。
pos[i]（1 <= i<= B->m）表示矩陣B的第i行的第一個非零元素在三元數組B的位置。
pos[i]（i == B->m + 1）是為了方便後面的計算而增加一個標記，類似於監視哨。
有了這個pos數組之後，我們再來計算公式1就可以提高效率了。
效率的提高表現在兩個方面：
1、我們看到代碼中（如下）雖然有兩層while循環，但是其實總共只執行了A->t次。
p=1;
while(p<=A->t)
{
……
while (p<=A->t&&A->data[p].row==crow)
{
……
p=p+1;
}
……
}
2、在第二層while循環中的for循環，由於pos數組的作用使得循環次數減少很多。
for(q=pos[brow];q<=pos[brow+1]-1;q++)
{
ccol=B->data[q].col;
ctemp[ccol]=ctemp[ccol] + A->data[p].v * B->data[q].v;
}
好了，以上就是這個稀疏矩陣乘法的大致過程，如有錯誤請指正。
*/
void MultS(TriTable *A,TriTable *B,TriTable *C)
{
int k,p,q,brow,crow,ccol;
int num[MAXSIZE],pos[MAXSIZE],ctemp[MAXSIZE];
if (A->n==B->m)
{
for(k=1;k<=B->m;k++)
num[k]=0;
for(k=1;k<=B->t;k++)
num[B->data[k].row]++;
pos[1]=1;
for(k=2;k<=B->m;k++) //這里將k<=B->t改為k<=B->m
pos[k]=pos[k-1]+num[k-1];
pos[1+B->m]=pos[B->m]+1; //這里將k<=B->t改為k<=B->m
C->m=A->m;
C->n=B->n;
C->t=0;
p=1;
while(p<=A->t)
{
crow=A->data[p].row;
for(k=1;k<=C->n;k++)
ctemp[k]=0;
while (p<=A->t&&A->data[p].row==crow)
{
brow=A->data[p].col;
for(q=pos[brow];q<=pos[brow+1]-1;q++)
{
ccol=B->data[q].col;
ctemp[ccol]=ctemp[ccol] + A->data[p].v * B->data[q].v;
}
p=p+1;
}
for(ccol=1;ccol<=B->n;ccol++)
if(ctemp[ccol]!=0)
{
C->t=C->t+1;
C->data[C->t].row=crow;
C->data[C->t].col=ccol;
C->data[C->t].v=ctemp[ccol];
}
}
}else
printf("these two arrat can't Multiplicate");
return;
}
6. 對稀疏矩陣進行壓縮存儲的目的是什麼
對稀疏矩陣進行壓縮存儲目的是節省存儲空間。
存儲矩陣的一般方法是採用二維數組，其優點是可以隨機地訪問每一個元素，因而能夠較容易地實現矩陣的各種運算。
但對於稀疏矩陣而言，若用二維數組來表示，會重復存儲了很多個0了，浪費空間，而且要花費時間來進行零元素的無效計算。所以必須考慮對稀疏矩陣進行壓縮存儲。
(6)稀疏矩陣演算法擴展閱讀
優點
稀疏矩陣的計算速度更快，因為MATLAB只對非零元素進行操作，這是稀疏矩陣的一個突出的優點。假設矩陣A,B中的矩陣一樣，計算2*A需要一百萬次的浮點運算，而計算2*B只需要2000次浮點運算。
因為MATLAB不能自動創建稀疏矩陣，所以要用特殊的命令來得到稀疏矩陣。算術和邏輯運算都適用於稀疏矩陣。對於一個用二維數組存儲的稀疏矩陣Amn，如果假設存儲每個數組元素需要L個位元組，那麼存儲整個矩陣需要m*n*L個位元組。
7. 稀疏矩陣
// algo5-1.cpp 實現演算法5.2的程序
// c1.h (程序名)
#include<string.h>
#include<ctype.h>數念
#include<malloc.h> // malloc()等
#include<limits.h> // INT_MAX等
#include<stdio.h> // EOF(=^Z或F6),NULL
#include<stdlib.h> // atoi()
#include<io.h> // eof()
#include<math.h> // floor(),ceil(),abs()
#include<process.h> // exit()
#include<iostream.h> // cout,cin
// 函數結果狀態代碼
#define TRUE 1
#define FALSE 0
#define OK 1
#define ERROR 0
#define INFEASIBLE -1
// #define OVERFLOW -2 因為在math.h中已定義OVERFLOW的值為3,故去掉此行
typedef int Status; // Status是函數的類型,其值是函數結果狀態代碼，如OK等
typedef int Boolean; // Boolean是布爾類型,其值是TRUE或FALSE
typedef int ElemType;
// c5-2.h 稀疏矩陣的三元組順序表存儲表示
#define MAXSIZE 100 // 非零元個數的最大值
struct Triple
{
int i,j; // 行下標,列下標
ElemType e; // 非零元素值
};
struct TSMatrix
{
Triple data[MAXSIZE+1]; // 非零元三元組表,data[0]未用
int mu,nu,tu; // 矩陣的行數、列數和非零元個數核哪
};
// bo5-2.cpp 三元組稀疏矩陣的基本操作,包括演算法5.1(9個)
Status CreateSMatrix(TSMatrix &M)
{ // 創建稀疏矩陣M
int i,m,n;
ElemType e;
Status k;
printf("請輸入矩陣的行數,列數,非零元素數：");
scanf("%d,%d,%d",&M.mu,&M.nu,&M.tu);
M.data[0].i=0; // 為以下比較順序做准備
for(i=1;i<=M.tu;i++)
{
do
{
printf("請按行序順序輸入第%d個非零元素所在的行(1～%d),列(1～%d),元素值：",i,M.mu,M.nu);
scanf("%d,%d,%d"改畢碼,&m,&n,&e);
k=0;
if(m<1||m>M.mu||n<1||n>M.nu) // 行或列超出范圍
k=1;
if(m<M.data[i-1].i||m==M.data[i-1].i&&n<=M.data[i-1].j) // 行或列的順序有錯
k=1;
}while(k);
M.data[i].i=m;
M.data[i].j=n;
M.data[i].e=e;
}
return OK;
}
void DestroySMatrix(TSMatrix &M)
{ // 銷毀稀疏矩陣M
M.mu=0;
M.nu=0;
M.tu=0;
}
void PrintSMatrix(TSMatrix M)
{ // 輸出稀疏矩陣M
int i;
printf("%d行%d列%d個非零元素。\n",M.mu,M.nu,M.tu);
printf("行 列 元素值\n");
for(i=1;i<=M.tu;i++)
printf("%2d%4d%8d\n",M.data[i].i,M.data[i].j,M.data[i].e);
}
Status CopySMatrix(TSMatrix M,TSMatrix &T)
{ // 由稀疏矩陣M復製得到T
T=M;
return OK;
}
int comp(int c1,int c2) // 另加
{ // AddSMatrix函數要用到
int i;
if(c1<c2)
i=1;
else if(c1==c2)
i=0;
else
i=-1;
return i;
}
Status AddSMatrix(TSMatrix M,TSMatrix N,TSMatrix &Q)
{ // 求稀疏矩陣的和Q=M+N
Triple *Mp,*Me,*Np,*Ne,*Qh,*Qe;
if(M.mu!=N.mu)
return ERROR;
if(M.nu!=N.nu)
return ERROR;
Q.mu=M.mu;
Q.nu=M.nu;
Mp=&M.data[1]; // Mp的初值指向矩陣M的非零元素首地址
Np=&N.data[1]; // Np的初值指向矩陣N的非零元素首地址
Me=&M.data[M.tu]; // Me指向矩陣M的非零元素尾地址
Ne=&N.data[N.tu]; // Ne指向矩陣N的非零元素尾地址
Qh=Qe=Q.data; // Qh、Qe的初值指向矩陣Q的非零元素首地址的前一地址
while(Mp<=Me&&Np<=Ne)
{
Qe++;
switch(comp(Mp->i,Np->i))
{
case 1: *Qe=*Mp;
Mp++;
break;
case 0: switch(comp(Mp->j,Np->j)) // M、N矩陣當前非零元素的行相等,繼續比較列
{
case 1: *Qe=*Mp;
Mp++;
break;
case 0: *Qe=*Mp;
Qe->e+=Np->e;
if(!Qe->e) // 元素值為0，不存入壓縮矩陣
Qe--;
Mp++;
Np++;
break;
case -1: *Qe=*Np;
Np++;
}
break;
case -1: *Qe=*Np;
Np++;
}
}
if(Mp>Me) // 矩陣M的元素全部處理完畢
while(Np<=Ne)
{
Qe++;
*Qe=*Np;
Np++;
}
if(Np>Ne) // 矩陣N的元素全部處理完畢
while(Mp<=Me)
{
Qe++;
*Qe=*Mp;
Mp++;
}
Q.tu=Qe-Qh; // 矩陣Q的非零元素個數
return OK;
}
Status SubtSMatrix(TSMatrix M,TSMatrix N,TSMatrix &Q)
{ // 求稀疏矩陣的差Q=M-N
int i;
for(i=1;i<=N.tu;i++)
N.data[i].e*=-1;
AddSMatrix(M,N,Q);
return OK;
}
Status MultSMatrix(TSMatrix M,TSMatrix N,TSMatrix &Q)
{ // 求稀疏矩陣的乘積Q=M*N
int i,j,h=M.mu,l=N.nu,Qn=0;
// h,l分別為矩陣Q的行、列值,Qn為矩陣Q的非零元素個數，初值為0
ElemType *Qe;
if(M.nu!=N.mu)
return ERROR;
Q.mu=M.mu;
Q.nu=N.nu;
Qe=(ElemType *)malloc(h*l*sizeof(ElemType)); // Qe為矩陣Q的臨時數組
// 矩陣Q的第i行j列的元素值存於*(Qe+(i-1)*l+j-1)中，初值為0
for(i=0;i<h*l;i++)
*(Qe+i)=0; // 賦初值0
for(i=1;i<=M.tu;i++) // 矩陣元素相乘，結果累加到Qe
for(j=1;j<=N.tu;j++)
if(M.data[i].j==N.data[j].i)
*(Qe+(M.data[i].i-1)*l+N.data[j].j-1)+=M.data[i].e*N.data[j].e;
for(i=1;i<=M.mu;i++)
for(j=1;j<=N.nu;j++)
if(*(Qe+(i-1)*l+j-1)!=0)
{
Qn++;
Q.data[Qn].e=*(Qe+(i-1)*l+j-1);
Q.data[Qn].i=i;
Q.data[Qn].j=j;
}
free(Qe);
Q.tu=Qn;
return OK;
}
Status TransposeSMatrix(TSMatrix M,TSMatrix &T)
{ // 求稀疏矩陣M的轉置矩陣T。演算法5.1
int p,q,col;
T.mu=M.nu;
T.nu=M.mu;
T.tu=M.tu;
if(T.tu)
{
q=1;
for(col=1;col<=M.nu;++col)
for(p=1;p<=M.tu;++p)
if(M.data[p].j==col)
{
T.data[q].i=M.data[p].j;
T.data[q].j=M.data[p].i;
T.data[q].e=M.data[p].e;
++q;
}
}
return OK;
}
Status FastTransposeSMatrix(TSMatrix M,TSMatrix &T)
{ // 快速求稀疏矩陣M的轉置矩陣T。演算法5.2
int p,q,t,col,*num,*cpot;
num=(int *)malloc((M.nu+1)*sizeof(int)); // 生成數組（[0]不用）
cpot=(int *)malloc((M.nu+1)*sizeof(int)); // 生成數組（[0]不用）
T.mu=M.nu;
T.nu=M.mu;
T.tu=M.tu;
if(T.tu)
{
for(col=1;col<=M.nu;++col)
num[col]=0; // 設初值
for(t=1;t<=M.tu;++t) // 求M中每一列含非零元素個數
++num[M.data[t].j];
cpot[1]=1;
for(col=2;col<=M.nu;++col) // 求第col列中第一個非零元在T.data中的序號
cpot[col]=cpot[col-1]+num[col-1];
for(p=1;p<=M.tu;++p)
{
col=M.data[p].j;
q=cpot[col];
T.data[q].i=M.data[p].j;
T.data[q].j=M.data[p].i;
T.data[q].e=M.data[p].e;
++cpot[col];
}
}
free(num);
free(cpot);
return OK;
}
void main()
{
TSMatrix A,B;
printf("創建矩陣A: ");
CreateSMatrix(A);
PrintSMatrix(A);
FastTransposeSMatrix(A,B);
printf("矩陣B(A的快速轉置): ");
PrintSMatrix(B);
DestroySMatrix(A);
DestroySMatrix(B);
}
轉置的，你參考下
8. 稀疏矩陣的壓縮存儲只需要存儲什麼
非零元素。
對於一個用二維數組存儲的稀疏矩陣Amn，如果假設存儲每個數組元素需要L個位元組，那麼存儲整個矩陣需要m*n*L個位元組。但是，這些存儲空間的大部分存放的是0元素，從而造成大量的空間浪費。為了節省存儲空間，可以只存儲其中的非0元素。
(8)稀疏矩陣演算法擴展閱讀
稀疏矩陣演算法的最大特點是通過只存儲和處理非零元素從而大幅度降低存儲空間需求以及計算復雜度，代價則是必須使用專門的稀疏矩陣壓縮存儲數據結構。稀疏矩陣演算法是典型的不規則演算法，計算訪存比很低，並且計算過程中的訪存軌跡與稀疏矩陣的稀疏結構相關。