oct_sparse.cpp

Go to the documentation of this file.
00001 #include <oct.h>
00002 #include "gs_oct.h"
00003 #include "grpc.h"
00004 
00005 int sparsematrix_octcsc_to_gscsr(
00006     int m, int n, int nnz, int datatype,
00007     double *a, double *ai, int *colind, int *rowptr,
00008     void *at, int *rowind, int *colptr) {
00009     int i, j, col, relpos;
00010     int *marker;
00011 
00012     marker = (int *) malloc(sizeof(int) * n);
00013     if (marker == NULL) return -1;
00014 
00015     /* Get counts of each column of A, and set up column pointers */
00016     for (i = 0; i < m; ++i)
00017         for (j = rowptr[i]; j < rowptr[i + 1]; ++j)
00018             ++marker[colind[j]];
00019 
00020     colptr[0] = 0;
00021 
00022     for (j = 0; j < n; ++j) {
00023         colptr[j + 1] = colptr[j] + marker[j];
00024         marker[j] = colptr[j];
00025     }
00026  
00027     /* Transfer the matrix into the compressed column storage. */
00028     for (i = 0; i < m; ++i) {
00029         for (j = rowptr[i]; j < rowptr[i + 1]; ++j) {
00030             col = colind[j];
00031             relpos = marker[col];
00032             rowind[relpos] = i;
00033 
00034             switch(datatype) {
00035                 case GS_INT:
00036                 ((int *) at)[relpos] = (int) ((double *) a)[j];
00037                 break;
00038                 case GS_FLOAT:
00039                     ((float *) at)[relpos] = (float) ((double *) a)[j];
00040                     break;
00041                 case GS_DOUBLE:
00042                     ((double *) at)[relpos] = (double) ((double *) a)[j];
00043                     break;
00044                 case GS_SCOMPLEX:
00045                     ((gs_scomplex *) at)[relpos].r = (float) ((double *) a)[j];
00046                     ((gs_scomplex *) at)[relpos].i = (float) ((double *) ai)[j];
00047                     break;
00048                 case GS_DCOMPLEX:
00049                     /* ((dcomplex *)at)[relpos] = ((dcomplex *)a)[j]; */
00050                     ((gs_dcomplex *) at)[relpos].r = (double) ((double *) a)[j];
00051                     ((gs_dcomplex *) at)[relpos].i = (double) ((double *) ai)[j];
00052                     break;
00053                 default:
00054                     return -1;
00055                     break;
00056             }
00057             ++marker[col];
00058         }
00059     }
00060 
00061     return 0;
00062 }
00063 
00064 
00065 int sparsematrix_gscsr_to_octcsc(
00066     int m, int n, int nnz, int datatype,
00067     void *a, int *colind, int *rowptr,
00068     double *at, double *ati, int *rowind, int *colptr) {
00069     int i, j, col, relpos;
00070     int *marker;
00071 
00072     marker = (int *)matlab_gs_calloc(n, sizeof(int));
00073     if (marker == NULL) return -1;
00074  
00075     /* Get counts of each column of A, and set up column pointers */
00076     for (i = 0; i < m; ++i)
00077         for (j = rowptr[i]; j < rowptr[i + 1]; ++j)
00078             ++marker[colind[j]];
00079 
00080     colptr[0] = 0;
00081 
00082     for (j = 0; j < n; ++j) {
00083         colptr[j + 1] = colptr[j] + marker[j];
00084         marker[j] = colptr[j];
00085     }
00086  
00087     /* Transfer the matrix into the compressed column storage. */
00088     for (i = 0; i < m; ++i) {
00089         for (j = rowptr[i]; j < rowptr[i + 1]; ++j) {
00090             col = colind[j];
00091             relpos = marker[col];
00092             rowind[relpos] = i;
00093       
00094             switch(datatype) {
00095                 case GS_INT:
00096                     ((double *) at)[relpos] = (double) ((int *) a)[j];
00097                     break;
00098                 case GS_FLOAT:
00099                     ((double *) at)[relpos] = (double) ((float *) a)[j];
00100                     break;
00101                 case GS_DOUBLE:
00102                     ((double *) at)[relpos] = ((double *) a)[j];
00103                     break;
00104                 case GS_SCOMPLEX:
00105                     ((double *) at)[relpos] = (float) (((gs_scomplex *) a)[j].r);
00106                     ((double *) ati)[relpos] = (float) (((gs_scomplex *) a)[j].i);
00107                     break;
00108                 case GS_DCOMPLEX:
00109                     ((double *) at)[relpos] = ((gs_dcomplex *) a)[j].r;
00110                     ((double *) ati)[relpos] = ((gs_dcomplex *) a)[j].i;
00111                     break;
00112                 default:
00113                     return -1;
00114                     break;
00115             }
00116             ++marker[col];
00117         }
00118     }
00119 
00120     matlab_gs_free(marker);
00121     return 0;
00122 }