12593348eSBarry Smith 235aab85fSBarry Smith #if !defined(__BAIJ_H) 335aab85fSBarry Smith #define __BAIJ_H 4da33ede1SBarry Smith #include "src/mat/matimpl.h" 5da33ede1SBarry Smith 62593348eSBarry Smith 72593348eSBarry Smith /* 835aab85fSBarry Smith MATSEQBAIJ format - Block compressed row storage. The i[] and j[] 9e24b481bSBarry Smith arrays start at 0. 102593348eSBarry Smith */ 112593348eSBarry Smith 126c6c5352SBarry Smith /* This header is shared by the SeqSBAIJ matrix */ 136c6c5352SBarry Smith #define SEQBAIJHEADER \ 146c6c5352SBarry Smith PetscTruth sorted; /* if true, rows are sorted by increasing columns */ \ 156c6c5352SBarry Smith PetscTruth roworiented; /* if true, row-oriented input, default */ \ 16690b6cddSBarry Smith PetscInt nonew; /* 1 don't add new nonzeros, -1 generate error on new */ \ 176c6c5352SBarry Smith PetscTruth singlemalloc; /* if true a, i, and j have been obtained with \ 186c6c5352SBarry Smith one big malloc */ \ 19521d7252SBarry Smith PetscInt bs2; /* square of block size */ \ 20690b6cddSBarry Smith PetscInt mbs,nbs; /* rows/bs, columns/bs */ \ 21690b6cddSBarry Smith PetscInt nz,maxnz; /* nonzeros, allocated nonzeros */ \ 22690b6cddSBarry Smith PetscInt *diag; /* pointers to diagonal elements */ \ 23690b6cddSBarry Smith PetscInt *i; /* pointer to beginning of each row */ \ 24690b6cddSBarry Smith PetscInt *imax; /* maximum space allocated for each row */ \ 25690b6cddSBarry Smith PetscInt *ilen; /* actual length of each row */ \ 26690b6cddSBarry Smith PetscInt *j; /* column values: j + i[k] - 1 is start of row k */ \ 276c6c5352SBarry Smith MatScalar *a; /* nonzero elements */ \ 286c6c5352SBarry Smith IS row,col,icol; /* index sets, used for reorderings */ \ 296c6c5352SBarry Smith PetscScalar *solve_work; /* work space used in MatSolve */ \ 30690b6cddSBarry Smith PetscInt reallocs; /* number of mallocs done during MatSetValues() \ 316c6c5352SBarry Smith as more values are set then were preallocated */ \ 326c6c5352SBarry Smith PetscScalar *mult_work; /* work array for matrix vector product*/ \ 336c6c5352SBarry Smith PetscScalar *saved_values; \ 346c6c5352SBarry Smith \ 356c6c5352SBarry Smith PetscTruth keepzeroedrows; /* if true, MatZeroRows() will not change nonzero structure */ \ 366c6c5352SBarry Smith \ 37690b6cddSBarry Smith PetscInt setvalueslen; /* only used for single precision */ \ 386c6c5352SBarry Smith MatScalar *setvaluescopy; /* area double precision values in MatSetValuesXXX() are copied \ 396c6c5352SBarry Smith before calling MatSetValuesXXX_SeqBAIJ_MatScalar() */ \ 406c6c5352SBarry Smith \ 416c6c5352SBarry Smith PetscTruth pivotinblocks; /* pivot inside factorization of each diagonal block */ \ 426c6c5352SBarry Smith \ 43690b6cddSBarry Smith PetscInt *xtoy,*xtoyB; /* map nonzero pattern of X into Y's, used by MatAXPY() */ \ 44a1d92eedSBarry Smith Mat XtoY; /* used by MatAXPY() */ \ 45b01c7715SBarry Smith PetscScalar *idiag; /* inverse of block diagonal */ \ 46*73e7a558SHong Zhang PetscTruth idiagvalid; /* if above has correct/current values */ \ 47*73e7a558SHong Zhang Mat_CompressedRow compressedrow; /* use compressed row format */ 486c6c5352SBarry Smith 496c6c5352SBarry Smith typedef struct { 506c6c5352SBarry Smith SEQBAIJHEADER 5135aab85fSBarry Smith } Mat_SeqBAIJ; 522593348eSBarry Smith 53dfbe8321SBarry Smith EXTERN PetscErrorCode MatILUFactorSymbolic_SeqBAIJ(Mat,IS,IS,MatFactorInfo*,Mat *); 54dfbe8321SBarry Smith EXTERN PetscErrorCode MatDuplicate_SeqBAIJ(Mat,MatDuplicateOption,Mat*); 55dfbe8321SBarry Smith EXTERN PetscErrorCode MatMissingDiagonal_SeqBAIJ(Mat); 56dfbe8321SBarry Smith EXTERN PetscErrorCode MatMarkDiagonal_SeqBAIJ(Mat); 57be3590efSBarry Smith 58dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorSymbolic_SeqBAIJ(Mat,IS,IS,MatFactorInfo*,Mat*); 59dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactor_SeqBAIJ(Mat,IS,IS,MatFactorInfo*); 60690b6cddSBarry Smith EXTERN PetscErrorCode MatIncreaseOverlap_SeqBAIJ(Mat,PetscInt,IS*,PetscInt); 61690b6cddSBarry Smith EXTERN PetscErrorCode MatGetSubMatrix_SeqBAIJ(Mat,IS,IS,PetscInt,MatReuse,Mat*); 62690b6cddSBarry Smith EXTERN PetscErrorCode MatGetSubMatrices_SeqBAIJ(Mat,PetscInt,const IS[],const IS[],MatReuse,Mat*[]); 63dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultTranspose_SeqBAIJ(Mat,Vec,Vec); 64dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultTransposeAdd_SeqBAIJ(Mat,Vec,Vec,Vec); 65dfbe8321SBarry Smith EXTERN PetscErrorCode MatScale_SeqBAIJ(const PetscScalar*,Mat); 66dfbe8321SBarry Smith EXTERN PetscErrorCode MatNorm_SeqBAIJ(Mat,NormType,PetscReal *); 67dfbe8321SBarry Smith EXTERN PetscErrorCode MatEqual_SeqBAIJ(Mat,Mat,PetscTruth*); 68dfbe8321SBarry Smith EXTERN PetscErrorCode MatGetDiagonal_SeqBAIJ(Mat,Vec); 69dfbe8321SBarry Smith EXTERN PetscErrorCode MatDiagonalScale_SeqBAIJ(Mat,Vec,Vec); 70dfbe8321SBarry Smith EXTERN PetscErrorCode MatGetInfo_SeqBAIJ(Mat,MatInfoType,MatInfo *); 71dfbe8321SBarry Smith EXTERN PetscErrorCode MatZeroEntries_SeqBAIJ(Mat); 7254138f6bSKris Buschelman 73dfbe8321SBarry Smith EXTERN PetscErrorCode MatSeqBAIJ_UpdateFactorNumeric_NaturalOrdering(Mat); 74dfbe8321SBarry Smith EXTERN PetscErrorCode MatSeqBAIJ_UpdateSolvers(Mat); 7554138f6bSKris Buschelman 76dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_Update(Mat,Vec,Vec); 77dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_1(Mat,Vec,Vec); 78dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_1_NaturalOrdering(Mat,Vec,Vec); 79dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_2(Mat,Vec,Vec); 80dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_2_NaturalOrdering(Mat,Vec,Vec); 81dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_3(Mat,Vec,Vec); 82dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_3_NaturalOrdering(Mat,Vec,Vec); 83dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_4(Mat,Vec,Vec); 84dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_4_NaturalOrdering(Mat,Vec,Vec); 8554138f6bSKris Buschelman #if defined(PETSC_HAVE_SSE) 86dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_4_SSE_Demotion(Mat,Vec,Vec); 87dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_4_NaturalOrdering_SSE_Demotion(Mat,Vec,Vec); 88dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_4_NaturalOrdering_SSE_Demotion_usj(Mat,Vec,Vec); 894110d75bSKris Buschelman #endif 90dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_5(Mat,Vec,Vec); 91dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_5_NaturalOrdering(Mat,Vec,Vec); 92dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_6(Mat,Vec,Vec); 93dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_6_NaturalOrdering(Mat,Vec,Vec); 94dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_7(Mat,Vec,Vec); 95dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_7_NaturalOrdering(Mat,Vec,Vec); 96dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolve_SeqBAIJ_N(Mat,Vec,Vec); 9754138f6bSKris Buschelman 98dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_Update(Mat,Vec,Vec); 99dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_1(Mat,Vec,Vec); 100dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_1_NaturalOrdering(Mat,Vec,Vec); 101dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_2(Mat,Vec,Vec); 102dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_2_NaturalOrdering(Mat,Vec,Vec); 103dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_3(Mat,Vec,Vec); 104dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_3_NaturalOrdering(Mat,Vec,Vec); 105dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_4(Mat,Vec,Vec); 106dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_4_NaturalOrdering(Mat,Vec,Vec); 107dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_5(Mat,Vec,Vec); 108dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_5_NaturalOrdering(Mat,Vec,Vec); 109dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_6(Mat,Vec,Vec); 110dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_6_NaturalOrdering(Mat,Vec,Vec); 111dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_7(Mat,Vec,Vec); 112dfbe8321SBarry Smith EXTERN PetscErrorCode MatSolveTranspose_SeqBAIJ_7_NaturalOrdering(Mat,Vec,Vec); 11354138f6bSKris Buschelman 114dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_1(Mat,Mat*); 115dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_2(Mat,Mat*); 116dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_2_NaturalOrdering(Mat,Mat*); 117dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_3(Mat,Mat*); 118dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_3_NaturalOrdering(Mat,Mat*); 119dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_4(Mat,Mat*); 120dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_4_NaturalOrdering(Mat,Mat*); 12154138f6bSKris Buschelman #if defined(PETSC_HAVE_SSE) 122dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_4_NaturalOrdering_SSE(Mat,Mat*); 123dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_4_NaturalOrdering_SSE_usj(Mat,Mat*); 12454138f6bSKris Buschelman #else 12554138f6bSKris Buschelman #endif 126dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_5(Mat,Mat*); 127dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_5_NaturalOrdering(Mat,Mat*); 128dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_6(Mat,Mat*); 129dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_6_NaturalOrdering(Mat,Mat*); 130dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_7(Mat,Mat*); 131dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_7_NaturalOrdering(Mat,Mat*); 132dfbe8321SBarry Smith EXTERN PetscErrorCode MatLUFactorNumeric_SeqBAIJ_N(Mat,Mat*); 13354138f6bSKris Buschelman 134dfbe8321SBarry Smith EXTERN PetscErrorCode MatMult_SeqBAIJ_1(Mat,Vec,Vec); 135dfbe8321SBarry Smith EXTERN PetscErrorCode MatMult_SeqBAIJ_2(Mat,Vec,Vec); 136dfbe8321SBarry Smith EXTERN PetscErrorCode MatMult_SeqBAIJ_3(Mat,Vec,Vec); 137dfbe8321SBarry Smith EXTERN PetscErrorCode MatMult_SeqBAIJ_4(Mat,Vec,Vec); 138dfbe8321SBarry Smith EXTERN PetscErrorCode MatMult_SeqBAIJ_5(Mat,Vec,Vec); 139dfbe8321SBarry Smith EXTERN PetscErrorCode MatMult_SeqBAIJ_6(Mat,Vec,Vec); 140dfbe8321SBarry Smith EXTERN PetscErrorCode MatMult_SeqBAIJ_7(Mat,Vec,Vec); 141dfbe8321SBarry Smith EXTERN PetscErrorCode MatMult_SeqBAIJ_N(Mat,Vec,Vec); 14254138f6bSKris Buschelman 143dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_1(Mat,Vec,Vec,Vec); 144dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_2(Mat,Vec,Vec,Vec); 145dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_3(Mat,Vec,Vec,Vec); 146dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_4(Mat,Vec,Vec,Vec); 147dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_5(Mat,Vec,Vec,Vec); 148dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_6(Mat,Vec,Vec,Vec); 149dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_7(Mat,Vec,Vec,Vec); 150dfbe8321SBarry Smith EXTERN PetscErrorCode MatMultAdd_SeqBAIJ_N(Mat,Vec,Vec,Vec); 151dfbe8321SBarry Smith EXTERN PetscErrorCode MatLoad_SeqBAIJ(PetscViewer,const MatType,Mat*); 1522593348eSBarry Smith 1532593348eSBarry Smith #endif 154