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