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