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