Actual source code: ex237.c

petsc-master 2020-08-25
Report Typos and Errors
  1: static char help[] = "Mini-app to benchmark matrix--matrix multiplication\n\n";

  3: /*
  4:   See the paper below for more information

  6:    "KSPHPDDM and PCHPDDM: Extending PETSc with Robust Overlapping Schwarz Preconditioners and Advanced Krylov Methods",
  7:    P. Jolivet, J. E. Roman, and S. Zampini (2020).
  8: */

 10:  #include <petsc.h>

 12: #if defined(PETSC_HAVE_MKL)
 13: #include <mkl.h>
 14: #define PetscStackCallMKLSparse(func, args) do {               \
 15:     sparse_status_t __ierr;                                    \
 16:     PetscStackPush(#func);                                     \
 17:     __func args;                                        \
 18:     PetscStackPop;                                             \
 19:     if (__ierr != SPARSE_STATUS_SUCCESS) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_LIB, "Error in %s(): error code %d", #func, (int)__ierr); \
 20:   } while (0)
 21: #else
 22: #define PetscStackCallMKLSparse(func, args) do {               \
 23:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "No MKL support"); \
 24:   } while (0)
 25: #endif

 27: int main(int argc, char** argv) {
 28:   Mat            A, C, D, E;
 29:   PetscInt       nbs = 10, ntype = 10, nN = 8, m, M, trial = 5;
 30:   PetscViewer    viewer;
 31:   PetscInt       bs[10], N[8];
 32:   char           *type[10];
 33:   PetscMPIInt    size;
 34:   PetscBool      flg, cuda, maij = PETSC_FALSE, check = PETSC_FALSE, trans = PETSC_FALSE, convert = PETSC_FALSE, mkl;
 35:   char           file[PETSC_MAX_PATH_LEN];

 38:   PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr;
 39:   if (ierr) return ierr;
 40:   MPI_Comm_size(PETSC_COMM_WORLD, &size);
 41:   if (size != 1) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only");
 42:   PetscOptionsGetString(NULL, NULL, "-f", file, PETSC_MAX_PATH_LEN, &flg);
 43:   if (!flg) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_USER_INPUT, "Must indicate binary file with the -f option");
 44:   PetscOptionsGetInt(NULL, NULL, "-trial", &trial, NULL);
 45:   PetscOptionsGetIntArray(NULL, NULL, "-bs", bs, &nbs, &flg);
 46:   if (!flg) {
 47:     nbs = 1;
 48:     bs[0] = 1;
 49:   }
 50:   PetscOptionsGetIntArray(NULL, NULL, "-N", N, &nN, &flg);
 51:   if (!flg) {
 52:     nN = 8;
 53:     N[0] = 1;  N[1] = 2;  N[2] = 4;  N[3] = 8;
 54:     N[4] = 16; N[5] = 32; N[6] = 64; N[7] = 128;
 55:   }
 56:   PetscOptionsGetStringArray(NULL, NULL, "-type", type, &ntype, &flg);
 57:   if (!flg) {
 58:     ntype = 1;
 59:     PetscStrallocpy(MATSEQAIJ, &type[0]);
 60:   }
 61:   PetscOptionsGetBool(NULL, NULL, "-check", &check, NULL);
 62:   PetscOptionsGetBool(NULL, NULL, "-trans", &trans, NULL);
 63:   PetscOptionsGetBool(NULL, NULL, "-convert_aij", &convert, NULL);
 64:   for (PetscInt j = 0; j < nbs; ++j) {
 65:     PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &viewer);
 66:     MatCreate(PETSC_COMM_WORLD, &A);
 67:     MatSetFromOptions(A);
 68:     MatLoad(A, viewer);
 69:     PetscViewerDestroy(&viewer);
 70:     PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJ, MATMPIAIJ, "");
 71:     if (!flg) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_USER_INPUT, "Must indicate a MatAIJ input matrix");
 72:     MatGetSize(A, &m, &M);
 73:     if (m == M) {
 74:       Mat oA;
 75:       MatTranspose(A, MAT_INITIAL_MATRIX, &oA);
 76:       MatAXPY(A, 1.0, oA, DIFFERENT_NONZERO_PATTERN);
 77:       MatDestroy(&oA);
 78:     }
 79:     MatGetLocalSize(A, &m, NULL);
 80:     MatGetSize(A, &M, NULL);
 81:     if (bs[j] > 1) {
 82:       Mat               T, Tt, B;
 83:       const PetscScalar *ptr;
 84:       PetscScalar       *val, *Aa;
 85:       const PetscInt    *Ai, *Aj;
 86:       PetscInt          An, i, k;
 87:       PetscBool         done;

 89:       MatCreateDense(PETSC_COMM_SELF, bs[j], bs[j], bs[j], bs[j], NULL, &T);
 90:       MatSetRandom(T, NULL);
 91:       MatTranspose(T, MAT_INITIAL_MATRIX, &Tt);
 92:       MatAXPY(T, 1.0, Tt, SAME_NONZERO_PATTERN);
 93:       MatDestroy(&Tt);
 94:       MatDenseGetArrayRead(T, &ptr);
 95:       MatGetRowIJ(A, 0, PETSC_FALSE, PETSC_FALSE, &An, &Ai, &Aj, &done);
 96:       if (!done || An != m) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Inconsistent sizes");
 97:       MatSeqAIJGetArray(A, &Aa);
 98:       MatCreate(PETSC_COMM_WORLD, &B);
 99:       MatSetType(B, MATSEQBAIJ);
100:       MatSetSizes(B, bs[j] * An, bs[j] * An, PETSC_DECIDE, PETSC_DECIDE);
101:       PetscMalloc1(Ai[An] * bs[j] * bs[j],&val);
102:       for (i = 0; i < Ai[An]; ++i)
103:         for (k = 0; k < bs[j] * bs[j]; ++k)
104:           val[i * bs[j] * bs[j] + k] = Aa[i] * ptr[k];
105:       MatSetOption(B, MAT_ROW_ORIENTED, PETSC_FALSE);
106:       MatSeqBAIJSetPreallocationCSR(B, bs[j], Ai, Aj, val);
107:       PetscFree(val);
108:       MatSeqAIJRestoreArray(A, &Aa);
109:       MatRestoreRowIJ(A, 0, PETSC_FALSE, PETSC_FALSE, &An, &Ai, &Aj, &done);
110:       MatDenseRestoreArrayRead(T, &ptr);
111:       MatDestroy(&T);
112:       MatDestroy(&A);
113:       A    = B;
114:     }
115:     /* reconvert back to SeqAIJ before converting to the desired type later */
116:     if (!convert) E = A;
117:     MatConvert(A, MATSEQAIJ, convert ? MAT_INITIAL_MATRIX : MAT_INPLACE_MATRIX, &E);
118:     MatSetOption(E, MAT_SYMMETRIC, PETSC_TRUE);
119:     for (PetscInt i = 0; i < ntype; ++i) {
120:       char        *tmp;
121:       PetscInt    *ia_ptr, *ja_ptr, k;
122:       PetscScalar *a_ptr;
123: #if defined(PETSC_HAVE_MKL)
124:       struct matrix_descr descr;
125:       sparse_matrix_t     spr;
126:       descr.type = SPARSE_MATRIX_TYPE_GENERAL;
127:       descr.diag = SPARSE_DIAG_NON_UNIT;
128: #endif
129:       if (convert) {
130:         MatDestroy(&A);
131:       }
132:       PetscStrstr(type[i], "mkl", &tmp);
133:       if (tmp) {
134:         size_t mlen, tlen;
135:         char base[256];

137:         mkl  = PETSC_TRUE;
138:         PetscStrlen(tmp, &mlen);
139:         PetscStrlen(type[i], &tlen);
140:         PetscStrncpy(base, type[i], tlen-mlen + 1);
141:         MatConvert(E, base, convert ? MAT_INITIAL_MATRIX : MAT_INPLACE_MATRIX, &A);
142:       } else {
143:         mkl  = PETSC_FALSE;
144:         PetscStrstr(type[i], "maij", &tmp);
145:         if (!tmp) {
146:           MatConvert(E, type[i], convert ? MAT_INITIAL_MATRIX : MAT_INPLACE_MATRIX, &A);
147:         } else {
148:           MatConvert(E, MATAIJ, convert ? MAT_INITIAL_MATRIX : MAT_INPLACE_MATRIX, &A);
149:           maij = PETSC_TRUE;
150:         }
151:       }
152:       PetscObjectTypeCompareAny((PetscObject)A, &cuda, MATSEQAIJCUSPARSE, MATMPIAIJCUSPARSE, "");
153:       if (mkl) {
154:         const PetscInt *Ai, *Aj;
155:         PetscInt       An;
156:         PetscBool      done;

158:         PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJ, MATSEQBAIJ, MATSEQSBAIJ, "");
159:         if (!flg) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_USER_INPUT, "Not implemented");
160:         PetscObjectTypeCompare((PetscObject)A, MATSEQAIJ, &flg);
161:         MatGetRowIJ(A, 0, PETSC_FALSE, flg ? PETSC_FALSE : PETSC_TRUE, &An, &Ai, &Aj, &done);
162:         if (!done) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Inconsistent sizes");
163:         PetscMalloc1(An + 1,&ia_ptr);
164:         PetscMalloc1(Ai[An],&ja_ptr);
165:         if (flg) { /* SeqAIJ */
166:           for (k = 0; k < An + 1; ++k) ia_ptr[k] = Ai[k];
167:           for (k = 0; k < Ai[An]; ++k) ja_ptr[k] = Aj[k];
168:           MatSeqAIJGetArray(A, &a_ptr);
169:           PetscStackCallMKLSparse(mkl_sparse_d_create_csr, (&spr, SPARSE_INDEX_BASE_ZERO, An, An, ia_ptr, ia_ptr + 1, ja_ptr, a_ptr));
170:         } else {
171:           PetscObjectTypeCompare((PetscObject)A, MATSEQBAIJ, &flg);
172:           if (flg) {
173:             for (k = 0; k < An + 1; ++k) ia_ptr[k] = Ai[k] + 1; /* Fortran indexing to maximize cases covered by _mm routines */
174:             for (k = 0; k < Ai[An]; ++k) ja_ptr[k] = Aj[k] + 1; /* Fortran indexing to maximize cases covered by _mm routines */
175:             MatSeqBAIJGetArray(A, &a_ptr);
176:             PetscStackCallMKLSparse(mkl_sparse_d_create_bsr, (&spr, SPARSE_INDEX_BASE_ONE, SPARSE_LAYOUT_COLUMN_MAJOR, An, An, bs[j], ia_ptr, ia_ptr + 1, ja_ptr, a_ptr));
177:           } else {
178:             PetscObjectTypeCompare((PetscObject)A, MATSEQSBAIJ, &flg);
179:             if (flg) {
180:               for (k = 0; k < An + 1; ++k) ia_ptr[k] = Ai[k] + 1; /* Fortran indexing to maximize cases covered by _mm routines */
181:               for (k = 0; k < Ai[An]; ++k) ja_ptr[k] = Aj[k] + 1; /* Fortran indexing to maximize cases covered by _mm routines */
182:               MatSeqSBAIJGetArray(A, &a_ptr);
183:               PetscStackCallMKLSparse(mkl_sparse_d_create_bsr, (&spr, SPARSE_INDEX_BASE_ONE, SPARSE_LAYOUT_COLUMN_MAJOR, An, An, bs[j], ia_ptr, ia_ptr + 1, ja_ptr, a_ptr));
184: #if defined(PETSC_HAVE_MKL)
185:               descr.type = SPARSE_MATRIX_TYPE_SYMMETRIC;
186:               descr.mode = SPARSE_FILL_MODE_UPPER;
187:               descr.diag = SPARSE_DIAG_NON_UNIT;
188: #endif
189:             }
190:           }
191:         }
192:         PetscObjectTypeCompare((PetscObject)A, MATSEQAIJ, &flg);
193:         MatRestoreRowIJ(A, 0, PETSC_FALSE, flg ? PETSC_FALSE : PETSC_TRUE, &An, &Ai, &Aj, &done);
194:       }

196:       MatViewFromOptions(A, NULL, "-A_view");

198:       for (k = 0; k < nN; ++k) {
199:         MatType       Atype, Ctype;
200:         PetscInt      AM, AN, CM, CN, t;
201:         PetscLogStage stage, tstage;
202:         char          stage_s[256];

204:         MatCreateDense(PETSC_COMM_WORLD, bs[j] * m, PETSC_DECIDE, bs[j] * M, N[k], NULL, &C);
205:         MatCreateDense(PETSC_COMM_WORLD, bs[j] * m, PETSC_DECIDE, bs[j] * M, N[k], NULL, &D);
206:         MatSetRandom(C, NULL);
207:         if (cuda) { /* convert to GPU if needed */
208:           MatConvert(C, MATDENSECUDA, MAT_INPLACE_MATRIX, &C);
209:           MatConvert(D, MATDENSECUDA, MAT_INPLACE_MATRIX, &D);
210:         }
211:         if (mkl) {
212:           if (N[k] > 1) PetscStackCallMKLSparse(mkl_sparse_set_mm_hint, (spr, SPARSE_OPERATION_NON_TRANSPOSE, descr, SPARSE_LAYOUT_COLUMN_MAJOR, N[k], 1 + trial));
213:           else          PetscStackCallMKLSparse(mkl_sparse_set_mv_hint, (spr, SPARSE_OPERATION_NON_TRANSPOSE, descr, 1 + trial));
214:           PetscStackCallMKLSparse(mkl_sparse_set_memory_hint, (spr, SPARSE_MEMORY_AGGRESSIVE));
215:           PetscStackCallMKLSparse(mkl_sparse_optimize, (spr));
216:         }
217:         MatGetType(A, &Atype);
218:         MatGetType(C, &Ctype);
219:         MatGetSize(A, &AM, &AN);
220:         MatGetSize(C, &CM, &CN);

222:         if (!maij || N[k] > 1) {
223:           PetscSNPrintf(stage_s, sizeof(stage_s), "type_%s-bs_%D-N_%02d", type[i], bs[j], (int)N[k]);
224:           PetscLogStageRegister(stage_s, &stage);
225:         }
226:         if (trans && N[k] > 1) {
227:           PetscSNPrintf(stage_s, sizeof(stage_s), "trans_type_%s-bs_%D-N_%02d", type[i], bs[j], (int)N[k]);
228:           PetscLogStageRegister(stage_s, &tstage);
229:         }

231:         /* A*B */
232:         if (N[k] > 1) {
233:           if (!maij) {
234:             MatProductCreateWithMat(A, C, NULL, D);
235:             MatProductSetType(D, MATPRODUCT_AB);
236:             MatProductSetFromOptions(D);
237:             MatProductSymbolic(D);
238:           }

240:           if (!mkl) {
241:             if (!maij) {
242:               MatProductNumeric(D);
243:               PetscPrintf(PETSC_COMM_WORLD, "Benchmarking MatProduct %s: with A %s %Dx%D and B %s %Dx%D\n", MatProductTypes[MATPRODUCT_AB], Atype, AM, AN, Ctype, CM, CN);
244:               PetscLogStagePush(stage);
245:               for (t = 0; t < trial; ++t) {
246:                 MatProductNumeric(D);
247:               }
248:               PetscLogStagePop();
249:             } else {
250:               Mat               E, Ct, Dt;
251:               Vec               cC, cD;
252:               const PetscScalar *c_ptr;
253:               PetscScalar       *d_ptr;
254:               MatCreateMAIJ(A, N[k], &E);
255:               MatDenseGetLocalMatrix(C, &Ct);
256:               MatDenseGetLocalMatrix(D, &Dt);
257:               MatTranspose(Ct, MAT_INPLACE_MATRIX, &Ct);
258:               MatTranspose(Dt, MAT_INPLACE_MATRIX, &Dt);
259:               MatDenseGetArrayRead(Ct, &c_ptr);
260:               MatDenseGetArrayWrite(Dt, &d_ptr);
261:               VecCreateMPIWithArray(PETSC_COMM_WORLD, 1, AM * N[k], PETSC_DECIDE, c_ptr, &cC);
262:               VecCreateMPIWithArray(PETSC_COMM_WORLD, 1, AM * N[k], PETSC_DECIDE, d_ptr, &cD);
263:               MatMult(E, cC, cD);
264:               PetscPrintf(PETSC_COMM_WORLD, "Benchmarking MatMult: with A %s %Dx%D and B %s %Dx%D\n", MATMAIJ, AM, AN, VECMPI, AM * N[k], 1);
265:               PetscLogStagePush(stage);
266:               for (t = 0; t < trial; ++t) {
267:                 MatMult(E, cC, cD);
268:               }
269:               PetscLogStagePop();
270:               VecDestroy(&cD);
271:               VecDestroy(&cC);
272:               MatDestroy(&E);
273:               MatDenseRestoreArrayWrite(Dt, &d_ptr);
274:               MatDenseRestoreArrayRead(Ct, &c_ptr);
275:               MatTranspose(Ct, MAT_INPLACE_MATRIX, &Ct);
276:               MatTranspose(Dt, MAT_INPLACE_MATRIX, &Dt);
277:             }
278:           } else {
279:             const PetscScalar *c_ptr;
280:             PetscScalar       *d_ptr;

282:             MatDenseGetArrayRead(C, &c_ptr);
283:             MatDenseGetArrayWrite(D, &d_ptr);
284:             PetscStackCallMKLSparse(mkl_sparse_d_mm,(SPARSE_OPERATION_NON_TRANSPOSE, 1.0, spr, descr, SPARSE_LAYOUT_COLUMN_MAJOR, c_ptr, CN, CM, 0.0, d_ptr, CM));
285:             PetscPrintf(PETSC_COMM_WORLD, "Benchmarking mkl_sparse_d_mm (COLUMN_MAJOR): with A %s %Dx%D and B %s %Dx%D\n", Atype, AM, AN, Ctype, CM, CN);
286:             PetscLogStagePush(stage);
287:             for (t = 0; t < trial; ++t) {
288:               PetscStackCallMKLSparse(mkl_sparse_d_mm, (SPARSE_OPERATION_NON_TRANSPOSE, 1.0, spr, descr, SPARSE_LAYOUT_COLUMN_MAJOR, c_ptr, CN, CM, 0.0, d_ptr, CM));
289:             }
290:             PetscLogStagePop();
291:             MatDenseRestoreArrayWrite(D, &d_ptr);
292:             MatDenseRestoreArrayRead(C, &c_ptr);
293:           }
294:         } else if (maij) {
295:           MatDestroy(&C);
296:           MatDestroy(&D);
297:           continue;
298:         } else if (!mkl) {
299:           Vec cC, cD;

301:           MatDenseGetColumnVecRead(C, 0, &cC);
302:           MatDenseGetColumnVecWrite(D, 0, &cD);
303:           MatMult(A, cC, cD);
304:           PetscPrintf(PETSC_COMM_WORLD, "Benchmarking MatMult: with A %s %Dx%D\n", Atype, AM, AN);
305:           PetscLogStagePush(stage);
306:           for (t = 0; t < trial; ++t) {
307:             MatMult(A, cC, cD);
308:           }
309:           PetscLogStagePop();
310:           MatDenseRestoreColumnVecRead(C, 0, &cC);
311:           MatDenseRestoreColumnVecWrite(D, 0, &cD);
312:         } else {
313:           const PetscScalar *c_ptr;
314:           PetscScalar       *d_ptr;

316:           MatDenseGetArrayRead(C, &c_ptr);
317:           MatDenseGetArrayWrite(D, &d_ptr);
318:           PetscPrintf(PETSC_COMM_WORLD, "Benchmarking mkl_sparse_d_mv: with A %s %Dx%D\n", Atype, AM, AN);
319:           PetscStackCallMKLSparse(mkl_sparse_d_mv, (SPARSE_OPERATION_NON_TRANSPOSE, 1.0, spr, descr, c_ptr, 0.0, d_ptr));
320:           PetscLogStagePush(stage);
321:           for (t = 0; t < trial; ++t) {
322:             PetscStackCallMKLSparse(mkl_sparse_d_mv, (SPARSE_OPERATION_NON_TRANSPOSE, 1.0, spr, descr, c_ptr, 0.0, d_ptr));
323:           }
324:           PetscLogStagePop();
325:           MatDenseRestoreArrayWrite(D, &d_ptr);
326:           MatDenseRestoreArrayRead(C, &c_ptr);
327:         }

329:         if (check) {
330:           MatMatMultEqual(A, C, D, 10, &flg);
331:           if (!flg) {
332:             MatType Dtype;

334:             MatGetType(D, &Dtype);
335:             PetscPrintf(PETSC_COMM_WORLD, "Error with A %s%s, C %s, D %s, Nk %D\n", Atype, mkl ? "mkl" : "", Ctype, Dtype, N[k]);
336:           }
337:         }

339:         /* MKL implementation seems buggy for ABt */
340:         /* A*Bt */
341:         if (!mkl && trans && N[k] > 1) {
342:           PetscObjectTypeCompareAny((PetscObject)A, &flg, MATSEQAIJ, MATMPIAIJ, "");
343:           if (flg) {
344:             MatTranspose(C, MAT_INPLACE_MATRIX, &C);
345:             MatGetType(C, &Ctype);
346:             if (!mkl) {
347:               MatProductCreateWithMat(A, C, NULL, D);
348:               MatProductSetType(D, MATPRODUCT_ABt);
349:               MatProductSetFromOptions(D);
350:               MatProductSymbolic(D);
351:               MatProductNumeric(D);
352:               PetscPrintf(PETSC_COMM_WORLD, "Benchmarking MatProduct %s: with A %s %Dx%D and Bt %s %Dx%D\n", MatProductTypes[MATPRODUCT_ABt], Atype, AM, AN, Ctype, CM, CN);
353:               PetscLogStagePush(tstage);
354:               for (t = 0; t < trial; ++t) {
355:                 MatProductNumeric(D);
356:               }
357:               PetscLogStagePop();
358:             } else {
359:               const PetscScalar *c_ptr;
360:               PetscScalar       *d_ptr;

362:               PetscStackCallMKLSparse(mkl_sparse_set_mm_hint, (spr, SPARSE_OPERATION_NON_TRANSPOSE, descr, SPARSE_LAYOUT_ROW_MAJOR, N[k], 1 + trial));
363:               PetscStackCallMKLSparse(mkl_sparse_optimize, (spr));
364:               MatDenseGetArrayRead(C, &c_ptr);
365:               MatDenseGetArrayWrite(D, &d_ptr);
366:               PetscPrintf(PETSC_COMM_WORLD, "Benchmarking mkl_sparse_d_mm (ROW_MAJOR): with A %s %Dx%D and B %s %Dx%D\n", Atype, AM, AN, Ctype, CM, CN);
367:               PetscStackCallMKLSparse(mkl_sparse_d_mm, (SPARSE_OPERATION_NON_TRANSPOSE, 1.0, spr, descr, SPARSE_LAYOUT_ROW_MAJOR, c_ptr, CN, CM, 0.0, d_ptr, CM));
368:               PetscLogStagePush(stage);
369:               for (t = 0; t < trial; ++t) {
370:                 PetscStackCallMKLSparse(mkl_sparse_d_mm, (SPARSE_OPERATION_NON_TRANSPOSE, 1.0, spr, descr, SPARSE_LAYOUT_ROW_MAJOR, c_ptr, CN, CM, 0.0, d_ptr, CM));
371:               }
372:               PetscLogStagePop();
373:               MatDenseRestoreArrayWrite(D, &d_ptr);
374:               MatDenseRestoreArrayRead(C, &c_ptr);
375:             }
376:           }
377:         }

379:         if (!mkl && trans && N[k] > 1 && flg && check) {
380:           MatMatTransposeMultEqual(A, C, D, 10, &flg);
381:           if (!flg) {
382:             MatType Dtype;
383:             MatGetType(D, &Dtype);
384:             PetscPrintf(PETSC_COMM_WORLD, "Error with A %s%s, C %s, D %s, Nk %D\n", Atype, mkl ? "mkl" : "", Ctype, Dtype, N[k]);
385:           }
386:         }
387:         MatDestroy(&C);
388:         MatDestroy(&D);
389:       }
390:       if (mkl) {
391:         PetscStackCallMKLSparse(mkl_sparse_destroy, (spr));
392:         PetscFree(ia_ptr);
393:         PetscFree(ja_ptr);
394:       }
395:       if (cuda && i != ntype - 1) {
396:         PetscPrintf(PETSC_COMM_WORLD, "AIJCUSPARSE must be last, otherwise MatConvert() to another MatType is too slow\n");
397:         break;
398:       }
399:     }
400:     if (E != A) {
401:       MatDestroy(&E);
402:     }
403:     MatDestroy(&A);
404:   }
405:   for (m = 0; m < ntype; ++m) {
406:     PetscFree(type[m]);
407:   }
408:   PetscFinalize();
409:   return 0;
410: }

412: /*TEST

414:    testset:
415:      nsize: 1
416:      requires: double !complex !define(PETSC_USE_64BIT_INDICES)
417:      filter: sed "/Benchmarking/d"
418:      args: -f ${wPETSC_DIR}/share/petsc/datafiles/matrices/spd-real-int32-float64 -bs 1,2,3 -N 1,2,18 -check -trans -convert_aij {{false true}shared output}
419:      test:
420:        suffix: basic
421:        args: -type aij,sbaij,baij
422:        output_file: output/ex237.out
423:      test:
424:        suffix: maij
425:        args: -type aij,maij
426:        output_file: output/ex237.out
427:      test:
428:        suffix: cuda
429:        requires: cuda
430:        args: -type aij,aijcusparse
431:        output_file: output/ex237.out
432:      test:
433:        suffix: mkl
434:        requires: mkl
435:        args: -type aij,aijmkl,baijmkl,sbaijmkl
436:        output_file: output/ex237.out

438: TEST*/