Actual source code: ex55.c
2: static char help[] = "Tests converting a matrix to another format with MatConvert().\n\n";
4: #include <petscmat.h>
5: /* Usage: mpiexec -n <np> ex55 -verbose <0 or 1> */
7: int main(int argc, char **args)
8: {
9: Mat C, A, B, D;
10: PetscInt i, j, ntypes, bs, mbs, m, block, d_nz = 6, o_nz = 3, col[3], row, verbose = 0;
11: PetscMPIInt size, rank;
12: MatType type[9];
13: char file[PETSC_MAX_PATH_LEN];
14: PetscViewer fd;
15: PetscBool equal, flg_loadmat, flg, issymmetric;
16: PetscScalar value[3];
19: PetscInitialize(&argc, &args, (char *)0, help);
20: PetscOptionsGetInt(NULL, NULL, "-verbose", &verbose, NULL);
21: PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg_loadmat);
22: MPI_Comm_size(PETSC_COMM_WORLD, &size);
23: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
25: PetscOptionsHasName(NULL, NULL, "-testseqaij", &flg);
26: if (flg) {
27: if (size == 1) {
28: type[0] = MATSEQAIJ;
29: } else {
30: type[0] = MATMPIAIJ;
31: }
32: } else {
33: type[0] = MATAIJ;
34: }
35: if (size == 1) {
36: ntypes = 3;
37: type[1] = MATSEQBAIJ;
38: type[2] = MATSEQSBAIJ;
39: } else {
40: ntypes = 3;
41: type[1] = MATMPIBAIJ;
42: type[2] = MATMPISBAIJ;
43: }
45: /* input matrix C */
46: if (flg_loadmat) {
47: /* Open binary file. */
48: PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd);
50: /* Load a baij matrix, then destroy the viewer. */
51: MatCreate(PETSC_COMM_WORLD, &C);
52: if (size == 1) {
53: MatSetType(C, MATSEQBAIJ);
54: } else {
55: MatSetType(C, MATMPIBAIJ);
56: }
57: MatSetFromOptions(C);
58: MatLoad(C, fd);
59: PetscViewerDestroy(&fd);
60: } else { /* Create a baij mat with bs>1 */
61: bs = 2;
62: mbs = 8;
63: PetscOptionsGetInt(NULL, NULL, "-mbs", &mbs, NULL);
64: PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL);
66: m = mbs * bs;
67: MatCreateBAIJ(PETSC_COMM_WORLD, bs, PETSC_DECIDE, PETSC_DECIDE, m, m, d_nz, NULL, o_nz, NULL, &C);
68: for (block = 0; block < mbs; block++) {
69: /* diagonal blocks */
70: value[0] = -1.0;
71: value[1] = 4.0;
72: value[2] = -1.0;
73: for (i = 1 + block * bs; i < bs - 1 + block * bs; i++) {
74: col[0] = i - 1;
75: col[1] = i;
76: col[2] = i + 1;
77: MatSetValues(C, 1, &i, 3, col, value, INSERT_VALUES);
78: }
79: i = bs - 1 + block * bs;
80: col[0] = bs - 2 + block * bs;
81: col[1] = bs - 1 + block * bs;
82: value[0] = -1.0;
83: value[1] = 4.0;
84: MatSetValues(C, 1, &i, 2, col, value, INSERT_VALUES);
86: i = 0 + block * bs;
87: col[0] = 0 + block * bs;
88: col[1] = 1 + block * bs;
89: value[0] = 4.0;
90: value[1] = -1.0;
91: MatSetValues(C, 1, &i, 2, col, value, INSERT_VALUES);
92: }
93: /* off-diagonal blocks */
94: value[0] = -1.0;
95: for (i = 0; i < (mbs - 1) * bs; i++) {
96: col[0] = i + bs;
97: MatSetValues(C, 1, &i, 1, col, value, INSERT_VALUES);
98: col[0] = i;
99: row = i + bs;
100: MatSetValues(C, 1, &row, 1, col, value, INSERT_VALUES);
101: }
102: MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY);
103: MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY);
104: }
105: {
106: /* Check the symmetry of C because it will be converted to a sbaij matrix */
107: Mat Ctrans;
108: MatTranspose(C, MAT_INITIAL_MATRIX, &Ctrans);
109: MatEqual(C, Ctrans, &flg);
110: /*
111: {
112: MatAXPY(C,1.,Ctrans,DIFFERENT_NONZERO_PATTERN);
113: flg = PETSC_TRUE;
114: }
115: */
116: MatSetOption(C, MAT_SYMMETRIC, flg);
117: MatDestroy(&Ctrans);
118: }
119: MatIsSymmetric(C, 0.0, &issymmetric);
120: MatViewFromOptions(C, NULL, "-view_mat");
122: /* convert C to other formats */
123: for (i = 0; i < ntypes; i++) {
124: PetscBool ismpisbaij, isseqsbaij;
126: PetscStrcmp(type[i], MATMPISBAIJ, &ismpisbaij);
127: PetscStrcmp(type[i], MATMPISBAIJ, &isseqsbaij);
128: if (!issymmetric && (ismpisbaij || isseqsbaij)) continue;
129: MatConvert(C, type[i], MAT_INITIAL_MATRIX, &A);
130: MatMultEqual(A, C, 10, &equal);
132: for (j = i + 1; j < ntypes; j++) {
133: PetscStrcmp(type[j], MATMPISBAIJ, &ismpisbaij);
134: PetscStrcmp(type[j], MATMPISBAIJ, &isseqsbaij);
135: if (!issymmetric && (ismpisbaij || isseqsbaij)) continue;
136: if (verbose > 0) PetscPrintf(PETSC_COMM_WORLD, " \n[%d] test conversion between %s and %s\n", rank, type[i], type[j]);
138: if (rank == 0 && verbose) printf("Convert %s A to %s B\n", type[i], type[j]);
139: MatConvert(A, type[j], MAT_INITIAL_MATRIX, &B);
140: /*
141: if (j == 2) {
142: PetscPrintf(PETSC_COMM_SELF," A: %s\n",type[i]);
143: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
144: PetscPrintf(PETSC_COMM_SELF," B: %s\n",type[j]);
145: MatView(B,PETSC_VIEWER_STDOUT_WORLD);
146: }
147: */
148: MatMultEqual(A, B, 10, &equal);
151: if (size == 1 || j != 2) { /* Matconvert from mpisbaij mat to other formats are not supported */
152: if (rank == 0 && verbose) printf("Convert %s B to %s D\n", type[j], type[i]);
153: MatConvert(B, type[i], MAT_INITIAL_MATRIX, &D);
154: MatMultEqual(B, D, 10, &equal);
157: MatDestroy(&D);
158: }
159: MatDestroy(&B);
160: MatDestroy(&D);
161: }
163: /* Test in-place convert */
164: if (size == 1) { /* size > 1 is not working yet! */
165: j = (i + 1) % ntypes;
166: PetscStrcmp(type[j], MATMPISBAIJ, &ismpisbaij);
167: PetscStrcmp(type[j], MATMPISBAIJ, &isseqsbaij);
168: if (!issymmetric && (ismpisbaij || isseqsbaij)) continue;
169: /* printf("[%d] i: %d, j: %d\n",rank,i,j); */
170: MatConvert(A, type[j], MAT_INPLACE_MATRIX, &A);
171: }
173: MatDestroy(&A);
174: }
176: /* test BAIJ to MATIS */
177: if (size > 1) {
178: MatType ctype;
180: MatGetType(C, &ctype);
181: MatConvert(C, MATIS, MAT_INITIAL_MATRIX, &A);
182: MatMultEqual(A, C, 10, &equal);
183: MatViewFromOptions(A, NULL, "-view_conv");
184: if (!equal) {
185: Mat C2;
187: MatConvert(A, ctype, MAT_INITIAL_MATRIX, &C2);
188: MatViewFromOptions(C2, NULL, "-view_conv_assembled");
189: MatAXPY(C2, -1., C, DIFFERENT_NONZERO_PATTERN);
190: MatChop(C2, PETSC_SMALL);
191: MatViewFromOptions(C2, NULL, "-view_err");
192: MatDestroy(&C2);
193: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error in conversion from BAIJ to MATIS");
194: }
195: MatConvert(C, MATIS, MAT_REUSE_MATRIX, &A);
196: MatMultEqual(A, C, 10, &equal);
197: MatViewFromOptions(A, NULL, "-view_conv");
198: if (!equal) {
199: Mat C2;
201: MatConvert(A, ctype, MAT_INITIAL_MATRIX, &C2);
202: MatViewFromOptions(C2, NULL, "-view_conv_assembled");
203: MatAXPY(C2, -1., C, DIFFERENT_NONZERO_PATTERN);
204: MatChop(C2, PETSC_SMALL);
205: MatViewFromOptions(C2, NULL, "-view_err");
206: MatDestroy(&C2);
207: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error in conversion from BAIJ to MATIS");
208: }
209: MatDestroy(&A);
210: MatDuplicate(C, MAT_COPY_VALUES, &A);
211: MatConvert(A, MATIS, MAT_INPLACE_MATRIX, &A);
212: MatViewFromOptions(A, NULL, "-view_conv");
213: MatMultEqual(A, C, 10, &equal);
214: if (!equal) {
215: Mat C2;
217: MatViewFromOptions(A, NULL, "-view_conv");
218: MatConvert(A, ctype, MAT_INITIAL_MATRIX, &C2);
219: MatViewFromOptions(C2, NULL, "-view_conv_assembled");
220: MatAXPY(C2, -1., C, DIFFERENT_NONZERO_PATTERN);
221: MatChop(C2, PETSC_SMALL);
222: MatViewFromOptions(C2, NULL, "-view_err");
223: MatDestroy(&C2);
224: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_PLIB, "Error in conversion from BAIJ to MATIS");
225: }
226: MatDestroy(&A);
227: }
228: MatDestroy(&C);
230: PetscFinalize();
231: return 0;
232: }
234: /*TEST
236: test:
238: test:
239: suffix: 2
240: nsize: 3
242: testset:
243: requires: parmetis
244: output_file: output/ex55_1.out
245: nsize: 3
246: args: -mat_is_disassemble_l2g_type nd -mat_partitioning_type parmetis
247: test:
248: suffix: matis_baij_parmetis_nd
249: test:
250: suffix: matis_aij_parmetis_nd
251: args: -testseqaij
252: test:
253: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
254: suffix: matis_poisson1_parmetis_nd
255: args: -f ${DATAFILESPATH}/matrices/poisson1
257: testset:
258: requires: ptscotch defined(PETSC_HAVE_SCOTCH_PARMETIS_V3_NODEND)
259: output_file: output/ex55_1.out
260: nsize: 4
261: args: -mat_is_disassemble_l2g_type nd -mat_partitioning_type ptscotch
262: test:
263: suffix: matis_baij_ptscotch_nd
264: test:
265: suffix: matis_aij_ptscotch_nd
266: args: -testseqaij
267: test:
268: requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
269: suffix: matis_poisson1_ptscotch_nd
270: args: -f ${DATAFILESPATH}/matrices/poisson1
272: TEST*/