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 */ 169371c9d4SSatish Balay PetscErrorCode MatNormDifference(Mat A, Mat B, PetscReal *norm) { 17c20d7725SJed Brown PetscFunctionBegin; 189566063dSJacob Faibussowitsch PetscCall(MatAXPY(B, -1.0, A, DIFFERENT_NONZERO_PATTERN)); 199566063dSJacob Faibussowitsch PetscCall(MatNorm(B, NORM_FROBENIUS, norm)); 20c20d7725SJed Brown PetscFunctionReturn(0); 21c20d7725SJed Brown } 22c20d7725SJed Brown 239371c9d4SSatish Balay int main(int argc, char **args) { 24c20d7725SJed Brown Mat A, A_save, B, C, P, C1, R; 25c20d7725SJed Brown PetscViewer viewer; 26c20d7725SJed Brown PetscMPIInt size, rank; 2720b3374bSStefano Zampini PetscInt i, j, *idxn, PM, PN = PETSC_DECIDE, rstart, rend; 28c20d7725SJed Brown PetscReal norm; 29c20d7725SJed Brown PetscRandom rdm; 3020b3374bSStefano Zampini char file[2][PETSC_MAX_PATH_LEN] = {"", ""}; 31c20d7725SJed Brown PetscScalar *a, rval, alpha; 32c20d7725SJed Brown PetscBool Test_MatMatMult = PETSC_TRUE, Test_MatTrMat = PETSC_TRUE, Test_MatMatTr = PETSC_TRUE; 3320b3374bSStefano Zampini PetscBool Test_MatPtAP = PETSC_TRUE, Test_MatRARt = PETSC_TRUE, flg, seqaij, flgA, flgB; 34c20d7725SJed Brown MatInfo info; 3520b3374bSStefano Zampini PetscInt nzp = 5; /* num of nonzeros in each row of P */ 36c20d7725SJed Brown MatType mattype; 3720b3374bSStefano Zampini const char *deft = MATAIJ; 3820b3374bSStefano Zampini char A_mattype[256], B_mattype[256]; 3920b3374bSStefano Zampini PetscInt mcheck = 10; 40c20d7725SJed Brown 41327415f7SBarry Smith PetscFunctionBeginUser; 429566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 439566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 449566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 45c20d7725SJed Brown 46c20d7725SJed Brown /* Load the matrices A_save and B */ 47d0609cedSBarry Smith PetscOptionsBegin(PETSC_COMM_WORLD, "", "", ""); 489566063dSJacob Faibussowitsch PetscCall(PetscOptionsBool("-test_rart", "Test MatRARt", "", Test_MatRARt, &Test_MatRARt, NULL)); 499566063dSJacob Faibussowitsch PetscCall(PetscOptionsInt("-PN", "Number of columns of P", "", PN, &PN, NULL)); 509566063dSJacob Faibussowitsch PetscCall(PetscOptionsInt("-mcheck", "Number of matmult checks", "", mcheck, &mcheck, NULL)); 519566063dSJacob Faibussowitsch PetscCall(PetscOptionsString("-fA", "Path for matrix A", "", file[0], file[0], sizeof(file[0]), &flg)); 5228b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_USER_INPUT, "Must indicate a file name for matrix A with the -fA option."); 539566063dSJacob Faibussowitsch PetscCall(PetscOptionsString("-fB", "Path for matrix B", "", file[1], file[1], sizeof(file[1]), &flg)); 549566063dSJacob Faibussowitsch PetscCall(PetscOptionsFList("-A_mat_type", "Matrix type", "MatSetType", MatList, deft, A_mattype, 256, &flgA)); 559566063dSJacob Faibussowitsch PetscCall(PetscOptionsFList("-B_mat_type", "Matrix type", "MatSetType", MatList, deft, B_mattype, 256, &flgB)); 56d0609cedSBarry Smith PetscOptionsEnd(); 57c20d7725SJed Brown 589566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file[0], FILE_MODE_READ, &viewer)); 599566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD, &A_save)); 609566063dSJacob Faibussowitsch PetscCall(MatLoad(A_save, viewer)); 619566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&viewer)); 62c20d7725SJed Brown 6320b3374bSStefano Zampini if (flg) { 649566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file[1], FILE_MODE_READ, &viewer)); 659566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD, &B)); 669566063dSJacob Faibussowitsch PetscCall(MatLoad(B, viewer)); 679566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&viewer)); 6820b3374bSStefano Zampini } else { 699566063dSJacob Faibussowitsch PetscCall(PetscObjectReference((PetscObject)A_save)); 7020b3374bSStefano Zampini B = A_save; 7120b3374bSStefano Zampini } 7220b3374bSStefano Zampini 739371c9d4SSatish Balay if (flgA) { PetscCall(MatConvert(A_save, A_mattype, MAT_INPLACE_MATRIX, &A_save)); } 749371c9d4SSatish Balay if (flgB) { PetscCall(MatConvert(B, B_mattype, MAT_INPLACE_MATRIX, &B)); } 759566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(A_save)); 769566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(B)); 77c20d7725SJed Brown 789566063dSJacob Faibussowitsch PetscCall(MatGetType(B, &mattype)); 79c20d7725SJed Brown 809566063dSJacob Faibussowitsch PetscCall(PetscMalloc(nzp * (sizeof(PetscInt) + sizeof(PetscScalar)), &idxn)); 81c20d7725SJed Brown a = (PetscScalar *)(idxn + nzp); 82c20d7725SJed Brown 839566063dSJacob Faibussowitsch PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rdm)); 849566063dSJacob Faibussowitsch PetscCall(PetscRandomSetFromOptions(rdm)); 85c20d7725SJed Brown 86c20d7725SJed Brown /* 1) MatMatMult() */ 87c20d7725SJed Brown /* ----------------*/ 88c20d7725SJed Brown if (Test_MatMatMult) { 899566063dSJacob Faibussowitsch PetscCall(MatDuplicate(A_save, MAT_COPY_VALUES, &A)); 90c20d7725SJed Brown 91c20d7725SJed Brown /* (1.1) Test developer API */ 929566063dSJacob Faibussowitsch PetscCall(MatProductCreate(A, B, NULL, &C)); 939566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C, "AB_")); 949566063dSJacob Faibussowitsch PetscCall(MatProductSetType(C, MATPRODUCT_AB)); 959566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(C, MATPRODUCTALGORITHMDEFAULT)); 969566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(C, PETSC_DEFAULT)); 979566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(C)); 98910fa13fSStefano Zampini /* we can inquire about MATOP_PRODUCTSYMBOLIC even if the destination matrix type has not been set yet */ 999566063dSJacob Faibussowitsch PetscCall(MatHasOperation(C, MATOP_PRODUCTSYMBOLIC, &flg)); 1009566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(C)); 1019566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 1029566063dSJacob Faibussowitsch PetscCall(MatMatMultEqual(A, B, C, mcheck, &flg)); 10328b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error in C=A*B"); 104c20d7725SJed Brown 105c20d7725SJed Brown /* Test reuse symbolic C */ 106c20d7725SJed Brown alpha = 0.9; 1079566063dSJacob Faibussowitsch PetscCall(MatScale(A, alpha)); 1089566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 109c20d7725SJed Brown 1109566063dSJacob Faibussowitsch PetscCall(MatMatMultEqual(A, B, C, mcheck, &flg)); 11128b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error in C=A*B"); 1129566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 113c20d7725SJed Brown 114c20d7725SJed Brown /* (1.2) Test user driver */ 1159566063dSJacob Faibussowitsch PetscCall(MatMatMult(A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 116c20d7725SJed Brown 117c20d7725SJed Brown /* Test MAT_REUSE_MATRIX - reuse symbolic C */ 118c20d7725SJed Brown alpha = 1.0; 119c20d7725SJed Brown for (i = 0; i < 2; i++) { 120c20d7725SJed Brown alpha -= 0.1; 1219566063dSJacob Faibussowitsch PetscCall(MatScale(A, alpha)); 1229566063dSJacob Faibussowitsch PetscCall(MatMatMult(A, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C)); 123c20d7725SJed Brown } 1249566063dSJacob Faibussowitsch PetscCall(MatMatMultEqual(A, B, C, mcheck, &flg)); 12528b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error: MatMatMult()"); 1269566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 1274417c5e8SHong Zhang 1284417c5e8SHong Zhang /* Test MatProductClear() */ 1299566063dSJacob Faibussowitsch PetscCall(MatProductClear(C)); 1309566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 131544a5e07SHong Zhang 132544a5e07SHong Zhang /* Test MatMatMult() for dense and aij matrices */ 1339566063dSJacob Faibussowitsch PetscCall(PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJ, MATMPIAIJ, "")); 13420b3374bSStefano Zampini if (flg) { 1359566063dSJacob Faibussowitsch PetscCall(MatConvert(A_save, MATDENSE, MAT_INITIAL_MATRIX, &A)); 1369566063dSJacob Faibussowitsch PetscCall(MatMatMult(A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 1379566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 1389566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 13920b3374bSStefano Zampini } 140c20d7725SJed Brown } 141c20d7725SJed Brown 142c20d7725SJed Brown /* Create P and R = P^T */ 143c20d7725SJed Brown /* --------------------- */ 1449566063dSJacob Faibussowitsch PetscCall(MatGetSize(B, &PM, NULL)); 14520b3374bSStefano Zampini if (PN < 0) PN = PM / 2; 1469566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD, &P)); 1479566063dSJacob Faibussowitsch PetscCall(MatSetSizes(P, PETSC_DECIDE, PETSC_DECIDE, PM, PN)); 1489566063dSJacob Faibussowitsch PetscCall(MatSetType(P, MATAIJ)); 1499566063dSJacob Faibussowitsch PetscCall(MatSeqAIJSetPreallocation(P, nzp, NULL)); 1509566063dSJacob Faibussowitsch PetscCall(MatMPIAIJSetPreallocation(P, nzp, NULL, nzp, NULL)); 1519566063dSJacob Faibussowitsch PetscCall(MatGetOwnershipRange(P, &rstart, &rend)); 1529371c9d4SSatish Balay for (i = 0; i < nzp; i++) { PetscCall(PetscRandomGetValue(rdm, &a[i])); } 153c20d7725SJed Brown for (i = rstart; i < rend; i++) { 154c20d7725SJed Brown for (j = 0; j < nzp; j++) { 1559566063dSJacob Faibussowitsch PetscCall(PetscRandomGetValue(rdm, &rval)); 156c20d7725SJed Brown idxn[j] = (PetscInt)(PetscRealPart(rval) * PN); 157c20d7725SJed Brown } 1589566063dSJacob Faibussowitsch PetscCall(MatSetValues(P, 1, &i, nzp, idxn, a, ADD_VALUES)); 159c20d7725SJed Brown } 1609566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(P, MAT_FINAL_ASSEMBLY)); 1619566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(P, MAT_FINAL_ASSEMBLY)); 162c20d7725SJed Brown 1639566063dSJacob Faibussowitsch PetscCall(MatTranspose(P, MAT_INITIAL_MATRIX, &R)); 1649566063dSJacob Faibussowitsch PetscCall(MatConvert(P, mattype, MAT_INPLACE_MATRIX, &P)); 1659566063dSJacob Faibussowitsch PetscCall(MatConvert(R, mattype, MAT_INPLACE_MATRIX, &R)); 1669566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(P)); 1679566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(R)); 168c20d7725SJed Brown 169c20d7725SJed Brown /* 2) MatTransposeMatMult() */ 170c20d7725SJed Brown /* ------------------------ */ 171c20d7725SJed Brown if (Test_MatTrMat) { 172c20d7725SJed Brown /* (2.1) Test developer driver C = P^T*B */ 1739566063dSJacob Faibussowitsch PetscCall(MatProductCreate(P, B, NULL, &C)); 1749566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C, "AtB_")); 1759566063dSJacob Faibussowitsch PetscCall(MatProductSetType(C, MATPRODUCT_AtB)); 1769566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(C, MATPRODUCTALGORITHMDEFAULT)); 1779566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(C, PETSC_DEFAULT)); 1789566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(C)); 1799566063dSJacob Faibussowitsch PetscCall(MatHasOperation(C, MATOP_PRODUCTSYMBOLIC, &flg)); 180263f2b91SStefano Zampini if (flg) { /* run tests if supported */ 1819566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(C)); /* equivalent to MatSetUp() */ 1829566063dSJacob Faibussowitsch PetscCall(MatSetOption(C, MAT_USE_INODES, PETSC_FALSE)); /* illustrate how to call MatSetOption() */ 1839566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 1849566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); /* test reuse symbolic C */ 18567b3012eSStefano Zampini 1869566063dSJacob Faibussowitsch PetscCall(MatTransposeMatMultEqual(P, B, C, mcheck, &flg)); 18728b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error: developer driver C = P^T*B"); 1889566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 189c20d7725SJed Brown 190c20d7725SJed Brown /* (2.2) Test user driver C = P^T*B */ 1919566063dSJacob Faibussowitsch PetscCall(MatTransposeMatMult(P, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 1929566063dSJacob Faibussowitsch PetscCall(MatTransposeMatMult(P, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C)); 1939566063dSJacob Faibussowitsch PetscCall(MatGetInfo(C, MAT_GLOBAL_SUM, &info)); 1949566063dSJacob Faibussowitsch PetscCall(MatProductClear(C)); 195c20d7725SJed Brown 196c20d7725SJed Brown /* Compare P^T*B and R*B */ 1979566063dSJacob Faibussowitsch PetscCall(MatMatMult(R, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C1)); 1989566063dSJacob Faibussowitsch PetscCall(MatNormDifference(C, C1, &norm)); 19908401ef6SPierre Jolivet PetscCheck(norm <= PETSC_SMALL, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatTransposeMatMult(): %g", (double)norm); 2009566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C1)); 201c20d7725SJed Brown 202c20d7725SJed Brown /* Test MatDuplicate() of C=P^T*B */ 2039566063dSJacob Faibussowitsch PetscCall(MatDuplicate(C, MAT_COPY_VALUES, &C1)); 2049566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C1)); 205263f2b91SStefano Zampini } else { 2069566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatTransposeMatMult not supported\n")); 207263f2b91SStefano Zampini } 2089566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 209c20d7725SJed Brown } 210c20d7725SJed Brown 21167b3012eSStefano Zampini /* 3) MatMatTransposeMult() */ 212c20d7725SJed Brown /* ------------------------ */ 213c20d7725SJed Brown if (Test_MatMatTr) { 214c20d7725SJed Brown /* C = B*R^T */ 2159566063dSJacob Faibussowitsch PetscCall(PetscObjectBaseTypeCompare((PetscObject)B, MATSEQAIJ, &seqaij)); 21620b3374bSStefano Zampini if (seqaij) { 2179566063dSJacob Faibussowitsch PetscCall(MatMatTransposeMult(B, R, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 2189566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C, "ABt_")); /* enable '-ABt_' for matrix C */ 2199566063dSJacob Faibussowitsch PetscCall(MatGetInfo(C, MAT_GLOBAL_SUM, &info)); 220c20d7725SJed Brown 221c20d7725SJed Brown /* Test MAT_REUSE_MATRIX - reuse symbolic C */ 2229566063dSJacob Faibussowitsch PetscCall(MatMatTransposeMult(B, R, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C)); 223c20d7725SJed Brown 224c20d7725SJed Brown /* Check */ 2259566063dSJacob Faibussowitsch PetscCall(MatMatMult(B, P, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C1)); 2269566063dSJacob Faibussowitsch PetscCall(MatNormDifference(C, C1, &norm)); 22708401ef6SPierre Jolivet PetscCheck(norm <= PETSC_SMALL, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatMatTransposeMult() %g", (double)norm); 2289566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C1)); 2299566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 230c20d7725SJed Brown } 231c20d7725SJed Brown } 232c20d7725SJed Brown 233c20d7725SJed Brown /* 4) Test MatPtAP() */ 234c20d7725SJed Brown /*-------------------*/ 235c20d7725SJed Brown if (Test_MatPtAP) { 2369566063dSJacob Faibussowitsch PetscCall(MatDuplicate(A_save, MAT_COPY_VALUES, &A)); 237c20d7725SJed Brown 238c20d7725SJed Brown /* (4.1) Test developer API */ 2399566063dSJacob Faibussowitsch PetscCall(MatProductCreate(A, P, NULL, &C)); 2409566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C, "PtAP_")); 2419566063dSJacob Faibussowitsch PetscCall(MatProductSetType(C, MATPRODUCT_PtAP)); 2429566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(C, MATPRODUCTALGORITHMDEFAULT)); 2439566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(C, PETSC_DEFAULT)); 2449566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(C)); 2459566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(C)); 2469566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 2479566063dSJacob Faibussowitsch PetscCall(MatPtAPMultEqual(A, P, C, mcheck, &flg)); 24828b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatProduct_PtAP"); 2499566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); /* reuse symbolic C */ 250c20d7725SJed Brown 2519566063dSJacob Faibussowitsch PetscCall(MatPtAPMultEqual(A, P, C, mcheck, &flg)); 25228b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatProduct_PtAP"); 2539566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 254c20d7725SJed Brown 255c20d7725SJed Brown /* (4.2) Test user driver */ 2569566063dSJacob Faibussowitsch PetscCall(MatPtAP(A, P, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 257c20d7725SJed Brown 258c20d7725SJed Brown /* Test MAT_REUSE_MATRIX - reuse symbolic C */ 259c20d7725SJed Brown alpha = 1.0; 260c20d7725SJed Brown for (i = 0; i < 2; i++) { 261c20d7725SJed Brown alpha -= 0.1; 2629566063dSJacob Faibussowitsch PetscCall(MatScale(A, alpha)); 2639566063dSJacob Faibussowitsch PetscCall(MatPtAP(A, P, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C)); 264c20d7725SJed Brown } 2659566063dSJacob Faibussowitsch PetscCall(MatPtAPMultEqual(A, P, C, mcheck, &flg)); 26628b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatPtAP"); 267c20d7725SJed Brown 268c20d7725SJed Brown /* 5) Test MatRARt() */ 269c20d7725SJed Brown /* ----------------- */ 270c20d7725SJed Brown if (Test_MatRARt) { 271c20d7725SJed Brown Mat RARt; 27267b3012eSStefano Zampini 27367b3012eSStefano Zampini /* (5.1) Test developer driver RARt = R*A*Rt */ 2749566063dSJacob Faibussowitsch PetscCall(MatProductCreate(A, R, NULL, &RARt)); 2759566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(RARt, "RARt_")); 2769566063dSJacob Faibussowitsch PetscCall(MatProductSetType(RARt, MATPRODUCT_RARt)); 2779566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(RARt, MATPRODUCTALGORITHMDEFAULT)); 2789566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(RARt, PETSC_DEFAULT)); 2799566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(RARt)); 2809566063dSJacob Faibussowitsch PetscCall(MatHasOperation(RARt, MATOP_PRODUCTSYMBOLIC, &flg)); 281263f2b91SStefano Zampini if (flg) { 2829566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(RARt)); /* equivalent to MatSetUp() */ 2839566063dSJacob Faibussowitsch PetscCall(MatSetOption(RARt, MAT_USE_INODES, PETSC_FALSE)); /* illustrate how to call MatSetOption() */ 2849566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(RARt)); 2859566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(RARt)); /* test reuse symbolic RARt */ 2869566063dSJacob Faibussowitsch PetscCall(MatDestroy(&RARt)); 28767b3012eSStefano Zampini 28867b3012eSStefano Zampini /* (2.2) Test user driver RARt = R*A*Rt */ 2899566063dSJacob Faibussowitsch PetscCall(MatRARt(A, R, MAT_INITIAL_MATRIX, 2.0, &RARt)); 2909566063dSJacob Faibussowitsch PetscCall(MatRARt(A, R, MAT_REUSE_MATRIX, 2.0, &RARt)); 29167b3012eSStefano Zampini 2929566063dSJacob Faibussowitsch PetscCall(MatNormDifference(C, RARt, &norm)); 29308401ef6SPierre Jolivet PetscCheck(norm <= PETSC_SMALL, PETSC_COMM_SELF, PETSC_ERR_PLIB, "|PtAP - RARt| = %g", (double)norm); 294263f2b91SStefano Zampini } else { 2959566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatRARt not supported\n")); 296263f2b91SStefano Zampini } 2979566063dSJacob Faibussowitsch PetscCall(MatDestroy(&RARt)); 298c20d7725SJed Brown } 299c20d7725SJed Brown 3009566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 3019566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 302c20d7725SJed Brown } 303c20d7725SJed Brown 304c20d7725SJed Brown /* Destroy objects */ 3059566063dSJacob Faibussowitsch PetscCall(PetscRandomDestroy(&rdm)); 3069566063dSJacob Faibussowitsch PetscCall(PetscFree(idxn)); 307c20d7725SJed Brown 3089566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A_save)); 3099566063dSJacob Faibussowitsch PetscCall(MatDestroy(&B)); 3109566063dSJacob Faibussowitsch PetscCall(MatDestroy(&P)); 3119566063dSJacob Faibussowitsch PetscCall(MatDestroy(&R)); 312c20d7725SJed Brown 3139566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 314b122ec5aSJacob Faibussowitsch return 0; 315c20d7725SJed Brown } 316c20d7725SJed Brown 317c20d7725SJed Brown /*TEST 318c20d7725SJed Brown test: 319c20d7725SJed Brown suffix: 1 320dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 321c20d7725SJed Brown args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium 322c20d7725SJed Brown output_file: output/ex62_1.out 323c20d7725SJed Brown 324c20d7725SJed Brown test: 325c20d7725SJed Brown suffix: 2_ab_scalable 326dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3273e662e0bSHong 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 328c20d7725SJed Brown output_file: output/ex62_1.out 329c20d7725SJed Brown 330c20d7725SJed Brown test: 331c20d7725SJed Brown suffix: 3_ab_scalable_fast 332dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3333e662e0bSHong Zhang args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm scalable_fast -matmatmult_via scalable_fast -matmattransmult_via color 334c20d7725SJed Brown output_file: output/ex62_1.out 335c20d7725SJed Brown 336c20d7725SJed Brown test: 337c20d7725SJed Brown suffix: 4_ab_heap 338dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3393e662e0bSHong 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 340c20d7725SJed Brown output_file: output/ex62_1.out 341c20d7725SJed Brown 342c20d7725SJed Brown test: 343c20d7725SJed Brown suffix: 5_ab_btheap 344dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3453e662e0bSHong Zhang args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm btheap -matmatmult_via btheap -matrart_via r*art 346c20d7725SJed Brown output_file: output/ex62_1.out 347c20d7725SJed Brown 348c20d7725SJed Brown test: 349c20d7725SJed Brown suffix: 6_ab_llcondensed 350dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3513e662e0bSHong Zhang args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm llcondensed -matmatmult_via llcondensed -matrart_via coloring_rart 352c20d7725SJed Brown output_file: output/ex62_1.out 353c20d7725SJed Brown 354c20d7725SJed Brown test: 355c20d7725SJed Brown suffix: 7_ab_rowmerge 356dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3573e662e0bSHong Zhang args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -AB_mat_product_algorithm rowmerge -matmatmult_via rowmerge 358c20d7725SJed Brown output_file: output/ex62_1.out 359c20d7725SJed Brown 360c20d7725SJed Brown test: 361c20d7725SJed Brown suffix: 8_ab_hypre 362dfd57a17SPierre Jolivet requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3633e662e0bSHong 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 364c20d7725SJed Brown output_file: output/ex62_1.out 365c20d7725SJed Brown 366c20d7725SJed Brown test: 367263f2b91SStefano Zampini suffix: hypre_medium 368263f2b91SStefano Zampini nsize: {{1 3}} 369263f2b91SStefano Zampini requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 370263f2b91SStefano Zampini args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -A_mat_type hypre -B_mat_type hypre -test_rart 0 371263f2b91SStefano Zampini output_file: output/ex62_hypre.out 372263f2b91SStefano Zampini 373263f2b91SStefano Zampini test: 374263f2b91SStefano Zampini suffix: hypre_tiny 375263f2b91SStefano Zampini nsize: {{1 3}} 376263f2b91SStefano Zampini requires: hypre !complex double !defined(PETSC_USE_64BIT_INDICES) 377263f2b91SStefano 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 378263f2b91SStefano Zampini output_file: output/ex62_hypre.out 379263f2b91SStefano Zampini 380263f2b91SStefano Zampini test: 381cec0a6c6SStefano Zampini suffix: 9_mkl 38215df37b9SStefano Zampini TODO: broken MatScale? 383b88df2e7SBarry Smith requires: mkl_sparse datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 384cec0a6c6SStefano Zampini args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -A_mat_type aijmkl -B_mat_type aijmkl 385cec0a6c6SStefano Zampini output_file: output/ex62_1.out 386cec0a6c6SStefano Zampini 387cec0a6c6SStefano Zampini test: 388c20d7725SJed Brown suffix: 10 389dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 390c20d7725SJed Brown nsize: 3 391c20d7725SJed Brown args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium 392c20d7725SJed Brown output_file: output/ex62_1.out 393c20d7725SJed Brown 394c20d7725SJed Brown test: 3953cea4f0aSStefano Zampini suffix: 10_backend 396dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3973cea4f0aSStefano Zampini nsize: 3 3983e662e0bSHong 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 3993cea4f0aSStefano Zampini output_file: output/ex62_1.out 4003cea4f0aSStefano Zampini 4013cea4f0aSStefano Zampini test: 402c20d7725SJed Brown suffix: 11_ab_scalable 403dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 404c20d7725SJed Brown nsize: 3 4053e662e0bSHong 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 406c20d7725SJed Brown output_file: output/ex62_1.out 407c20d7725SJed Brown 408c20d7725SJed Brown test: 409c20d7725SJed Brown suffix: 12_ab_seqmpi 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 seqmpi -matmatmult_via seqmpi -AtB_mat_product_algorithm at*b -mattransposematmult_via at*b 413c20d7725SJed Brown output_file: output/ex62_1.out 414c20d7725SJed Brown 415c20d7725SJed Brown test: 416c20d7725SJed Brown suffix: 13_ab_hypre 417dfd57a17SPierre Jolivet requires: hypre 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 hypre -matmatmult_via hypre -PtAP_mat_product_algorithm hypre -matptap_via hypre 420c20d7725SJed Brown output_file: output/ex62_1.out 421c20d7725SJed Brown 42220b3374bSStefano Zampini test: 42320b3374bSStefano Zampini suffix: 14_seqaij 424dfd57a17SPierre Jolivet requires: !complex double !defined(PETSC_USE_64BIT_INDICES) 42520b3374bSStefano Zampini args: -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system 42620b3374bSStefano Zampini output_file: output/ex62_1.out 42720b3374bSStefano Zampini 42820b3374bSStefano Zampini test: 42920b3374bSStefano Zampini suffix: 14_seqaijcusparse 430dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4311a2c6b5cSJunchao 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 43220b3374bSStefano Zampini output_file: output/ex62_1.out 43320b3374bSStefano Zampini 43420b3374bSStefano Zampini test: 43565e4b4d4SStefano Zampini suffix: 14_seqaijcusparse_cpu 436dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4373e662e0bSHong 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 43865e4b4d4SStefano Zampini output_file: output/ex62_1.out 43965e4b4d4SStefano Zampini 44065e4b4d4SStefano Zampini test: 4415cb69cabSStefano Zampini suffix: 14_mpiaijcusparse_seq 4425cb69cabSStefano Zampini nsize: 1 443dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4441a2c6b5cSJunchao 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 4455cb69cabSStefano Zampini output_file: output/ex62_1.out 4465cb69cabSStefano Zampini 4475cb69cabSStefano Zampini test: 44865e4b4d4SStefano Zampini suffix: 14_mpiaijcusparse_seq_cpu 44965e4b4d4SStefano Zampini nsize: 1 450dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4514005f4d6SPierre 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 45265e4b4d4SStefano Zampini output_file: output/ex62_1.out 45365e4b4d4SStefano Zampini 45465e4b4d4SStefano Zampini test: 4555cb69cabSStefano Zampini suffix: 14_mpiaijcusparse 4565cb69cabSStefano Zampini nsize: 3 457dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4581a2c6b5cSJunchao 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 4595cb69cabSStefano Zampini output_file: output/ex62_1.out 4605cb69cabSStefano Zampini 4615cb69cabSStefano Zampini test: 46265e4b4d4SStefano Zampini suffix: 14_mpiaijcusparse_cpu 46365e4b4d4SStefano Zampini nsize: 3 464dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4654005f4d6SPierre 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 46665e4b4d4SStefano Zampini output_file: output/ex62_1.out 46765e4b4d4SStefano Zampini 46865e4b4d4SStefano Zampini test: 469076ba34aSJunchao Zhang nsize: {{1 3}} 470076ba34aSJunchao Zhang suffix: 14_aijkokkos 471*dcfd994dSJunchao Zhang requires: kokkos_kernels !complex double !defined(PETSC_USE_64BIT_INDICES) 47220b3374bSStefano 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 47320b3374bSStefano Zampini output_file: output/ex62_1.out 47420b3374bSStefano Zampini 4755cb69cabSStefano Zampini # these tests use matrices with many zero rows 47620b3374bSStefano Zampini test: 47720b3374bSStefano Zampini suffix: 15_seqaijcusparse 478dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath 4791a2c6b5cSJunchao Zhang args: -A_mat_type aijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith 48020b3374bSStefano Zampini output_file: output/ex62_1.out 4815cb69cabSStefano Zampini 4825cb69cabSStefano Zampini test: 4835cb69cabSStefano Zampini suffix: 15_mpiaijcusparse_seq 4845cb69cabSStefano Zampini nsize: 1 485dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath 4861a2c6b5cSJunchao Zhang args: -A_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith 4875cb69cabSStefano Zampini output_file: output/ex62_1.out 4885cb69cabSStefano Zampini 4895cb69cabSStefano Zampini test: 4905cb69cabSStefano Zampini nsize: 3 4915cb69cabSStefano Zampini suffix: 15_mpiaijcusparse 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 496c20d7725SJed Brown TEST*/ 497