1 2 static char help[] = "Tests converting a matrix to another format with MatConvert().\n\n"; 3 4 #include <petscmat.h> 5 /* Usage: mpiexec -n <np> ex55 -verbose <0 or 1> */ 6 7 int main(int argc,char **args) 8 { 9 Mat C,A,B,D; 10 PetscErrorCode ierr; 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]; 18 19 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 20 ierr = PetscOptionsGetInt(NULL,NULL,"-verbose",&verbose,NULL);CHKERRQ(ierr); 21 ierr = PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg_loadmat);CHKERRQ(ierr); 22 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr); 23 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 24 25 ierr = PetscOptionsHasName(NULL,NULL,"-testseqaij",&flg);CHKERRQ(ierr); 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 } 44 45 /* input matrix C */ 46 if (flg_loadmat) { 47 /* Open binary file. */ 48 ierr = PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd);CHKERRQ(ierr); 49 50 /* Load a baij matrix, then destroy the viewer. */ 51 ierr = MatCreate(PETSC_COMM_WORLD,&C);CHKERRQ(ierr); 52 if (size == 1) { 53 ierr = MatSetType(C,MATSEQBAIJ);CHKERRQ(ierr); 54 } else { 55 ierr = MatSetType(C,MATMPIBAIJ);CHKERRQ(ierr); 56 } 57 ierr = MatSetFromOptions(C);CHKERRQ(ierr); 58 ierr = MatLoad(C,fd);CHKERRQ(ierr); 59 ierr = PetscViewerDestroy(&fd);CHKERRQ(ierr); 60 } else { /* Create a baij mat with bs>1 */ 61 bs = 2; mbs=8; 62 ierr = PetscOptionsGetInt(NULL,NULL,"-mbs",&mbs,NULL);CHKERRQ(ierr); 63 ierr = PetscOptionsGetInt(NULL,NULL,"-bs",&bs,NULL);CHKERRQ(ierr); 64 if (bs <= 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG," bs must be >1 in this case"); 65 m = mbs*bs; 66 ierr = MatCreateBAIJ(PETSC_COMM_WORLD,bs,PETSC_DECIDE,PETSC_DECIDE,m,m,d_nz,NULL,o_nz,NULL,&C);CHKERRQ(ierr); 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 ierr = MatSetValues(C,1,&i,3,col,value,INSERT_VALUES);CHKERRQ(ierr); 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 ierr = MatSetValues(C,1,&i,2,col,value,INSERT_VALUES);CHKERRQ(ierr); 77 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 ierr = MatSetValues(C,1,&i,2,col,value,INSERT_VALUES);CHKERRQ(ierr); 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 ierr = MatSetValues(C,1,&i,1,col,value,INSERT_VALUES);CHKERRQ(ierr); 87 col[0]=i; row=i+bs; 88 ierr = MatSetValues(C,1,&row,1,col,value,INSERT_VALUES);CHKERRQ(ierr); 89 } 90 ierr = MatAssemblyBegin(C,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 91 ierr = MatAssemblyEnd(C,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 92 } 93 { 94 /* Check the symmetry of C because it will be converted to a sbaij matrix */ 95 Mat Ctrans; 96 ierr = MatTranspose(C,MAT_INITIAL_MATRIX,&Ctrans);CHKERRQ(ierr); 97 ierr = MatEqual(C,Ctrans,&flg);CHKERRQ(ierr); 98 /* 99 { 100 ierr = MatAXPY(C,1.,Ctrans,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr); 101 flg = PETSC_TRUE; 102 } 103 */ 104 ierr = MatSetOption(C,MAT_SYMMETRIC,flg);CHKERRQ(ierr); 105 ierr = MatDestroy(&Ctrans);CHKERRQ(ierr); 106 } 107 ierr = MatIsSymmetric(C,0.0,&issymmetric);CHKERRQ(ierr); 108 ierr = MatViewFromOptions(C,NULL,"-view_mat");CHKERRQ(ierr); 109 110 /* convert C to other formats */ 111 for (i=0; i<ntypes; i++) { 112 PetscBool ismpisbaij,isseqsbaij; 113 114 ierr = PetscStrcmp(type[i],MATMPISBAIJ,&ismpisbaij);CHKERRQ(ierr); 115 ierr = PetscStrcmp(type[i],MATMPISBAIJ,&isseqsbaij);CHKERRQ(ierr); 116 if (!issymmetric && (ismpisbaij || isseqsbaij)) continue; 117 ierr = MatConvert(C,type[i],MAT_INITIAL_MATRIX,&A);CHKERRQ(ierr); 118 ierr = MatMultEqual(A,C,10,&equal);CHKERRQ(ierr); 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 ierr = PetscStrcmp(type[j],MATMPISBAIJ,&ismpisbaij);CHKERRQ(ierr); 122 ierr = PetscStrcmp(type[j],MATMPISBAIJ,&isseqsbaij);CHKERRQ(ierr); 123 if (!issymmetric && (ismpisbaij || isseqsbaij)) continue; 124 if (verbose>0) { 125 ierr = PetscPrintf(PETSC_COMM_WORLD," \n[%d] test conversion between %s and %s\n",rank,type[i],type[j]);CHKERRQ(ierr); 126 } 127 128 if (!rank && verbose) printf("Convert %s A to %s B\n",type[i],type[j]); 129 ierr = MatConvert(A,type[j],MAT_INITIAL_MATRIX,&B);CHKERRQ(ierr); 130 /* 131 if (j == 2) { 132 ierr = PetscPrintf(PETSC_COMM_SELF," A: %s\n",type[i]);CHKERRQ(ierr); 133 ierr = MatView(A,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 134 ierr = PetscPrintf(PETSC_COMM_SELF," B: %s\n",type[j]);CHKERRQ(ierr); 135 ierr = MatView(B,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 136 } 137 */ 138 ierr = MatMultEqual(A,B,10,&equal);CHKERRQ(ierr); 139 if (!equal) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in conversion from %s to %s",type[i],type[j]); 140 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 ierr = MatConvert(B,type[i],MAT_INITIAL_MATRIX,&D);CHKERRQ(ierr); 144 ierr = MatMultEqual(B,D,10,&equal);CHKERRQ(ierr); 145 if (!equal) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"Error in conversion from %s to %s",type[j],type[i]); 146 147 ierr = MatDestroy(&D);CHKERRQ(ierr); 148 } 149 ierr = MatDestroy(&B);CHKERRQ(ierr); 150 ierr = MatDestroy(&D);CHKERRQ(ierr); 151 } 152 153 /* Test in-place convert */ 154 if (size == 1) { /* size > 1 is not working yet! */ 155 j = (i+1)%ntypes; 156 ierr = PetscStrcmp(type[j],MATMPISBAIJ,&ismpisbaij);CHKERRQ(ierr); 157 ierr = PetscStrcmp(type[j],MATMPISBAIJ,&isseqsbaij);CHKERRQ(ierr); 158 if (!issymmetric && (ismpisbaij || isseqsbaij)) continue; 159 /* printf("[%d] i: %d, j: %d\n",rank,i,j); */ 160 ierr = MatConvert(A,type[j],MAT_INPLACE_MATRIX,&A);CHKERRQ(ierr); 161 } 162 163 ierr = MatDestroy(&A);CHKERRQ(ierr); 164 } 165 166 /* test BAIJ to MATIS */ 167 if (size > 1) { 168 MatType ctype; 169 170 ierr = MatGetType(C,&ctype);CHKERRQ(ierr); 171 ierr = MatConvert(C,MATIS,MAT_INITIAL_MATRIX,&A);CHKERRQ(ierr); 172 ierr = MatMultEqual(A,C,10,&equal);CHKERRQ(ierr); 173 ierr = MatViewFromOptions(A,NULL,"-view_conv");CHKERRQ(ierr); 174 if (!equal) { 175 Mat C2; 176 177 ierr = MatConvert(A,ctype,MAT_INITIAL_MATRIX,&C2);CHKERRQ(ierr); 178 ierr = MatViewFromOptions(C2,NULL,"-view_conv_assembled");CHKERRQ(ierr); 179 ierr = MatAXPY(C2,-1.,C,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr); 180 ierr = MatChop(C2,PETSC_SMALL);CHKERRQ(ierr); 181 ierr = MatViewFromOptions(C2,NULL,"-view_err");CHKERRQ(ierr); 182 ierr = MatDestroy(&C2);CHKERRQ(ierr); 183 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in conversion from BAIJ to MATIS"); 184 } 185 ierr = MatConvert(C,MATIS,MAT_REUSE_MATRIX,&A);CHKERRQ(ierr); 186 ierr = MatMultEqual(A,C,10,&equal);CHKERRQ(ierr); 187 ierr = MatViewFromOptions(A,NULL,"-view_conv");CHKERRQ(ierr); 188 if (!equal) { 189 Mat C2; 190 191 ierr = MatConvert(A,ctype,MAT_INITIAL_MATRIX,&C2);CHKERRQ(ierr); 192 ierr = MatViewFromOptions(C2,NULL,"-view_conv_assembled");CHKERRQ(ierr); 193 ierr = MatAXPY(C2,-1.,C,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr); 194 ierr = MatChop(C2,PETSC_SMALL);CHKERRQ(ierr); 195 ierr = MatViewFromOptions(C2,NULL,"-view_err");CHKERRQ(ierr); 196 ierr = MatDestroy(&C2);CHKERRQ(ierr); 197 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in conversion from BAIJ to MATIS"); 198 } 199 ierr = MatDestroy(&A);CHKERRQ(ierr); 200 ierr = MatDuplicate(C,MAT_COPY_VALUES,&A);CHKERRQ(ierr); 201 ierr = MatConvert(A,MATIS,MAT_INPLACE_MATRIX,&A);CHKERRQ(ierr); 202 ierr = MatViewFromOptions(A,NULL,"-view_conv");CHKERRQ(ierr); 203 ierr = MatMultEqual(A,C,10,&equal);CHKERRQ(ierr); 204 if (!equal) { 205 Mat C2; 206 207 ierr = MatViewFromOptions(A,NULL,"-view_conv");CHKERRQ(ierr); 208 ierr = MatConvert(A,ctype,MAT_INITIAL_MATRIX,&C2);CHKERRQ(ierr); 209 ierr = MatViewFromOptions(C2,NULL,"-view_conv_assembled");CHKERRQ(ierr); 210 ierr = MatAXPY(C2,-1.,C,DIFFERENT_NONZERO_PATTERN);CHKERRQ(ierr); 211 ierr = MatChop(C2,PETSC_SMALL);CHKERRQ(ierr); 212 ierr = MatViewFromOptions(C2,NULL,"-view_err");CHKERRQ(ierr); 213 ierr = MatDestroy(&C2);CHKERRQ(ierr); 214 SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Error in conversion from BAIJ to MATIS"); 215 } 216 ierr = MatDestroy(&A);CHKERRQ(ierr); 217 } 218 ierr = MatDestroy(&C);CHKERRQ(ierr); 219 220 ierr = PetscFinalize(); 221 return ierr; 222 } 223 224 /*TEST 225 226 test: 227 228 test: 229 suffix: 2 230 nsize: 3 231 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 246 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 261 262 TEST*/ 263