171063593SHong Zhang /* 26d0b6147SHong Zhang Defines matrix-matrix-matrix product routines for SeqAIJ matrices 371063593SHong Zhang D = A * B * C 471063593SHong Zhang */ 571063593SHong Zhang #include <../src/mat/impls/aij/seq/aij.h> /*I "petscmat.h" I*/ 66d0b6147SHong Zhang 76d0b6147SHong Zhang #undef __FUNCT__ 86d0b6147SHong Zhang #define __FUNCT__ "MatDestroy_SeqAIJ_MatMatMatMult" 96d0b6147SHong Zhang PetscErrorCode MatDestroy_SeqAIJ_MatMatMatMult(Mat A) 106d0b6147SHong Zhang { 116d0b6147SHong Zhang Mat_SeqAIJ *a = (Mat_SeqAIJ*)A->data; 126d0b6147SHong Zhang Mat_MatMatMatMult *matmatmatmult=a->matmatmatmult; 136d0b6147SHong Zhang PetscErrorCode ierr; 146d0b6147SHong Zhang 156d0b6147SHong Zhang PetscFunctionBegin; 166d0b6147SHong Zhang ierr = MatDestroy(&matmatmatmult->BC);CHKERRQ(ierr); 176d0b6147SHong Zhang ierr = matmatmatmult->destroy(A);CHKERRQ(ierr); 186d0b6147SHong Zhang ierr = PetscFree(matmatmatmult);CHKERRQ(ierr); 196d0b6147SHong Zhang PetscFunctionReturn(0); 206d0b6147SHong Zhang } 2171063593SHong Zhang 2271063593SHong Zhang #undef __FUNCT__ 2371063593SHong Zhang #define __FUNCT__ "MatMatMatMult_SeqAIJ_SeqAIJ_SeqAIJ" 2471063593SHong Zhang PetscErrorCode MatMatMatMult_SeqAIJ_SeqAIJ_SeqAIJ(Mat A,Mat B,Mat C,MatReuse scall,PetscReal fill,Mat *D) 2571063593SHong Zhang { 266d0b6147SHong Zhang PetscErrorCode ierr; 2771063593SHong Zhang 2871063593SHong Zhang PetscFunctionBegin; 2971063593SHong Zhang if (scall == MAT_INITIAL_MATRIX) { 306d0b6147SHong Zhang ierr = PetscLogEventBegin(MAT_MatMatMultSymbolic,A,B,C,0);CHKERRQ(ierr); 316d0b6147SHong Zhang ierr = MatMatMatMultSymbolic_SeqAIJ_SeqAIJ_SeqAIJ(A,B,C,fill,D);CHKERRQ(ierr); 326d0b6147SHong Zhang ierr = PetscLogEventEnd(MAT_MatMatMultSymbolic,A,B,C,0);CHKERRQ(ierr); 3371063593SHong Zhang } 346d0b6147SHong Zhang ierr = PetscLogEventBegin(MAT_MatMatMultNumeric,A,B,C,0);CHKERRQ(ierr); 356d0b6147SHong Zhang ierr = MatMatMatMultNumeric_SeqAIJ_SeqAIJ_SeqAIJ(A,B,C,*D);CHKERRQ(ierr); 366d0b6147SHong Zhang ierr = PetscLogEventEnd(MAT_MatMatMultNumeric,A,B,C,0);CHKERRQ(ierr); 3771063593SHong Zhang PetscFunctionReturn(0); 3871063593SHong Zhang } 3971063593SHong Zhang 4071063593SHong Zhang #undef __FUNCT__ 4171063593SHong Zhang #define __FUNCT__ "MatMatMatMultSymbolic_SeqAIJ_SeqAIJ_SeqAIJ" 4271063593SHong Zhang PetscErrorCode MatMatMatMultSymbolic_SeqAIJ_SeqAIJ_SeqAIJ(Mat A,Mat B,Mat C,PetscReal fill,Mat *D) 4371063593SHong Zhang { 4471063593SHong Zhang PetscErrorCode ierr; 456d0b6147SHong Zhang Mat BC; 466d0b6147SHong Zhang Mat_MatMatMatMult *matmatmatmult; 476d0b6147SHong Zhang Mat_SeqAIJ *d; 486d0b6147SHong Zhang PetscBool scalable=PETSC_TRUE; 4971063593SHong Zhang 5071063593SHong Zhang PetscFunctionBegin; 516d0b6147SHong Zhang ierr = PetscObjectOptionsBegin((PetscObject)B);CHKERRQ(ierr); 520298fd71SBarry Smith ierr = PetscOptionsBool("-matmatmatmult_scalable","Use a scalable but slower D=A*B*C","",scalable,&scalable,NULL);CHKERRQ(ierr); 536d0b6147SHong Zhang ierr = PetscOptionsEnd();CHKERRQ(ierr); 546d0b6147SHong Zhang if (scalable) { 556d0b6147SHong Zhang ierr = MatMatMultSymbolic_SeqAIJ_SeqAIJ_Scalable(B,C,fill,&BC);CHKERRQ(ierr); 566d0b6147SHong Zhang ierr = MatMatMultSymbolic_SeqAIJ_SeqAIJ_Scalable(A,BC,fill,D);CHKERRQ(ierr); 576d0b6147SHong Zhang } else { 586d0b6147SHong Zhang ierr = MatMatMultSymbolic_SeqAIJ_SeqAIJ(B,C,fill,&BC);CHKERRQ(ierr); 596d0b6147SHong Zhang ierr = MatMatMultSymbolic_SeqAIJ_SeqAIJ(A,BC,fill,D);CHKERRQ(ierr); 606d0b6147SHong Zhang } 616d0b6147SHong Zhang 626d0b6147SHong Zhang /* create struct Mat_MatMatMatMult and attached it to *D */ 63*b00a9115SJed Brown ierr = PetscNew(&matmatmatmult);CHKERRQ(ierr); 642205254eSKarl Rupp 656d0b6147SHong Zhang matmatmatmult->BC = BC; 666d0b6147SHong Zhang matmatmatmult->destroy = (*D)->ops->destroy; 676d0b6147SHong Zhang d = (Mat_SeqAIJ*)(*D)->data; 686d0b6147SHong Zhang d->matmatmatmult = matmatmatmult; 696d0b6147SHong Zhang 706d0b6147SHong Zhang (*D)->ops->matmatmultnumeric = MatMatMatMultNumeric_SeqAIJ_SeqAIJ_SeqAIJ; 716d0b6147SHong Zhang (*D)->ops->destroy = MatDestroy_SeqAIJ_MatMatMatMult; 7271063593SHong Zhang PetscFunctionReturn(0); 7371063593SHong Zhang } 7471063593SHong Zhang 7571063593SHong Zhang #undef __FUNCT__ 7671063593SHong Zhang #define __FUNCT__ "MatMatMatMultNumeric_SeqAIJ_SeqAIJ_SeqAIJ" 7771063593SHong Zhang PetscErrorCode MatMatMatMultNumeric_SeqAIJ_SeqAIJ_SeqAIJ(Mat A,Mat B,Mat C,Mat D) 7871063593SHong Zhang { 7971063593SHong Zhang PetscErrorCode ierr; 806d0b6147SHong Zhang Mat_SeqAIJ *d =(Mat_SeqAIJ*)D->data; 816d0b6147SHong Zhang Mat_MatMatMatMult *matmatmatmult=d->matmatmatmult; 826d0b6147SHong Zhang Mat BC = matmatmatmult->BC; 8371063593SHong Zhang 8471063593SHong Zhang PetscFunctionBegin; 856d0b6147SHong Zhang ierr = (BC->ops->matmultnumeric)(B,C,BC);CHKERRQ(ierr); 866d0b6147SHong Zhang ierr = (D->ops->matmultnumeric)(A,BC,D);CHKERRQ(ierr); 8771063593SHong Zhang PetscFunctionReturn(0); 8871063593SHong Zhang } 89