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