xref: /petsc/src/mat/tests/ex62.c (revision 28b400f66ebc7ae0049166a2294dfcd3df27e64b)
1c20d7725SJed Brown 
2c20d7725SJed Brown static char help[] = "Test Matrix products for AIJ matrices\n\
3c20d7725SJed Brown Input arguments are:\n\
4c20d7725SJed Brown   -fA <input_file> -fB <input_file> -fC <input_file>: file to load\n\n";
5c20d7725SJed Brown /* Example of usage:
6c20d7725SJed Brown    ./ex62 -fA <A_binary> -fB <B_binary>
7c20d7725SJed Brown    mpiexec -n 3 ./ex62 -fA medium -fB medium
8c20d7725SJed Brown */
9c20d7725SJed Brown 
10c20d7725SJed Brown #include <petscmat.h>
11c20d7725SJed Brown 
12c20d7725SJed Brown /*
13c20d7725SJed Brown      B = A - B
14c20d7725SJed Brown      norm = norm(B)
15c20d7725SJed Brown */
16c20d7725SJed Brown PetscErrorCode MatNormDifference(Mat A,Mat B,PetscReal *norm)
17c20d7725SJed Brown {
18c20d7725SJed Brown   PetscFunctionBegin;
195f80ce2aSJacob Faibussowitsch   CHKERRQ(MatAXPY(B,-1.0,A,DIFFERENT_NONZERO_PATTERN));
205f80ce2aSJacob Faibussowitsch   CHKERRQ(MatNorm(B,NORM_FROBENIUS,norm));
21c20d7725SJed Brown   PetscFunctionReturn(0);
22c20d7725SJed Brown }
23c20d7725SJed Brown 
24c20d7725SJed Brown int main(int argc,char **args)
25c20d7725SJed Brown {
26c20d7725SJed Brown   Mat            A,A_save,B,C,P,C1,R;
27c20d7725SJed Brown   PetscViewer    viewer;
28c20d7725SJed Brown   PetscErrorCode ierr;
29c20d7725SJed Brown   PetscMPIInt    size,rank;
3020b3374bSStefano Zampini   PetscInt       i,j,*idxn,PM,PN = PETSC_DECIDE,rstart,rend;
31c20d7725SJed Brown   PetscReal      norm;
32c20d7725SJed Brown   PetscRandom    rdm;
3320b3374bSStefano Zampini   char           file[2][PETSC_MAX_PATH_LEN] = { "", ""};
34c20d7725SJed Brown   PetscScalar    *a,rval,alpha;
35c20d7725SJed Brown   PetscBool      Test_MatMatMult=PETSC_TRUE,Test_MatTrMat=PETSC_TRUE,Test_MatMatTr=PETSC_TRUE;
3620b3374bSStefano Zampini   PetscBool      Test_MatPtAP=PETSC_TRUE,Test_MatRARt=PETSC_TRUE,flg,seqaij,flgA,flgB;
37c20d7725SJed Brown   MatInfo        info;
3820b3374bSStefano Zampini   PetscInt       nzp  = 5; /* num of nonzeros in each row of P */
39c20d7725SJed Brown   MatType        mattype;
4020b3374bSStefano Zampini   const char     *deft = MATAIJ;
4120b3374bSStefano Zampini   char           A_mattype[256], B_mattype[256];
4220b3374bSStefano Zampini   PetscInt       mcheck = 10;
43c20d7725SJed Brown 
44c20d7725SJed Brown   ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
455f80ce2aSJacob Faibussowitsch   CHKERRMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size));
465f80ce2aSJacob Faibussowitsch   CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank));
47c20d7725SJed Brown 
48c20d7725SJed Brown   /*  Load the matrices A_save and B */
4920b3374bSStefano Zampini   ierr = PetscOptionsBegin(PETSC_COMM_WORLD,"","","");CHKERRQ(ierr);
505f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsBool("-test_rart","Test MatRARt","",Test_MatRARt,&Test_MatRARt,NULL));
515f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsInt("-PN","Number of columns of P","",PN,&PN,NULL));
525f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsInt("-mcheck","Number of matmult checks","",mcheck,&mcheck,NULL));
535f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsString("-fA","Path for matrix A","",file[0],file[0],sizeof(file[0]),&flg));
54*28b400f6SJacob Faibussowitsch   PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Must indicate a file name for matrix A with the -fA option.");
555f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsString("-fB","Path for matrix B","",file[1],file[1],sizeof(file[1]),&flg));
565f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsFList("-A_mat_type","Matrix type","MatSetType",MatList,deft,A_mattype,256,&flgA));
575f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscOptionsFList("-B_mat_type","Matrix type","MatSetType",MatList,deft,B_mattype,256,&flgB));
5820b3374bSStefano Zampini   ierr = PetscOptionsEnd();CHKERRQ(ierr);
59c20d7725SJed Brown 
605f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[0],FILE_MODE_READ,&viewer));
615f80ce2aSJacob Faibussowitsch   CHKERRQ(MatCreate(PETSC_COMM_WORLD,&A_save));
625f80ce2aSJacob Faibussowitsch   CHKERRQ(MatLoad(A_save,viewer));
635f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscViewerDestroy(&viewer));
64c20d7725SJed Brown 
6520b3374bSStefano Zampini   if (flg) {
665f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[1],FILE_MODE_READ,&viewer));
675f80ce2aSJacob Faibussowitsch     CHKERRQ(MatCreate(PETSC_COMM_WORLD,&B));
685f80ce2aSJacob Faibussowitsch     CHKERRQ(MatLoad(B,viewer));
695f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscViewerDestroy(&viewer));
7020b3374bSStefano Zampini   } else {
715f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscObjectReference((PetscObject)A_save));
7220b3374bSStefano Zampini     B = A_save;
7320b3374bSStefano Zampini   }
7420b3374bSStefano Zampini 
7520b3374bSStefano Zampini   if (flgA) {
765f80ce2aSJacob Faibussowitsch     CHKERRQ(MatConvert(A_save,A_mattype,MAT_INPLACE_MATRIX,&A_save));
7720b3374bSStefano Zampini   }
7820b3374bSStefano Zampini   if (flgB) {
795f80ce2aSJacob Faibussowitsch     CHKERRQ(MatConvert(B,B_mattype,MAT_INPLACE_MATRIX,&B));
8020b3374bSStefano Zampini   }
815f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetFromOptions(A_save));
825f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetFromOptions(B));
83c20d7725SJed Brown 
845f80ce2aSJacob Faibussowitsch   CHKERRQ(MatGetType(B,&mattype));
85c20d7725SJed Brown 
865f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscMalloc(nzp*(sizeof(PetscInt)+sizeof(PetscScalar)),&idxn));
87c20d7725SJed Brown   a    = (PetscScalar*)(idxn + nzp);
88c20d7725SJed Brown 
895f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscRandomCreate(PETSC_COMM_WORLD,&rdm));
905f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscRandomSetFromOptions(rdm));
91c20d7725SJed Brown 
92c20d7725SJed Brown   /* 1) MatMatMult() */
93c20d7725SJed Brown   /* ----------------*/
94c20d7725SJed Brown   if (Test_MatMatMult) {
955f80ce2aSJacob Faibussowitsch     CHKERRQ(MatDuplicate(A_save,MAT_COPY_VALUES,&A));
96c20d7725SJed Brown 
97c20d7725SJed Brown     /* (1.1) Test developer API */
985f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductCreate(A,B,NULL,&C));
995f80ce2aSJacob Faibussowitsch     CHKERRQ(MatSetOptionsPrefix(C,"AB_"));
1005f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetType(C,MATPRODUCT_AB));
1015f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetAlgorithm(C,MATPRODUCTALGORITHMDEFAULT));
1025f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetFill(C,PETSC_DEFAULT));
1035f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetFromOptions(C));
104910fa13fSStefano Zampini     /* we can inquire about MATOP_PRODUCTSYMBOLIC even if the destination matrix type has not been set yet */
1055f80ce2aSJacob Faibussowitsch     CHKERRQ(MatHasOperation(C,MATOP_PRODUCTSYMBOLIC,&flg));
1065f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSymbolic(C));
1075f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductNumeric(C));
1085f80ce2aSJacob Faibussowitsch     CHKERRQ(MatMatMultEqual(A,B,C,mcheck,&flg));
109*28b400f6SJacob Faibussowitsch     PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in C=A*B");
110c20d7725SJed Brown 
111c20d7725SJed Brown     /* Test reuse symbolic C */
112c20d7725SJed Brown     alpha = 0.9;
1135f80ce2aSJacob Faibussowitsch     CHKERRQ(MatScale(A,alpha));
1145f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductNumeric(C));
115c20d7725SJed Brown 
1165f80ce2aSJacob Faibussowitsch     CHKERRQ(MatMatMultEqual(A,B,C,mcheck,&flg));
117*28b400f6SJacob Faibussowitsch     PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in C=A*B");
1185f80ce2aSJacob Faibussowitsch     CHKERRQ(MatDestroy(&C));
119c20d7725SJed Brown 
120c20d7725SJed Brown     /* (1.2) Test user driver */
1215f80ce2aSJacob Faibussowitsch     CHKERRQ(MatMatMult(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C));
122c20d7725SJed Brown 
123c20d7725SJed Brown     /* Test MAT_REUSE_MATRIX - reuse symbolic C */
124c20d7725SJed Brown     alpha = 1.0;
125c20d7725SJed Brown     for (i=0; i<2; i++) {
126c20d7725SJed Brown       alpha -= 0.1;
1275f80ce2aSJacob Faibussowitsch       CHKERRQ(MatScale(A,alpha));
1285f80ce2aSJacob Faibussowitsch       CHKERRQ(MatMatMult(A,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C));
129c20d7725SJed Brown     }
1305f80ce2aSJacob Faibussowitsch     CHKERRQ(MatMatMultEqual(A,B,C,mcheck,&flg));
131*28b400f6SJacob Faibussowitsch     PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error: MatMatMult()");
1325f80ce2aSJacob Faibussowitsch     CHKERRQ(MatDestroy(&A));
1334417c5e8SHong Zhang 
1344417c5e8SHong Zhang     /* Test MatProductClear() */
1355f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductClear(C));
1365f80ce2aSJacob Faibussowitsch     CHKERRQ(MatDestroy(&C));
137544a5e07SHong Zhang 
138544a5e07SHong Zhang     /* Test MatMatMult() for dense and aij matrices */
1395f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscObjectTypeCompareAny((PetscObject)A,&flg,MATSEQAIJ,MATMPIAIJ,""));
14020b3374bSStefano Zampini     if (flg) {
1415f80ce2aSJacob Faibussowitsch       CHKERRQ(MatConvert(A_save,MATDENSE,MAT_INITIAL_MATRIX,&A));
1425f80ce2aSJacob Faibussowitsch       CHKERRQ(MatMatMult(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C));
1435f80ce2aSJacob Faibussowitsch       CHKERRQ(MatDestroy(&C));
1445f80ce2aSJacob Faibussowitsch       CHKERRQ(MatDestroy(&A));
14520b3374bSStefano Zampini     }
146c20d7725SJed Brown   }
147c20d7725SJed Brown 
148c20d7725SJed Brown   /* Create P and R = P^T  */
149c20d7725SJed Brown   /* --------------------- */
1505f80ce2aSJacob Faibussowitsch   CHKERRQ(MatGetSize(B,&PM,NULL));
15120b3374bSStefano Zampini   if (PN < 0) PN = PM/2;
1525f80ce2aSJacob Faibussowitsch   CHKERRQ(MatCreate(PETSC_COMM_WORLD,&P));
1535f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetSizes(P,PETSC_DECIDE,PETSC_DECIDE,PM,PN));
1545f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetType(P,MATAIJ));
1555f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSeqAIJSetPreallocation(P,nzp,NULL));
1565f80ce2aSJacob Faibussowitsch   CHKERRQ(MatMPIAIJSetPreallocation(P,nzp,NULL,nzp,NULL));
1575f80ce2aSJacob Faibussowitsch   CHKERRQ(MatGetOwnershipRange(P,&rstart,&rend));
158c20d7725SJed Brown   for (i=0; i<nzp; i++) {
1595f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscRandomGetValue(rdm,&a[i]));
160c20d7725SJed Brown   }
161c20d7725SJed Brown   for (i=rstart; i<rend; i++) {
162c20d7725SJed Brown     for (j=0; j<nzp; j++) {
1635f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscRandomGetValue(rdm,&rval));
164c20d7725SJed Brown       idxn[j] = (PetscInt)(PetscRealPart(rval)*PN);
165c20d7725SJed Brown     }
1665f80ce2aSJacob Faibussowitsch     CHKERRQ(MatSetValues(P,1,&i,nzp,idxn,a,ADD_VALUES));
167c20d7725SJed Brown   }
1685f80ce2aSJacob Faibussowitsch   CHKERRQ(MatAssemblyBegin(P,MAT_FINAL_ASSEMBLY));
1695f80ce2aSJacob Faibussowitsch   CHKERRQ(MatAssemblyEnd(P,MAT_FINAL_ASSEMBLY));
170c20d7725SJed Brown 
1715f80ce2aSJacob Faibussowitsch   CHKERRQ(MatTranspose(P,MAT_INITIAL_MATRIX,&R));
1725f80ce2aSJacob Faibussowitsch   CHKERRQ(MatConvert(P,mattype,MAT_INPLACE_MATRIX,&P));
1735f80ce2aSJacob Faibussowitsch   CHKERRQ(MatConvert(R,mattype,MAT_INPLACE_MATRIX,&R));
1745f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetFromOptions(P));
1755f80ce2aSJacob Faibussowitsch   CHKERRQ(MatSetFromOptions(R));
176c20d7725SJed Brown 
177c20d7725SJed Brown   /* 2) MatTransposeMatMult() */
178c20d7725SJed Brown   /* ------------------------ */
179c20d7725SJed Brown   if (Test_MatTrMat) {
180c20d7725SJed Brown     /* (2.1) Test developer driver C = P^T*B */
1815f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductCreate(P,B,NULL,&C));
1825f80ce2aSJacob Faibussowitsch     CHKERRQ(MatSetOptionsPrefix(C,"AtB_"));
1835f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetType(C,MATPRODUCT_AtB));
1845f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetAlgorithm(C,MATPRODUCTALGORITHMDEFAULT));
1855f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetFill(C,PETSC_DEFAULT));
1865f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetFromOptions(C));
1875f80ce2aSJacob Faibussowitsch     CHKERRQ(MatHasOperation(C,MATOP_PRODUCTSYMBOLIC,&flg));
188263f2b91SStefano Zampini     if (flg) { /* run tests if supported */
1895f80ce2aSJacob Faibussowitsch       CHKERRQ(MatProductSymbolic(C)); /* equivalent to MatSetUp() */
1905f80ce2aSJacob Faibussowitsch       CHKERRQ(MatSetOption(C,MAT_USE_INODES,PETSC_FALSE)); /* illustrate how to call MatSetOption() */
1915f80ce2aSJacob Faibussowitsch       CHKERRQ(MatProductNumeric(C));
1925f80ce2aSJacob Faibussowitsch       CHKERRQ(MatProductNumeric(C)); /* test reuse symbolic C */
19367b3012eSStefano Zampini 
1945f80ce2aSJacob Faibussowitsch       CHKERRQ(MatTransposeMatMultEqual(P,B,C,mcheck,&flg));
195*28b400f6SJacob Faibussowitsch       PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error: developer driver C = P^T*B");
1965f80ce2aSJacob Faibussowitsch       CHKERRQ(MatDestroy(&C));
197c20d7725SJed Brown 
198c20d7725SJed Brown       /* (2.2) Test user driver C = P^T*B */
1995f80ce2aSJacob Faibussowitsch       CHKERRQ(MatTransposeMatMult(P,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C));
2005f80ce2aSJacob Faibussowitsch       CHKERRQ(MatTransposeMatMult(P,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C));
2015f80ce2aSJacob Faibussowitsch       CHKERRQ(MatGetInfo(C,MAT_GLOBAL_SUM,&info));
2025f80ce2aSJacob Faibussowitsch       CHKERRQ(MatProductClear(C));
203c20d7725SJed Brown 
204c20d7725SJed Brown       /* Compare P^T*B and R*B */
2055f80ce2aSJacob Faibussowitsch       CHKERRQ(MatMatMult(R,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C1));
2065f80ce2aSJacob Faibussowitsch       CHKERRQ(MatNormDifference(C,C1,&norm));
2072c71b3e2SJacob Faibussowitsch       PetscCheckFalse(norm > PETSC_SMALL,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatTransposeMatMult(): %g",(double)norm);
2085f80ce2aSJacob Faibussowitsch       CHKERRQ(MatDestroy(&C1));
209c20d7725SJed Brown 
210c20d7725SJed Brown       /* Test MatDuplicate() of C=P^T*B */
2115f80ce2aSJacob Faibussowitsch       CHKERRQ(MatDuplicate(C,MAT_COPY_VALUES,&C1));
2125f80ce2aSJacob Faibussowitsch       CHKERRQ(MatDestroy(&C1));
213263f2b91SStefano Zampini     } else {
2145f80ce2aSJacob Faibussowitsch       CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"MatTransposeMatMult not supported\n"));
215263f2b91SStefano Zampini     }
2165f80ce2aSJacob Faibussowitsch     CHKERRQ(MatDestroy(&C));
217c20d7725SJed Brown   }
218c20d7725SJed Brown 
21967b3012eSStefano Zampini   /* 3) MatMatTransposeMult() */
220c20d7725SJed Brown   /* ------------------------ */
221c20d7725SJed Brown   if (Test_MatMatTr) {
222c20d7725SJed Brown     /* C = B*R^T */
2235f80ce2aSJacob Faibussowitsch     CHKERRQ(PetscObjectBaseTypeCompare((PetscObject)B,MATSEQAIJ,&seqaij));
22420b3374bSStefano Zampini     if (seqaij) {
2255f80ce2aSJacob Faibussowitsch       CHKERRQ(MatMatTransposeMult(B,R,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C));
2265f80ce2aSJacob Faibussowitsch       CHKERRQ(MatSetOptionsPrefix(C,"ABt_")); /* enable '-ABt_' for matrix C */
2275f80ce2aSJacob Faibussowitsch       CHKERRQ(MatGetInfo(C,MAT_GLOBAL_SUM,&info));
228c20d7725SJed Brown 
229c20d7725SJed Brown       /* Test MAT_REUSE_MATRIX - reuse symbolic C */
2305f80ce2aSJacob Faibussowitsch       CHKERRQ(MatMatTransposeMult(B,R,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C));
231c20d7725SJed Brown 
232c20d7725SJed Brown       /* Check */
2335f80ce2aSJacob Faibussowitsch       CHKERRQ(MatMatMult(B,P,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C1));
2345f80ce2aSJacob Faibussowitsch       CHKERRQ(MatNormDifference(C,C1,&norm));
2352c71b3e2SJacob Faibussowitsch       PetscCheckFalse(norm > PETSC_SMALL,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatMatTransposeMult() %g",(double)norm);
2365f80ce2aSJacob Faibussowitsch       CHKERRQ(MatDestroy(&C1));
2375f80ce2aSJacob Faibussowitsch       CHKERRQ(MatDestroy(&C));
238c20d7725SJed Brown     }
239c20d7725SJed Brown   }
240c20d7725SJed Brown 
241c20d7725SJed Brown   /* 4) Test MatPtAP() */
242c20d7725SJed Brown   /*-------------------*/
243c20d7725SJed Brown   if (Test_MatPtAP) {
2445f80ce2aSJacob Faibussowitsch     CHKERRQ(MatDuplicate(A_save,MAT_COPY_VALUES,&A));
245c20d7725SJed Brown 
246c20d7725SJed Brown     /* (4.1) Test developer API */
2475f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductCreate(A,P,NULL,&C));
2485f80ce2aSJacob Faibussowitsch     CHKERRQ(MatSetOptionsPrefix(C,"PtAP_"));
2495f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetType(C,MATPRODUCT_PtAP));
2505f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetAlgorithm(C,MATPRODUCTALGORITHMDEFAULT));
2515f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetFill(C,PETSC_DEFAULT));
2525f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSetFromOptions(C));
2535f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductSymbolic(C));
2545f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductNumeric(C));
2555f80ce2aSJacob Faibussowitsch     CHKERRQ(MatPtAPMultEqual(A,P,C,mcheck,&flg));
256*28b400f6SJacob Faibussowitsch     PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatProduct_PtAP");
2575f80ce2aSJacob Faibussowitsch     CHKERRQ(MatProductNumeric(C)); /* reuse symbolic C */
258c20d7725SJed Brown 
2595f80ce2aSJacob Faibussowitsch     CHKERRQ(MatPtAPMultEqual(A,P,C,mcheck,&flg));
260*28b400f6SJacob Faibussowitsch     PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatProduct_PtAP");
2615f80ce2aSJacob Faibussowitsch     CHKERRQ(MatDestroy(&C));
262c20d7725SJed Brown 
263c20d7725SJed Brown     /* (4.2) Test user driver */
2645f80ce2aSJacob Faibussowitsch     CHKERRQ(MatPtAP(A,P,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C));
265c20d7725SJed Brown 
266c20d7725SJed Brown     /* Test MAT_REUSE_MATRIX - reuse symbolic C */
267c20d7725SJed Brown     alpha = 1.0;
268c20d7725SJed Brown     for (i=0; i<2; i++) {
269c20d7725SJed Brown       alpha -= 0.1;
2705f80ce2aSJacob Faibussowitsch       CHKERRQ(MatScale(A,alpha));
2715f80ce2aSJacob Faibussowitsch       CHKERRQ(MatPtAP(A,P,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C));
272c20d7725SJed Brown     }
2735f80ce2aSJacob Faibussowitsch     CHKERRQ(MatPtAPMultEqual(A,P,C,mcheck,&flg));
274*28b400f6SJacob Faibussowitsch     PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatPtAP");
275c20d7725SJed Brown 
276c20d7725SJed Brown     /* 5) Test MatRARt() */
277c20d7725SJed Brown     /* ----------------- */
278c20d7725SJed Brown     if (Test_MatRARt) {
279c20d7725SJed Brown       Mat RARt;
28067b3012eSStefano Zampini 
28167b3012eSStefano Zampini       /* (5.1) Test developer driver RARt = R*A*Rt */
2825f80ce2aSJacob Faibussowitsch       CHKERRQ(MatProductCreate(A,R,NULL,&RARt));
2835f80ce2aSJacob Faibussowitsch       CHKERRQ(MatSetOptionsPrefix(RARt,"RARt_"));
2845f80ce2aSJacob Faibussowitsch       CHKERRQ(MatProductSetType(RARt,MATPRODUCT_RARt));
2855f80ce2aSJacob Faibussowitsch       CHKERRQ(MatProductSetAlgorithm(RARt,MATPRODUCTALGORITHMDEFAULT));
2865f80ce2aSJacob Faibussowitsch       CHKERRQ(MatProductSetFill(RARt,PETSC_DEFAULT));
2875f80ce2aSJacob Faibussowitsch       CHKERRQ(MatProductSetFromOptions(RARt));
2885f80ce2aSJacob Faibussowitsch       CHKERRQ(MatHasOperation(RARt,MATOP_PRODUCTSYMBOLIC,&flg));
289263f2b91SStefano Zampini       if (flg) {
2905f80ce2aSJacob Faibussowitsch         CHKERRQ(MatProductSymbolic(RARt)); /* equivalent to MatSetUp() */
2915f80ce2aSJacob Faibussowitsch         CHKERRQ(MatSetOption(RARt,MAT_USE_INODES,PETSC_FALSE)); /* illustrate how to call MatSetOption() */
2925f80ce2aSJacob Faibussowitsch         CHKERRQ(MatProductNumeric(RARt));
2935f80ce2aSJacob Faibussowitsch         CHKERRQ(MatProductNumeric(RARt)); /* test reuse symbolic RARt */
2945f80ce2aSJacob Faibussowitsch         CHKERRQ(MatDestroy(&RARt));
29567b3012eSStefano Zampini 
29667b3012eSStefano Zampini         /* (2.2) Test user driver RARt = R*A*Rt */
2975f80ce2aSJacob Faibussowitsch         CHKERRQ(MatRARt(A,R,MAT_INITIAL_MATRIX,2.0,&RARt));
2985f80ce2aSJacob Faibussowitsch         CHKERRQ(MatRARt(A,R,MAT_REUSE_MATRIX,2.0,&RARt));
29967b3012eSStefano Zampini 
3005f80ce2aSJacob Faibussowitsch         CHKERRQ(MatNormDifference(C,RARt,&norm));
3012c71b3e2SJacob Faibussowitsch         PetscCheckFalse(norm > PETSC_SMALL,PETSC_COMM_SELF,PETSC_ERR_PLIB,"|PtAP - RARt| = %g",(double)norm);
302263f2b91SStefano Zampini       } else {
3035f80ce2aSJacob Faibussowitsch         CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"MatRARt not supported\n"));
304263f2b91SStefano Zampini       }
3055f80ce2aSJacob Faibussowitsch       CHKERRQ(MatDestroy(&RARt));
306c20d7725SJed Brown     }
307c20d7725SJed Brown 
3085f80ce2aSJacob Faibussowitsch     CHKERRQ(MatDestroy(&A));
3095f80ce2aSJacob Faibussowitsch     CHKERRQ(MatDestroy(&C));
310c20d7725SJed Brown   }
311c20d7725SJed Brown 
312c20d7725SJed Brown   /* Destroy objects */
3135f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscRandomDestroy(&rdm));
3145f80ce2aSJacob Faibussowitsch   CHKERRQ(PetscFree(idxn));
315c20d7725SJed Brown 
3165f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&A_save));
3175f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&B));
3185f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&P));
3195f80ce2aSJacob Faibussowitsch   CHKERRQ(MatDestroy(&R));
320c20d7725SJed Brown 
321263f2b91SStefano Zampini   ierr = PetscFinalize();
322c20d7725SJed Brown   return ierr;
323c20d7725SJed Brown }
324c20d7725SJed Brown 
325c20d7725SJed Brown /*TEST
326c20d7725SJed Brown    test:
327c20d7725SJed Brown      suffix: 1
328dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
329c20d7725SJed Brown      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium
330c20d7725SJed Brown      output_file: output/ex62_1.out
331c20d7725SJed Brown 
332c20d7725SJed Brown    test:
333c20d7725SJed Brown      suffix: 2_ab_scalable
334dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
3353e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm scalable -matmatmult_via scalable -AtB_mat_product_algorithm outerproduct -mattransposematmult_via outerproduct
336c20d7725SJed Brown      output_file: output/ex62_1.out
337c20d7725SJed Brown 
338c20d7725SJed Brown    test:
339c20d7725SJed Brown      suffix: 3_ab_scalable_fast
340dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
3413e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm scalable_fast -matmatmult_via scalable_fast -matmattransmult_via color
342c20d7725SJed Brown      output_file: output/ex62_1.out
343c20d7725SJed Brown 
344c20d7725SJed Brown    test:
345c20d7725SJed Brown      suffix: 4_ab_heap
346dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
3473e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm heap -matmatmult_via heap -PtAP_mat_product_algorithm rap -matptap_via rap
348c20d7725SJed Brown      output_file: output/ex62_1.out
349c20d7725SJed Brown 
350c20d7725SJed Brown    test:
351c20d7725SJed Brown      suffix: 5_ab_btheap
352dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
3533e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm btheap -matmatmult_via btheap -matrart_via r*art
354c20d7725SJed Brown      output_file: output/ex62_1.out
355c20d7725SJed Brown 
356c20d7725SJed Brown    test:
357c20d7725SJed Brown      suffix: 6_ab_llcondensed
358dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
3593e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm llcondensed -matmatmult_via llcondensed -matrart_via coloring_rart
360c20d7725SJed Brown      output_file: output/ex62_1.out
361c20d7725SJed Brown 
362c20d7725SJed Brown    test:
363c20d7725SJed Brown      suffix: 7_ab_rowmerge
364dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
3653e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm rowmerge -matmatmult_via rowmerge
366c20d7725SJed Brown      output_file: output/ex62_1.out
367c20d7725SJed Brown 
368c20d7725SJed Brown    test:
369c20d7725SJed Brown      suffix: 8_ab_hypre
370dfd57a17SPierre Jolivet      requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
3713e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm hypre -matmatmult_via hypre -PtAP_mat_product_algorithm hypre -matptap_via hypre
372c20d7725SJed Brown      output_file: output/ex62_1.out
373c20d7725SJed Brown 
374c20d7725SJed Brown    test:
375263f2b91SStefano Zampini      suffix: hypre_medium
376263f2b91SStefano Zampini      nsize: {{1 3}}
377263f2b91SStefano Zampini      requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
378263f2b91SStefano Zampini      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -A_mat_type hypre -B_mat_type hypre -test_rart 0
379263f2b91SStefano Zampini      output_file: output/ex62_hypre.out
380263f2b91SStefano Zampini 
381263f2b91SStefano Zampini    test:
382263f2b91SStefano Zampini      suffix: hypre_tiny
383263f2b91SStefano Zampini      nsize: {{1 3}}
384263f2b91SStefano Zampini      requires: hypre !complex double !defined(PETSC_USE_64BIT_INDICES)
385263f2b91SStefano Zampini      args: -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -A_mat_type hypre -B_mat_type hypre -test_rart 0
386263f2b91SStefano Zampini      output_file: output/ex62_hypre.out
387263f2b91SStefano Zampini 
388263f2b91SStefano Zampini    test:
389cec0a6c6SStefano Zampini      suffix: 9_mkl
39015df37b9SStefano Zampini      TODO: broken MatScale?
391b88df2e7SBarry Smith      requires: mkl_sparse datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
392cec0a6c6SStefano Zampini      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -A_mat_type aijmkl -B_mat_type aijmkl
393cec0a6c6SStefano Zampini      output_file: output/ex62_1.out
394cec0a6c6SStefano Zampini 
395cec0a6c6SStefano Zampini    test:
396c20d7725SJed Brown      suffix: 10
397dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
398c20d7725SJed Brown      nsize: 3
399c20d7725SJed Brown      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium
400c20d7725SJed Brown      output_file: output/ex62_1.out
401c20d7725SJed Brown 
402c20d7725SJed Brown    test:
4033cea4f0aSStefano Zampini      suffix: 10_backend
404dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
4053cea4f0aSStefano Zampini      nsize: 3
4063e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm backend -matmatmult_via backend -AtB_mat_product_algorithm backend -mattransposematmult_via backend -PtAP_mat_product_algorithm backend -matptap_via backend
4073cea4f0aSStefano Zampini      output_file: output/ex62_1.out
4083cea4f0aSStefano Zampini 
4093cea4f0aSStefano Zampini    test:
410c20d7725SJed Brown      suffix: 11_ab_scalable
411dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
412c20d7725SJed Brown      nsize: 3
4133e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm scalable -matmatmult_via scalable -AtB_mat_product_algorithm scalable -mattransposematmult_via scalable
414c20d7725SJed Brown      output_file: output/ex62_1.out
415c20d7725SJed Brown 
416c20d7725SJed Brown    test:
417c20d7725SJed Brown      suffix: 12_ab_seqmpi
418dfd57a17SPierre Jolivet      requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
419c20d7725SJed Brown      nsize: 3
4203e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm seqmpi -matmatmult_via seqmpi -AtB_mat_product_algorithm at*b -mattransposematmult_via at*b
421c20d7725SJed Brown      output_file: output/ex62_1.out
422c20d7725SJed Brown 
423c20d7725SJed Brown    test:
424c20d7725SJed Brown      suffix: 13_ab_hypre
425dfd57a17SPierre Jolivet      requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
426c20d7725SJed Brown      nsize: 3
4273e662e0bSHong Zhang      args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm hypre -matmatmult_via hypre -PtAP_mat_product_algorithm hypre -matptap_via hypre
428c20d7725SJed Brown      output_file: output/ex62_1.out
429c20d7725SJed Brown 
43020b3374bSStefano Zampini    test:
43120b3374bSStefano Zampini      suffix: 14_seqaij
432dfd57a17SPierre Jolivet      requires: !complex double !defined(PETSC_USE_64BIT_INDICES)
43320b3374bSStefano Zampini      args: -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system
43420b3374bSStefano Zampini      output_file: output/ex62_1.out
43520b3374bSStefano Zampini 
43620b3374bSStefano Zampini    test:
43720b3374bSStefano Zampini      suffix: 14_seqaijcusparse
438dfd57a17SPierre Jolivet      requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES)
4391a2c6b5cSJunchao Zhang      args: -A_mat_type aijcusparse -B_mat_type aijcusparse -mat_form_explicit_transpose -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system
44020b3374bSStefano Zampini      output_file: output/ex62_1.out
44120b3374bSStefano Zampini 
44220b3374bSStefano Zampini    test:
44365e4b4d4SStefano Zampini      suffix: 14_seqaijcusparse_cpu
444dfd57a17SPierre Jolivet      requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES)
4453e662e0bSHong Zhang      args: -A_mat_type aijcusparse -B_mat_type aijcusparse -mat_form_explicit_transpose -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -AB_mat_product_algorithm_backend_cpu -matmatmult_backend_cpu -PtAP_mat_product_algorithm_backend_cpu -matptap_backend_cpu -RARt_mat_product_algorithm_backend_cpu -matrart_backend_cpu
44665e4b4d4SStefano Zampini      output_file: output/ex62_1.out
44765e4b4d4SStefano Zampini 
44865e4b4d4SStefano Zampini    test:
4495cb69cabSStefano Zampini      suffix: 14_mpiaijcusparse_seq
4505cb69cabSStefano Zampini      nsize: 1
451dfd57a17SPierre Jolivet      requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES)
4521a2c6b5cSJunchao Zhang      args: -A_mat_type mpiaijcusparse -B_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system
4535cb69cabSStefano Zampini      output_file: output/ex62_1.out
4545cb69cabSStefano Zampini 
4555cb69cabSStefano Zampini    test:
45665e4b4d4SStefano Zampini      suffix: 14_mpiaijcusparse_seq_cpu
45765e4b4d4SStefano Zampini      nsize: 1
458dfd57a17SPierre Jolivet      requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES)
4593e662e0bSHong Zhang      args: -A_mat_type mpiaijcusparse -B_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -AB_mat_product_algorithm_backend_cpu -matmatmult_backend_cpu -PtAP_mat_product_algorithm_backend_cpu -matptap_backend_cpu
46065e4b4d4SStefano Zampini      output_file: output/ex62_1.out
46165e4b4d4SStefano Zampini 
46265e4b4d4SStefano Zampini    test:
4635cb69cabSStefano Zampini      suffix: 14_mpiaijcusparse
4645cb69cabSStefano Zampini      nsize: 3
465dfd57a17SPierre Jolivet      requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES)
4661a2c6b5cSJunchao Zhang      args: -A_mat_type mpiaijcusparse -B_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system
4675cb69cabSStefano Zampini      output_file: output/ex62_1.out
4685cb69cabSStefano Zampini 
4695cb69cabSStefano Zampini    test:
47065e4b4d4SStefano Zampini      suffix: 14_mpiaijcusparse_cpu
47165e4b4d4SStefano Zampini      nsize: 3
472dfd57a17SPierre Jolivet      requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES)
4733e662e0bSHong Zhang      args: -A_mat_type mpiaijcusparse -B_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -AB_mat_product_algorithm_backend_cpu -matmatmult_backend_cpu -PtAP_mat_product_algorithm_backend_cpu -matptap_backend_cpu
47465e4b4d4SStefano Zampini      output_file: output/ex62_1.out
47565e4b4d4SStefano Zampini 
47665e4b4d4SStefano Zampini    test:
477076ba34aSJunchao Zhang      nsize: {{1 3}}
478076ba34aSJunchao Zhang      suffix: 14_aijkokkos
4793078479eSJunchao Zhang      requires: !sycl kokkos_kernels !complex double !defined(PETSC_USE_64BIT_INDICES)
48020b3374bSStefano Zampini      args: -A_mat_type aijkokkos -B_mat_type aijkokkos -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system
48120b3374bSStefano Zampini      output_file: output/ex62_1.out
48220b3374bSStefano Zampini 
4835cb69cabSStefano Zampini    # these tests use matrices with many zero rows
48420b3374bSStefano Zampini    test:
48520b3374bSStefano Zampini      suffix: 15_seqaijcusparse
486dfd57a17SPierre Jolivet      requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath
4871a2c6b5cSJunchao Zhang      args: -A_mat_type aijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith
48820b3374bSStefano Zampini      output_file: output/ex62_1.out
4895cb69cabSStefano Zampini 
4905cb69cabSStefano Zampini    test:
4915cb69cabSStefano Zampini      suffix: 15_mpiaijcusparse_seq
4925cb69cabSStefano Zampini      nsize: 1
493dfd57a17SPierre Jolivet      requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath
4941a2c6b5cSJunchao Zhang      args: -A_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith
4955cb69cabSStefano Zampini      output_file: output/ex62_1.out
4965cb69cabSStefano Zampini 
4975cb69cabSStefano Zampini    test:
4985cb69cabSStefano Zampini      nsize: 3
4995cb69cabSStefano Zampini      suffix: 15_mpiaijcusparse
500dfd57a17SPierre Jolivet      requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath
5011a2c6b5cSJunchao Zhang      args: -A_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith
5025cb69cabSStefano Zampini      output_file: output/ex62_1.out
5035cb69cabSStefano Zampini 
504c20d7725SJed Brown TEST*/
505