Actual source code: ex10.c

petsc-master 2017-05-24
Report Typos and Errors

  2: static char help[] = "Reads a PETSc matrix and vector from a file and solves a linear system.\n\
  3: This version first preloads and solves a small system, then loads \n\
  4: another (larger) system and solves it as well.  This example illustrates\n\
  5: preloading of instructions with the smaller system so that more accurate\n\
  6: performance monitoring can be done with the larger one (that actually\n\
  7: is the system of interest).  See the 'Performance Hints' chapter of the\n\
  8: users manual for a discussion of preloading.  Input parameters include\n\
  9:   -f0 <input_file> : first file to load (small system)\n\
 10:   -f1 <input_file> : second file to load (larger system)\n\n\
 11:   -nearnulldim <0> : number of vectors in the near-null space immediately following matrix\n\n\
 12:   -trans  : solve transpose system instead\n\n";
 13: /*
 14:   This code can be used to test PETSc interface to other packages.\n\
 15:   Examples of command line options:       \n\
 16:    ./ex10 -f0 <datafile> -ksp_type preonly  \n\
 17:         -help -ksp_view                  \n\
 18:         -num_numfac <num_numfac> -num_rhs <num_rhs> \n\
 19:         -ksp_type preonly -pc_type lu -pc_factor_mat_solver_package superlu or superlu_dist or mumps \n\
 20:         -ksp_type preonly -pc_type cholesky -pc_factor_mat_solver_package mumps \n\
 21:    mpiexec -n <np> ./ex10 -f0 <datafile> -ksp_type cg -pc_type asm -pc_asm_type basic -sub_pc_type icc -mat_type sbaij
 22:  \n\n";
 23: */
 24: /*T
 25:    Concepts: KSP^solving a linear system
 26:    Processors: n
 27:    requires: !complex 
 28: T*/





 34: /*
 35:   Include "petscksp.h" so that we can use KSP solvers.  Note that this file
 36:   automatically includes:
 37:      petscsys.h       - base PETSc routines   petscvec.h - vectors
 38:      petscmat.h - matrices
 39:      petscis.h     - index sets            petscksp.h - Krylov subspace methods
 40:      petscviewer.h - viewers               petscpc.h  - preconditioners
 41: */
 42:  #include <petscksp.h>

 44: int main(int argc,char **args)
 45: {
 46:   KSP            ksp;             /* linear solver context */
 47:   Mat            A;               /* matrix */
 48:   Vec            x,b,u;           /* approx solution, RHS, exact solution */
 49:   PetscViewer    fd;              /* viewer */
 50:   char           file[4][PETSC_MAX_PATH_LEN];     /* input file name */
 51:   PetscBool      table     =PETSC_FALSE,flg,trans=PETSC_FALSE,initialguess = PETSC_FALSE;
 52:   PetscBool      outputSoln=PETSC_FALSE;
 54:   PetscInt       its,num_numfac,m,n,M,nearnulldim = 0;
 55:   PetscReal      norm;
 56:   PetscBool      preload=PETSC_TRUE,isSymmetric,cknorm=PETSC_FALSE,initialguessfile = PETSC_FALSE;
 57:   PetscMPIInt    rank;
 58:   char           initialguessfilename[PETSC_MAX_PATH_LEN];

 60:   PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
 61:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 62:   PetscOptionsGetBool(NULL,NULL,"-table",&table,NULL);
 63:   PetscOptionsGetBool(NULL,NULL,"-trans",&trans,NULL);
 64:   PetscOptionsGetBool(NULL,NULL,"-initialguess",&initialguess,NULL);
 65:   PetscOptionsGetBool(NULL,NULL,"-output_solution",&outputSoln,NULL);
 66:   PetscOptionsGetString(NULL,NULL,"-initialguessfilename",initialguessfilename,PETSC_MAX_PATH_LEN,&initialguessfile);
 67:   PetscOptionsGetInt(NULL,NULL,"-nearnulldim",&nearnulldim,NULL);

 69:   /*
 70:      Determine files from which we read the two linear systems
 71:      (matrix and right-hand-side vector).
 72:   */
 73:   PetscOptionsGetString(NULL,NULL,"-f",file[0],PETSC_MAX_PATH_LEN,&flg);
 74:   if (flg) {
 75:     PetscStrcpy(file[1],file[0]);
 76:     preload = PETSC_FALSE;
 77:   } else {
 78:     PetscOptionsGetString(NULL,NULL,"-f0",file[0],PETSC_MAX_PATH_LEN,&flg);
 79:     if (!flg) SETERRQ(PETSC_COMM_WORLD,1,"Must indicate binary file with the -f0 or -f option");
 80:     PetscOptionsGetString(NULL,NULL,"-f1",file[1],PETSC_MAX_PATH_LEN,&flg);
 81:     if (!flg) preload = PETSC_FALSE;   /* don't bother with second system */
 82:   }

 84:   /* -----------------------------------------------------------
 85:                   Beginning of linear solver loop
 86:      ----------------------------------------------------------- */
 87:   /*
 88:      Loop through the linear solve 2 times.
 89:       - The intention here is to preload and solve a small system;
 90:         then load another (larger) system and solve it as well.
 91:         This process preloads the instructions with the smaller
 92:         system so that more accurate performance monitoring (via
 93:         -log_view) can be done with the larger one (that actually
 94:         is the system of interest).
 95:   */
 96:   PetscPreLoadBegin(preload,"Load system");

 98:   /* - - - - - - - - - - - New Stage - - - - - - - - - - - - -
 99:                          Load system
100:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

102:   /*
103:      Open binary file.  Note that we use FILE_MODE_READ to indicate
104:      reading from this file.
105:   */
106:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[PetscPreLoadIt],FILE_MODE_READ,&fd);

108:   /*
109:      Load the matrix and vector; then destroy the viewer.
110:   */
111:   MatCreate(PETSC_COMM_WORLD,&A);
112:   MatSetFromOptions(A);
113:   MatLoad(A,fd);
114:   if (nearnulldim) {
115:     MatNullSpace nullsp;
116:     Vec          *nullvecs;
117:     PetscInt     i;
118:     PetscMalloc1(nearnulldim,&nullvecs);
119:     for (i=0; i<nearnulldim; i++) {
120:       VecCreate(PETSC_COMM_WORLD,&nullvecs[i]);
121:       VecLoad(nullvecs[i],fd);
122:     }
123:     MatNullSpaceCreate(PETSC_COMM_WORLD,PETSC_FALSE,nearnulldim,nullvecs,&nullsp);
124:     MatSetNearNullSpace(A,nullsp);
125:     for (i=0; i<nearnulldim; i++) {VecDestroy(&nullvecs[i]);}
126:     PetscFree(nullvecs);
127:     MatNullSpaceDestroy(&nullsp);
128:   }

130:   flg  = PETSC_FALSE;
131:   PetscOptionsGetString(NULL,NULL,"-rhs",file[2],PETSC_MAX_PATH_LEN,&flg);
132:   VecCreate(PETSC_COMM_WORLD,&b);
133:   if (flg) {   /* rhs is stored in a separate file */
134:     if (file[2][0] == '0' || file[2][0] == 0) {
135:       PetscInt    m;
136:       PetscScalar one = 1.0;
137:       PetscInfo(0,"Using vector of ones for RHS\n");
138:       MatGetLocalSize(A,&m,NULL);
139:       VecSetSizes(b,m,PETSC_DECIDE);
140:       VecSetFromOptions(b);
141:       VecSet(b,one);
142:     } else {
143:       PetscViewerDestroy(&fd);
144:       PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[2],FILE_MODE_READ,&fd);
145:       VecSetFromOptions(b);
146:       VecLoad(b,fd);
147:     }
148:   } else {   /* rhs is stored in the same file as matrix */
149:     VecSetFromOptions(b);
150:     VecLoad(b,fd);
151:   }
152:   PetscViewerDestroy(&fd);

154:   /* Make A singular for testing zero-pivot of ilu factorization        */
155:   /* Example: ./ex10 -f0 <datafile> -test_zeropivot -set_row_zero -pc_factor_shift_type <shift_type> */
156:   flg  = PETSC_FALSE;
157:   PetscOptionsGetBool(NULL,NULL, "-test_zeropivot", &flg,NULL);
158:   if (flg) {
159:     PetscInt          row=0;
160:     PetscBool         flg1=PETSC_FALSE;
161:     PetscScalar       zero=0.0;

163:     PetscOptionsGetBool(NULL,NULL, "-set_row_zero", &flg1,NULL);
164:     if (flg1) {   /* set a row as zeros */
165:       MatSetOption(A,MAT_KEEP_NONZERO_PATTERN,PETSC_TRUE);
166:       MatZeroRows(A,1,&row,0.0,NULL,NULL);
167:     } else if (!rank) {
168:       /* set (row,row) entry as zero */
169:       MatSetValues(A,1,&row,1,&row,&zero,INSERT_VALUES);
170:     }
171:     MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
172:     MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
173:   }

175:   /* Check whether A is symmetric */
176:   flg  = PETSC_FALSE;
177:   PetscOptionsGetBool(NULL,NULL, "-check_symmetry", &flg,NULL);
178:   if (flg) {
179:     Mat Atrans;
180:     MatTranspose(A, MAT_INITIAL_MATRIX,&Atrans);
181:     MatEqual(A, Atrans, &isSymmetric);
182:     if (isSymmetric) {
183:       MatSetOption(A,MAT_SYMMETRIC,PETSC_TRUE);
184:     } else {
185:       PetscPrintf(PETSC_COMM_WORLD,"Warning: A is non-symmetric \n");
186:     }
187:     MatDestroy(&Atrans);
188:   }

190:   /*
191:      If the loaded matrix is larger than the vector (due to being padded
192:      to match the block size of the system), then create a new padded vector.
193:   */

195:   MatGetLocalSize(A,&m,&n);
196:   /*  if (m != n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ, "This example is not intended for rectangular matrices (%d, %d)", m, n);*/
197:   MatGetSize(A,&M,NULL);
198:   VecGetSize(b,&m);
199:   if (M != m) {   /* Create a new vector b by padding the old one */
200:     PetscInt    j,mvec,start,end,indx;
201:     Vec         tmp;
202:     PetscScalar *bold;

204:     VecCreate(PETSC_COMM_WORLD,&tmp);
205:     VecSetSizes(tmp,n,PETSC_DECIDE);
206:     VecSetFromOptions(tmp);
207:     VecGetOwnershipRange(b,&start,&end);
208:     VecGetLocalSize(b,&mvec);
209:     VecGetArray(b,&bold);
210:     for (j=0; j<mvec; j++) {
211:       indx = start+j;
212:       VecSetValues(tmp,1,&indx,bold+j,INSERT_VALUES);
213:     }
214:     VecRestoreArray(b,&bold);
215:     VecDestroy(&b);
216:     VecAssemblyBegin(tmp);
217:     VecAssemblyEnd(tmp);
218:     b    = tmp;
219:   }

221:   MatCreateVecs(A,&x,NULL);
222:   VecDuplicate(b,&u);
223:   if (initialguessfile) {
224:     PetscViewer viewer2;
225:     PetscViewerBinaryOpen(PETSC_COMM_WORLD,initialguessfilename,FILE_MODE_READ,&viewer2);
226:     VecLoad(x,viewer2);
227:     PetscViewerDestroy(&viewer2);
228:     initialguess = PETSC_TRUE;
229:   } else if (initialguess) {
230:     VecSet(x,1.0);
231:   } else {
232:     VecSet(x,0.0);
233:   }


236:   /* Check scaling in A */
237:   flg  = PETSC_FALSE;
238:   PetscOptionsGetBool(NULL,NULL, "-check_scaling", &flg,NULL);
239:   if (flg) {
240:     Vec       max, min;
241:     PetscInt  idx;
242:     PetscReal val;

244:     VecDuplicate(x, &max);
245:     VecDuplicate(x, &min);
246:     MatGetRowMaxAbs(A, max, NULL);
247:     MatGetRowMinAbs(A, min, NULL);
248:     {
249:       PetscViewer viewer;

251:       PetscViewerASCIIOpen(PETSC_COMM_WORLD, "max.data", &viewer);
252:       VecView(max, viewer);
253:       PetscViewerDestroy(&viewer);
254:       PetscViewerASCIIOpen(PETSC_COMM_WORLD, "min.data", &viewer);
255:       VecView(min, viewer);
256:       PetscViewerDestroy(&viewer);
257:     }
258:     VecView(max, PETSC_VIEWER_DRAW_WORLD);
259:     VecMax(max, &idx, &val);
260:     PetscPrintf(PETSC_COMM_WORLD, "Largest max row element %g at row %D\n", (double)val, idx);
261:     VecView(min, PETSC_VIEWER_DRAW_WORLD);
262:     VecMin(min, &idx, &val);
263:     PetscPrintf(PETSC_COMM_WORLD, "Smallest min row element %g at row %D\n", (double)val, idx);
264:     VecMin(max, &idx, &val);
265:     PetscPrintf(PETSC_COMM_WORLD, "Smallest max row element %g at row %D\n", (double)val, idx);
266:     VecPointwiseDivide(max, max, min);
267:     VecMax(max, &idx, &val);
268:     PetscPrintf(PETSC_COMM_WORLD, "Largest row ratio %g at row %D\n", (double)val, idx);
269:     VecView(max, PETSC_VIEWER_DRAW_WORLD);
270:     VecDestroy(&max);
271:     VecDestroy(&min);
272:   }

274:   /*  MatView(A,PETSC_VIEWER_STDOUT_WORLD); */
275:   /* - - - - - - - - - - - New Stage - - - - - - - - - - - - -
276:                     Setup solve for system
277:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
278:   /*
279:      Conclude profiling last stage; begin profiling next stage.
280:   */
281:   PetscPreLoadStage("KSPSetUpSolve");

283:   /*
284:      Create linear solver; set operators; set runtime options.
285:   */
286:   KSPCreate(PETSC_COMM_WORLD,&ksp);
287:   KSPSetInitialGuessNonzero(ksp,initialguess);
288:   num_numfac = 1;
289:   PetscOptionsGetInt(NULL,NULL,"-num_numfac",&num_numfac,NULL);
290:   while (num_numfac--) {
291:     PetscBool lsqr;
292:     char      str[32];
293:     PetscOptionsGetString(NULL,NULL,"-ksp_type",str,32,&lsqr);
294:     if (lsqr) {
295:       PetscStrcmp("lsqr",str,&lsqr);
296:     }
297:     if (lsqr) {
298:       Mat BtB;
299:       MatTransposeMatMult(A,A,MAT_INITIAL_MATRIX,4,&BtB);
300:       KSPSetOperators(ksp,A,BtB);
301:       MatDestroy(&BtB);
302:     } else {
303:       KSPSetOperators(ksp,A,A);
304:     }
305:     KSPSetFromOptions(ksp);

307:     /*
308:      Here we explicitly call KSPSetUp() and KSPSetUpOnBlocks() to
309:      enable more precise profiling of setting up the preconditioner.
310:      These calls are optional, since both will be called within
311:      KSPSolve() if they haven't been called already.
312:     */
313:     KSPSetUp(ksp);
314:     KSPSetUpOnBlocks(ksp);

316:     /* - - - - - - - - - - - New Stage - - - - - - - - - - - - -
317:                          Solve system
318:       - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

320:     /*
321:      Solve linear system; 
322:     */
323:     if (trans) {
324:       KSPSolveTranspose(ksp,b,x);
325:       KSPGetIterationNumber(ksp,&its);
326:     } else {
327:       PetscInt num_rhs=1;
328:       PetscOptionsGetInt(NULL,NULL,"-num_rhs",&num_rhs,NULL);
329:       cknorm = PETSC_FALSE;
330:       PetscOptionsGetBool(NULL,NULL,"-cknorm",&cknorm,NULL);
331:       while (num_rhs--) {
332:         if (num_rhs == 1) VecSet(x,0.0);
333:         KSPSolve(ksp,b,x);
334:       }
335:       KSPGetIterationNumber(ksp,&its);
336:       if (cknorm) {     /* Check error for each rhs */
337:         if (trans) {
338:           MatMultTranspose(A,x,u);
339:         } else {
340:           MatMult(A,x,u);
341:         }
342:         VecAXPY(u,-1.0,b);
343:         VecNorm(u,NORM_2,&norm);
344:         PetscPrintf(PETSC_COMM_WORLD,"  Number of iterations = %3D\n",its);
345:         if (!PetscIsNanScalar(norm)) {
346:           if (norm < 1.e-12) {
347:             PetscPrintf(PETSC_COMM_WORLD,"  Residual norm < 1.e-12\n");
348:           } else {
349:             PetscPrintf(PETSC_COMM_WORLD,"  Residual norm %g\n",(double)norm);
350:           }
351:         }
352:       }
353:     }   /* while (num_rhs--) */

355:     /* - - - - - - - - - - - New Stage - - - - - - - - - - - - -
356:           Check error, print output, free data structures.
357:      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */

359:     /*
360:        Check error
361:     */
362:     if (trans) {
363:       MatMultTranspose(A,x,u);
364:     } else {
365:       MatMult(A,x,u);
366:     }
367:     VecAXPY(u,-1.0,b);
368:     VecNorm(u,NORM_2,&norm);
369:     /*
370:      Write output (optinally using table for solver details).
371:       - PetscPrintf() handles output for multiprocessor jobs
372:         by printing from only one processor in the communicator.
373:       - KSPView() prints information about the linear solver.
374:     */
375:     if (table) {
376:       char        *matrixname,kspinfo[120];
377:       PetscViewer viewer;

379:       /*
380:        Open a string viewer; then write info to it.
381:       */
382:       PetscViewerStringOpen(PETSC_COMM_WORLD,kspinfo,120,&viewer);
383:       KSPView(ksp,viewer);
384:       PetscStrrchr(file[PetscPreLoadIt],'/',&matrixname);
385:       PetscPrintf(PETSC_COMM_WORLD,"%-8.8s %3D %2.0e %s \n",matrixname,its,norm,kspinfo);

387:       /*
388:         Destroy the viewer
389:       */
390:       PetscViewerDestroy(&viewer);
391:     } else {
392:       PetscPrintf(PETSC_COMM_WORLD,"Number of iterations = %3D\n",its);
393:       if (!PetscIsNanScalar(norm)) {
394:         if (norm < 1.e-12 && !PetscIsNanScalar((PetscScalar)norm)) {
395:           PetscPrintf(PETSC_COMM_WORLD,"  Residual norm < 1.e-12\n");
396:         } else {
397:           PetscPrintf(PETSC_COMM_WORLD,"Residual norm %g\n",(double)norm);
398:         }
399:       }
400:     }
401:     PetscOptionsGetString(NULL,NULL,"-solution",file[3],PETSC_MAX_PATH_LEN,&flg);
402:     if (flg) {
403:       PetscViewer viewer;
404:       Vec         xstar;
405:       PetscReal   norm;

407:       PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[3],FILE_MODE_READ,&viewer);
408:       VecCreate(PETSC_COMM_WORLD,&xstar);
409:       VecLoad(xstar,viewer);
410:       VecAXPY(xstar, -1.0, x);
411:       VecNorm(xstar, NORM_2, &norm);
412:       PetscPrintf(PETSC_COMM_WORLD, "Error norm %g\n", (double)norm);
413:       VecDestroy(&xstar);
414:       PetscViewerDestroy(&viewer);
415:     }
416:     if (outputSoln) {
417:       PetscViewer viewer;

419:       PetscViewerBinaryOpen(PETSC_COMM_WORLD,"solution.petsc",FILE_MODE_WRITE,&viewer);
420:       VecView(x, viewer);
421:       PetscViewerDestroy(&viewer);
422:     }

424:     flg  = PETSC_FALSE;
425:     PetscOptionsGetBool(NULL,NULL, "-ksp_reason", &flg,NULL);
426:     if (flg) {
427:       KSPConvergedReason reason;
428:       KSPGetConvergedReason(ksp,&reason);
429:       PetscPrintf(PETSC_COMM_WORLD,"KSPConvergedReason: %D\n", reason);
430:     }

432:   }   /* while (num_numfac--) */

434:   /*
435:      Free work space.  All PETSc objects should be destroyed when they
436:      are no longer needed.
437:   */
438:   MatDestroy(&A); VecDestroy(&b);
439:   VecDestroy(&u); VecDestroy(&x);
440:   KSPDestroy(&ksp);
441:   PetscPreLoadEnd();
442:   /* -----------------------------------------------------------
443:                       End of linear solver loop
444:      ----------------------------------------------------------- */

446:   PetscFinalize();
447:   return ierr;
448: }


451: /*TEST
452:    
453:    testset:
454:       suffix: 1
455:       nsize: 2
456:       args: -f0 ${wPETSC_DIR}/share/petsc/datafiles/matrices/spd-real-int@PETSC_INDEX_SIZE@-float@PETSC_SCALAR_SIZE@

458:    testset:
459:       suffix: 1a
460:       args: -f0 ${wPETSC_DIR}/share/petsc/datafiles/matrices/spd-real-int@PETSC_INDEX_SIZE@-float@PETSC_SCALAR_SIZE@

462:    testset:
463:       nsize: 2
464:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
465:       args: -f0 ${DATAFILESPATH}/matrices/medium
466:       args:  -ksp_type bicg
467:       test:
468:          suffix: 2
469:       test:
470:          suffix: 3
471:          args: -pc_type asm 
472:       
473:    testset:
474:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
475:       args: -f0 ${DATAFILESPATH}/matrices/medium
476:       args: -ksp_type bicg
477:       test:
478:          suffix: 4
479:          args: -pc_type lu 
480:       test:
481:          suffix: 5
482:    
483:    testset:
484:       suffix: 6
485:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
486:       args: -f0 ${DATAFILESPATH}/matrices/fem1
487:       args: -pc_factor_levels 2 -pc_factor_fill 1.73 -ksp_gmres_cgs_refinement_type refine_always
488:    
489:    testset:
490:       suffix: 7
491:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
492:       args: -f0 ${DATAFILESPATH}/matrices/medium 
493:       args: -viewer_binary_skip_info -mat_type seqbaij 
494:       args: -matload_block_size {{2 3 4 5 6 7 8}separate output}
495:       args: -ksp_max_it 100 -ksp_gmres_cgs_refinement_type refine_always 
496:       args: -ksp_rtol 1.0e-15 -ksp_monitor_short
497:       test:
498:          suffix: a
499:       test:
500:          suffix: b
501:          args: -pc_factor_mat_ordering_type nd
502:       test:
503:          suffix: c
504:          args: -pc_factor_levels 1

506:    testset:
507:       suffix: 7_d
508:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
509:       args: -f0 ${DATAFILESPATH}/matrices/medium 
510:       args: -viewer_binary_skip_info -mat_type seqbaij 
511:       args: -matload_block_size {{2 3 4 5 6 7 8}shared output}
512:       args: -ksp_type preonly -pc_type lu
513:    
514:    testset:
515:       suffix: 8
516:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
517:       args: -f0 ${DATAFILESPATH}/matrices/medium  
518:       args: -ksp_diagonal_scale -pc_type eisenstat -ksp_monitor_short -ksp_diagonal_scale_fix -ksp_gmres_cgs_refinement_type refine_always -mat_no_inode
519:    
520:    testset:
521:       suffix: 9
522:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
523:       args: -f0 ${DATAFILESPATH}/matrices/medium 
524:       args: -viewer_binary_skip_info  -matload_block_size {{1 2 3 4 5 6 7}separate output} -ksp_max_it 100 -ksp_gmres_cgs_refinement_type refine_always -ksp_rtol 1.0e-15 -ksp_monitor_short
525:       test:
526:          suffix: a
527:          args: -mat_type seqbaij
528:       test:
529:          suffix: b
530:          args: -mat_type seqbaij -trans
531:       test:
532:          suffix: c
533:          nsize: 2
534:          args: -mat_type mpibaij 
535:       test:
536:          suffix: d
537:          nsize: 2
538:          args: -mat_type mpibaij -trans
539:       test:
540:          suffix: e
541:          nsize: 3
542:          args: -mat_type mpibaij 
543:       test:
544:          suffix: f
545:          nsize: 3
546:          args: -mat_type mpibaij -trans
547:    

549:    testset:
550:       suffix: 10
551:       nsize: 2
552:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
553:       args: -ksp_type fgmres -pc_type ksp -f0 ${DATAFILESPATH}/matrices/medium -ksp_fgmres_modifypcksp -ksp_monitor_short
554:    
555:    testset:
556:       TODO: Need to determine goal of this test
557:       suffix: 11
558:       nsize: 2
559:       args: -f0 http://ftp.mcs.anl.gov/pub/petsc/matrices/testmatrix.gz
560:    
561:    testset:
562:       suffix: 12
563:       requires: matlab
564:       args: -pc_type lu -pc_factor_mat_solver_package matlab -f0 ${DATAFILESPATH}/matrices/arco1
565:    
566:    testset:
567:       suffix: 13
568:       requires: lusol
569:       args: -f0 ${DATAFILESPATH}/matrices/arco1
570:       args: -mat_type lusol -pc_type lu
571:    
572:    testset:
573:       nsize: 3
574:       args: -f0 ${DATAFILESPATH}/matrices/medium
575:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
576:       test:
577:          suffix: 14
578:          requires: spai 
579:          args: -pc_type spai 
580:       test:
581:          suffix: 15
582:          requires: hypre
583:          args: -pc_type hypre -pc_hypre_type pilut 
584:       test:
585:          suffix: 16
586:          requires: hypre
587:          args: -pc_type hypre -pc_hypre_type parasails
588:       test:
589:          suffix: 17
590:          requires: hypre
591:          args: -pc_type hypre -pc_hypre_type boomeramg 
592:    
593:    testset:
594:       suffix: 19
595:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
596:       args: -f0 ${DATAFILESPATH}/matrices/poisson1
597:       args: -ksp_type cg -pc_type icc 
598:       args: -pc_factor_levels {{0 2 4}separate output}
599:       test:
600:       test:
601:          args: -mat_type seqsbaij
602:    

604:    testset:
605:       suffix: ILU
606:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
607:       args: -f0 ${DATAFILESPATH}/matrices/small
608:       args: -pc_factor_levels 1
609:       test:
610:       test:
611:          # This is tested against regular ILU (used to be denoted ILUBAIJ)
612:          args: -mat_type baij
613:    
614:    testset:
615:       suffix: aijcusparse
616:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) cusparse
617:       args: -f0 ${DATAFILESPATH}/matrices/medium -ksp_monitor_short -ksp_view -mat_view ascii::ascii_info -mat_type aijcusparse -pc_factor_mat_solver_package cusparse -pc_type ilu
618:    
619:    testset:
620:       TODO: No output file. Need to determine if deprecated
621:       suffix: asm_viennacl
622:       nsize: 2
623:       requires: viennacl
624:       args: -pc_type asm -pc_asm_sub_mat_type aijviennacl -f0 ${wPETSC_DIR}/share/petsc/datafiles/matrices/spd-real-int${PETSC_INDEX_SIZE}-float${PETSC_SCALAR_SIZE}
625:    
626:    testset:
627:       nsize: 2
628:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) hypre
629:       args: -f0 ${DATAFILESPATH}/matrices/poisson2.gz -ksp_monitor_short -ksp_rtol 1.E-9 -pc_type hypre -pc_hypre_type boomeramg 
630:       test:
631:          suffix: boomeramg_euclid
632:          args: -pc_hypre_boomeramg_smooth_type Euclid -pc_hypre_boomeramg_smooth_num_levels 2 -pc_hypre_boomeramg_eu_level 1 -pc_hypre_boomeramg_eu_droptolerance 0.01 
633:          TODO: Need to determine if deprecated
634:       test:
635:          suffix: boomeramg_euclid_bj
636:          args: -pc_hypre_boomeramg_smooth_type Euclid -pc_hypre_boomeramg_smooth_num_levels 2 -pc_hypre_boomeramg_eu_level 1 -pc_hypre_boomeramg_eu_droptolerance 0.01 -pc_hypre_boomeramg_eu_bj 
637:          TODO: Need to determine if deprecated
638:       test:
639:          suffix: boomeramg_parasails
640:          args: -pc_hypre_boomeramg_smooth_type ParaSails -pc_hypre_boomeramg_smooth_num_levels 2 
641:       test:
642:          suffix: boomeramg_pilut
643:          args: -pc_hypre_boomeramg_smooth_type Pilut -pc_hypre_boomeramg_smooth_num_levels 2 
644:       test:
645:          suffix: boomeramg_schwarz
646:          args: -pc_hypre_boomeramg_smooth_type Schwarz-smoothers 
647:       
648:    testset:
649:       suffix: cg_singlereduction
650:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
651:       args: -f0 ${DATAFILESPATH}/matrices/small
652:       args: -mat_type mpisbaij -ksp_type cg -pc_type eisenstat -ksp_monitor_short -ksp_converged_reason
653:       test:
654:       test:
655:          args: -ksp_cg_single_reduction
656:       
657:    testset:
658:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
659:       args: -f0 ${DATAFILESPATH}/matrices/poisson2.gz
660:       args: -ksp_monitor_short -pc_type icc
661:       test:
662:          suffix: cr
663:          args: -ksp_type cr
664:       test:
665:          suffix: lcd
666:          args: -ksp_type lcd
667:    
668:    testset:
669:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
670:       args: -f0 ${DATAFILESPATH}/matrices/small
671:       args: -ksp_monitor_short -ksp_view -mat_view ascii::ascii_info 
672:       test:
673:          suffix: seqaijcrl
674:          args: -mat_type seqaijcrl
675:       test:
676:          suffix: seqaijperm
677:          args: -mat_type seqaijperm

679:    testset:
680:       nsize: 2
681:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
682:       args: -f0 ${DATAFILESPATH}/matrices/small
683:       args: -ksp_monitor_short -ksp_view 
684:       # Different output files
685:       test:
686:          suffix: mpiaijcrl
687:          args: -mat_type mpiaijcrl
688:       test:
689:          suffix: mpiaijperm
690:          args: -mat_type mpiaijperm
691:    
692:    testset:
693:       nsize: 8
694:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
695:       args: -ksp_monitor_short -ksp_view
696:       test:
697:          suffix: xxt
698:          args: -f0 ${DATAFILESPATH}/matrices/poisson1 -check_symmetry -ksp_type cg -pc_type tfs
699:       test:
700:          suffix: xyt
701:          args: -f0 ${DATAFILESPATH}/matrices/arco1 -ksp_type gmres -pc_type tfs

703:    testset:
704:       # The output file here is the same as mumps
705:       suffix: mumps_cholesky
706:       output_file: output/ex10_mumps.out
707:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) mumps
708:       args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type cholesky -pc_factor_mat_solver_package mumps -num_numfac 2 -num_rhs 2
709:       nsize: {{1 2}}
710:       test:
711:          args: -mat_type sbaij -mat_ignore_lower_triangular
712:       test:
713:          args: -mat_type aij
714:       test:
715:          args: -mat_type aij -matload_spd
716:       
717:    testset:
718:       # The output file here is the same as mumps
719:       suffix: mumps_lu
720:       output_file: output/ex10_mumps.out
721:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) mumps
722:       args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -pc_factor_mat_solver_package mumps -num_numfac 2 -num_rhs 2
723:       test:
724:          args: -mat_type seqaij
725:       test:
726:          nsize: 2
727:          args: -mat_type mpiaij
728:       test:
729:          args: -mat_type seqbaij -matload_block_size 2
730:       test:
731:          nsize: 2
732:          args: -mat_type mpibaij -matload_block_size 2
733:       test:
734:          args: -mat_type aij -mat_mumps_icntl_7 5
735:          TODO: Need to determine if deprecated
736:       test:
737:          nsize: 2
738:          args: -mat_type mpiaij -mat_mumps_icntl_28 2 -mat_mumps_icntl_29 2
739:          TODO: Need to determine if deprecated
740:       
741:    testset:
742:       # The output file here is the same as mumps
743:       suffix: mumps_redundant
744:       output_file: output/ex10_mumps_redundant.out
745:       nsize: 8
746:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) mumps
747:       args: -f0 ${DATAFILESPATH}/matrices/medium -ksp_type preonly -pc_type redundant -pc_redundant_number {{8 7 6 5 4 3 2 1}} -redundant_pc_factor_mat_solver_package mumps -num_numfac 2 -num_rhs 2
748:    
749:    testset:
750:       suffix: pastix_cholesky
751:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) pastix
752:       output_file: output/ex10_mumps.out
753:       nsize: {{1 2}}
754:       args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_factor_mat_solver_package pastix -num_numfac 2 -num_rhs 2 -pc_type cholesky -mat_type sbaij -mat_ignore_lower_triangular
755:       
756:    testset:
757:       suffix: pastix_lu
758:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) pastix
759:       args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -pc_factor_mat_solver_package pastix -num_numfac 2 -num_rhs 2
760:       output_file: output/ex10_mumps.out
761:       test:
762:          args: -mat_type seqaij
763:       test:
764:          nsize: 2
765:          args: -mat_type mpiaij
766:    
767:    testset:
768:       suffix: pastix_redundant
769:       output_file: output/ex10_mumps_redundant.out
770:       nsize: 8
771:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) pastix
772:       args: -f0 ${DATAFILESPATH}/matrices/medium -ksp_type preonly -pc_type redundant -pc_redundant_number {{8 7 6 5 4 3 2 1}} -redundant_pc_factor_mat_solver_package pastix -num_numfac 2 -num_rhs 2
773:    
774:      
775:    testset:
776:       suffix: superlu_dist_lu
777:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) superlu_dist
778:       output_file: output/ex10_mumps.out
779:       args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -pc_factor_mat_solver_package superlu_dist -num_numfac 2 -num_rhs 2
780:       nsize: {{1 2}}
781:       
782:    testset:
783:       suffix: superlu_dist_redundant
784:       nsize: 8
785:       output_file: output/ex10_mumps_redundant.out
786:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) superlu
787:       args: -f0 ${DATAFILESPATH}/matrices/medium -ksp_type preonly -pc_type redundant -pc_redundant_number {{8 7 6 5 4 3 2 1}} -redundant_pc_factor_mat_solver_package superlu_dist -num_numfac 2 -num_rhs 2
788:    
789:    testset:
790:       suffix: superlu_lu
791:       output_file: output/ex10_mumps.out
792:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) superlu
793:       args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -pc_factor_mat_solver_package superlu -num_numfac 2 -num_rhs 2
794:    
795:    testset:
796:       suffix: umfpack
797:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) suitesparse
798:       args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -mat_type seqaij -pc_factor_mat_solver_package umfpack -num_numfac 2 -num_rhs 2
799:    
800:      
801:    testset:
802:       suffix: zeropivot
803:       requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) mumps
804:       args: -f0 ${DATAFILESPATH}/matrices/small -test_zeropivot -ksp_converged_reason -ksp_type fgmres -pc_type ksp
805:       test:
806:          nsize: 3
807:          args: -ksp_pc_type bjacobi
808:       test:
809:          nsize: 2
810:          args: -ksp_ksp_type cg -ksp_pc_type bjacobi -ksp_pc_bjacobi_blocks 1
811:       #test:
812:          #nsize: 3
813:          #args: -ksp_ksp_converged_reason -ksp_pc_type bjacobi -ksp_sub_ksp_converged_reason
814:          #TODO: Need to determine if deprecated
815: TEST*/