18965ea79SLois Curfman McInnes 270f55243SBarry Smith #include "src/mat/impls/dense/seq/dense.h" 38965ea79SLois Curfman McInnes 4052cd425SLois Curfman McInnes /* Data stuctures for basic parallel dense matrix */ 5052cd425SLois Curfman McInnes 6052cd425SLois Curfman McInnes /* Structure to hold the information for factorization of a dense matrix */ 726b0ef5dSLois Curfman McInnes /* Most of this info is used in the pipe send/recv routines */ 8052cd425SLois Curfman McInnes typedef struct { 913f74950SBarry Smith PetscInt nlnr; /* number of local rows downstream */ 1013f74950SBarry Smith PetscInt nrend; /* rend for downstream processor */ 1113f74950SBarry Smith PetscInt nbr,pnbr; /* Down and upstream neighbors */ 1213f74950SBarry Smith PetscInt *tag; /* message tags */ 1313f74950SBarry Smith PetscInt currow; /* current row number */ 1413f74950SBarry Smith PetscInt phase; /* phase (used to indicate tag) */ 1513f74950SBarry Smith PetscInt up; /* Are we moving up or down in row number? */ 1613f74950SBarry Smith PetscInt use_bcast; /* Are we broadcasting max length? */ 1713f74950SBarry Smith PetscInt nsend; /* number of sends */ 1813f74950SBarry Smith PetscInt nrecv; /* number of receives */ 1926b0ef5dSLois Curfman McInnes 2026b0ef5dSLois Curfman McInnes /* data initially in matrix context */ 2113f74950SBarry Smith PetscInt k; /* Blocking factor (unused as yet) */ 2213f74950SBarry Smith PetscInt k2; /* Blocking factor for solves */ 23ea709b57SSatish Balay PetscScalar *temp; 2413f74950SBarry Smith PetscInt nlptr; 2513f74950SBarry Smith PetscInt *lrows; 2613f74950SBarry Smith PetscInt *nlrows; 2713f74950SBarry Smith PetscInt *pivots; 2826b0ef5dSLois Curfman McInnes } FactorCtx; 29052cd425SLois Curfman McInnes 30052cd425SLois Curfman McInnes #define PIPEPHASE (ctx->phase == 0) 3167e560aaSBarry Smith 328965ea79SLois Curfman McInnes typedef struct { 3313f74950SBarry Smith PetscInt nvec; /* this is the n size for the vector one multiplies with */ 34bad5519cSBarry Smith Mat A; /* local submatrix */ 3513f74950SBarry Smith PetscMPIInt size; /* size of communicator */ 3613f74950SBarry Smith PetscMPIInt rank; /* rank of proc in communicator */ 378965ea79SLois Curfman McInnes /* The following variables are used for matrix assembly */ 38273d9f13SBarry Smith PetscTruth donotstash; /* Flag indicationg if values should be stashed */ 398965ea79SLois Curfman McInnes MPI_Request *send_waits; /* array of send requests */ 408965ea79SLois Curfman McInnes MPI_Request *recv_waits; /* array of receive requests */ 4113f74950SBarry Smith PetscInt nsends,nrecvs; /* numbers of sends and receives */ 42ea709b57SSatish Balay PetscScalar *svalues,*rvalues; /* sending and receiving data */ 4313f74950SBarry Smith PetscInt rmax; /* maximum message length */ 448965ea79SLois Curfman McInnes 458965ea79SLois Curfman McInnes /* The following variables are used for matrix-vector products */ 468965ea79SLois Curfman McInnes 478965ea79SLois Curfman McInnes Vec lvec; /* local vector */ 48bad5519cSBarry Smith VecScatter Mvctx; /* scatter context for vector */ 49052cd425SLois Curfman McInnes 50273d9f13SBarry Smith PetscTruth roworiented; /* if true, row oriented input (default) */ 51622d7880SLois Curfman McInnes FactorCtx *factor; /* factorization context */ 528965ea79SLois Curfman McInnes } Mat_MPIDense; 5397304618SKris Buschelman 54f69a0ea3SMatthew Knepley EXTERN PetscErrorCode MatLoad_MPIDense(PetscViewer, MatType,Mat*); 55dfbe8321SBarry Smith EXTERN PetscErrorCode MatSetUpMultiply_MPIDense(Mat); 5613f74950SBarry Smith EXTERN PetscErrorCode MatGetSubMatrices_MPIDense(Mat,PetscInt,const IS[],const IS[],MatReuse,Mat *[]); 576e4ee0c6SHong Zhang EXTERN PetscErrorCode MatEqual_MPIDense(Mat,Mat,PetscTruth*); 58*4ae313f4SHong Zhang EXTERN PetscErrorCode MatMatMultSymbolic_MPIDense_MPIDense(Mat,Mat,PetscReal,Mat*); 59*4ae313f4SHong Zhang EXTERN PetscErrorCode MatMatMult_MPIAIJ_MPIDense(Mat,Mat,MatReuse,PetscReal,Mat*); 60*4ae313f4SHong Zhang EXTERN PetscErrorCode MatMatMultSymbolic_MPIAIJ_MPIDense(Mat,Mat,PetscReal,Mat*); 61*4ae313f4SHong Zhang EXTERN PetscErrorCode MatMatMultNumeric_MPIAIJ_MPIDense(Mat,Mat,Mat); 62