Actual source code: sbaij.h
4: #include include/private/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: SEQAIJHEADER(PetscScalar);
14: SEQBAIJHEADER;
15: PetscInt *inew; /* pointer to beginning of each row of reordered matrix */
16: PetscInt *jnew; /* column values: jnew + i[k] is start of row k */
17: MatScalar *anew; /* nonzero diagonal and superdiagonal elements of reordered matrix */
18: PetscScalar *solves_work; /* work space used in MatSolves */
19: PetscInt solves_work_n;/* size of solves_work */
20: PetscScalar *relax_work;
21: PetscInt *a2anew; /* map used for symm permutation */
22: PetscTruth permute; /* if true, a non-trivial permutation is used for factorization */
23: PetscTruth ignore_ltriangular; /* if true, ignore the lower triangular values inserted by users */
24: PetscTruth getrow_utriangular; /* if true, MatGetRow_SeqSBAIJ() is enabled to get the upper part of the row */
25: } Mat_SeqSBAIJ;
28: EXTERN PetscErrorCode MatSeqSBAIJSetPreallocation_SeqSBAIJ(Mat,PetscInt,PetscInt,PetscInt*);
30: EXTERN PetscErrorCode MatCholeskyFactorSymbolic_SeqSBAIJ(Mat,IS,MatFactorInfo*,Mat*);
31: EXTERN PetscErrorCode MatCholeskyFactor_SeqSBAIJ(Mat,IS,MatFactorInfo*);
32: EXTERN PetscErrorCode MatICCFactorSymbolic_SeqSBAIJ(Mat,IS,MatFactorInfo*,Mat *);
33: EXTERN PetscErrorCode MatDuplicate_SeqSBAIJ(Mat,MatDuplicateOption,Mat*);
34: EXTERN PetscErrorCode MatMarkDiagonal_SeqSBAIJ(Mat);
35: EXTERN PetscErrorCode MatIncreaseOverlap_SeqSBAIJ(Mat,PetscInt,IS[],PetscInt);
36: EXTERN PetscErrorCode MatGetSubMatrix_SeqSBAIJ(Mat,IS,IS,PetscInt,MatReuse,Mat*);
37: EXTERN PetscErrorCode MatGetSubMatrices_SeqSBAIJ(Mat,PetscInt,const IS[],const IS[],MatReuse,Mat*[]);
38: EXTERN PetscErrorCode MatScale_SeqSBAIJ(Mat,PetscScalar);
39: EXTERN PetscErrorCode MatNorm_SeqSBAIJ(Mat,NormType,PetscReal *);
40: EXTERN PetscErrorCode MatEqual_SeqSBAIJ(Mat,Mat,PetscTruth*);
41: EXTERN PetscErrorCode MatGetDiagonal_SeqSBAIJ(Mat,Vec);
42: EXTERN PetscErrorCode MatDiagonalScale_SeqSBAIJ(Mat,Vec,Vec);
43: EXTERN PetscErrorCode MatGetInfo_SeqSBAIJ(Mat,MatInfoType,MatInfo *);
44: EXTERN PetscErrorCode MatZeroEntries_SeqSBAIJ(Mat);
45: EXTERN PetscErrorCode MatGetRowMax_SeqSBAIJ(Mat,Vec);
46: EXTERN PetscErrorCode MatGetInertia_SeqSBAIJ(Mat,PetscInt*,PetscInt*,PetscInt*);
48: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_1_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
49: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_1_NaturalOrdering(Mat,Vec,Vec);
50: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_1_NaturalOrdering(Mat,Vec,Vec);
51: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_1_NaturalOrdering(Mat,Vec,Vec);
52: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_1(Mat,Vec,Vec);
53: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_1(Mat,Vec,Vec);
55: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_2_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
56: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_2_NaturalOrdering(Mat,Vec,Vec);
57: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_2_NaturalOrdering(Mat,Vec,Vec);
58: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_2_NaturalOrdering(Mat,Vec,Vec);
60: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_3_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
61: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_3_NaturalOrdering(Mat,Vec,Vec);
62: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_3_NaturalOrdering(Mat,Vec,Vec);
63: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_3_NaturalOrdering(Mat,Vec,Vec);
65: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_4_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
66: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_4_NaturalOrdering(Mat,Vec,Vec);
67: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_4_NaturalOrdering(Mat,Vec,Vec);
68: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_4_NaturalOrdering(Mat,Vec,Vec);
70: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_5_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
71: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_5_NaturalOrdering(Mat,Vec,Vec);
72: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_5_NaturalOrdering(Mat,Vec,Vec);
73: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_5_NaturalOrdering(Mat,Vec,Vec);
75: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_6_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
76: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_6_NaturalOrdering(Mat,Vec,Vec);
77: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_6_NaturalOrdering(Mat,Vec,Vec);
78: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_6_NaturalOrdering(Mat,Vec,Vec);
80: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_7_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
81: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_7_NaturalOrdering(Mat,Vec,Vec);
82: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_7_NaturalOrdering(Mat,Vec,Vec);
83: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_7_NaturalOrdering(Mat,Vec,Vec);
85: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_N_NaturalOrdering(Mat,MatFactorInfo*,Mat*);
86: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_N_NaturalOrdering(Mat,Vec,Vec);
87: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_N_NaturalOrdering(Mat,Vec,Vec);
88: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_N_NaturalOrdering(Mat,Vec,Vec);
89: EXTERN PetscErrorCode MatForwardSolve_SeqSBAIJ_N(Mat,Vec,Vec);
90: EXTERN PetscErrorCode MatBackwardSolve_SeqSBAIJ_N(Mat,Vec,Vec);
92: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_N(Mat,MatFactorInfo*,Mat*);
93: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_1(Mat,MatFactorInfo*,Mat*);
94: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_2(Mat,MatFactorInfo*,Mat*);
95: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_3(Mat,MatFactorInfo*,Mat*);
96: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_4(Mat,MatFactorInfo*,Mat*);
97: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_5(Mat,MatFactorInfo*,Mat*);
98: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_6(Mat,MatFactorInfo*,Mat*);
99: EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqSBAIJ_7(Mat,MatFactorInfo*,Mat*);
101: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_N(Mat,Vec,Vec);
102: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_1(Mat,Vec,Vec);
103: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_2(Mat,Vec,Vec);
104: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_3(Mat,Vec,Vec);
105: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_4(Mat,Vec,Vec);
106: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_5(Mat,Vec,Vec);
107: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_6(Mat,Vec,Vec);
108: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_7(Mat,Vec,Vec);
110: EXTERN PetscErrorCode MatSolves_SeqSBAIJ_1(Mat,Vecs,Vecs);
112: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_1_NaturalOrdering(Mat,Vec,Vec);
113: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_2_NaturalOrdering(Mat,Vec,Vec);
114: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_3_NaturalOrdering(Mat,Vec,Vec);
115: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_4_NaturalOrdering(Mat,Vec,Vec);
116: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_5_NaturalOrdering(Mat,Vec,Vec);
117: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_6_NaturalOrdering(Mat,Vec,Vec);
118: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_7_NaturalOrdering(Mat,Vec,Vec);
119: EXTERN PetscErrorCode MatSolve_SeqSBAIJ_N_NaturalOrdering(Mat,Vec,Vec);
121: EXTERN PetscErrorCode MatMult_SeqSBAIJ_1(Mat,Vec,Vec);
122: EXTERN PetscErrorCode MatMult_SeqSBAIJ_2(Mat,Vec,Vec);
123: EXTERN PetscErrorCode MatMult_SeqSBAIJ_3(Mat,Vec,Vec);
124: EXTERN PetscErrorCode MatMult_SeqSBAIJ_4(Mat,Vec,Vec);
125: EXTERN PetscErrorCode MatMult_SeqSBAIJ_5(Mat,Vec,Vec);
126: EXTERN PetscErrorCode MatMult_SeqSBAIJ_6(Mat,Vec,Vec);
127: EXTERN PetscErrorCode MatMult_SeqSBAIJ_7(Mat,Vec,Vec);
128: EXTERN PetscErrorCode MatMult_SeqSBAIJ_N(Mat,Vec,Vec);
130: EXTERN PetscErrorCode MatMultAdd_SeqSBAIJ_1(Mat,Vec,Vec,Vec);
131: EXTERN PetscErrorCode MatMultAdd_SeqSBAIJ_2(Mat,Vec,Vec,Vec);
132: EXTERN PetscErrorCode MatMultAdd_SeqSBAIJ_3(Mat,Vec,Vec,Vec);
133: EXTERN PetscErrorCode MatMultAdd_SeqSBAIJ_4(Mat,Vec,Vec,Vec);
134: EXTERN PetscErrorCode MatMultAdd_SeqSBAIJ_5(Mat,Vec,Vec,Vec);
135: EXTERN PetscErrorCode MatMultAdd_SeqSBAIJ_6(Mat,Vec,Vec,Vec);
136: EXTERN PetscErrorCode MatMultAdd_SeqSBAIJ_7(Mat,Vec,Vec,Vec);
137: EXTERN PetscErrorCode MatMultAdd_SeqSBAIJ_N(Mat,Vec,Vec,Vec);
139: EXTERN PetscErrorCode MatRelax_SeqSBAIJ(Mat,Vec,PetscReal,MatSORType,PetscReal,PetscInt,PetscInt,Vec);
140: EXTERN PetscErrorCode MatLoad_SeqSBAIJ(PetscViewer, MatType,Mat*);
142: #endif