Actual source code: sbaij.h
4: #include src/mat/matimpl.h
5: #include src/mat/impls/baij/seq/baij.h
7: /*
8: MATSEQSBAIJ format - Block compressed row storage. The i[] and j[]
9: arrays start at 0.
10: */
12: typedef struct {
13: SEQBAIJHEADER
14: PetscInt *inew; /* pointer to beginning of each row of reordered matrix */
15: PetscInt *jnew; /* column values: jnew + i[k] is start of row k */
16: MatScalar *anew; /* nonzero diagonal and superdiagonal elements of reordered matrix */
17: PetscScalar *solves_work; /* work space used in MatSolves */
18: PetscInt solves_work_n;/* size of solves_work */
19: PetscInt *a2anew; /* map used for symm permutation */
20: PetscTruth permute; /* if true, a non-trivial permutation is used for factorization */
21: PetscTruth ignore_ltriangular; /* if true, ignore the lower triangular values inserted by users */
22: PetscTruth getrow_utriangular; /* if true, MatGetRow_SeqSBAIJ() is enabled to get the upper part of the row */
23: } Mat_SeqSBAIJ;
26: EXTERN PetscErrorCode MatSeqSBAIJSetPreallocation_SeqSBAIJ(Mat,PetscInt,PetscInt,PetscInt*);
28: EXTERN PetscErrorCode MatICCFactorSymbolic_SeqSBAIJ(Mat,IS,MatFactorInfo*,Mat *);
29: EXTERN PetscErrorCode MatDuplicate_SeqSBAIJ(Mat,MatDuplicateOption,Mat*);
30: EXTERN PetscErrorCode MatMarkDiagonal_SeqSBAIJ(Mat);
32: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_1_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
33: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_1_NaturalOrdering(Mat,Vec,Vec);
35: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_2_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
36: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_2_NaturalOrdering(Mat,Vec,Vec);
38: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_3_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
39: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_3_NaturalOrdering(Mat,Vec,Vec);
41: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_4_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
42: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_4_NaturalOrdering(Mat,Vec,Vec);
44: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_5_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
45: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_5_NaturalOrdering(Mat,Vec,Vec);
47: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_6_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
48: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_6_NaturalOrdering(Mat,Vec,Vec);
50: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_7_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
51: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_7_NaturalOrdering(Mat,Vec,Vec);
53: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_N_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
54: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_N_NaturalOrdering(Mat,Vec,Vec);
56: EXTERN PetscErrorCode MatRelax_SeqSBAIJ(Mat,Vec,PetscReal,MatSORType,PetscReal,PetscInt,PetscInt,Vec);
57: EXTERN PetscErrorCode MatLoad_SeqSBAIJ(PetscViewer, MatType,Mat*);
59: #endif