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