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