1 2 static char help[] = "Tests the use of MatZeroRows() for uniprocessor matrices.\n\n"; 3 4 #include <petscmat.h> 5 6 int main(int argc, char **args) 7 { 8 Mat C; 9 PetscInt i, j, m = 5, n = 5, Ii, J; 10 PetscScalar v, five = 5.0; 11 IS isrow; 12 PetscBool keepnonzeropattern; 13 14 PetscFunctionBeginUser; 15 PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 16 /* create the matrix for the five point stencil, YET AGAIN*/ 17 PetscCall(MatCreate(PETSC_COMM_SELF, &C)); 18 PetscCall(MatSetSizes(C, PETSC_DECIDE, PETSC_DECIDE, m * n, m * n)); 19 PetscCall(MatSetFromOptions(C)); 20 PetscCall(MatSetUp(C)); 21 for (i = 0; i < m; i++) { 22 for (j = 0; j < n; j++) { 23 v = -1.0; 24 Ii = j + n * i; 25 if (i > 0) { 26 J = Ii - n; 27 PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 28 } 29 if (i < m - 1) { 30 J = Ii + n; 31 PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 32 } 33 if (j > 0) { 34 J = Ii - 1; 35 PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 36 } 37 if (j < n - 1) { 38 J = Ii + 1; 39 PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 40 } 41 v = 4.0; 42 PetscCall(MatSetValues(C, 1, &Ii, 1, &Ii, &v, INSERT_VALUES)); 43 } 44 } 45 PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY)); 46 PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY)); 47 48 PetscCall(ISCreateStride(PETSC_COMM_SELF, (m * n) / 2, 0, 2, &isrow)); 49 50 PetscCall(PetscOptionsHasName(NULL, NULL, "-keep_nonzero_pattern", &keepnonzeropattern)); 51 if (keepnonzeropattern) PetscCall(MatSetOption(C, MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE)); 52 53 PetscCall(MatZeroRowsIS(C, isrow, five, 0, 0)); 54 55 PetscCall(MatView(C, PETSC_VIEWER_STDOUT_SELF)); 56 57 PetscCall(ISDestroy(&isrow)); 58 PetscCall(MatDestroy(&C)); 59 PetscCall(PetscFinalize()); 60 return 0; 61 } 62 63 /*TEST 64 65 test: 66 67 test: 68 suffix: 2 69 args: -mat_type seqbaij -mat_block_size 5 70 71 test: 72 suffix: 3 73 args: -keep_nonzero_pattern 74 75 test: 76 suffix: 4 77 args: -keep_nonzero_pattern -mat_type seqbaij -mat_block_size 5 78 79 TEST*/ 80