1 #ifndef MPIDENSE_H 2 #define MPIDENSE_H 3 4 #include <../src/mat/impls/dense/seq/dense.h> 5 #include <petscsf.h> 6 7 /* Data structures for basic parallel dense matrix */ 8 9 typedef struct { /* used by MatMatMultxxx_MPIDense_MPIDense() */ 10 Mat Ae, Be, Ce; /* matrix in Elemental format */ 11 } Mat_MatMultDense; 12 13 typedef struct { /* used by MatTransposeMatMultXXX_MPIDense_MPIDense() */ 14 PetscScalar *sendbuf; 15 Mat atb; 16 PetscMPIInt *recvcounts; 17 PetscMPIInt tag; 18 } Mat_TransMatMultDense; 19 20 typedef struct { /* used by MatMatTransposeMultxxx_MPIDense_MPIDense() */ 21 PetscScalar *buf[2]; 22 PetscMPIInt tag; 23 PetscMPIInt *recvcounts; 24 PetscMPIInt *recvdispls; 25 PetscInt alg; /* algorithm used */ 26 } Mat_MatTransMultDense; 27 28 typedef struct { 29 Mat A; /* local submatrix */ 30 31 /* The following variables are used for matrix assembly */ 32 PetscBool donotstash; /* Flag indicating if values should be stashed */ 33 MPI_Request *send_waits; /* array of send requests */ 34 MPI_Request *recv_waits; /* array of receive requests */ 35 PetscInt nsends, nrecvs; /* numbers of sends and receives */ 36 PetscScalar *svalues, *rvalues; /* sending and receiving data */ 37 PetscInt rmax; /* maximum message length */ 38 39 /* The following variables are used for matrix-vector products */ 40 Vec lvec; /* local vector */ 41 PetscSF Mvctx; /* for mat-mult communications */ 42 PetscBool roworiented; /* if true, row oriented input (default) */ 43 44 /* Support for MatDenseGetColumnVec and MatDenseGetSubMatrix */ 45 Mat cmat; /* matrix representation of a given subset of columns */ 46 Vec cvec; /* vector representation of a given column */ 47 const PetscScalar *ptrinuse; /* holds array to be restored (just a placeholder) */ 48 PetscInt vecinuse; /* if cvec is in use (col = vecinuse-1) */ 49 PetscInt matinuse; /* if cmat is in use (cbegin = matinuse-1) */ 50 } Mat_MPIDense; 51 52 PETSC_INTERN PetscErrorCode MatSetUpMultiply_MPIDense(Mat); 53 PETSC_INTERN PetscErrorCode MatCreateSubMatrices_MPIDense(Mat, PetscInt, const IS[], const IS[], MatReuse, Mat *[]); 54 PETSC_INTERN PetscErrorCode MatProductSetFromOptions_MPIAIJ_MPIDense(Mat); 55 56 PETSC_INTERN PetscErrorCode MatProductSetFromOptions_MPIDense_MPIAIJ(Mat); 57 PETSC_INTERN PetscErrorCode MatProductSetFromOptions_MPIDense(Mat); 58 59 #if defined(PETSC_HAVE_ELEMENTAL) 60 PETSC_INTERN PetscErrorCode MatMatMultSymbolic_Elemental(Mat, Mat, PetscReal, Mat); 61 PETSC_INTERN PetscErrorCode MatMatMultNumeric_Elemental(Mat, Mat, Mat); 62 #endif 63 PETSC_INTERN PetscErrorCode MatConvert_SeqDense_SeqAIJ(Mat A, MatType newtype, MatReuse reuse, Mat *newmat); 64 65 PETSC_INTERN PetscErrorCode MatShift_MPIDense(Mat, PetscScalar); 66 PETSC_INTERN PetscErrorCode MatDenseGetColumnVecWrite_MPIDense(Mat, PetscInt, Vec *); 67 PETSC_INTERN PetscErrorCode MatDenseGetColumnVecRead_MPIDense(Mat, PetscInt, Vec *); 68 PETSC_INTERN PetscErrorCode MatDenseGetColumnVec_MPIDense(Mat, PetscInt, Vec *); 69 PETSC_INTERN PetscErrorCode MatDenseRestoreColumnVecWrite_MPIDense(Mat, PetscInt, Vec *); 70 PETSC_INTERN PetscErrorCode MatDenseRestoreColumnVecRead_MPIDense(Mat, PetscInt, Vec *); 71 PETSC_INTERN PetscErrorCode MatDenseRestoreColumnVec_MPIDense(Mat, PetscInt, Vec *); 72 73 PETSC_INTERN PetscErrorCode MatCreate_MPIDense(Mat); 74 PETSC_INTERN PetscErrorCode MatGetDiagonal_MPIDense(Mat, Vec); 75 76 #if PetscDefined(HAVE_CUDA) 77 PETSC_INTERN PetscErrorCode MatConvert_MPIDense_MPIDenseCUDA(Mat, MatType, MatReuse, Mat *); 78 #endif 79 80 #if PetscDefined(HAVE_HIP) 81 PETSC_INTERN PetscErrorCode MatConvert_MPIDense_MPIDenseHIP(Mat, MatType, MatReuse, Mat *); 82 #endif 83 84 #endif 85