1 2 static char help[] = "Tests converting a SBAIJ matrix to BAIJ format with MatConvert. Modified from ex55.c\n\n"; 3 4 #include <petscmat.h> 5 /* Usage: ./ex141 -mat_view */ 6 7 int main(int argc, char **args) 8 { 9 Mat C, B; 10 PetscInt i, bs = 2, mbs, m, block, d_nz = 6, col[3]; 11 PetscMPIInt size; 12 char file[PETSC_MAX_PATH_LEN]; 13 PetscViewer fd; 14 PetscBool equal, loadmat; 15 PetscScalar value[4]; 16 PetscInt ridx[2], cidx[2]; 17 18 PetscFunctionBeginUser; 19 PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 20 PetscCall(PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &loadmat)); 21 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 22 PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!"); 23 24 /* input matrix C */ 25 if (loadmat) { 26 /* Open binary file. Load a sbaij matrix, then destroy the viewer. */ 27 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd)); 28 PetscCall(MatCreate(PETSC_COMM_WORLD, &C)); 29 PetscCall(MatSetType(C, MATSEQSBAIJ)); 30 PetscCall(MatLoad(C, fd)); 31 PetscCall(PetscViewerDestroy(&fd)); 32 } else { /* Create a sbaij mat with bs>1 */ 33 mbs = 8; 34 PetscCall(PetscOptionsGetInt(NULL, NULL, "-mbs", &mbs, NULL)); 35 m = mbs * bs; 36 PetscCall(MatCreate(PETSC_COMM_WORLD, &C)); 37 PetscCall(MatSetSizes(C, PETSC_DECIDE, PETSC_DECIDE, m, m)); 38 PetscCall(MatSetType(C, MATSBAIJ)); 39 PetscCall(MatSetFromOptions(C)); 40 PetscCall(MatSeqSBAIJSetPreallocation(C, bs, d_nz, NULL)); 41 PetscCall(MatSetUp(C)); 42 PetscCall(MatSetOption(C, MAT_IGNORE_LOWER_TRIANGULAR, PETSC_TRUE)); 43 44 for (block = 0; block < mbs; block++) { 45 /* diagonal blocks */ 46 value[0] = -1.0; 47 value[1] = 4.0; 48 value[2] = -1.0; 49 for (i = 1 + block * bs; i < bs - 1 + block * bs; i++) { 50 col[0] = i - 1; 51 col[1] = i; 52 col[2] = i + 1; 53 PetscCall(MatSetValues(C, 1, &i, 3, col, value, INSERT_VALUES)); 54 } 55 i = bs - 1 + block * bs; 56 col[0] = bs - 2 + block * bs; 57 col[1] = bs - 1 + block * bs; 58 59 value[0] = -1.0; 60 value[1] = 4.0; 61 PetscCall(MatSetValues(C, 1, &i, 2, col, value, INSERT_VALUES)); 62 63 i = 0 + block * bs; 64 col[0] = 0 + block * bs; 65 col[1] = 1 + block * bs; 66 67 value[0] = 4.0; 68 value[1] = -1.0; 69 PetscCall(MatSetValues(C, 1, &i, 2, col, value, INSERT_VALUES)); 70 } 71 /* off-diagonal blocks */ 72 value[0] = -1.0; 73 value[1] = -0.1; 74 value[2] = 0.0; 75 value[3] = -1.0; /* row-oriented */ 76 for (block = 0; block < mbs - 1; block++) { 77 for (i = 0; i < bs; i++) { 78 ridx[i] = block * bs + i; 79 cidx[i] = (block + 1) * bs + i; 80 } 81 PetscCall(MatSetValues(C, bs, ridx, bs, cidx, value, INSERT_VALUES)); 82 } 83 PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY)); 84 PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY)); 85 } 86 87 /* convert C to BAIJ format */ 88 PetscCall(MatConvert(C, MATSEQBAIJ, MAT_INITIAL_MATRIX, &B)); 89 PetscCall(MatMultEqual(B, C, 10, &equal)); 90 PetscCheck(equal, PETSC_COMM_SELF, PETSC_ERR_PLIB, "MatConvert fails!"); 91 92 PetscCall(MatDestroy(&B)); 93 PetscCall(MatDestroy(&C)); 94 PetscCall(PetscFinalize()); 95 return 0; 96 } 97 98 /*TEST 99 100 test: 101 output_file: output/ex141.out 102 103 TEST*/ 104