1*f996eeb8SHong Zhang /* 2*f996eeb8SHong Zhang Defines matrix-matrix-matrix product routines for MPIAIJ matrices 3*f996eeb8SHong Zhang D = A * B * C 4*f996eeb8SHong Zhang */ 5*f996eeb8SHong Zhang #include <../src/mat/impls/aij/mpi/mpiaij.h> /*I "petscmat.h" I*/ 6*f996eeb8SHong Zhang 7*f996eeb8SHong Zhang #undef __FUNCT__ 8*f996eeb8SHong Zhang #define __FUNCT__ "MatDestroy_MPIAIJ_MatMatMatMult" 9*f996eeb8SHong Zhang PetscErrorCode MatDestroy_MPIAIJ_MatMatMatMult(Mat A) 10*f996eeb8SHong Zhang { 11*f996eeb8SHong Zhang Mat_MPIAIJ *a = (Mat_MPIAIJ*)A->data; 12*f996eeb8SHong Zhang Mat_MatMatMatMult *matmatmatmult=a->matmatmatmult; 13*f996eeb8SHong Zhang PetscErrorCode ierr; 14*f996eeb8SHong Zhang 15*f996eeb8SHong Zhang PetscFunctionBegin; 16*f996eeb8SHong Zhang ierr = MatDestroy(&matmatmatmult->BC);CHKERRQ(ierr); 17*f996eeb8SHong Zhang ierr = matmatmatmult->destroy(A);CHKERRQ(ierr); 18*f996eeb8SHong Zhang ierr = PetscFree(matmatmatmult);CHKERRQ(ierr); 19*f996eeb8SHong Zhang PetscFunctionReturn(0); 20*f996eeb8SHong Zhang } 21*f996eeb8SHong Zhang 22*f996eeb8SHong Zhang #undef __FUNCT__ 23*f996eeb8SHong Zhang #define __FUNCT__ "MatMatMatMult_MPIAIJ_MPIAIJ_MPIAIJ" 24*f996eeb8SHong Zhang PetscErrorCode MatMatMatMult_MPIAIJ_MPIAIJ_MPIAIJ(Mat A,Mat B,Mat C,MatReuse scall,PetscReal fill,Mat *D) 25*f996eeb8SHong Zhang { 26*f996eeb8SHong Zhang PetscErrorCode ierr; 27*f996eeb8SHong Zhang 28*f996eeb8SHong Zhang PetscFunctionBegin; 29*f996eeb8SHong Zhang if (scall == MAT_INITIAL_MATRIX){ 30*f996eeb8SHong Zhang ierr = PetscLogEventBegin(MAT_MatMatMultSymbolic,A,B,C,0);CHKERRQ(ierr); 31*f996eeb8SHong Zhang ierr = MatMatMatMultSymbolic_MPIAIJ_MPIAIJ_MPIAIJ(A,B,C,fill,D);CHKERRQ(ierr); 32*f996eeb8SHong Zhang ierr = PetscLogEventEnd(MAT_MatMatMultSymbolic,A,B,C,0);CHKERRQ(ierr); 33*f996eeb8SHong Zhang } 34*f996eeb8SHong Zhang ierr = PetscLogEventBegin(MAT_MatMatMultNumeric,A,B,C,0);CHKERRQ(ierr); 35*f996eeb8SHong Zhang ierr = MatMatMatMultNumeric_MPIAIJ_MPIAIJ_MPIAIJ(A,B,C,*D);CHKERRQ(ierr); 36*f996eeb8SHong Zhang ierr = PetscLogEventEnd(MAT_MatMatMultNumeric,A,B,C,0);CHKERRQ(ierr); 37*f996eeb8SHong Zhang PetscFunctionReturn(0); 38*f996eeb8SHong Zhang } 39*f996eeb8SHong Zhang 40*f996eeb8SHong Zhang #undef __FUNCT__ 41*f996eeb8SHong Zhang #define __FUNCT__ "MatMatMatMultSymbolic_MPIAIJ_MPIAIJ_MPIAIJ" 42*f996eeb8SHong Zhang PetscErrorCode MatMatMatMultSymbolic_MPIAIJ_MPIAIJ_MPIAIJ(Mat A,Mat B,Mat C,PetscReal fill,Mat *D) 43*f996eeb8SHong Zhang { 44*f996eeb8SHong Zhang PetscErrorCode ierr; 45*f996eeb8SHong Zhang Mat BC; 46*f996eeb8SHong Zhang Mat_MatMatMatMult *matmatmatmult; 47*f996eeb8SHong Zhang Mat_MPIAIJ *d; 48*f996eeb8SHong Zhang PetscBool scalable=PETSC_TRUE; 49*f996eeb8SHong Zhang 50*f996eeb8SHong Zhang PetscFunctionBegin; 51*f996eeb8SHong Zhang ierr = PetscObjectOptionsBegin((PetscObject)B);CHKERRQ(ierr); 52*f996eeb8SHong Zhang ierr = PetscOptionsBool("-matmatmatmult_scalable","Use a scalable but slower D=A*B*C","",scalable,&scalable,PETSC_NULL);CHKERRQ(ierr); 53*f996eeb8SHong Zhang ierr = PetscOptionsEnd();CHKERRQ(ierr); 54*f996eeb8SHong Zhang if (scalable){ 55*f996eeb8SHong Zhang ierr = MatMatMultSymbolic_MPIAIJ_MPIAIJ_Scalable(B,C,fill,&BC);CHKERRQ(ierr); 56*f996eeb8SHong Zhang ierr = MatMatMultSymbolic_MPIAIJ_MPIAIJ_Scalable(A,BC,fill,D);CHKERRQ(ierr); 57*f996eeb8SHong Zhang } else { 58*f996eeb8SHong Zhang ierr = MatMatMultSymbolic_MPIAIJ_MPIAIJ(B,C,fill,&BC);CHKERRQ(ierr); 59*f996eeb8SHong Zhang ierr = MatMatMultSymbolic_MPIAIJ_MPIAIJ(A,BC,fill,D);CHKERRQ(ierr); 60*f996eeb8SHong Zhang } 61*f996eeb8SHong Zhang 62*f996eeb8SHong Zhang /* create struct Mat_MatMatMatMult and attached it to *D */ 63*f996eeb8SHong Zhang ierr = PetscNew(Mat_MatMatMatMult,&matmatmatmult);CHKERRQ(ierr); 64*f996eeb8SHong Zhang matmatmatmult->BC = BC; 65*f996eeb8SHong Zhang matmatmatmult->destroy = (*D)->ops->destroy; 66*f996eeb8SHong Zhang d = (Mat_MPIAIJ*)(*D)->data; 67*f996eeb8SHong Zhang d->matmatmatmult = matmatmatmult; 68*f996eeb8SHong Zhang 69*f996eeb8SHong Zhang (*D)->ops->matmatmultnumeric = MatMatMatMultNumeric_MPIAIJ_MPIAIJ_MPIAIJ; 70*f996eeb8SHong Zhang (*D)->ops->destroy = MatDestroy_MPIAIJ_MatMatMatMult; 71*f996eeb8SHong Zhang PetscFunctionReturn(0); 72*f996eeb8SHong Zhang } 73*f996eeb8SHong Zhang 74*f996eeb8SHong Zhang #undef __FUNCT__ 75*f996eeb8SHong Zhang #define __FUNCT__ "MatMatMatMultNumeric_MPIAIJ_MPIAIJ_MPIAIJ" 76*f996eeb8SHong Zhang PetscErrorCode MatMatMatMultNumeric_MPIAIJ_MPIAIJ_MPIAIJ(Mat A,Mat B,Mat C,Mat D) 77*f996eeb8SHong Zhang { 78*f996eeb8SHong Zhang PetscErrorCode ierr; 79*f996eeb8SHong Zhang Mat_MPIAIJ *d=(Mat_MPIAIJ*)D->data; 80*f996eeb8SHong Zhang Mat_MatMatMatMult *matmatmatmult=d->matmatmatmult; 81*f996eeb8SHong Zhang Mat BC= matmatmatmult->BC; 82*f996eeb8SHong Zhang 83*f996eeb8SHong Zhang PetscFunctionBegin; 84*f996eeb8SHong Zhang ierr = (BC->ops->matmultnumeric)(B,C,BC);CHKERRQ(ierr); 85*f996eeb8SHong Zhang ierr = (D->ops->matmultnumeric)(A,BC,D);CHKERRQ(ierr); 86*f996eeb8SHong Zhang PetscFunctionReturn(0); 87*f996eeb8SHong Zhang } 88