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 */ 16*d71ae5a4SJacob Faibussowitsch PetscErrorCode MatNormDifference(Mat A, Mat B, PetscReal *norm) 17*d71ae5a4SJacob Faibussowitsch { 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 24*d71ae5a4SJacob Faibussowitsch int main(int argc, char **args) 25*d71ae5a4SJacob Faibussowitsch { 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 43327415f7SBarry Smith PetscFunctionBeginUser; 449566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 459566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 469566063dSJacob Faibussowitsch PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 47c20d7725SJed Brown 48c20d7725SJed Brown /* Load the matrices A_save and B */ 49d0609cedSBarry Smith PetscOptionsBegin(PETSC_COMM_WORLD, "", "", ""); 509566063dSJacob Faibussowitsch PetscCall(PetscOptionsBool("-test_rart", "Test MatRARt", "", Test_MatRARt, &Test_MatRARt, NULL)); 519566063dSJacob Faibussowitsch PetscCall(PetscOptionsInt("-PN", "Number of columns of P", "", PN, &PN, NULL)); 529566063dSJacob Faibussowitsch PetscCall(PetscOptionsInt("-mcheck", "Number of matmult checks", "", mcheck, &mcheck, NULL)); 539566063dSJacob Faibussowitsch PetscCall(PetscOptionsString("-fA", "Path for matrix A", "", file[0], file[0], sizeof(file[0]), &flg)); 5428b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_USER_INPUT, "Must indicate a file name for matrix A with the -fA option."); 559566063dSJacob Faibussowitsch PetscCall(PetscOptionsString("-fB", "Path for matrix B", "", file[1], file[1], sizeof(file[1]), &flg)); 569566063dSJacob Faibussowitsch PetscCall(PetscOptionsFList("-A_mat_type", "Matrix type", "MatSetType", MatList, deft, A_mattype, 256, &flgA)); 579566063dSJacob Faibussowitsch PetscCall(PetscOptionsFList("-B_mat_type", "Matrix type", "MatSetType", MatList, deft, B_mattype, 256, &flgB)); 58d0609cedSBarry Smith PetscOptionsEnd(); 59c20d7725SJed Brown 609566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file[0], FILE_MODE_READ, &viewer)); 619566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD, &A_save)); 629566063dSJacob Faibussowitsch PetscCall(MatLoad(A_save, viewer)); 639566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&viewer)); 64c20d7725SJed Brown 6520b3374bSStefano Zampini if (flg) { 669566063dSJacob Faibussowitsch PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file[1], FILE_MODE_READ, &viewer)); 679566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD, &B)); 689566063dSJacob Faibussowitsch PetscCall(MatLoad(B, viewer)); 699566063dSJacob Faibussowitsch PetscCall(PetscViewerDestroy(&viewer)); 7020b3374bSStefano Zampini } else { 719566063dSJacob Faibussowitsch PetscCall(PetscObjectReference((PetscObject)A_save)); 7220b3374bSStefano Zampini B = A_save; 7320b3374bSStefano Zampini } 7420b3374bSStefano Zampini 7548a46eb9SPierre Jolivet if (flgA) PetscCall(MatConvert(A_save, A_mattype, MAT_INPLACE_MATRIX, &A_save)); 7648a46eb9SPierre Jolivet if (flgB) PetscCall(MatConvert(B, B_mattype, MAT_INPLACE_MATRIX, &B)); 779566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(A_save)); 789566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(B)); 79c20d7725SJed Brown 809566063dSJacob Faibussowitsch PetscCall(MatGetType(B, &mattype)); 81c20d7725SJed Brown 829566063dSJacob Faibussowitsch PetscCall(PetscMalloc(nzp * (sizeof(PetscInt) + sizeof(PetscScalar)), &idxn)); 83c20d7725SJed Brown a = (PetscScalar *)(idxn + nzp); 84c20d7725SJed Brown 859566063dSJacob Faibussowitsch PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rdm)); 869566063dSJacob Faibussowitsch PetscCall(PetscRandomSetFromOptions(rdm)); 87c20d7725SJed Brown 88c20d7725SJed Brown /* 1) MatMatMult() */ 89c20d7725SJed Brown /* ----------------*/ 90c20d7725SJed Brown if (Test_MatMatMult) { 919566063dSJacob Faibussowitsch PetscCall(MatDuplicate(A_save, MAT_COPY_VALUES, &A)); 92c20d7725SJed Brown 93c20d7725SJed Brown /* (1.1) Test developer API */ 949566063dSJacob Faibussowitsch PetscCall(MatProductCreate(A, B, NULL, &C)); 959566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C, "AB_")); 969566063dSJacob Faibussowitsch PetscCall(MatProductSetType(C, MATPRODUCT_AB)); 979566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(C, MATPRODUCTALGORITHMDEFAULT)); 989566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(C, PETSC_DEFAULT)); 999566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(C)); 100910fa13fSStefano Zampini /* we can inquire about MATOP_PRODUCTSYMBOLIC even if the destination matrix type has not been set yet */ 1019566063dSJacob Faibussowitsch PetscCall(MatHasOperation(C, MATOP_PRODUCTSYMBOLIC, &flg)); 1029566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(C)); 1039566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 1049566063dSJacob Faibussowitsch PetscCall(MatMatMultEqual(A, B, C, mcheck, &flg)); 10528b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error in C=A*B"); 106c20d7725SJed Brown 107c20d7725SJed Brown /* Test reuse symbolic C */ 108c20d7725SJed Brown alpha = 0.9; 1099566063dSJacob Faibussowitsch PetscCall(MatScale(A, alpha)); 1109566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 111c20d7725SJed Brown 1129566063dSJacob Faibussowitsch PetscCall(MatMatMultEqual(A, B, C, mcheck, &flg)); 11328b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error in C=A*B"); 1149566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 115c20d7725SJed Brown 116c20d7725SJed Brown /* (1.2) Test user driver */ 1179566063dSJacob Faibussowitsch PetscCall(MatMatMult(A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 118c20d7725SJed Brown 119c20d7725SJed Brown /* Test MAT_REUSE_MATRIX - reuse symbolic C */ 120c20d7725SJed Brown alpha = 1.0; 121c20d7725SJed Brown for (i = 0; i < 2; i++) { 122c20d7725SJed Brown alpha -= 0.1; 1239566063dSJacob Faibussowitsch PetscCall(MatScale(A, alpha)); 1249566063dSJacob Faibussowitsch PetscCall(MatMatMult(A, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C)); 125c20d7725SJed Brown } 1269566063dSJacob Faibussowitsch PetscCall(MatMatMultEqual(A, B, C, mcheck, &flg)); 12728b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error: MatMatMult()"); 1289566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 1294417c5e8SHong Zhang 1304417c5e8SHong Zhang /* Test MatProductClear() */ 1319566063dSJacob Faibussowitsch PetscCall(MatProductClear(C)); 1329566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 133544a5e07SHong Zhang 134544a5e07SHong Zhang /* Test MatMatMult() for dense and aij matrices */ 1359566063dSJacob Faibussowitsch PetscCall(PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJ, MATMPIAIJ, "")); 13620b3374bSStefano Zampini if (flg) { 1379566063dSJacob Faibussowitsch PetscCall(MatConvert(A_save, MATDENSE, MAT_INITIAL_MATRIX, &A)); 1389566063dSJacob Faibussowitsch PetscCall(MatMatMult(A, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 1399566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 1409566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 14120b3374bSStefano Zampini } 142c20d7725SJed Brown } 143c20d7725SJed Brown 144c20d7725SJed Brown /* Create P and R = P^T */ 145c20d7725SJed Brown /* --------------------- */ 1469566063dSJacob Faibussowitsch PetscCall(MatGetSize(B, &PM, NULL)); 14720b3374bSStefano Zampini if (PN < 0) PN = PM / 2; 1489566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_WORLD, &P)); 1499566063dSJacob Faibussowitsch PetscCall(MatSetSizes(P, PETSC_DECIDE, PETSC_DECIDE, PM, PN)); 1509566063dSJacob Faibussowitsch PetscCall(MatSetType(P, MATAIJ)); 1519566063dSJacob Faibussowitsch PetscCall(MatSeqAIJSetPreallocation(P, nzp, NULL)); 1529566063dSJacob Faibussowitsch PetscCall(MatMPIAIJSetPreallocation(P, nzp, NULL, nzp, NULL)); 1539566063dSJacob Faibussowitsch PetscCall(MatGetOwnershipRange(P, &rstart, &rend)); 15448a46eb9SPierre Jolivet for (i = 0; i < nzp; i++) PetscCall(PetscRandomGetValue(rdm, &a[i])); 155c20d7725SJed Brown for (i = rstart; i < rend; i++) { 156c20d7725SJed Brown for (j = 0; j < nzp; j++) { 1579566063dSJacob Faibussowitsch PetscCall(PetscRandomGetValue(rdm, &rval)); 158c20d7725SJed Brown idxn[j] = (PetscInt)(PetscRealPart(rval) * PN); 159c20d7725SJed Brown } 1609566063dSJacob Faibussowitsch PetscCall(MatSetValues(P, 1, &i, nzp, idxn, a, ADD_VALUES)); 161c20d7725SJed Brown } 1629566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(P, MAT_FINAL_ASSEMBLY)); 1639566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(P, MAT_FINAL_ASSEMBLY)); 164c20d7725SJed Brown 1659566063dSJacob Faibussowitsch PetscCall(MatTranspose(P, MAT_INITIAL_MATRIX, &R)); 1669566063dSJacob Faibussowitsch PetscCall(MatConvert(P, mattype, MAT_INPLACE_MATRIX, &P)); 1679566063dSJacob Faibussowitsch PetscCall(MatConvert(R, mattype, MAT_INPLACE_MATRIX, &R)); 1689566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(P)); 1699566063dSJacob Faibussowitsch PetscCall(MatSetFromOptions(R)); 170c20d7725SJed Brown 171c20d7725SJed Brown /* 2) MatTransposeMatMult() */ 172c20d7725SJed Brown /* ------------------------ */ 173c20d7725SJed Brown if (Test_MatTrMat) { 174c20d7725SJed Brown /* (2.1) Test developer driver C = P^T*B */ 1759566063dSJacob Faibussowitsch PetscCall(MatProductCreate(P, B, NULL, &C)); 1769566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C, "AtB_")); 1779566063dSJacob Faibussowitsch PetscCall(MatProductSetType(C, MATPRODUCT_AtB)); 1789566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(C, MATPRODUCTALGORITHMDEFAULT)); 1799566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(C, PETSC_DEFAULT)); 1809566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(C)); 1819566063dSJacob Faibussowitsch PetscCall(MatHasOperation(C, MATOP_PRODUCTSYMBOLIC, &flg)); 182263f2b91SStefano Zampini if (flg) { /* run tests if supported */ 1839566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(C)); /* equivalent to MatSetUp() */ 1849566063dSJacob Faibussowitsch PetscCall(MatSetOption(C, MAT_USE_INODES, PETSC_FALSE)); /* illustrate how to call MatSetOption() */ 1859566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 1869566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); /* test reuse symbolic C */ 18767b3012eSStefano Zampini 1889566063dSJacob Faibussowitsch PetscCall(MatTransposeMatMultEqual(P, B, C, mcheck, &flg)); 18928b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error: developer driver C = P^T*B"); 1909566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 191c20d7725SJed Brown 192c20d7725SJed Brown /* (2.2) Test user driver C = P^T*B */ 1939566063dSJacob Faibussowitsch PetscCall(MatTransposeMatMult(P, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 1949566063dSJacob Faibussowitsch PetscCall(MatTransposeMatMult(P, B, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C)); 1959566063dSJacob Faibussowitsch PetscCall(MatGetInfo(C, MAT_GLOBAL_SUM, &info)); 1969566063dSJacob Faibussowitsch PetscCall(MatProductClear(C)); 197c20d7725SJed Brown 198c20d7725SJed Brown /* Compare P^T*B and R*B */ 1999566063dSJacob Faibussowitsch PetscCall(MatMatMult(R, B, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C1)); 2009566063dSJacob Faibussowitsch PetscCall(MatNormDifference(C, C1, &norm)); 20108401ef6SPierre Jolivet PetscCheck(norm <= PETSC_SMALL, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatTransposeMatMult(): %g", (double)norm); 2029566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C1)); 203c20d7725SJed Brown 204c20d7725SJed Brown /* Test MatDuplicate() of C=P^T*B */ 2059566063dSJacob Faibussowitsch PetscCall(MatDuplicate(C, MAT_COPY_VALUES, &C1)); 2069566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C1)); 207263f2b91SStefano Zampini } else { 2089566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatTransposeMatMult not supported\n")); 209263f2b91SStefano Zampini } 2109566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 211c20d7725SJed Brown } 212c20d7725SJed Brown 21367b3012eSStefano Zampini /* 3) MatMatTransposeMult() */ 214c20d7725SJed Brown /* ------------------------ */ 215c20d7725SJed Brown if (Test_MatMatTr) { 216c20d7725SJed Brown /* C = B*R^T */ 2179566063dSJacob Faibussowitsch PetscCall(PetscObjectBaseTypeCompare((PetscObject)B, MATSEQAIJ, &seqaij)); 21820b3374bSStefano Zampini if (seqaij) { 2199566063dSJacob Faibussowitsch PetscCall(MatMatTransposeMult(B, R, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 2209566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C, "ABt_")); /* enable '-ABt_' for matrix C */ 2219566063dSJacob Faibussowitsch PetscCall(MatGetInfo(C, MAT_GLOBAL_SUM, &info)); 222c20d7725SJed Brown 223c20d7725SJed Brown /* Test MAT_REUSE_MATRIX - reuse symbolic C */ 2249566063dSJacob Faibussowitsch PetscCall(MatMatTransposeMult(B, R, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C)); 225c20d7725SJed Brown 226c20d7725SJed Brown /* Check */ 2279566063dSJacob Faibussowitsch PetscCall(MatMatMult(B, P, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C1)); 2289566063dSJacob Faibussowitsch PetscCall(MatNormDifference(C, C1, &norm)); 22908401ef6SPierre Jolivet PetscCheck(norm <= PETSC_SMALL, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatMatTransposeMult() %g", (double)norm); 2309566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C1)); 2319566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 232c20d7725SJed Brown } 233c20d7725SJed Brown } 234c20d7725SJed Brown 235c20d7725SJed Brown /* 4) Test MatPtAP() */ 236c20d7725SJed Brown /*-------------------*/ 237c20d7725SJed Brown if (Test_MatPtAP) { 2389566063dSJacob Faibussowitsch PetscCall(MatDuplicate(A_save, MAT_COPY_VALUES, &A)); 239c20d7725SJed Brown 240c20d7725SJed Brown /* (4.1) Test developer API */ 2419566063dSJacob Faibussowitsch PetscCall(MatProductCreate(A, P, NULL, &C)); 2429566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(C, "PtAP_")); 2439566063dSJacob Faibussowitsch PetscCall(MatProductSetType(C, MATPRODUCT_PtAP)); 2449566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(C, MATPRODUCTALGORITHMDEFAULT)); 2459566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(C, PETSC_DEFAULT)); 2469566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(C)); 2479566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(C)); 2489566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); 2499566063dSJacob Faibussowitsch PetscCall(MatPtAPMultEqual(A, P, C, mcheck, &flg)); 25028b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatProduct_PtAP"); 2519566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(C)); /* reuse symbolic C */ 252c20d7725SJed Brown 2539566063dSJacob Faibussowitsch PetscCall(MatPtAPMultEqual(A, P, C, mcheck, &flg)); 25428b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatProduct_PtAP"); 2559566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 256c20d7725SJed Brown 257c20d7725SJed Brown /* (4.2) Test user driver */ 2589566063dSJacob Faibussowitsch PetscCall(MatPtAP(A, P, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 259c20d7725SJed Brown 260c20d7725SJed Brown /* Test MAT_REUSE_MATRIX - reuse symbolic C */ 261c20d7725SJed Brown alpha = 1.0; 262c20d7725SJed Brown for (i = 0; i < 2; i++) { 263c20d7725SJed Brown alpha -= 0.1; 2649566063dSJacob Faibussowitsch PetscCall(MatScale(A, alpha)); 2659566063dSJacob Faibussowitsch PetscCall(MatPtAP(A, P, MAT_REUSE_MATRIX, PETSC_DEFAULT, &C)); 266c20d7725SJed Brown } 2679566063dSJacob Faibussowitsch PetscCall(MatPtAPMultEqual(A, P, C, mcheck, &flg)); 26828b400f6SJacob Faibussowitsch PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in MatPtAP"); 269c20d7725SJed Brown 270c20d7725SJed Brown /* 5) Test MatRARt() */ 271c20d7725SJed Brown /* ----------------- */ 272c20d7725SJed Brown if (Test_MatRARt) { 273c20d7725SJed Brown Mat RARt; 27467b3012eSStefano Zampini 27567b3012eSStefano Zampini /* (5.1) Test developer driver RARt = R*A*Rt */ 2769566063dSJacob Faibussowitsch PetscCall(MatProductCreate(A, R, NULL, &RARt)); 2779566063dSJacob Faibussowitsch PetscCall(MatSetOptionsPrefix(RARt, "RARt_")); 2789566063dSJacob Faibussowitsch PetscCall(MatProductSetType(RARt, MATPRODUCT_RARt)); 2799566063dSJacob Faibussowitsch PetscCall(MatProductSetAlgorithm(RARt, MATPRODUCTALGORITHMDEFAULT)); 2809566063dSJacob Faibussowitsch PetscCall(MatProductSetFill(RARt, PETSC_DEFAULT)); 2819566063dSJacob Faibussowitsch PetscCall(MatProductSetFromOptions(RARt)); 2829566063dSJacob Faibussowitsch PetscCall(MatHasOperation(RARt, MATOP_PRODUCTSYMBOLIC, &flg)); 283263f2b91SStefano Zampini if (flg) { 2849566063dSJacob Faibussowitsch PetscCall(MatProductSymbolic(RARt)); /* equivalent to MatSetUp() */ 2859566063dSJacob Faibussowitsch PetscCall(MatSetOption(RARt, MAT_USE_INODES, PETSC_FALSE)); /* illustrate how to call MatSetOption() */ 2869566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(RARt)); 2879566063dSJacob Faibussowitsch PetscCall(MatProductNumeric(RARt)); /* test reuse symbolic RARt */ 2889566063dSJacob Faibussowitsch PetscCall(MatDestroy(&RARt)); 28967b3012eSStefano Zampini 29067b3012eSStefano Zampini /* (2.2) Test user driver RARt = R*A*Rt */ 2919566063dSJacob Faibussowitsch PetscCall(MatRARt(A, R, MAT_INITIAL_MATRIX, 2.0, &RARt)); 2929566063dSJacob Faibussowitsch PetscCall(MatRARt(A, R, MAT_REUSE_MATRIX, 2.0, &RARt)); 29367b3012eSStefano Zampini 2949566063dSJacob Faibussowitsch PetscCall(MatNormDifference(C, RARt, &norm)); 29508401ef6SPierre Jolivet PetscCheck(norm <= PETSC_SMALL, PETSC_COMM_SELF, PETSC_ERR_PLIB, "|PtAP - RARt| = %g", (double)norm); 296263f2b91SStefano Zampini } else { 2979566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatRARt not supported\n")); 298263f2b91SStefano Zampini } 2999566063dSJacob Faibussowitsch PetscCall(MatDestroy(&RARt)); 300c20d7725SJed Brown } 301c20d7725SJed Brown 3029566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A)); 3039566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 304c20d7725SJed Brown } 305c20d7725SJed Brown 306c20d7725SJed Brown /* Destroy objects */ 3079566063dSJacob Faibussowitsch PetscCall(PetscRandomDestroy(&rdm)); 3089566063dSJacob Faibussowitsch PetscCall(PetscFree(idxn)); 309c20d7725SJed Brown 3109566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A_save)); 3119566063dSJacob Faibussowitsch PetscCall(MatDestroy(&B)); 3129566063dSJacob Faibussowitsch PetscCall(MatDestroy(&P)); 3139566063dSJacob Faibussowitsch PetscCall(MatDestroy(&R)); 314c20d7725SJed Brown 3159566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 316b122ec5aSJacob Faibussowitsch return 0; 317c20d7725SJed Brown } 318c20d7725SJed Brown 319c20d7725SJed Brown /*TEST 320c20d7725SJed Brown test: 321c20d7725SJed Brown suffix: 1 322dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 323c20d7725SJed Brown args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium 324c20d7725SJed Brown output_file: output/ex62_1.out 325c20d7725SJed Brown 326c20d7725SJed Brown test: 327c20d7725SJed Brown suffix: 2_ab_scalable 328dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3293e662e0bSHong 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 330c20d7725SJed Brown output_file: output/ex62_1.out 331c20d7725SJed Brown 332c20d7725SJed Brown test: 333c20d7725SJed Brown suffix: 3_ab_scalable_fast 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_fast -matmatmult_via scalable_fast -matmattransmult_via color 336c20d7725SJed Brown output_file: output/ex62_1.out 337c20d7725SJed Brown 338c20d7725SJed Brown test: 339c20d7725SJed Brown suffix: 4_ab_heap 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 heap -matmatmult_via heap -PtAP_mat_product_algorithm rap -matptap_via rap 342c20d7725SJed Brown output_file: output/ex62_1.out 343c20d7725SJed Brown 344c20d7725SJed Brown test: 345c20d7725SJed Brown suffix: 5_ab_btheap 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 btheap -matmatmult_via btheap -matrart_via r*art 348c20d7725SJed Brown output_file: output/ex62_1.out 349c20d7725SJed Brown 350c20d7725SJed Brown test: 351c20d7725SJed Brown suffix: 6_ab_llcondensed 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 llcondensed -matmatmult_via llcondensed -matrart_via coloring_rart 354c20d7725SJed Brown output_file: output/ex62_1.out 355c20d7725SJed Brown 356c20d7725SJed Brown test: 357c20d7725SJed Brown suffix: 7_ab_rowmerge 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 rowmerge -matmatmult_via rowmerge 360c20d7725SJed Brown output_file: output/ex62_1.out 361c20d7725SJed Brown 362c20d7725SJed Brown test: 363c20d7725SJed Brown suffix: 8_ab_hypre 364dfd57a17SPierre Jolivet requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3653e662e0bSHong 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 366c20d7725SJed Brown output_file: output/ex62_1.out 367c20d7725SJed Brown 368c20d7725SJed Brown test: 369263f2b91SStefano Zampini suffix: hypre_medium 370263f2b91SStefano Zampini nsize: {{1 3}} 371263f2b91SStefano Zampini requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 372263f2b91SStefano Zampini args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -A_mat_type hypre -B_mat_type hypre -test_rart 0 373263f2b91SStefano Zampini output_file: output/ex62_hypre.out 374263f2b91SStefano Zampini 375263f2b91SStefano Zampini test: 376263f2b91SStefano Zampini suffix: hypre_tiny 377263f2b91SStefano Zampini nsize: {{1 3}} 378263f2b91SStefano Zampini requires: hypre !complex double !defined(PETSC_USE_64BIT_INDICES) 379263f2b91SStefano 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 380263f2b91SStefano Zampini output_file: output/ex62_hypre.out 381263f2b91SStefano Zampini 382263f2b91SStefano Zampini test: 383cec0a6c6SStefano Zampini suffix: 9_mkl 38415df37b9SStefano Zampini TODO: broken MatScale? 385b88df2e7SBarry Smith requires: mkl_sparse datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 386cec0a6c6SStefano Zampini args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium -A_mat_type aijmkl -B_mat_type aijmkl 387cec0a6c6SStefano Zampini output_file: output/ex62_1.out 388cec0a6c6SStefano Zampini 389cec0a6c6SStefano Zampini test: 390c20d7725SJed Brown suffix: 10 391dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 392c20d7725SJed Brown nsize: 3 393c20d7725SJed Brown args: -fA ${DATAFILESPATH}/matrices/medium -fB ${DATAFILESPATH}/matrices/medium 394c20d7725SJed Brown output_file: output/ex62_1.out 395c20d7725SJed Brown 396c20d7725SJed Brown test: 3973cea4f0aSStefano Zampini suffix: 10_backend 398dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 3993cea4f0aSStefano Zampini nsize: 3 4003e662e0bSHong 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 4013cea4f0aSStefano Zampini output_file: output/ex62_1.out 4023cea4f0aSStefano Zampini 4033cea4f0aSStefano Zampini test: 404c20d7725SJed Brown suffix: 11_ab_scalable 405dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 406c20d7725SJed Brown nsize: 3 4073e662e0bSHong 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 408c20d7725SJed Brown output_file: output/ex62_1.out 409c20d7725SJed Brown 410c20d7725SJed Brown test: 411c20d7725SJed Brown suffix: 12_ab_seqmpi 412dfd57a17SPierre Jolivet requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 413c20d7725SJed Brown nsize: 3 4143e662e0bSHong 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 415c20d7725SJed Brown output_file: output/ex62_1.out 416c20d7725SJed Brown 417c20d7725SJed Brown test: 418c20d7725SJed Brown suffix: 13_ab_hypre 419dfd57a17SPierre Jolivet requires: hypre datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES) 420c20d7725SJed Brown nsize: 3 4213e662e0bSHong 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 422c20d7725SJed Brown output_file: output/ex62_1.out 423c20d7725SJed Brown 42420b3374bSStefano Zampini test: 42520b3374bSStefano Zampini suffix: 14_seqaij 426dfd57a17SPierre Jolivet requires: !complex double !defined(PETSC_USE_64BIT_INDICES) 42720b3374bSStefano Zampini args: -fA ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system -fB ${wPETSC_DIR}/share/petsc/datafiles/matrices/tiny_system 42820b3374bSStefano Zampini output_file: output/ex62_1.out 42920b3374bSStefano Zampini 43020b3374bSStefano Zampini test: 43120b3374bSStefano Zampini suffix: 14_seqaijcusparse 432dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4331a2c6b5cSJunchao 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 43420b3374bSStefano Zampini output_file: output/ex62_1.out 43520b3374bSStefano Zampini 43620b3374bSStefano Zampini test: 43765e4b4d4SStefano Zampini suffix: 14_seqaijcusparse_cpu 438dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4393e662e0bSHong 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 44065e4b4d4SStefano Zampini output_file: output/ex62_1.out 44165e4b4d4SStefano Zampini 44265e4b4d4SStefano Zampini test: 4435cb69cabSStefano Zampini suffix: 14_mpiaijcusparse_seq 4445cb69cabSStefano Zampini nsize: 1 445dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4461a2c6b5cSJunchao 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 4475cb69cabSStefano Zampini output_file: output/ex62_1.out 4485cb69cabSStefano Zampini 4495cb69cabSStefano Zampini test: 45065e4b4d4SStefano Zampini suffix: 14_mpiaijcusparse_seq_cpu 45165e4b4d4SStefano Zampini nsize: 1 452dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4534005f4d6SPierre 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 45465e4b4d4SStefano Zampini output_file: output/ex62_1.out 45565e4b4d4SStefano Zampini 45665e4b4d4SStefano Zampini test: 4575cb69cabSStefano Zampini suffix: 14_mpiaijcusparse 4585cb69cabSStefano Zampini nsize: 3 459dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4601a2c6b5cSJunchao 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 4615cb69cabSStefano Zampini output_file: output/ex62_1.out 4625cb69cabSStefano Zampini 4635cb69cabSStefano Zampini test: 46465e4b4d4SStefano Zampini suffix: 14_mpiaijcusparse_cpu 46565e4b4d4SStefano Zampini nsize: 3 466dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) 4674005f4d6SPierre 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 46865e4b4d4SStefano Zampini output_file: output/ex62_1.out 46965e4b4d4SStefano Zampini 47065e4b4d4SStefano Zampini test: 471076ba34aSJunchao Zhang nsize: {{1 3}} 472076ba34aSJunchao Zhang suffix: 14_aijkokkos 473dcfd994dSJunchao Zhang requires: kokkos_kernels !complex double !defined(PETSC_USE_64BIT_INDICES) 47420b3374bSStefano 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 47520b3374bSStefano Zampini output_file: output/ex62_1.out 47620b3374bSStefano Zampini 4775cb69cabSStefano Zampini # these tests use matrices with many zero rows 47820b3374bSStefano Zampini test: 47920b3374bSStefano Zampini suffix: 15_seqaijcusparse 480dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath 4811a2c6b5cSJunchao Zhang args: -A_mat_type aijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith 48220b3374bSStefano Zampini output_file: output/ex62_1.out 4835cb69cabSStefano Zampini 4845cb69cabSStefano Zampini test: 4855cb69cabSStefano Zampini suffix: 15_mpiaijcusparse_seq 4865cb69cabSStefano Zampini nsize: 1 487dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath 4881a2c6b5cSJunchao Zhang args: -A_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith 4895cb69cabSStefano Zampini output_file: output/ex62_1.out 4905cb69cabSStefano Zampini 4915cb69cabSStefano Zampini test: 4925cb69cabSStefano Zampini nsize: 3 4935cb69cabSStefano Zampini suffix: 15_mpiaijcusparse 494dfd57a17SPierre Jolivet requires: cuda !complex double !defined(PETSC_USE_64BIT_INDICES) datafilespath 4951a2c6b5cSJunchao Zhang args: -A_mat_type mpiaijcusparse -mat_form_explicit_transpose -fA ${DATAFILESPATH}/matrices/matmatmult/A4.BGriffith 4965cb69cabSStefano Zampini output_file: output/ex62_1.out 4975cb69cabSStefano Zampini 498c20d7725SJed Brown TEST*/ 499