Actual source code: ex141.c
2: static char help[] = "Tests converting a SBAIJ matrix to BAIJ format with MatConvert. Modified from ex55.c\n\n";
4: #include <petscmat.h>
5: /* Usage: ./ex141 -mat_view */
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];
19: PetscInitialize(&argc, &args, (char *)0, help);
20: PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &loadmat);
21: MPI_Comm_size(PETSC_COMM_WORLD, &size);
24: /* input matrix C */
25: if (loadmat) {
26: /* Open binary file. Load a sbaij matrix, then destroy the viewer. */
27: PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd);
28: MatCreate(PETSC_COMM_WORLD, &C);
29: MatSetType(C, MATSEQSBAIJ);
30: MatLoad(C, fd);
31: PetscViewerDestroy(&fd);
32: } else { /* Create a sbaij mat with bs>1 */
33: mbs = 8;
34: PetscOptionsGetInt(NULL, NULL, "-mbs", &mbs, NULL);
35: m = mbs * bs;
36: MatCreate(PETSC_COMM_WORLD, &C);
37: MatSetSizes(C, PETSC_DECIDE, PETSC_DECIDE, m, m);
38: MatSetType(C, MATSBAIJ);
39: MatSetFromOptions(C);
40: MatSeqSBAIJSetPreallocation(C, bs, d_nz, NULL);
41: MatSetUp(C);
42: MatSetOption(C, MAT_IGNORE_LOWER_TRIANGULAR, PETSC_TRUE);
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: 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;
59: value[0] = -1.0;
60: value[1] = 4.0;
61: MatSetValues(C, 1, &i, 2, col, value, INSERT_VALUES);
63: i = 0 + block * bs;
64: col[0] = 0 + block * bs;
65: col[1] = 1 + block * bs;
67: value[0] = 4.0;
68: value[1] = -1.0;
69: 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: MatSetValues(C, bs, ridx, bs, cidx, value, INSERT_VALUES);
82: }
83: MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY);
84: MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY);
85: }
87: /* convert C to BAIJ format */
88: MatConvert(C, MATSEQBAIJ, MAT_INITIAL_MATRIX, &B);
89: MatMultEqual(B, C, 10, &equal);
92: MatDestroy(&B);
93: MatDestroy(&C);
94: PetscFinalize();
95: return 0;
96: }
98: /*TEST
100: test:
101: output_file: output/ex141.out
103: TEST*/