1*b2f1ab58SBarry Smith /* 2*b2f1ab58SBarry Smith Define type spbas_matrix: sparse matrices using pointers 3*b2f1ab58SBarry Smith 4*b2f1ab58SBarry Smith Global matrix information 5*b2f1ab58SBarry Smith nrows, ncols: dimensions 6*b2f1ab58SBarry Smith nnz : number of nonzeros (in whole matrix) 7*b2f1ab58SBarry Smith col_idx_type: storage scheme for column numbers 8*b2f1ab58SBarry Smith SPBAS_COLUMN_NUMBERS: 9*b2f1ab58SBarry Smith array icol contains column indices: 10*b2f1ab58SBarry Smith A(i,icol[i][j]) = values[i][j] 11*b2f1ab58SBarry Smith SPBAS_DIAGONAL_OFFSETS: 12*b2f1ab58SBarry Smith array icol contains diagonal offsets: 13*b2f1ab58SBarry Smith A(i,i+icol[i][j]) = values[i][j] 14*b2f1ab58SBarry Smith SPBAS_OFFSET_ARRAY: 15*b2f1ab58SBarry Smith array icol contains offsets wrt array 16*b2f1ab58SBarry Smith icol0: 17*b2f1ab58SBarry Smith A(i,icol0[i]+icol[i][j]) = values[i][j] 18*b2f1ab58SBarry Smith 19*b2f1ab58SBarry Smith Information about each row 20*b2f1ab58SBarry Smith row_nnz : number of nonzeros for each row 21*b2f1ab58SBarry Smith icol0 : column index offset (when needed, otherwise NULL) 22*b2f1ab58SBarry Smith icols : array of diagonal offsets for each row, as descibed 23*b2f1ab58SBarry Smith for col_idx_type, above 24*b2f1ab58SBarry Smith values : array of matrix entries for each row 25*b2f1ab58SBarry Smith when values == NULL, this matrix is really 26*b2f1ab58SBarry Smith a sparseness pattern, not a matrix 27*b2f1ab58SBarry Smith 28*b2f1ab58SBarry Smith The other fields describe the way in which the data are stored 29*b2f1ab58SBarry Smith in memory. 30*b2f1ab58SBarry Smith 31*b2f1ab58SBarry Smith block_data : The pointers icols[i] all point to places in a 32*b2f1ab58SBarry Smith single allocated array. Only for icols[0] was 33*b2f1ab58SBarry Smith malloc called. Freeing icols[0] will free 34*b2f1ab58SBarry Smith all other icols=arrays as well. 35*b2f1ab58SBarry Smith Same for arrays values[i] 36*b2f1ab58SBarry Smith */ 37*b2f1ab58SBarry Smith 38*b2f1ab58SBarry Smith #define SPBAS_COLUMN_NUMBERS (0) 39*b2f1ab58SBarry Smith #define SPBAS_DIAGONAL_OFFSETS (1) 40*b2f1ab58SBarry Smith #define SPBAS_OFFSET_ARRAY (2) 41*b2f1ab58SBarry Smith 42*b2f1ab58SBarry Smith #define NEGATIVE_DIAGONAL (-42) 43*b2f1ab58SBarry Smith 44*b2f1ab58SBarry Smith typedef struct { 45*b2f1ab58SBarry Smith PetscInt nrows; 46*b2f1ab58SBarry Smith PetscInt ncols; 47*b2f1ab58SBarry Smith PetscInt nnz; 48*b2f1ab58SBarry Smith PetscInt col_idx_type; 49*b2f1ab58SBarry Smith 50*b2f1ab58SBarry Smith PetscInt *row_nnz; 51*b2f1ab58SBarry Smith PetscInt *icol0; 52*b2f1ab58SBarry Smith PetscInt **icols; 53*b2f1ab58SBarry Smith PetscScalar **values; 54*b2f1ab58SBarry Smith 55*b2f1ab58SBarry Smith PetscTruth block_data; 56*b2f1ab58SBarry Smith PetscInt n_alloc_icol; 57*b2f1ab58SBarry Smith PetscInt n_alloc_val; 58*b2f1ab58SBarry Smith PetscInt *alloc_icol; 59*b2f1ab58SBarry Smith PetscScalar *alloc_val; 60*b2f1ab58SBarry Smith } spbas_matrix ; 61*b2f1ab58SBarry Smith 62*b2f1ab58SBarry Smith 63*b2f1ab58SBarry Smith /* 64*b2f1ab58SBarry Smith spbas_compress_pattern: 65*b2f1ab58SBarry Smith calculate a compressed sparseness pattern for a sparseness pattern 66*b2f1ab58SBarry Smith given in compressed row storage. The compressed sparseness pattern may 67*b2f1ab58SBarry Smith require (much) less memory. 68*b2f1ab58SBarry Smith 69*b2f1ab58SBarry Smith spbas_memory_requirement: 70*b2f1ab58SBarry Smith Calculate the number of bytes needed to store tha matrix 71*b2f1ab58SBarry Smith 72*b2f1ab58SBarry Smith spbas_incomplete_cholesky: 73*b2f1ab58SBarry Smith Incomplete Cholesky decomposition 74*b2f1ab58SBarry Smith 75*b2f1ab58SBarry Smith spbas_delete: 76*b2f1ab58SBarry Smith de-allocate the arrays owned by this matrix 77*b2f1ab58SBarry Smith 78*b2f1ab58SBarry Smith spbas_matrix_to_crs: 79*b2f1ab58SBarry Smith Convert an spbas_matrix to compessed row storage 80*b2f1ab58SBarry Smith 81*b2f1ab58SBarry Smith spbas_dump: 82*b2f1ab58SBarry Smith Print the matrix in i,j,val-format 83*b2f1ab58SBarry Smith 84*b2f1ab58SBarry Smith spbas_transpose: 85*b2f1ab58SBarry Smith Return the transpose of a matrix 86*b2f1ab58SBarry Smith 87*b2f1ab58SBarry Smith spbas_pattern_only: 88*b2f1ab58SBarry Smith Return the sparseness pattern (matrix without values) of a 89*b2f1ab58SBarry Smith compressed row storage 90*b2f1ab58SBarry Smith */ 91*b2f1ab58SBarry Smith PetscErrorCode spbas_compress_pattern(PetscInt*,PetscInt*,PetscInt,PetscInt,PetscInt,spbas_matrix*,PetscScalar*); 92*b2f1ab58SBarry Smith long int spbas_memory_requirement( spbas_matrix); 93*b2f1ab58SBarry Smith PetscErrorCode spbas_delete(spbas_matrix); 94*b2f1ab58SBarry Smith PetscErrorCode spbas_incomplete_cholesky(Mat,const PetscInt*,const PetscInt*,spbas_matrix,PetscScalar,PetscScalar,spbas_matrix*); 95*b2f1ab58SBarry Smith PetscErrorCode spbas_matrix_to_crs(spbas_matrix, MatScalar **,PetscInt **,PetscInt**); 96*b2f1ab58SBarry Smith PetscErrorCode spbas_dump(const char*,spbas_matrix); 97*b2f1ab58SBarry Smith PetscErrorCode spbas_transpose(spbas_matrix,spbas_matrix*); 98*b2f1ab58SBarry Smith PetscErrorCode spbas_apply_reordering(spbas_matrix*, const PetscInt *, const PetscInt*); 99*b2f1ab58SBarry Smith PetscErrorCode spbas_pattern_only(PetscInt, PetscInt, PetscInt *, PetscInt*, spbas_matrix *); 100*b2f1ab58SBarry Smith PetscErrorCode spbas_power (spbas_matrix, PetscInt, spbas_matrix *); 101*b2f1ab58SBarry Smith PetscErrorCode spbas_keep_upper(spbas_matrix *); 102*b2f1ab58SBarry Smith 103