Actual source code: ex55.c

petsc-3.14.6 2021-03-30
Report Typos and Errors

  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,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; 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*/