xref: /petsc/src/mat/impls/aij/seq/matmatmatmult.c (revision 6d0b6147cb0f49956a99006e1af308da93e18f1c)
171063593SHong Zhang /*
2*6d0b6147SHong 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*/
6*6d0b6147SHong Zhang 
7*6d0b6147SHong Zhang #undef __FUNCT__
8*6d0b6147SHong Zhang #define __FUNCT__ "MatDestroy_SeqAIJ_MatMatMatMult"
9*6d0b6147SHong Zhang PetscErrorCode MatDestroy_SeqAIJ_MatMatMatMult(Mat A)
10*6d0b6147SHong Zhang {
11*6d0b6147SHong Zhang   Mat_SeqAIJ         *a = (Mat_SeqAIJ*)A->data;
12*6d0b6147SHong Zhang   Mat_MatMatMatMult  *matmatmatmult=a->matmatmatmult;
13*6d0b6147SHong Zhang   PetscErrorCode     ierr;
14*6d0b6147SHong Zhang 
15*6d0b6147SHong Zhang   PetscFunctionBegin;
16*6d0b6147SHong Zhang   ierr = MatDestroy(&matmatmatmult->BC);CHKERRQ(ierr);
17*6d0b6147SHong Zhang   ierr = matmatmatmult->destroy(A);CHKERRQ(ierr);
18*6d0b6147SHong Zhang   ierr = PetscFree(matmatmatmult);CHKERRQ(ierr);
19*6d0b6147SHong Zhang   PetscFunctionReturn(0);
20*6d0b6147SHong 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 {
26*6d0b6147SHong Zhang   PetscErrorCode ierr;
2771063593SHong Zhang 
2871063593SHong Zhang   PetscFunctionBegin;
2971063593SHong Zhang   if (scall == MAT_INITIAL_MATRIX){
30*6d0b6147SHong Zhang     ierr = PetscLogEventBegin(MAT_MatMatMultSymbolic,A,B,C,0);CHKERRQ(ierr);
31*6d0b6147SHong Zhang     ierr = MatMatMatMultSymbolic_SeqAIJ_SeqAIJ_SeqAIJ(A,B,C,fill,D);CHKERRQ(ierr);
32*6d0b6147SHong Zhang     ierr = PetscLogEventEnd(MAT_MatMatMultSymbolic,A,B,C,0);CHKERRQ(ierr);
3371063593SHong Zhang   }
34*6d0b6147SHong Zhang   ierr = PetscLogEventBegin(MAT_MatMatMultNumeric,A,B,C,0);CHKERRQ(ierr);
35*6d0b6147SHong Zhang   ierr = MatMatMatMultNumeric_SeqAIJ_SeqAIJ_SeqAIJ(A,B,C,*D);CHKERRQ(ierr);
36*6d0b6147SHong 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;
45*6d0b6147SHong Zhang   Mat                BC;
46*6d0b6147SHong Zhang   Mat_MatMatMatMult  *matmatmatmult;
47*6d0b6147SHong Zhang   Mat_SeqAIJ         *d;
48*6d0b6147SHong Zhang   PetscBool          scalable=PETSC_TRUE;
4971063593SHong Zhang 
5071063593SHong Zhang   PetscFunctionBegin;
51*6d0b6147SHong Zhang   ierr = PetscObjectOptionsBegin((PetscObject)B);CHKERRQ(ierr);
52*6d0b6147SHong Zhang     ierr = PetscOptionsBool("-matmatmatmult_scalable","Use a scalable but slower D=A*B*C","",scalable,&scalable,PETSC_NULL);CHKERRQ(ierr);
53*6d0b6147SHong Zhang   ierr = PetscOptionsEnd();CHKERRQ(ierr);
54*6d0b6147SHong Zhang   if (scalable){
55*6d0b6147SHong Zhang     ierr = MatMatMultSymbolic_SeqAIJ_SeqAIJ_Scalable(B,C,fill,&BC);CHKERRQ(ierr);
56*6d0b6147SHong Zhang     ierr = MatMatMultSymbolic_SeqAIJ_SeqAIJ_Scalable(A,BC,fill,D);CHKERRQ(ierr);
57*6d0b6147SHong Zhang   } else {
58*6d0b6147SHong Zhang     ierr = MatMatMultSymbolic_SeqAIJ_SeqAIJ(B,C,fill,&BC);CHKERRQ(ierr);
59*6d0b6147SHong Zhang     ierr = MatMatMultSymbolic_SeqAIJ_SeqAIJ(A,BC,fill,D);CHKERRQ(ierr);
60*6d0b6147SHong Zhang   }
61*6d0b6147SHong Zhang 
62*6d0b6147SHong Zhang   /* create struct Mat_MatMatMatMult and attached it to *D */
63*6d0b6147SHong Zhang   ierr = PetscNew(Mat_MatMatMatMult,&matmatmatmult);CHKERRQ(ierr);
64*6d0b6147SHong Zhang   matmatmatmult->BC      = BC;
65*6d0b6147SHong Zhang   matmatmatmult->destroy = (*D)->ops->destroy;
66*6d0b6147SHong Zhang   d                      = (Mat_SeqAIJ*)(*D)->data;
67*6d0b6147SHong Zhang   d->matmatmatmult       = matmatmatmult;
68*6d0b6147SHong Zhang 
69*6d0b6147SHong Zhang   (*D)->ops->matmatmultnumeric = MatMatMatMultNumeric_SeqAIJ_SeqAIJ_SeqAIJ;
70*6d0b6147SHong Zhang   (*D)->ops->destroy           = MatDestroy_SeqAIJ_MatMatMatMult;
7171063593SHong Zhang   PetscFunctionReturn(0);
7271063593SHong Zhang }
7371063593SHong Zhang 
7471063593SHong Zhang #undef __FUNCT__
7571063593SHong Zhang #define __FUNCT__ "MatMatMatMultNumeric_SeqAIJ_SeqAIJ_SeqAIJ"
7671063593SHong Zhang PetscErrorCode MatMatMatMultNumeric_SeqAIJ_SeqAIJ_SeqAIJ(Mat A,Mat B,Mat C,Mat D)
7771063593SHong Zhang {
7871063593SHong Zhang   PetscErrorCode    ierr;
79*6d0b6147SHong Zhang   Mat_SeqAIJ        *d=(Mat_SeqAIJ*)D->data;
80*6d0b6147SHong Zhang   Mat_MatMatMatMult *matmatmatmult=d->matmatmatmult;
81*6d0b6147SHong Zhang   Mat               BC= matmatmatmult->BC;
8271063593SHong Zhang 
8371063593SHong Zhang   PetscFunctionBegin;
84*6d0b6147SHong Zhang   ierr = (BC->ops->matmultnumeric)(B,C,BC);CHKERRQ(ierr);
85*6d0b6147SHong Zhang   ierr = (D->ops->matmultnumeric)(A,BC,D);CHKERRQ(ierr);
8671063593SHong Zhang   PetscFunctionReturn(0);
8771063593SHong Zhang }
88