1 2 #if !defined(__BAIJ_H) 3 #define __BAIJ_H 4 #include "src/mat/matimpl.h" 5 6 7 /* 8 MATSEQBAIJ format - Block compressed row storage. The i[] and j[] 9 arrays start at 0. 10 */ 11 12 /* This header is shared by the SeqSBAIJ matrix */ 13 #define SEQBAIJHEADER \ 14 PetscTruth sorted; /* if true, rows are sorted by increasing columns */ \ 15 PetscTruth roworiented; /* if true, row-oriented input, default */ \ 16 PetscInt nonew; /* 1 don't add new nonzeros, -1 generate error on new */ \ 17 PetscTruth singlemalloc; /* if true a, i, and j have been obtained with \ 18 one big malloc */ \ 19 PetscInt bs2; /* square of block size */ \ 20 PetscInt mbs,nbs; /* rows/bs, columns/bs */ \ 21 PetscInt nz,maxnz; /* nonzeros, allocated nonzeros */ \ 22 PetscInt *diag; /* pointers to diagonal elements */ \ 23 PetscInt *i; /* pointer to beginning of each row */ \ 24 PetscInt *imax; /* maximum space allocated for each row */ \ 25 PetscInt *ilen; /* actual length of each row */ \ 26 PetscInt *j; /* column values: j + i[k] - 1 is start of row k */ \ 27 MatScalar *a; /* nonzero elements */ \ 28 IS row,col,icol; /* index sets, used for reorderings */ \ 29 PetscScalar *solve_work; /* work space used in MatSolve */ \ 30 PetscInt reallocs; /* number of mallocs done during MatSetValues() \ 31 as more values are set then were preallocated */ \ 32 PetscScalar *mult_work; /* work array for matrix vector product*/ \ 33 PetscScalar *saved_values; \ 34 \ 35 PetscTruth keepzeroedrows; /* if true, MatZeroRows() will not change nonzero structure */ \ 36 Mat sbaijMat; /* mat in sbaij format */ \ 37 \ 38 PetscInt setvalueslen; /* only used for single precision */ \ 39 MatScalar *setvaluescopy; /* area double precision values in MatSetValuesXXX() are copied \ 40 before calling MatSetValuesXXX_SeqBAIJ_MatScalar() */ \ 41 \ 42 PetscTruth pivotinblocks; /* pivot inside factorization of each diagonal block */ \ 43 \ 44 PetscInt *xtoy,*xtoyB; /* map nonzero pattern of X into Y's, used by MatAXPY() */ \ 45 Mat XtoY; /* used by MatAXPY() */ \ 46 PetscScalar *idiag; /* inverse of block diagonal */ \ 47 PetscTruth idiagvalid; /* if above has correct/current values */ \ 48 Mat_CompressedRow compressedrow; /* use compressed row format */ \ 49 PetscTruth freedata; /* free the i,j,a data when the matrix is destroyed; true by default */ 50 51 typedef struct { 52 SEQBAIJHEADER 53 } Mat_SeqBAIJ; 54 55 EXTERN_C_BEGIN 56 EXTERN PetscErrorCode MatSeqBAIJSetPreallocation_SeqBAIJ(Mat,PetscInt,PetscInt,PetscInt*); 57 EXTERN_C_END 58 EXTERN PetscErrorCode MatILUFactorSymbolic_SeqBAIJ(Mat,IS,IS,MatFactorInfo*,Mat *); 59 EXTERN PetscErrorCode MatICCFactorSymbolic_SeqBAIJ(Mat,IS,MatFactorInfo*,Mat *); 60 EXTERN PetscErrorCode MatCholeskyFactorSymbolic_SeqBAIJ(Mat,IS,MatFactorInfo*,Mat*); 61 EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqBAIJ_N(Mat,MatFactorInfo*,Mat*); 62 EXTERN PetscErrorCode MatCholeskyFactorNumeric_SeqBAIJ_N_NaturalOrdering(Mat,MatFactorInfo*,Mat*); 63 EXTERN PetscErrorCode MatDuplicate_SeqBAIJ(Mat,MatDuplicateOption,Mat*); 64 EXTERN PetscErrorCode MatMissingDiagonal_SeqBAIJ(Mat); 65 EXTERN PetscErrorCode MatMarkDiagonal_SeqBAIJ(Mat); 66 67 EXTERN PetscErrorCode MatLUFactorSymbolic_SeqBAIJ(Mat,IS,IS,MatFactorInfo*,Mat*); 68 EXTERN PetscErrorCode MatLUFactor_SeqBAIJ(Mat,IS,IS,MatFactorInfo*); 69 EXTERN PetscErrorCode MatIncreaseOverlap_SeqBAIJ(Mat,PetscInt,IS*,PetscInt); 70 EXTERN PetscErrorCode MatGetSubMatrix_SeqBAIJ(Mat,IS,IS,PetscInt,MatReuse,Mat*); 71 EXTERN PetscErrorCode MatGetSubMatrices_SeqBAIJ(Mat,PetscInt,const IS[],const IS[],MatReuse,Mat*[]); 72 EXTERN PetscErrorCode MatMultTranspose_SeqBAIJ(Mat,Vec,Vec); 73 EXTERN PetscErrorCode MatMultTransposeAdd_SeqBAIJ(Mat,Vec,Vec,Vec); 74 EXTERN PetscErrorCode MatScale_SeqBAIJ(Mat,PetscScalar); 75 EXTERN PetscErrorCode MatNorm_SeqBAIJ(Mat,NormType,PetscReal *); 76 EXTERN PetscErrorCode MatEqual_SeqBAIJ(Mat,Mat,PetscTruth*); 77 EXTERN PetscErrorCode MatGetDiagonal_SeqBAIJ(Mat,Vec); 78 EXTERN PetscErrorCode MatDiagonalScale_SeqBAIJ(Mat,Vec,Vec); 79 EXTERN PetscErrorCode MatGetInfo_SeqBAIJ(Mat,MatInfoType,MatInfo *); 80 EXTERN PetscErrorCode MatZeroEntries_SeqBAIJ(Mat); 81 82 EXTERN PetscErrorCode MatSeqBAIJ_UpdateFactorNumeric_NaturalOrdering(Mat); 83 EXTERN PetscErrorCode MatSeqBAIJ_UpdateSolvers(Mat); 84 85 EXTERN PetscErrorCode MatSolve_SeqBAIJ_Update(Mat,Vec,Vec); 86 EXTERN PetscErrorCode MatSolve_SeqBAIJ_1(Mat,Vec,Vec); 87 EXTERN PetscErrorCode MatSolve_SeqBAIJ_1_NaturalOrdering(Mat,Vec,Vec); 88 EXTERN PetscErrorCode MatSolve_SeqBAIJ_2(Mat,Vec,Vec); 89 EXTERN PetscErrorCode MatSolve_SeqBAIJ_2_NaturalOrdering(Mat,Vec,Vec); 90 EXTERN PetscErrorCode MatSolve_SeqBAIJ_3(Mat,Vec,Vec); 91 EXTERN PetscErrorCode MatSolve_SeqBAIJ_3_NaturalOrdering(Mat,Vec,Vec); 92 EXTERN PetscErrorCode MatSolve_SeqBAIJ_4(Mat,Vec,Vec); 93 EXTERN PetscErrorCode MatSolve_SeqBAIJ_4_NaturalOrdering(Mat,Vec,Vec); 94 #if defined(PETSC_HAVE_SSE) 95 EXTERN PetscErrorCode MatSolve_SeqBAIJ_4_SSE_Demotion(Mat,Vec,Vec); 96 EXTERN PetscErrorCode MatSolve_SeqBAIJ_4_NaturalOrdering_SSE_Demotion(Mat,Vec,Vec); 97 EXTERN PetscErrorCode MatSolve_SeqBAIJ_4_NaturalOrdering_SSE_Demotion_usj(Mat,Vec,Vec); 98 #endif 99 EXTERN PetscErrorCode MatSolve_SeqBAIJ_5(Mat,Vec,Vec); 100 EXTERN PetscErrorCode MatSolve_SeqBAIJ_5_NaturalOrdering(Mat,Vec,Vec); 101 EXTERN PetscErrorCode MatSolve_SeqBAIJ_6(Mat,Vec,Vec); 102 EXTERN PetscErrorCode MatSolve_SeqBAIJ_6_NaturalOrdering(Mat,Vec,Vec); 103 EXTERN PetscErrorCode MatSolve_SeqBAIJ_7(Mat,Vec,Vec); 104 EXTERN PetscErrorCode MatSolve_SeqBAIJ_7_NaturalOrdering(Mat,Vec,Vec); 105 EXTERN PetscErrorCode MatSolve_SeqBAIJ_N(Mat,Vec,Vec); 106 107 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_Update(Mat,Vec,Vec); 108 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_1(Mat,Vec,Vec); 109 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_1_NaturalOrdering(Mat,Vec,Vec); 110 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_2(Mat,Vec,Vec); 111 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_2_NaturalOrdering(Mat,Vec,Vec); 112 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_3(Mat,Vec,Vec); 113 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_3_NaturalOrdering(Mat,Vec,Vec); 114 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_4(Mat,Vec,Vec); 115 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_4_NaturalOrdering(Mat,Vec,Vec); 116 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_5(Mat,Vec,Vec); 117 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_5_NaturalOrdering(Mat,Vec,Vec); 118 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_6(Mat,Vec,Vec); 119 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_6_NaturalOrdering(Mat,Vec,Vec); 120 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_7(Mat,Vec,Vec); 121 EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_7_NaturalOrdering(Mat,Vec,Vec); 122 123 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_1(Mat,MatFactorInfo*,Mat*); 124 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_2(Mat,MatFactorInfo*,Mat*); 125 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_2_NaturalOrdering(Mat,MatFactorInfo*,Mat*); 126 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_3(Mat,MatFactorInfo*,Mat*); 127 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_3_NaturalOrdering(Mat,MatFactorInfo*,Mat*); 128 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_4(Mat,MatFactorInfo*,Mat*); 129 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_4_NaturalOrdering(Mat,MatFactorInfo*,Mat*); 130 #if defined(PETSC_HAVE_SSE) 131 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_4_NaturalOrdering_SSE(Mat,MatFactorInfo*,Mat*); 132 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_4_NaturalOrdering_SSE_usj(Mat,MatFactorInfo*,Mat*); 133 #else 134 #endif 135 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_5(Mat,MatFactorInfo*,Mat*); 136 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_5_NaturalOrdering(Mat,MatFactorInfo*,Mat*); 137 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_6(Mat,MatFactorInfo*,Mat*); 138 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_6_NaturalOrdering(Mat,MatFactorInfo*,Mat*); 139 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_7(Mat,MatFactorInfo*,Mat*); 140 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_7_NaturalOrdering(Mat,MatFactorInfo*,Mat*); 141 EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_N(Mat,MatFactorInfo*,Mat*); 142 143 EXTERN PetscErrorCode MatMult_SeqBAIJ_1(Mat,Vec,Vec); 144 EXTERN PetscErrorCode MatMult_SeqBAIJ_2(Mat,Vec,Vec); 145 EXTERN PetscErrorCode MatMult_SeqBAIJ_3(Mat,Vec,Vec); 146 EXTERN PetscErrorCode MatMult_SeqBAIJ_4(Mat,Vec,Vec); 147 EXTERN PetscErrorCode MatMult_SeqBAIJ_5(Mat,Vec,Vec); 148 EXTERN PetscErrorCode MatMult_SeqBAIJ_6(Mat,Vec,Vec); 149 EXTERN PetscErrorCode MatMult_SeqBAIJ_7(Mat,Vec,Vec); 150 EXTERN PetscErrorCode MatMult_SeqBAIJ_N(Mat,Vec,Vec); 151 152 EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_1(Mat,Vec,Vec,Vec); 153 EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_2(Mat,Vec,Vec,Vec); 154 EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_3(Mat,Vec,Vec,Vec); 155 EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_4(Mat,Vec,Vec,Vec); 156 EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_5(Mat,Vec,Vec,Vec); 157 EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_6(Mat,Vec,Vec,Vec); 158 EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_7(Mat,Vec,Vec,Vec); 159 EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_N(Mat,Vec,Vec,Vec); 160 EXTERN PetscErrorCode MatLoad_SeqBAIJ(PetscViewer, MatType,Mat*); 161 162 #endif 163