1 static char help[] = "Tests converting a parallel AIJ formatted matrix to the parallel Row format.\n\ 2 This also tests MatGetRow() and MatRestoreRow() for the parallel case.\n\n"; 3 4 #include <petscmat.h> 5 6 int main(int argc, char **args) 7 { 8 Mat C, A; 9 PetscInt i, j, m = 3, n = 2, Ii, J, rstart, rend, nz; 10 PetscMPIInt rank, size; 11 const PetscInt *idx; 12 PetscScalar v; 13 const PetscScalar *values; 14 15 PetscFunctionBeginUser; 16 PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 17 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 18 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 19 n = 2 * size; 20 21 /* create the matrix for the five point stencil, YET AGAIN*/ 22 PetscCall(MatCreate(PETSC_COMM_WORLD, &C)); 23 PetscCall(MatSetSizes(C, PETSC_DECIDE, PETSC_DECIDE, m * n, m * n)); 24 PetscCall(MatSetFromOptions(C)); 25 PetscCall(MatMPIAIJSetPreallocation(C, 5, NULL, 5, NULL)); 26 PetscCall(MatSeqAIJSetPreallocation(C, 5, NULL)); 27 for (i = 0; i < m; i++) { 28 for (j = 2 * rank; j < 2 * rank + 2; j++) { 29 v = -1.0; 30 Ii = j + n * i; 31 if (i > 0) { 32 J = Ii - n; 33 PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 34 } 35 if (i < m - 1) { 36 J = Ii + n; 37 PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 38 } 39 if (j > 0) { 40 J = Ii - 1; 41 PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 42 } 43 if (j < n - 1) { 44 J = Ii + 1; 45 PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 46 } 47 v = 4.0; 48 PetscCall(MatSetValues(C, 1, &Ii, 1, &Ii, &v, INSERT_VALUES)); 49 } 50 } 51 PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY)); 52 PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY)); 53 PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_INFO)); 54 PetscCall(MatView(C, PETSC_VIEWER_STDOUT_WORLD)); 55 PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); 56 PetscCall(MatView(C, PETSC_VIEWER_STDOUT_WORLD)); 57 58 PetscCall(MatGetOwnershipRange(C, &rstart, &rend)); 59 for (i = rstart; i < rend; i++) { 60 PetscCall(MatGetRow(C, i, &nz, &idx, &values)); 61 PetscCall(PetscSynchronizedFPrintf(PETSC_COMM_WORLD, stdout, "[%d] get row %" PetscInt_FMT ": ", rank, i)); 62 for (j = 0; j < nz; j++) PetscCall(PetscSynchronizedFPrintf(PETSC_COMM_WORLD, stdout, "%" PetscInt_FMT " %g ", idx[j], (double)PetscRealPart(values[j]))); 63 PetscCall(PetscSynchronizedFPrintf(PETSC_COMM_WORLD, stdout, "\n")); 64 PetscCall(MatRestoreRow(C, i, &nz, &idx, &values)); 65 } 66 PetscCall(PetscSynchronizedFlush(PETSC_COMM_WORLD, stdout)); 67 68 PetscCall(MatConvert(C, MATSAME, MAT_INITIAL_MATRIX, &A)); 69 PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_INFO)); 70 PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 71 PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); 72 PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 73 74 PetscCall(MatDestroy(&A)); 75 PetscCall(MatDestroy(&C)); 76 PetscCall(PetscFinalize()); 77 return 0; 78 } 79 80 /*TEST 81 82 test: 83 args: -mat_type seqaij 84 85 TEST*/ 86