1 /* 2 Defines matrix-matrix-matrix product routines for MPIAIJ matrices 3 D = A * B * C 4 */ 5 #include <../src/mat/impls/aij/mpi/mpiaij.h> /*I "petscmat.h" I*/ 6 7 #if defined(PETSC_HAVE_HYPRE) 8 PETSC_INTERN PetscErrorCode MatTransposeMatMatMultSymbolic_AIJ_AIJ_AIJ_wHYPRE(Mat,Mat,Mat,PetscReal,Mat*); 9 PETSC_INTERN PetscErrorCode MatTransposeMatMatMultNumeric_AIJ_AIJ_AIJ_wHYPRE(Mat,Mat,Mat,Mat); 10 11 PETSC_INTERN PetscErrorCode MatMatMatMult_Transpose_AIJ_AIJ(Mat R,Mat A,Mat P,MatReuse scall, PetscReal fill, Mat *RAP) 12 { 13 Mat Rt; 14 PetscBool flg; 15 PetscErrorCode ierr; 16 17 PetscFunctionBegin; 18 ierr = MatTransposeGetMat(R,&Rt);CHKERRQ(ierr); 19 ierr = PetscObjectTypeCompareAny((PetscObject)Rt,&flg,MATSEQAIJ,MATMPIAIJ,NULL);CHKERRQ(ierr); 20 if (!flg) SETERRQ1(PetscObjectComm((PetscObject)Rt),PETSC_ERR_SUP,"Not for matrix type %s",((PetscObject)Rt)->type_name); 21 if (scall == MAT_INITIAL_MATRIX) { 22 ierr = PetscLogEventBegin(MAT_MatMatMultSymbolic,R,A,P,0);CHKERRQ(ierr); 23 ierr = MatTransposeMatMatMultSymbolic_AIJ_AIJ_AIJ_wHYPRE(Rt,A,P,fill,RAP);CHKERRQ(ierr); 24 ierr = PetscLogEventEnd(MAT_MatMatMultSymbolic,R,A,P,0);CHKERRQ(ierr); 25 } 26 ierr = PetscLogEventBegin(MAT_MatMatMultNumeric,R,A,P,0);CHKERRQ(ierr); 27 ierr = MatTransposeMatMatMultNumeric_AIJ_AIJ_AIJ_wHYPRE(Rt,A,P,*RAP);CHKERRQ(ierr); 28 ierr = PetscLogEventEnd(MAT_MatMatMultNumeric,R,A,P,0);CHKERRQ(ierr); 29 PetscFunctionReturn(0); 30 } 31 #endif 32 33 PetscErrorCode MatDestroy_MPIAIJ_MatMatMatMult(Mat A) 34 { 35 Mat_MPIAIJ *a = (Mat_MPIAIJ*)A->data; 36 Mat_MatMatMatMult *matmatmatmult=a->matmatmatmult; 37 PetscErrorCode ierr; 38 39 PetscFunctionBegin; 40 ierr = MatDestroy(&matmatmatmult->BC);CHKERRQ(ierr); 41 ierr = matmatmatmult->destroy(A);CHKERRQ(ierr); 42 ierr = PetscFree(matmatmatmult);CHKERRQ(ierr); 43 PetscFunctionReturn(0); 44 } 45 46 PETSC_INTERN PetscErrorCode MatMatMatMult_MPIAIJ_MPIAIJ_MPIAIJ(Mat A,Mat B,Mat C,MatReuse scall,PetscReal fill,Mat *D) 47 { 48 PetscErrorCode ierr; 49 50 PetscFunctionBegin; 51 if (scall == MAT_INITIAL_MATRIX) { 52 ierr = PetscLogEventBegin(MAT_MatMatMultSymbolic,A,B,C,0);CHKERRQ(ierr); 53 ierr = MatMatMatMultSymbolic_MPIAIJ_MPIAIJ_MPIAIJ(A,B,C,fill,D);CHKERRQ(ierr); 54 ierr = PetscLogEventEnd(MAT_MatMatMultSymbolic,A,B,C,0);CHKERRQ(ierr); 55 } 56 ierr = PetscLogEventBegin(MAT_MatMatMultNumeric,A,B,C,0);CHKERRQ(ierr); 57 ierr = ((*D)->ops->matmatmultnumeric)(A,B,C,*D);CHKERRQ(ierr); 58 ierr = PetscLogEventEnd(MAT_MatMatMultNumeric,A,B,C,0);CHKERRQ(ierr); 59 PetscFunctionReturn(0); 60 } 61 62 PetscErrorCode MatMatMatMultSymbolic_MPIAIJ_MPIAIJ_MPIAIJ(Mat A,Mat B,Mat C,PetscReal fill,Mat *D) 63 { 64 PetscErrorCode ierr; 65 Mat BC; 66 Mat_MatMatMatMult *matmatmatmult; 67 Mat_MPIAIJ *d; 68 PetscBool scalable=PETSC_TRUE; 69 70 PetscFunctionBegin; 71 ierr = PetscObjectOptionsBegin((PetscObject)B);CHKERRQ(ierr); 72 ierr = PetscOptionsBool("-matmatmatmult_scalable","Use a scalable but slower D=A*B*C","",scalable,&scalable,NULL);CHKERRQ(ierr); 73 ierr = PetscOptionsEnd();CHKERRQ(ierr); 74 if (scalable) { 75 ierr = MatMatMultSymbolic_MPIAIJ_MPIAIJ(B,C,fill,&BC);CHKERRQ(ierr); 76 ierr = MatMatMultSymbolic_MPIAIJ_MPIAIJ(A,BC,fill,D);CHKERRQ(ierr); 77 } else { 78 ierr = MatMatMultSymbolic_MPIAIJ_MPIAIJ_nonscalable(B,C,fill,&BC);CHKERRQ(ierr); 79 ierr = MatMatMultSymbolic_MPIAIJ_MPIAIJ_nonscalable(A,BC,fill,D);CHKERRQ(ierr); 80 } 81 82 /* create struct Mat_MatMatMatMult and attached it to *D */ 83 ierr = PetscNew(&matmatmatmult);CHKERRQ(ierr); 84 85 matmatmatmult->BC = BC; 86 matmatmatmult->destroy = (*D)->ops->destroy; 87 d = (Mat_MPIAIJ*)(*D)->data; 88 d->matmatmatmult = matmatmatmult; 89 90 (*D)->ops->matmatmultnumeric = MatMatMatMultNumeric_MPIAIJ_MPIAIJ_MPIAIJ; 91 (*D)->ops->destroy = MatDestroy_MPIAIJ_MatMatMatMult; 92 PetscFunctionReturn(0); 93 } 94 95 PetscErrorCode MatMatMatMultNumeric_MPIAIJ_MPIAIJ_MPIAIJ(Mat A,Mat B,Mat C,Mat D) 96 { 97 PetscErrorCode ierr; 98 Mat_MPIAIJ *d = (Mat_MPIAIJ*)D->data; 99 Mat_MatMatMatMult *matmatmatmult = d->matmatmatmult; 100 Mat BC = matmatmatmult->BC; 101 102 PetscFunctionBegin; 103 ierr = (BC->ops->matmultnumeric)(B,C,BC);CHKERRQ(ierr); 104 ierr = (D->ops->matmultnumeric)(A,BC,D);CHKERRQ(ierr); 105 PetscFunctionReturn(0); 106 } 107 108 PetscErrorCode MatRARt_MPIAIJ_MPIAIJ(Mat A,Mat R,MatReuse scall,PetscReal fill,Mat *C) 109 { 110 PetscErrorCode ierr; 111 Mat Rt; 112 113 PetscFunctionBegin; 114 ierr = MatTranspose(R,MAT_INITIAL_MATRIX,&Rt);CHKERRQ(ierr); 115 if (scall == MAT_INITIAL_MATRIX) { 116 ierr = MatMatMatMultSymbolic_MPIAIJ_MPIAIJ_MPIAIJ(R,A,Rt,fill,C);CHKERRQ(ierr); 117 } 118 ierr = MatMatMatMultNumeric_MPIAIJ_MPIAIJ_MPIAIJ(R,A,Rt,*C);CHKERRQ(ierr); 119 ierr = MatDestroy(&Rt);CHKERRQ(ierr); 120 PetscFunctionReturn(0); 121 } 122