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