1 2 static char help[] = "Tests MatTranspose() and MatEqual() for MPIAIJ matrices.\n\n"; 3 4 #include <petscmat.h> 5 6 int main(int argc, char **argv) 7 { 8 Mat A, B; 9 PetscInt m = 7, n, i, rstart, rend, cols[3]; 10 PetscScalar v[3]; 11 PetscBool equal; 12 const char *eq[2]; 13 14 PetscFunctionBeginUser; 15 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 16 PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL)); 17 n = m; 18 19 /* ------- Assemble matrix, --------- */ 20 21 PetscCall(MatCreateAIJ(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, m, n, 0, 0, 0, 0, &A)); 22 PetscCall(MatSetOption(A, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE)); 23 PetscCall(MatGetOwnershipRange(A, &rstart, &rend)); 24 if (!rstart) { 25 cols[0] = 0; 26 cols[1] = 1; 27 v[0] = 2.0; 28 v[1] = -1.0; 29 PetscCall(MatSetValues(A, 1, &rstart, 2, cols, v, INSERT_VALUES)); 30 rstart++; 31 } 32 if (rend == m) { 33 rend--; 34 cols[0] = rend - 1; 35 cols[1] = rend; 36 v[0] = -1.0; 37 v[1] = 2.0; 38 PetscCall(MatSetValues(A, 1, &rend, 2, cols, v, INSERT_VALUES)); 39 } 40 v[0] = -1.0; 41 v[1] = 2.0; 42 v[2] = -1.0; 43 for (i = rstart; i < rend; i++) { 44 cols[0] = i - 1; 45 cols[1] = i; 46 cols[2] = i + 1; 47 PetscCall(MatSetValues(A, 1, &i, 3, cols, v, INSERT_VALUES)); 48 } 49 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 50 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 51 52 PetscCall(MatTranspose(A, MAT_INITIAL_MATRIX, &B)); 53 54 PetscCall(MatEqual(A, B, &equal)); 55 56 eq[0] = "not equal"; 57 eq[1] = "equal"; 58 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Matrices are %s\n", eq[equal])); 59 60 PetscCall(MatTranspose(A, MAT_REUSE_MATRIX, &B)); 61 PetscCall(MatEqual(A, B, &equal)); 62 if (!equal) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "MatTranspose with MAT_REUSE_MATRIX failed")); 63 64 /* Free data structures */ 65 PetscCall(MatDestroy(&A)); 66 PetscCall(MatDestroy(&B)); 67 68 PetscCall(PetscFinalize()); 69 return 0; 70 } 71 72 /*TEST 73 74 test: 75 76 test: 77 suffix: 2 78 nsize: 2 79 output_file: output/ex58_1.out 80 81 TEST*/ 82