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; 199566063dSJacob Faibussowitsch PetscCall(MatAXPY(B,-1.0,A,DIFFERENT_NONZERO_PATTERN)); 209566063dSJacob Faibussowitsch PetscCall(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 PetscMPIInt size,rank; 2920b3374bSStefano Zampini PetscInt i,j,*idxn,PM,PN = PETSC_DECIDE,rstart,rend; 30c20d7725SJed Brown PetscReal norm; 31c20d7725SJed Brown PetscRandom rdm; 3220b3374bSStefano Zampini char file[2][PETSC_MAX_PATH_LEN] = { "", ""}; 33c20d7725SJed Brown PetscScalar *a,rval,alpha; 34c20d7725SJed Brown PetscBool Test_MatMatMult=PETSC_TRUE,Test_MatTrMat=PETSC_TRUE,Test_MatMatTr=PETSC_TRUE; 3520b3374bSStefano Zampini PetscBool Test_MatPtAP=PETSC_TRUE,Test_MatRARt=PETSC_TRUE,flg,seqaij,flgA,flgB; 36c20d7725SJed Brown MatInfo info; 3720b3374bSStefano Zampini PetscInt nzp = 5; /* num of nonzeros in each row of P */ 38c20d7725SJed Brown MatType mattype; 3920b3374bSStefano Zampini const char *deft = MATAIJ; 4020b3374bSStefano Zampini char A_mattype[256], B_mattype[256]; 4120b3374bSStefano Zampini PetscInt mcheck = 10; 42c20d7725SJed Brown 439566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc,&args,(char*)0,help)); 449566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 459566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 46c20d7725SJed Brown 47c20d7725SJed Brown /* Load the matrices A_save and B */ 48*d0609cedSBarry Smith PetscOptionsBegin(PETSC_COMM_WORLD,"","",""); 499566063dSJacob Faibussowitsch PetscCall(PetscOptionsBool("-test_rart","Test MatRARt","",Test_MatRARt,&Test_MatRARt,NULL)); 509566063dSJacob Faibussowitsch PetscCall(PetscOptionsInt("-PN","Number of columns of P","",PN,&PN,NULL)); 519566063dSJacob Faibussowitsch PetscCall(PetscOptionsInt("-mcheck","Number of matmult checks","",mcheck,&mcheck,NULL)); 529566063dSJacob Faibussowitsch PetscCall(PetscOptionsString("-fA","Path for matrix A","",file[0],file[0],sizeof(file[0]),&flg)); 5328b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_USER_INPUT,"Must indicate a file name for matrix A with the -fA option."); 549566063dSJacob Faibussowitsch PetscCall(PetscOptionsString("-fB","Path for matrix B","",file[1],file[1],sizeof(file[1]),&flg)); 559566063dSJacob Faibussowitsch PetscCall(PetscOptionsFList("-A_mat_type","Matrix type","MatSetType",MatList,deft,A_mattype,256,&flgA)); 569566063dSJacob Faibussowitsch PetscCall(PetscOptionsFList("-B_mat_type","Matrix type","MatSetType",MatList,deft,B_mattype,256,&flgB)); 57*d0609cedSBarry Smith PetscOptionsEnd(); 58c20d7725SJed Brown 599566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[0],FILE_MODE_READ,&viewer)); 609566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD,&A_save)); 619566063dSJacob Faibussowitsch PetscCall(MatLoad(A_save,viewer)); 629566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&viewer)); 63c20d7725SJed Brown 6420b3374bSStefano Zampini if (flg) { 659566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[1],FILE_MODE_READ,&viewer)); 669566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD,&B)); 679566063dSJacob Faibussowitsch PetscCall(MatLoad(B,viewer)); 689566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&viewer)); 6920b3374bSStefano Zampini } else { 709566063dSJacob Faibussowitsch PetscCall(PetscObjectReference((PetscObject)A_save)); 7120b3374bSStefano Zampini B = A_save; 7220b3374bSStefano Zampini } 7320b3374bSStefano Zampini 7420b3374bSStefano Zampini if (flgA) { 759566063dSJacob Faibussowitsch PetscCall(MatConvert(A_save,A_mattype,MAT_INPLACE_MATRIX,&A_save)); 7620b3374bSStefano Zampini } 7720b3374bSStefano Zampini if (flgB) { 789566063dSJacob Faibussowitsch PetscCall(MatConvert(B,B_mattype,MAT_INPLACE_MATRIX,&B)); 7920b3374bSStefano Zampini } 809566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(A_save)); 819566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(B)); 82c20d7725SJed Brown 839566063dSJacob Faibussowitsch PetscCall(MatGetType(B,&mattype)); 84c20d7725SJed Brown 859566063dSJacob Faibussowitsch PetscCall(PetscMalloc(nzp*(sizeof(PetscInt)+sizeof(PetscScalar)),&idxn)); 86c20d7725SJed Brown a = (PetscScalar*)(idxn + nzp); 87c20d7725SJed Brown 889566063dSJacob Faibussowitsch PetscCall(PetscRandomCreate(PETSC_COMM_WORLD,&rdm)); 899566063dSJacob Faibussowitsch PetscCall(PetscRandomSetFromOptions(rdm)); 90c20d7725SJed Brown 91c20d7725SJed Brown /* 1) MatMatMult() */ 92c20d7725SJed Brown /* ----------------*/ 93c20d7725SJed Brown if (Test_MatMatMult) { 949566063dSJacob Faibussowitsch PetscCall(MatDuplicate(A_save,MAT_COPY_VALUES,&A)); 95c20d7725SJed Brown 96c20d7725SJed Brown /* (1.1) Test developer API */ 979566063dSJacob Faibussowitsch PetscCall(MatProductCreate(A,B,NULL,&C)); 989566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C,"AB_")); 999566063dSJacob Faibussowitsch PetscCall(MatProductSetType(C,MATPRODUCT_AB)); 1009566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(C,MATPRODUCTALGORITHMDEFAULT)); 1019566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(C,PETSC_DEFAULT)); 1029566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(C)); 103910fa13fSStefano Zampini /* we can inquire about MATOP_PRODUCTSYMBOLIC even if the destination matrix type has not been set yet */ 1049566063dSJacob Faibussowitsch PetscCall(MatHasOperation(C,MATOP_PRODUCTSYMBOLIC,&flg)); 1059566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(C)); 1069566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 1079566063dSJacob Faibussowitsch PetscCall(MatMatMultEqual(A,B,C,mcheck,&flg)); 10828b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in C=A*B"); 109c20d7725SJed Brown 110c20d7725SJed Brown /* Test reuse symbolic C */ 111c20d7725SJed Brown alpha = 0.9; 1129566063dSJacob Faibussowitsch PetscCall(MatScale(A,alpha)); 1139566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 114c20d7725SJed Brown 1159566063dSJacob Faibussowitsch PetscCall(MatMatMultEqual(A,B,C,mcheck,&flg)); 11628b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in C=A*B"); 1179566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 118c20d7725SJed Brown 119c20d7725SJed Brown /* (1.2) Test user driver */ 1209566063dSJacob Faibussowitsch PetscCall(MatMatMult(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C)); 121c20d7725SJed Brown 122c20d7725SJed Brown /* Test MAT_REUSE_MATRIX - reuse symbolic C */ 123c20d7725SJed Brown alpha = 1.0; 124c20d7725SJed Brown for (i=0; i<2; i++) { 125c20d7725SJed Brown alpha -= 0.1; 1269566063dSJacob Faibussowitsch PetscCall(MatScale(A,alpha)); 1279566063dSJacob Faibussowitsch PetscCall(MatMatMult(A,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C)); 128c20d7725SJed Brown } 1299566063dSJacob Faibussowitsch PetscCall(MatMatMultEqual(A,B,C,mcheck,&flg)); 13028b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error: MatMatMult()"); 1319566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 1324417c5e8SHong Zhang 1334417c5e8SHong Zhang /* Test MatProductClear() */ 1349566063dSJacob Faibussowitsch PetscCall(MatProductClear(C)); 1359566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 136544a5e07SHong Zhang 137544a5e07SHong Zhang /* Test MatMatMult() for dense and aij matrices */ 1389566063dSJacob Faibussowitsch PetscCall(PetscObjectTypeCompareAny((PetscObject)A,&flg,MATSEQAIJ,MATMPIAIJ,"")); 13920b3374bSStefano Zampini if (flg) { 1409566063dSJacob Faibussowitsch PetscCall(MatConvert(A_save,MATDENSE,MAT_INITIAL_MATRIX,&A)); 1419566063dSJacob Faibussowitsch PetscCall(MatMatMult(A,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C)); 1429566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 1439566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 14420b3374bSStefano Zampini } 145c20d7725SJed Brown } 146c20d7725SJed Brown 147c20d7725SJed Brown /* Create P and R = P^T */ 148c20d7725SJed Brown /* --------------------- */ 1499566063dSJacob Faibussowitsch PetscCall(MatGetSize(B,&PM,NULL)); 15020b3374bSStefano Zampini if (PN < 0) PN = PM/2; 1519566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD,&P)); 1529566063dSJacob Faibussowitsch PetscCall(MatSetSizes(P,PETSC_DECIDE,PETSC_DECIDE,PM,PN)); 1539566063dSJacob Faibussowitsch PetscCall(MatSetType(P,MATAIJ)); 1549566063dSJacob Faibussowitsch PetscCall(MatSeqAIJSetPreallocation(P,nzp,NULL)); 1559566063dSJacob Faibussowitsch PetscCall(MatMPIAIJSetPreallocation(P,nzp,NULL,nzp,NULL)); 1569566063dSJacob Faibussowitsch PetscCall(MatGetOwnershipRange(P,&rstart,&rend)); 157c20d7725SJed Brown for (i=0; i<nzp; i++) { 1589566063dSJacob Faibussowitsch PetscCall(PetscRandomGetValue(rdm,&a[i])); 159c20d7725SJed Brown } 160c20d7725SJed Brown for (i=rstart; i<rend; i++) { 161c20d7725SJed Brown for (j=0; j<nzp; j++) { 1629566063dSJacob Faibussowitsch PetscCall(PetscRandomGetValue(rdm,&rval)); 163c20d7725SJed Brown idxn[j] = (PetscInt)(PetscRealPart(rval)*PN); 164c20d7725SJed Brown } 1659566063dSJacob Faibussowitsch PetscCall(MatSetValues(P,1,&i,nzp,idxn,a,ADD_VALUES)); 166c20d7725SJed Brown } 1679566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(P,MAT_FINAL_ASSEMBLY)); 1689566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(P,MAT_FINAL_ASSEMBLY)); 169c20d7725SJed Brown 1709566063dSJacob Faibussowitsch PetscCall(MatTranspose(P,MAT_INITIAL_MATRIX,&R)); 1719566063dSJacob Faibussowitsch PetscCall(MatConvert(P,mattype,MAT_INPLACE_MATRIX,&P)); 1729566063dSJacob Faibussowitsch PetscCall(MatConvert(R,mattype,MAT_INPLACE_MATRIX,&R)); 1739566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(P)); 1749566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(R)); 175c20d7725SJed Brown 176c20d7725SJed Brown /* 2) MatTransposeMatMult() */ 177c20d7725SJed Brown /* ------------------------ */ 178c20d7725SJed Brown if (Test_MatTrMat) { 179c20d7725SJed Brown /* (2.1) Test developer driver C = P^T*B */ 1809566063dSJacob Faibussowitsch PetscCall(MatProductCreate(P,B,NULL,&C)); 1819566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C,"AtB_")); 1829566063dSJacob Faibussowitsch PetscCall(MatProductSetType(C,MATPRODUCT_AtB)); 1839566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(C,MATPRODUCTALGORITHMDEFAULT)); 1849566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(C,PETSC_DEFAULT)); 1859566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(C)); 1869566063dSJacob Faibussowitsch PetscCall(MatHasOperation(C,MATOP_PRODUCTSYMBOLIC,&flg)); 187263f2b91SStefano Zampini if (flg) { /* run tests if supported */ 1889566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(C)); /* equivalent to MatSetUp() */ 1899566063dSJacob Faibussowitsch PetscCall(MatSetOption(C,MAT_USE_INODES,PETSC_FALSE)); /* illustrate how to call MatSetOption() */ 1909566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 1919566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); /* test reuse symbolic C */ 19267b3012eSStefano Zampini 1939566063dSJacob Faibussowitsch PetscCall(MatTransposeMatMultEqual(P,B,C,mcheck,&flg)); 19428b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error: developer driver C = P^T*B"); 1959566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 196c20d7725SJed Brown 197c20d7725SJed Brown /* (2.2) Test user driver C = P^T*B */ 1989566063dSJacob Faibussowitsch PetscCall(MatTransposeMatMult(P,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C)); 1999566063dSJacob Faibussowitsch PetscCall(MatTransposeMatMult(P,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C)); 2009566063dSJacob Faibussowitsch PetscCall(MatGetInfo(C,MAT_GLOBAL_SUM,&info)); 2019566063dSJacob Faibussowitsch PetscCall(MatProductClear(C)); 202c20d7725SJed Brown 203c20d7725SJed Brown /* Compare P^T*B and R*B */ 2049566063dSJacob Faibussowitsch PetscCall(MatMatMult(R,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C1)); 2059566063dSJacob Faibussowitsch PetscCall(MatNormDifference(C,C1,&norm)); 20608401ef6SPierre Jolivet PetscCheck(norm <= PETSC_SMALL,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatTransposeMatMult(): %g",(double)norm); 2079566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C1)); 208c20d7725SJed Brown 209c20d7725SJed Brown /* Test MatDuplicate() of C=P^T*B */ 2109566063dSJacob Faibussowitsch PetscCall(MatDuplicate(C,MAT_COPY_VALUES,&C1)); 2119566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C1)); 212263f2b91SStefano Zampini } else { 2139566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD,"MatTransposeMatMult not supported\n")); 214263f2b91SStefano Zampini } 2159566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 216c20d7725SJed Brown } 217c20d7725SJed Brown 21867b3012eSStefano Zampini /* 3) MatMatTransposeMult() */ 219c20d7725SJed Brown /* ------------------------ */ 220c20d7725SJed Brown if (Test_MatMatTr) { 221c20d7725SJed Brown /* C = B*R^T */ 2229566063dSJacob Faibussowitsch PetscCall(PetscObjectBaseTypeCompare((PetscObject)B,MATSEQAIJ,&seqaij)); 22320b3374bSStefano Zampini if (seqaij) { 2249566063dSJacob Faibussowitsch PetscCall(MatMatTransposeMult(B,R,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C)); 2259566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C,"ABt_")); /* enable '-ABt_' for matrix C */ 2269566063dSJacob Faibussowitsch PetscCall(MatGetInfo(C,MAT_GLOBAL_SUM,&info)); 227c20d7725SJed Brown 228c20d7725SJed Brown /* Test MAT_REUSE_MATRIX - reuse symbolic C */ 2299566063dSJacob Faibussowitsch PetscCall(MatMatTransposeMult(B,R,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C)); 230c20d7725SJed Brown 231c20d7725SJed Brown /* Check */ 2329566063dSJacob Faibussowitsch PetscCall(MatMatMult(B,P,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C1)); 2339566063dSJacob Faibussowitsch PetscCall(MatNormDifference(C,C1,&norm)); 23408401ef6SPierre Jolivet PetscCheck(norm <= PETSC_SMALL,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatMatTransposeMult() %g",(double)norm); 2359566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C1)); 2369566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 237c20d7725SJed Brown } 238c20d7725SJed Brown } 239c20d7725SJed Brown 240c20d7725SJed Brown /* 4) Test MatPtAP() */ 241c20d7725SJed Brown /*-------------------*/ 242c20d7725SJed Brown if (Test_MatPtAP) { 2439566063dSJacob Faibussowitsch PetscCall(MatDuplicate(A_save,MAT_COPY_VALUES,&A)); 244c20d7725SJed Brown 245c20d7725SJed Brown /* (4.1) Test developer API */ 2469566063dSJacob Faibussowitsch PetscCall(MatProductCreate(A,P,NULL,&C)); 2479566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C,"PtAP_")); 2489566063dSJacob Faibussowitsch PetscCall(MatProductSetType(C,MATPRODUCT_PtAP)); 2499566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(C,MATPRODUCTALGORITHMDEFAULT)); 2509566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(C,PETSC_DEFAULT)); 2519566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(C)); 2529566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(C)); 2539566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 2549566063dSJacob Faibussowitsch PetscCall(MatPtAPMultEqual(A,P,C,mcheck,&flg)); 25528b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatProduct_PtAP"); 2569566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); /* reuse symbolic C */ 257c20d7725SJed Brown 2589566063dSJacob Faibussowitsch PetscCall(MatPtAPMultEqual(A,P,C,mcheck,&flg)); 25928b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatProduct_PtAP"); 2609566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 261c20d7725SJed Brown 262c20d7725SJed Brown /* (4.2) Test user driver */ 2639566063dSJacob Faibussowitsch PetscCall(MatPtAP(A,P,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&C)); 264c20d7725SJed Brown 265c20d7725SJed Brown /* Test MAT_REUSE_MATRIX - reuse symbolic C */ 266c20d7725SJed Brown alpha = 1.0; 267c20d7725SJed Brown for (i=0; i<2; i++) { 268c20d7725SJed Brown alpha -= 0.1; 2699566063dSJacob Faibussowitsch PetscCall(MatScale(A,alpha)); 2709566063dSJacob Faibussowitsch PetscCall(MatPtAP(A,P,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C)); 271c20d7725SJed Brown } 2729566063dSJacob Faibussowitsch PetscCall(MatPtAPMultEqual(A,P,C,mcheck,&flg)); 27328b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_WORLD,PETSC_ERR_PLIB,"Error in MatPtAP"); 274c20d7725SJed Brown 275c20d7725SJed Brown /* 5) Test MatRARt() */ 276c20d7725SJed Brown /* ----------------- */ 277c20d7725SJed Brown if (Test_MatRARt) { 278c20d7725SJed Brown Mat RARt; 27967b3012eSStefano Zampini 28067b3012eSStefano Zampini /* (5.1) Test developer driver RARt = R*A*Rt */ 2819566063dSJacob Faibussowitsch PetscCall(MatProductCreate(A,R,NULL,&RARt)); 2829566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(RARt,"RARt_")); 2839566063dSJacob Faibussowitsch PetscCall(MatProductSetType(RARt,MATPRODUCT_RARt)); 2849566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(RARt,MATPRODUCTALGORITHMDEFAULT)); 2859566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(RARt,PETSC_DEFAULT)); 2869566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(RARt)); 2879566063dSJacob Faibussowitsch PetscCall(MatHasOperation(RARt,MATOP_PRODUCTSYMBOLIC,&flg)); 288263f2b91SStefano Zampini if (flg) { 2899566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(RARt)); /* equivalent to MatSetUp() */ 2909566063dSJacob Faibussowitsch PetscCall(MatSetOption(RARt,MAT_USE_INODES,PETSC_FALSE)); /* illustrate how to call MatSetOption() */ 2919566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(RARt)); 2929566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(RARt)); /* test reuse symbolic RARt */ 2939566063dSJacob Faibussowitsch PetscCall(MatDestroy(&RARt)); 29467b3012eSStefano Zampini 29567b3012eSStefano Zampini /* (2.2) Test user driver RARt = R*A*Rt */ 2969566063dSJacob Faibussowitsch PetscCall(MatRARt(A,R,MAT_INITIAL_MATRIX,2.0,&RARt)); 2979566063dSJacob Faibussowitsch PetscCall(MatRARt(A,R,MAT_REUSE_MATRIX,2.0,&RARt)); 29867b3012eSStefano Zampini 2999566063dSJacob Faibussowitsch PetscCall(MatNormDifference(C,RARt,&norm)); 30008401ef6SPierre Jolivet PetscCheck(norm <= PETSC_SMALL,PETSC_COMM_SELF,PETSC_ERR_PLIB,"|PtAP - RARt| = %g",(double)norm); 301263f2b91SStefano Zampini } else { 3029566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD,"MatRARt not supported\n")); 303263f2b91SStefano Zampini } 3049566063dSJacob Faibussowitsch PetscCall(MatDestroy(&RARt)); 305c20d7725SJed Brown } 306c20d7725SJed Brown 3079566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 3089566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 309c20d7725SJed Brown } 310c20d7725SJed Brown 311c20d7725SJed Brown /* Destroy objects */ 3129566063dSJacob Faibussowitsch PetscCall(PetscRandomDestroy(&rdm)); 3139566063dSJacob Faibussowitsch PetscCall(PetscFree(idxn)); 314c20d7725SJed Brown 3159566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A_save)); 3169566063dSJacob Faibussowitsch PetscCall(MatDestroy(&B)); 3179566063dSJacob Faibussowitsch PetscCall(MatDestroy(&P)); 3189566063dSJacob Faibussowitsch PetscCall(MatDestroy(&R)); 319c20d7725SJed Brown 3209566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 321b122ec5aSJacob Faibussowitsch return 0; 322c20d7725SJed Brown } 323c20d7725SJed Brown 324c20d7725SJed Brown /*TEST 325c20d7725SJed Brown test: 326c20d7725SJed Brown suffix: 1 327dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 328c20d7725SJed Brown args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium 329c20d7725SJed Brown output_file: output/ex62_1.out 330c20d7725SJed Brown 331c20d7725SJed Brown test: 332c20d7725SJed Brown suffix: 2_ab_scalable 333dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3343e662e0bSHong 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 335c20d7725SJed Brown output_file: output/ex62_1.out 336c20d7725SJed Brown 337c20d7725SJed Brown test: 338c20d7725SJed Brown suffix: 3_ab_scalable_fast 339dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3403e662e0bSHong Zhang args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm scalable_fast -matmatmult_via scalable_fast -matmattransmult_via color 341c20d7725SJed Brown output_file: output/ex62_1.out 342c20d7725SJed Brown 343c20d7725SJed Brown test: 344c20d7725SJed Brown suffix: 4_ab_heap 345dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3463e662e0bSHong 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 347c20d7725SJed Brown output_file: output/ex62_1.out 348c20d7725SJed Brown 349c20d7725SJed Brown test: 350c20d7725SJed Brown suffix: 5_ab_btheap 351dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3523e662e0bSHong Zhang args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm btheap -matmatmult_via btheap -matrart_via r*art 353c20d7725SJed Brown output_file: output/ex62_1.out 354c20d7725SJed Brown 355c20d7725SJed Brown test: 356c20d7725SJed Brown suffix: 6_ab_llcondensed 357dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3583e662e0bSHong Zhang args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm llcondensed -matmatmult_via llcondensed -matrart_via coloring_rart 359c20d7725SJed Brown output_file: output/ex62_1.out 360c20d7725SJed Brown 361c20d7725SJed Brown test: 362c20d7725SJed Brown suffix: 7_ab_rowmerge 363dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3643e662e0bSHong Zhang args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm rowmerge -matmatmult_via rowmerge 365c20d7725SJed Brown output_file: output/ex62_1.out 366c20d7725SJed Brown 367c20d7725SJed Brown test: 368c20d7725SJed Brown suffix: 8_ab_hypre 369dfd57a17SPierre Jolivet requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3703e662e0bSHong 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 371c20d7725SJed Brown output_file: output/ex62_1.out 372c20d7725SJed Brown 373c20d7725SJed Brown test: 374263f2b91SStefano Zampini suffix: hypre_medium 375263f2b91SStefano Zampini nsize: {{1 3}} 376263f2b91SStefano Zampini requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 377263f2b91SStefano Zampini args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -A_mat_type hypre -B_mat_type hypre -test_rart 0 378263f2b91SStefano Zampini output_file: output/ex62_hypre.out 379263f2b91SStefano Zampini 380263f2b91SStefano Zampini test: 381263f2b91SStefano Zampini suffix: hypre_tiny 382263f2b91SStefano Zampini nsize: {{1 3}} 383263f2b91SStefano Zampini requires: hypre !complex double !defined(PETSC_USE_64BIT_INDICES) 384263f2b91SStefano 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 385263f2b91SStefano Zampini output_file: output/ex62_hypre.out 386263f2b91SStefano Zampini 387263f2b91SStefano Zampini test: 388cec0a6c6SStefano Zampini suffix: 9_mkl 38915df37b9SStefano Zampini TODO: broken MatScale? 390b88df2e7SBarry Smith requires: mkl_sparse datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 391cec0a6c6SStefano Zampini args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -A_mat_type aijmkl -B_mat_type aijmkl 392cec0a6c6SStefano Zampini output_file: output/ex62_1.out 393cec0a6c6SStefano Zampini 394cec0a6c6SStefano Zampini test: 395c20d7725SJed Brown suffix: 10 396dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 397c20d7725SJed Brown nsize: 3 398c20d7725SJed Brown args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium 399c20d7725SJed Brown output_file: output/ex62_1.out 400c20d7725SJed Brown 401c20d7725SJed Brown test: 4023cea4f0aSStefano Zampini suffix: 10_backend 403dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 4043cea4f0aSStefano Zampini nsize: 3 4053e662e0bSHong 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 4063cea4f0aSStefano Zampini output_file: output/ex62_1.out 4073cea4f0aSStefano Zampini 4083cea4f0aSStefano Zampini test: 409c20d7725SJed Brown suffix: 11_ab_scalable 410dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 411c20d7725SJed Brown nsize: 3 4123e662e0bSHong 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 413c20d7725SJed Brown output_file: output/ex62_1.out 414c20d7725SJed Brown 415c20d7725SJed Brown test: 416c20d7725SJed Brown suffix: 12_ab_seqmpi 417dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 418c20d7725SJed Brown nsize: 3 4193e662e0bSHong 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 420c20d7725SJed Brown output_file: output/ex62_1.out 421c20d7725SJed Brown 422c20d7725SJed Brown test: 423c20d7725SJed Brown suffix: 13_ab_hypre 424dfd57a17SPierre Jolivet requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 425c20d7725SJed Brown nsize: 3 4263e662e0bSHong 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 427c20d7725SJed Brown output_file: output/ex62_1.out 428c20d7725SJed Brown 42920b3374bSStefano Zampini test: 43020b3374bSStefano Zampini suffix: 14_seqaij 431dfd57a17SPierre Jolivet requires: !complex double !defined(PETSC_USE_64BIT_INDICES) 43220b3374bSStefano Zampini args: -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system 43320b3374bSStefano Zampini output_file: output/ex62_1.out 43420b3374bSStefano Zampini 43520b3374bSStefano Zampini test: 43620b3374bSStefano Zampini suffix: 14_seqaijcusparse 437dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4381a2c6b5cSJunchao 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 43920b3374bSStefano Zampini output_file: output/ex62_1.out 44020b3374bSStefano Zampini 44120b3374bSStefano Zampini test: 44265e4b4d4SStefano Zampini suffix: 14_seqaijcusparse_cpu 443dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4443e662e0bSHong 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 44565e4b4d4SStefano Zampini output_file: output/ex62_1.out 44665e4b4d4SStefano Zampini 44765e4b4d4SStefano Zampini test: 4485cb69cabSStefano Zampini suffix: 14_mpiaijcusparse_seq 4495cb69cabSStefano Zampini nsize: 1 450dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4511a2c6b5cSJunchao 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 4525cb69cabSStefano Zampini output_file: output/ex62_1.out 4535cb69cabSStefano Zampini 4545cb69cabSStefano Zampini test: 45565e4b4d4SStefano Zampini suffix: 14_mpiaijcusparse_seq_cpu 45665e4b4d4SStefano Zampini nsize: 1 457dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4584005f4d6SPierre Jolivet 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 -test_rart 0 45965e4b4d4SStefano Zampini output_file: output/ex62_1.out 46065e4b4d4SStefano Zampini 46165e4b4d4SStefano Zampini test: 4625cb69cabSStefano Zampini suffix: 14_mpiaijcusparse 4635cb69cabSStefano Zampini nsize: 3 464dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4651a2c6b5cSJunchao 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 4665cb69cabSStefano Zampini output_file: output/ex62_1.out 4675cb69cabSStefano Zampini 4685cb69cabSStefano Zampini test: 46965e4b4d4SStefano Zampini suffix: 14_mpiaijcusparse_cpu 47065e4b4d4SStefano Zampini nsize: 3 471dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4724005f4d6SPierre Jolivet 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 -test_rart 0 47365e4b4d4SStefano Zampini output_file: output/ex62_1.out 47465e4b4d4SStefano Zampini 47565e4b4d4SStefano Zampini test: 476076ba34aSJunchao Zhang nsize: {{1 3}} 477076ba34aSJunchao Zhang suffix: 14_aijkokkos 4783078479eSJunchao Zhang requires: !sycl kokkos_kernels !complex double !defined(PETSC_USE_64BIT_INDICES) 47920b3374bSStefano 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 48020b3374bSStefano Zampini output_file: output/ex62_1.out 48120b3374bSStefano Zampini 4825cb69cabSStefano Zampini # these tests use matrices with many zero rows 48320b3374bSStefano Zampini test: 48420b3374bSStefano Zampini suffix: 15_seqaijcusparse 485dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath 4861a2c6b5cSJunchao Zhang args: -A_mat_type aijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith 48720b3374bSStefano Zampini output_file: output/ex62_1.out 4885cb69cabSStefano Zampini 4895cb69cabSStefano Zampini test: 4905cb69cabSStefano Zampini suffix: 15_mpiaijcusparse_seq 4915cb69cabSStefano Zampini nsize: 1 492dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath 4931a2c6b5cSJunchao Zhang args: -A_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith 4945cb69cabSStefano Zampini output_file: output/ex62_1.out 4955cb69cabSStefano Zampini 4965cb69cabSStefano Zampini test: 4975cb69cabSStefano Zampini nsize: 3 4985cb69cabSStefano Zampini suffix: 15_mpiaijcusparse 499dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath 5001a2c6b5cSJunchao Zhang args: -A_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith 5015cb69cabSStefano Zampini output_file: output/ex62_1.out 5025cb69cabSStefano Zampini 503c20d7725SJed Brown TEST*/ 504