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