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