1 2 static char help[] = "Test saving SeqSBAIJ matrix that is missing diagonal entries."; 3 4 #include <petscmat.h> 5 6 int main(int argc, char **args) 7 { 8 Mat A; 9 PetscInt bs = 3, m = 4, i, j, val = 10, row[2], col[3], rstart; 10 PetscMPIInt size; 11 PetscScalar x[6][9]; 12 13 PetscFunctionBeginUser; 14 PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 15 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 16 PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "Test is only for sequential"); 17 PetscCall(MatCreateSeqSBAIJ(PETSC_COMM_SELF, bs, m * bs, m * bs, 1, NULL, &A)); 18 PetscCall(MatSetOption(A, MAT_IGNORE_LOWER_TRIANGULAR, PETSC_TRUE)); 19 rstart = 0; 20 21 row[0] = rstart + 0; 22 row[1] = rstart + 2; 23 col[0] = rstart + 0; 24 col[1] = rstart + 1; 25 col[2] = rstart + 3; 26 for (i = 0; i < 6; i++) { 27 for (j = 0; j < 9; j++) x[i][j] = (PetscScalar)val++; 28 } 29 PetscCall(MatSetValuesBlocked(A, 2, row, 3, col, &x[0][0], INSERT_VALUES)); 30 PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 31 PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 32 PetscCall(MatView(A, PETSC_VIEWER_BINARY_WORLD)); 33 34 PetscCall(MatDestroy(&A)); 35 PetscCall(PetscFinalize()); 36 return 0; 37 } 38