Actual source code: ex55.c
petsc-3.13.6 2020-09-29
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;
11: PetscInt i,j,ntypes,bs,mbs,m,block,d_nz=6, o_nz=3,col[3],row,verbose=0;
12: PetscMPIInt size,rank;
13: MatType type[9];
14: char file[PETSC_MAX_PATH_LEN];
15: PetscViewer fd;
16: PetscBool equal,flg_loadmat,flg,issymmetric;
17: PetscScalar value[3];
19: PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
20: PetscOptionsGetInt(NULL,NULL,"-verbose",&verbose,NULL);
21: PetscOptionsGetString(NULL,NULL,"-f",file,PETSC_MAX_PATH_LEN,&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; mbs=8;
62: PetscOptionsGetInt(NULL,NULL,"-mbs",&mbs,NULL);
63: PetscOptionsGetInt(NULL,NULL,"-bs",&bs,NULL);
64: if (bs <= 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG," bs must be >1 in this case");
65: m = mbs*bs;
66: MatCreateBAIJ(PETSC_COMM_WORLD,bs,PETSC_DECIDE,PETSC_DECIDE,m,m,d_nz,NULL,o_nz,NULL,&C);
67: for (block=0; block<mbs; block++) {
68: /* diagonal blocks */
69: value[0] = -1.0; value[1] = 4.0; value[2] = -1.0;
70: for (i=1+block*bs; i<bs-1+block*bs; i++) {
71: col[0] = i-1; col[1] = i; col[2] = i+1;
72: MatSetValues(C,1,&i,3,col,value,INSERT_VALUES);
73: }
74: i = bs - 1+block*bs; col[0] = bs - 2+block*bs; col[1] = bs - 1+block*bs;
75: value[0]=-1.0; value[1]=4.0;
76: MatSetValues(C,1,&i,2,col,value,INSERT_VALUES);
78: i = 0+block*bs; col[0] = 0+block*bs; col[1] = 1+block*bs;
79: value[0]=4.0; value[1] = -1.0;
80: MatSetValues(C,1,&i,2,col,value,INSERT_VALUES);
81: }
82: /* off-diagonal blocks */
83: value[0]=-1.0;
84: for (i=0; i<(mbs-1)*bs; i++) {
85: col[0]=i+bs;
86: MatSetValues(C,1,&i,1,col,value,INSERT_VALUES);
87: col[0]=i; row=i+bs;
88: MatSetValues(C,1,&row,1,col,value,INSERT_VALUES);
89: }
90: MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);
91: MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);
92: }
93: {
94: /* Check the symmetry of C because it will be converted to a sbaij matrix */
95: Mat Ctrans;
96: MatTranspose(C,MAT_INITIAL_MATRIX,&Ctrans);
97: MatEqual(C,Ctrans,&flg);
98: /*
99: {
100: MatAXPY(C,1.,Ctrans,DIFFERENT_NONZERO_PATTERN);
101: flg = PETSC_TRUE;
102: }
103: */
104: MatSetOption(C,MAT_SYMMETRIC,flg);
105: MatDestroy(&Ctrans);
106: }
107: MatIsSymmetric(C,0.0,&issymmetric);
108: MatViewFromOptions(C,NULL,"-view_mat");
110: /* convert C to other formats */
111: for (i=0; i<ntypes; i++) {
112: PetscBool ismpisbaij,isseqsbaij;
114: PetscStrcmp(type[i],MATMPISBAIJ,&ismpisbaij);
115: PetscStrcmp(type[i],MATMPISBAIJ,&isseqsbaij);
116: if (!issymmetric && (ismpisbaij || isseqsbaij)) continue;
117: MatConvert(C,type[i],MAT_INITIAL_MATRIX,&A);
118: MatMultEqual(A,C,10,&equal);
119: if (!equal) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in conversion from BAIJ to %s",type[i]);
120: for (j=i+1; j<ntypes; j++) {
121: PetscStrcmp(type[j],MATMPISBAIJ,&ismpisbaij);
122: PetscStrcmp(type[j],MATMPISBAIJ,&isseqsbaij);
123: if (!issymmetric && (ismpisbaij || isseqsbaij)) continue;
124: if (verbose>0) {
125: PetscPrintf(PETSC_COMM_WORLD," \n[%d] test conversion between %s and %s\n",rank,type[i],type[j]);
126: }
128: if (!rank && verbose) printf("Convert %s A to %s B\n",type[i],type[j]);
129: MatConvert(A,type[j],MAT_INITIAL_MATRIX,&B);
130: /*
131: if (j == 2) {
132: PetscPrintf(PETSC_COMM_SELF," A: %s\n",type[i]);
133: MatView(A,PETSC_VIEWER_STDOUT_WORLD);
134: PetscPrintf(PETSC_COMM_SELF," B: %s\n",type[j]);
135: MatView(B,PETSC_VIEWER_STDOUT_WORLD);
136: }
137: */
138: MatMultEqual(A,B,10,&equal);
139: if (!equal) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in conversion from %s to %s",type[i],type[j]);
141: if (size == 1 || j != 2) { /* Matconvert from mpisbaij mat to other formats are not supported */
142: if (!rank && verbose) printf("Convert %s B to %s D\n",type[j],type[i]);
143: MatConvert(B,type[i],MAT_INITIAL_MATRIX,&D);
144: MatMultEqual(B,D,10,&equal);
145: if (!equal) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in conversion from %s to %s",type[j],type[i]);
147: MatDestroy(&D);
148: }
149: MatDestroy(&B);
150: MatDestroy(&D);
151: }
153: /* Test in-place convert */
154: if (size == 1) { /* size > 1 is not working yet! */
155: j = (i+1)%ntypes;
156: PetscStrcmp(type[j],MATMPISBAIJ,&ismpisbaij);
157: PetscStrcmp(type[j],MATMPISBAIJ,&isseqsbaij);
158: if (!issymmetric && (ismpisbaij || isseqsbaij)) continue;
159: /* printf("[%d] i: %d, j: %d\n",rank,i,j); */
160: MatConvert(A,type[j],MAT_INPLACE_MATRIX,&A);
161: }
163: MatDestroy(&A);
164: }
166: /* test BAIJ to MATIS */
167: if (size > 1) {
168: MatType ctype;
170: MatGetType(C,&ctype);
171: MatConvert(C,MATIS,MAT_INITIAL_MATRIX,&A);
172: MatMultEqual(A,C,10,&equal);
173: MatViewFromOptions(A,NULL,"-view_conv");
174: if (!equal) {
175: Mat C2;
177: MatConvert(A,ctype,MAT_INITIAL_MATRIX,&C2);
178: MatViewFromOptions(C2,NULL,"-view_conv_assembled");
179: MatAXPY(C2,-1.,C,DIFFERENT_NONZERO_PATTERN);
180: MatChop(C2,PETSC_SMALL);
181: MatViewFromOptions(C2,NULL,"-view_err");
182: MatDestroy(&C2);
183: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in conversion from BAIJ to MATIS");
184: }
185: MatConvert(C,MATIS,MAT_REUSE_MATRIX,&A);
186: MatMultEqual(A,C,10,&equal);
187: MatViewFromOptions(A,NULL,"-view_conv");
188: if (!equal) {
189: Mat C2;
191: MatConvert(A,ctype,MAT_INITIAL_MATRIX,&C2);
192: MatViewFromOptions(C2,NULL,"-view_conv_assembled");
193: MatAXPY(C2,-1.,C,DIFFERENT_NONZERO_PATTERN);
194: MatChop(C2,PETSC_SMALL);
195: MatViewFromOptions(C2,NULL,"-view_err");
196: MatDestroy(&C2);
197: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in conversion from BAIJ to MATIS");
198: }
199: MatDestroy(&A);
200: MatDuplicate(C,MAT_COPY_VALUES,&A);
201: MatConvert(A,MATIS,MAT_INPLACE_MATRIX,&A);
202: MatViewFromOptions(A,NULL,"-view_conv");
203: MatMultEqual(A,C,10,&equal);
204: if (!equal) {
205: Mat C2;
207: MatViewFromOptions(A,NULL,"-view_conv");
208: MatConvert(A,ctype,MAT_INITIAL_MATRIX,&C2);
209: MatViewFromOptions(C2,NULL,"-view_conv_assembled");
210: MatAXPY(C2,-1.,C,DIFFERENT_NONZERO_PATTERN);
211: MatChop(C2,PETSC_SMALL);
212: MatViewFromOptions(C2,NULL,"-view_err");
213: MatDestroy(&C2);
214: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in conversion from BAIJ to MATIS");
215: }
216: MatDestroy(&A);
217: }
218: MatDestroy(&C);
220: PetscFinalize();
221: return ierr;
222: }
224: /*TEST
226: test:
228: test:
229: suffix: 2
230: nsize: 3
232: testset:
233: requires: parmetis
234: output_file: output/ex55_1.out
235: nsize: 3
236: args: -mat_is_disassemble_l2g_type nd -mat_partitioning_type parmetis
237: test:
238: suffix: matis_baij_parmetis_nd
239: test:
240: suffix: matis_aij_parmetis_nd
241: args: -testseqaij
242: test:
243: requires: datafilespath !complex double !define(PETSC_USE_64BIT_INDICES)
244: suffix: matis_poisson1_parmetis_nd
245: args: -f ${DATAFILESPATH}/matrices/poisson1
247: testset:
248: requires: ptscotch define(PETSC_HAVE_SCOTCH_PARMETIS_V3_NODEND)
249: output_file: output/ex55_1.out
250: nsize: 4
251: args: -mat_is_disassemble_l2g_type nd -mat_partitioning_type ptscotch
252: test:
253: suffix: matis_baij_ptscotch_nd
254: test:
255: suffix: matis_aij_ptscotch_nd
256: args: -testseqaij
257: test:
258: requires: datafilespath !complex double !define(PETSC_USE_64BIT_INDICES)
259: suffix: matis_poisson1_ptscotch_nd
260: args: -f ${DATAFILESPATH}/matrices/poisson1
262: TEST*/