1 2 static char help[] = "Tests MatSeq(B)AIJSetColumnIndices().\n\n"; 3 4 #include <petscmat.h> 5 6 /* 7 Generate the following matrix: 8 9 1 0 3 10 1 2 3 11 0 0 3 12 */ 13 int main(int argc, char **args) 14 { 15 Mat A; 16 PetscScalar v; 17 PetscInt i, j, rowlens[] = {2, 3, 1}, cols[] = {0, 2, 0, 1, 2, 2}; 18 PetscBool flg; 19 20 PetscFunctionBeginUser; 21 PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 22 PetscCall(PetscOptionsHasName(NULL, NULL, "-baij", &flg)); 23 if (flg) { 24 PetscCall(MatCreateSeqBAIJ(PETSC_COMM_WORLD, 1, 3, 3, 0, rowlens, &A)); 25 PetscCall(MatSeqBAIJSetColumnIndices(A, cols)); 26 } else { 27 PetscCall(MatCreateSeqAIJ(PETSC_COMM_WORLD, 3, 3, 0, rowlens, &A)); 28 PetscCall(MatSeqAIJSetColumnIndices(A, cols)); 29 } 30 31 i = 0; 32 j = 0; 33 v = 1.0; 34 PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES)); 35 i = 0; 36 j = 2; 37 v = 3.0; 38 PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES)); 39 40 i = 1; 41 j = 0; 42 v = 1.0; 43 PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES)); 44 i = 1; 45 j = 1; 46 v = 2.0; 47 PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES)); 48 i = 1; 49 j = 2; 50 v = 3.0; 51 PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES)); 52 53 i = 2; 54 j = 2; 55 v = 3.0; 56 PetscCall(MatSetValues(A, 1, &i, 1, &j, &v, INSERT_VALUES)); 57 58 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 59 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 60 PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 61 62 PetscCall(MatDestroy(&A)); 63 PetscCall(PetscFinalize()); 64 return 0; 65 } 66 67 /*TEST 68 69 test: 70 71 test: 72 suffix: 2 73 args: -baij 74 75 TEST*/ 76