1 static char help[] = "Tests MatCreateSubMatrices().\n\n"; 2 3 #include <petscmat.h> 4 5 int main(int argc, char **args) 6 { 7 Mat A, B, *Bsub; 8 PetscInt i, j, m = 6, n = 6, N = 36, Ii, J; 9 PetscScalar v; 10 IS isrow; 11 12 PetscFunctionBeginUser; 13 PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 14 PetscCall(MatCreateSeqAIJ(PETSC_COMM_WORLD, N, N, 5, NULL, &A)); 15 for (i = 0; i < m; i++) { 16 for (j = 0; j < n; j++) { 17 v = -1.0; 18 Ii = j + n * i; 19 if (i > 0) { 20 J = Ii - n; 21 PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 22 } 23 if (i < m - 1) { 24 J = Ii + n; 25 PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 26 } 27 if (j > 0) { 28 J = Ii - 1; 29 PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 30 } 31 if (j < n - 1) { 32 J = Ii + 1; 33 PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES)); 34 } 35 v = 4.0; 36 PetscCall(MatSetValues(A, 1, &Ii, 1, &Ii, &v, INSERT_VALUES)); 37 } 38 } 39 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 40 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 41 PetscCall(MatView(A, PETSC_VIEWER_STDOUT_SELF)); 42 43 /* take the first diagonal block */ 44 PetscCall(ISCreateStride(PETSC_COMM_WORLD, m, 0, 1, &isrow)); 45 PetscCall(MatCreateSubMatrices(A, 1, &isrow, &isrow, MAT_INITIAL_MATRIX, &Bsub)); 46 B = *Bsub; 47 PetscCall(ISDestroy(&isrow)); 48 PetscCall(MatView(B, PETSC_VIEWER_STDOUT_SELF)); 49 PetscCall(MatDestroySubMatrices(1, &Bsub)); 50 51 /* take a strided block */ 52 PetscCall(ISCreateStride(PETSC_COMM_WORLD, m, 0, 2, &isrow)); 53 PetscCall(MatCreateSubMatrices(A, 1, &isrow, &isrow, MAT_INITIAL_MATRIX, &Bsub)); 54 B = *Bsub; 55 PetscCall(ISDestroy(&isrow)); 56 PetscCall(MatView(B, PETSC_VIEWER_STDOUT_SELF)); 57 PetscCall(MatDestroySubMatrices(1, &Bsub)); 58 59 /* take the last block */ 60 PetscCall(ISCreateStride(PETSC_COMM_WORLD, m, N - m - 1, 1, &isrow)); 61 PetscCall(MatCreateSubMatrices(A, 1, &isrow, &isrow, MAT_INITIAL_MATRIX, &Bsub)); 62 B = *Bsub; 63 PetscCall(ISDestroy(&isrow)); 64 PetscCall(MatView(B, PETSC_VIEWER_STDOUT_SELF)); 65 66 PetscCall(MatDestroySubMatrices(1, &Bsub)); 67 PetscCall(MatDestroy(&A)); 68 69 PetscCall(PetscFinalize()); 70 return 0; 71 } 72 73 /*TEST 74 75 test: 76 77 TEST*/ 78