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 *);