1 2 #include <../src/mat/impls/dense/seq/dense.h> 3 4 /* Data stuctures for basic parallel dense matrix */ 5 6 typedef struct { /* used by MatTransposeMatMult_MPIDense_MPIDense() */ 7 PetscScalar *sendbuf,*recvbuf; 8 PetscInt *recvcounts; 9 Mat C_loc; 10 PetscErrorCode (*destroy)(Mat); 11 } Mat_TransMatMultDense; 12 13 typedef struct { 14 PetscInt nvec; /* this is the n size for the vector one multiplies with */ 15 Mat A; /* local submatrix */ 16 PetscMPIInt size; /* size of communicator */ 17 PetscMPIInt rank; /* rank of proc in communicator */ 18 /* The following variables are used for matrix assembly */ 19 PetscBool donotstash; /* Flag indicationg if values should be stashed */ 20 MPI_Request *send_waits; /* array of send requests */ 21 MPI_Request *recv_waits; /* array of receive requests */ 22 PetscInt nsends,nrecvs; /* numbers of sends and receives */ 23 PetscScalar *svalues,*rvalues; /* sending and receiving data */ 24 PetscInt rmax; /* maximum message length */ 25 26 /* The following variables are used for matrix-vector products */ 27 Vec lvec; /* local vector */ 28 VecScatter Mvctx; /* scatter context for vector */ 29 PetscBool roworiented; /* if true, row oriented input (default) */ 30 31 Mat_MatTransMatMult *atb; /* used by MatTransposeMatMult_MPIAIJ_MPIDense */ 32 Mat_TransMatMultDense *atbdense; /* used by MatTransposeMatMult_MPIDense_MPIDense */ 33 } Mat_MPIDense; 34 35 PETSC_INTERN PetscErrorCode MatLoad_MPIDense(Mat,PetscViewer); 36 PETSC_INTERN PetscErrorCode MatSetUpMultiply_MPIDense(Mat); 37 PETSC_INTERN PetscErrorCode MatGetSubMatrices_MPIDense(Mat,PetscInt,const IS[],const IS[],MatReuse,Mat *[]); 38 PETSC_INTERN PetscErrorCode MatEqual_MPIDense(Mat,Mat,PetscBool*); 39 PETSC_INTERN PetscErrorCode MatMatMultSymbolic_MPIDense_MPIDense(Mat,Mat,PetscReal,Mat*); 40 PETSC_INTERN PetscErrorCode MatMatMult_MPIAIJ_MPIDense(Mat,Mat,MatReuse,PetscReal,Mat*); 41 PETSC_INTERN PetscErrorCode MatMatMultSymbolic_MPIAIJ_MPIDense(Mat,Mat,PetscReal,Mat*); 42 PETSC_INTERN PetscErrorCode MatMatMultNumeric_MPIAIJ_MPIDense(Mat,Mat,Mat); 43 PETSC_INTERN PetscErrorCode MatTransposeMatMult_MPIDense_MPIDense(Mat,Mat,MatReuse,PetscReal,Mat*); 44 PETSC_INTERN PetscErrorCode MatTransposeMatMultSymbolic_MPIDense_MPIDense(Mat,Mat,PetscReal,Mat*); 45 PETSC_INTERN PetscErrorCode MatTransposeMatMultNumeric_MPIDense_MPIDense(Mat,Mat,Mat); 46