1 static char help[] = "Tests MatPermute() in parallel.\n\n"; 2 /* Results: 3 Sequential: 4 - seqaij: correct permutation 5 - seqbaij: permutation not supported for this MATTYPE 6 - seqsbaij: permutation not supported for this MATTYPE 7 Parallel: 8 - mpiaij: correct permutation 9 - mpibaij: correct permutation 10 - mpisbaij: permutation not supported for this MATTYPE 11 */ 12 13 #include <petscmat.h> 14 15 int main(int argc, char **argv) { 16 const struct { 17 PetscInt i, j; 18 PetscScalar v; 19 } entries[] = { 20 {0, 3, 1.}, 21 {1, 2, 2.}, 22 {2, 1, 3.}, 23 {2, 5, 4.}, 24 {3, 0, 5.}, 25 {3, 6, 6.}, 26 {4, 1, 7.}, 27 {4, 4, 8.} 28 }; 29 const PetscInt ixrow[5] = {4, 2, 1, 0, 3}, ixcol[7] = {5, 3, 6, 1, 2, 0, 4}; 30 Mat A, B; 31 PetscInt i, rstart, rend, cstart, cend; 32 IS isrow, iscol; 33 PetscViewer viewer; 34 PetscBool view_sparse; 35 36 PetscFunctionBeginUser; 37 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 38 /* ------- Assemble matrix, --------- */ 39 PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 40 PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, 5, 7)); 41 PetscCall(MatSetFromOptions(A)); 42 PetscCall(MatSetUp(A)); 43 PetscCall(MatGetOwnershipRange(A, &rstart, &rend)); 44 PetscCall(MatGetOwnershipRangeColumn(A, &cstart, &cend)); 45 46 for (i = 0; i < (PetscInt)PETSC_STATIC_ARRAY_LENGTH(entries); i++) { PetscCall(MatSetValue(A, entries[i].i, entries[i].j, entries[i].v, INSERT_VALUES)); } 47 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 48 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 49 50 /* ------ Prepare index sets ------ */ 51 PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, rend - rstart, ixrow + rstart, PETSC_USE_POINTER, &isrow)); 52 PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, cend - cstart, ixcol + cstart, PETSC_USE_POINTER, &iscol)); 53 PetscCall(ISSetPermutation(isrow)); 54 PetscCall(ISSetPermutation(iscol)); 55 56 PetscCall(PetscViewerASCIIGetStdout(PETSC_COMM_WORLD, &viewer)); 57 view_sparse = PETSC_FALSE; 58 PetscCall(PetscOptionsGetBool(NULL, NULL, "-view_sparse", &view_sparse, NULL)); 59 if (!view_sparse) { PetscCall(PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_DENSE)); } 60 PetscCall(PetscViewerASCIIPrintf(viewer, "Original matrix\n")); 61 PetscCall(MatView(A, viewer)); 62 63 PetscCall(MatPermute(A, isrow, iscol, &B)); 64 PetscCall(PetscViewerASCIIPrintf(viewer, "Permuted matrix\n")); 65 PetscCall(MatView(B, viewer)); 66 67 if (!view_sparse) { PetscCall(PetscViewerPopFormat(viewer)); } 68 PetscCall(PetscViewerASCIIPrintf(viewer, "Row permutation\n")); 69 PetscCall(ISView(isrow, viewer)); 70 PetscCall(PetscViewerASCIIPrintf(viewer, "Column permutation\n")); 71 PetscCall(ISView(iscol, viewer)); 72 73 /* Free data structures */ 74 PetscCall(ISDestroy(&isrow)); 75 PetscCall(ISDestroy(&iscol)); 76 PetscCall(MatDestroy(&A)); 77 PetscCall(MatDestroy(&B)); 78 79 PetscCall(PetscFinalize()); 80 return 0; 81 } 82 83 /*TEST 84 85 build: 86 requires: !complex 87 88 test: 89 args: -view_sparse 90 91 test: 92 suffix: 2 93 nsize: 2 94 args: -view_sparse 95 96 test: 97 suffix: 2b 98 nsize: 2 99 args: -mat_type baij -view_sparse 100 101 test: 102 suffix: 3 103 nsize: 3 104 args: -view_sparse 105 106 test: 107 suffix: 3b 108 nsize: 3 109 args: -mat_type baij -view_sparse 110 111 test: 112 suffix: dense 113 args: -mat_type dense 114 115 TEST*/ 116