Actual source code: bdiag.h

  1: /* $Id: bdiag.h,v 1.31 2001/08/07 03:02:53 balay Exp $ */

 5:  #include src/mat/matimpl.h


  8: /*
  9:    Mat_SeqBDiag (MATSEQBDIAG) - block-diagonal format, where each diagonal
 10:    element consists of a square block of size bs  x bs.  Dense storage
 11:    within each block is in column-major order.  The diagonals are the
 12:    full length of the matrix.  As a special case, blocks of size bs=1
 13:    (scalars) are supported as well.
 14: */

 16: typedef struct {
 17:   int          mblock,nblock;    /* block rows and columns */
 18:   int          nonew;            /* if true, no new nonzeros allowed in matrix */
 19:   int          nonew_diag;       /* if true, no new diagonals allowed in matrix */
 20:   int          nz,maxnz;         /* nonzeros, allocated nonzeros */
 21:   int          nd;               /* number of block diagonals */
 22:   int          mainbd;           /* the number of the main block diagonal */
 23:   int          bs;               /* Each diagonal element is an bs x bs matrix */
 24:   int          *diag;            /* value of (row-col)/bs for each diagonal */
 25:   int          *bdlen;           /* block-length of each diagonal */
 26:   int          ndim;             /* diagonals come from an ndim pde (if 0, ignore) */
 27:   int          ndims[3];         /* sizes of the mesh if ndim > 0 */
 28:   PetscTruth   user_alloc;       /* true if the user provided the diagonals */
 29:   int          *colloc;          /* holds column locations if using MatGetRow */
 30:   PetscScalar  **diagv;          /* The actual diagonals */
 31:   PetscScalar  *dvalue;          /* Used to hold a row if MatGetRow is used */
 32:   int          *pivot;           /* pivots for LU factorization (temporary loc) */
 33:   PetscTruth   roworiented;      /* inputs to MatSetValue() are row oriented (default = 1) */
 34:   int          reallocs;         /* number of allocations during MatSetValues */
 35: } Mat_SeqBDiag;

 37: EXTERN int MatNorm_SeqBDiag_Columns(Mat,PetscReal*,int);
 38: EXTERN int MatMult_SeqBDiag_1(Mat,Vec,Vec);
 39: EXTERN int MatMult_SeqBDiag_2(Mat,Vec,Vec);
 40: EXTERN int MatMult_SeqBDiag_3(Mat,Vec,Vec);
 41: EXTERN int MatMult_SeqBDiag_4(Mat,Vec,Vec);
 42: EXTERN int MatMult_SeqBDiag_5(Mat,Vec,Vec);
 43: EXTERN int MatMult_SeqBDiag_N(Mat A,Vec,Vec);
 44: EXTERN int MatMultAdd_SeqBDiag_1(Mat,Vec,Vec,Vec);
 45: EXTERN int MatMultAdd_SeqBDiag_2(Mat,Vec,Vec,Vec);
 46: EXTERN int MatMultAdd_SeqBDiag_3(Mat,Vec,Vec,Vec);
 47: EXTERN int MatMultAdd_SeqBDiag_4(Mat,Vec,Vec,Vec);
 48: EXTERN int MatMultAdd_SeqBDiag_5(Mat,Vec,Vec,Vec);
 49: EXTERN int MatMultAdd_SeqBDiag_N(Mat A,Vec,Vec,Vec);
 50: EXTERN int MatMultTranspose_SeqBDiag_1(Mat,Vec,Vec);
 51: EXTERN int MatMultTranspose_SeqBDiag_N(Mat A,Vec,Vec);
 52: EXTERN int MatMultTransposeAdd_SeqBDiag_1(Mat,Vec,Vec,Vec);
 53: EXTERN int MatMultTransposeAdd_SeqBDiag_N(Mat A,Vec,Vec,Vec);
 54: EXTERN int MatSetValues_SeqBDiag_1(Mat,int,const int [],int,const int [],const PetscScalar [],InsertMode);
 55: EXTERN int MatSetValues_SeqBDiag_N(Mat,int,const int [],int,const int [],const PetscScalar [],InsertMode);
 56: EXTERN int MatGetValues_SeqBDiag_1(Mat,int,const int [],int,const int [],PetscScalar []);
 57: EXTERN int MatGetValues_SeqBDiag_N(Mat,int,const int [],int,const int [],PetscScalar []);
 58: EXTERN int MatRelax_SeqBDiag_1(Mat,Vec,PetscReal,MatSORType,PetscReal,int,int,Vec);
 59: EXTERN int MatRelax_SeqBDiag_N(Mat,Vec,PetscReal,MatSORType,PetscReal,int,int,Vec);
 60: EXTERN int MatView_SeqBDiag(Mat,PetscViewer);
 61: EXTERN int MatGetInfo_SeqBDiag(Mat,MatInfoType,MatInfo*);
 62: EXTERN int MatGetRow_SeqBDiag(Mat,int,int *,int **,PetscScalar **);
 63: EXTERN int MatRestoreRow_SeqBDiag(Mat,int,int *,int **,PetscScalar **);
 64: EXTERN int MatTranspose_SeqBDiag(Mat,Mat *);
 65: EXTERN int MatNorm_SeqBDiag(Mat,NormType,PetscReal *);
 66: EXTERN int MatLUFactorSymbolic_SeqBDiag(Mat,IS,IS,MatFactorInfo*,Mat*);
 67: EXTERN int MatILUFactorSymbolic_SeqBDiag(Mat,IS,IS,MatFactorInfo*,Mat*);
 68: EXTERN int MatILUFactor_SeqBDiag(Mat,IS,IS,MatFactorInfo*);
 69: EXTERN int MatLUFactorNumeric_SeqBDiag_N(Mat,Mat*);
 70: EXTERN int MatLUFactorNumeric_SeqBDiag_1(Mat,Mat*);
 71: EXTERN int MatSolve_SeqBDiag_1(Mat,Vec,Vec);
 72: EXTERN int MatSolve_SeqBDiag_2(Mat,Vec,Vec);
 73: EXTERN int MatSolve_SeqBDiag_3(Mat,Vec,Vec);
 74: EXTERN int MatSolve_SeqBDiag_4(Mat,Vec,Vec);
 75: EXTERN int MatSolve_SeqBDiag_5(Mat,Vec,Vec);
 76: EXTERN int MatSolve_SeqBDiag_N(Mat,Vec,Vec);
 77: EXTERN int MatLoad_SeqBDiag(PetscViewer,const MatType,Mat*);

 79: #endif