2014-10-26 34 views
1

我有一个存储在compressed row storage(CRS)中的相对较大的(例如5000行×8000列)和稀疏矩阵。我正在尝试获取其compressed column storage(CCS)表格。如何将压缩行存储转换为稀疏矩阵的压缩列存储?

这样做是否已经有一个标准算法?一种选择可以是从CRS重建整个矩阵(4000万条目),然后使用简单的算法来获得其CCS。但是,这种时间复杂性非常糟糕,我计划在更大的矩阵上使用这种算法。有关如何做到这一点的任何其他想法?

回答

2

可能效率不高的数值方法代码,但我想出了这似乎工作:

#include <stdio.h> 
#include <string.h> 

#define COLS 6 
#define SIZE(a) (sizeof(a)/sizeof(*(a))) 

int main() { 
    float f[] = {10,-2, 3, 9, 3, 7, 8, 7, 3, 8, 7, 5, 8, 9, 9,13, 4, 2, 1}; 
    int c[] = { 0, 4, 0, 1, 5, 1, 2, 3, 0, 2, 3, 4, 1, 3, 4, 5, 1, 4, 5}; 
    int r[] = { 0, 2,  5,  8,   12,   16,  19}; 
    float nf[SIZE(f)]; 
    int nc[COLS+1] = {0}; 
    int nr[SIZE(f)]; 
    int nn[COLS+1]; 

    int rr[SIZE(f)]; 
    for (int k = 0, i = 0; i < SIZE(r); i++) 
    for (int j = 0; j < r[i+1] - r[i]; j++) 
     rr[k++] = i; 

    for (int i = 0; i < SIZE(f); i++) 
    nc[c[i]+1]++; 
    for (int i = 1; i <= COLS; i++) 
    nc[i] += nc[i-1]; 
    memcpy(nn, nc, sizeof(nc)); 

    for (int i = 0; i < SIZE(f); i++) { 
    int x = nn[c[i]]++; 
    nf[x] = f[i]; 
    nr[x] = rr[i]; 
    } 

    for (int i = 0; i < SIZE(nf); i++) printf("%2.0f ", nf[i]); 
    putchar('\n'); 
    for (int i = 0; i < SIZE(nr); i++) printf("%2d ", nr[i]); 
    putchar('\n'); 
    for (int i = 0; i < SIZE(nc); i++) printf("%2d ", nc[i]); 
    putchar('\n'); 

    return 0; 
} 
1

似乎有一种类似于标准方法的东西,因为在数值配方中描述了一种算法。我会在这里引用代码给你的想法,而更多的细节你可以参考2.7章节。的第三版。

NRsparseMat NRsparseMat::transpose() const { 
    Int i,j,k,index,m=nrows,n=ncols; 
    NRsparseMat at(n,m,nvals); //Initialized to zero. 

    //First find the column lengths for AT , i.e. the row lengths of A. 
    VecInt count(m,0); //Temporary counters for each row of A. 
    for (i=0;i<n;i++) 
     for (j=col_ptr[i];j<col_ptr[i+1];j++) { 
      k=row_ind[j]; 
      count[k]++; 
     } 

    for (j=0;j<m;j++) //Now set at.col_ptr. 0th entry stays 0. 
     at.col_ptr[j+1]=at.col_ptr[j]+count[j]; 

    for(j=0;j<m;j++) //Reset counters to zero. 
     count[j]=0; 

    for (i=0;i<n;i++) //Main loop. 
     for (j=col_ptr[i];j<col_ptr[i+1];j++) { 
      k=row_ind[j]; 
      index=at.col_ptr[k]+count[k]; //Element’s position in column of AT . 
      at.row_ind[index]=i; 
      at.val[index]=val[j]; 
      count[k]++; //Increment counter for next element in that column. 
     } 
    return at; 
} 

对于我个人用的,我通常从数值方法去除它重新编写代码的具体类型定义(如IntVecInt),重命名,重新格式化等