Actual source code: ex72.c
petsc-3.9.4 2018-09-11
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: ./ex72 -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_type superlu or superlu_dist or mumps \n\
20: -ksp_type preonly -pc_type cholesky -pc_factor_mat_solver_type mumps \n\
21: mpiexec -n <np> ./ex72 -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: T*/
33: /*
34: Include "petscksp.h" so that we can use KSP solvers. Note that this file
35: automatically includes:
36: petscsys.h - base PETSc routines petscvec.h - vectors
37: petscmat.h - matrices
38: petscis.h - index sets petscksp.h - Krylov subspace methods
39: petscviewer.h - viewers petscpc.h - preconditioners
40: */
41: #include <petscksp.h>
43: int main(int argc,char **args)
44: {
45: KSP ksp; /* linear solver context */
46: Mat A; /* matrix */
47: Vec x,b,u; /* approx solution, RHS, exact solution */
48: PetscViewer viewer; /* viewer */
49: char file[4][PETSC_MAX_PATH_LEN]; /* input file name */
50: PetscBool table =PETSC_FALSE,flg,trans=PETSC_FALSE,initialguess = PETSC_FALSE;
51: PetscBool outputSoln=PETSC_FALSE,constantnullspace = PETSC_FALSE;
53: PetscInt its,num_numfac,m,n,M,nearnulldim = 0;
54: PetscReal norm;
55: PetscBool preload=PETSC_TRUE,isSymmetric,cknorm=PETSC_FALSE,initialguessfile = PETSC_FALSE;
56: PetscMPIInt rank;
57: char initialguessfilename[PETSC_MAX_PATH_LEN];
59: PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
60: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
61: PetscOptionsGetBool(NULL,NULL,"-table",&table,NULL);
62: PetscOptionsGetBool(NULL,NULL,"-constantnullspace",&constantnullspace,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,&viewer);
108: /*
109: Load the matrix and vector; then destroy the viewer.
110: */
111: MatCreate(PETSC_COMM_WORLD,&A);
112: MatSetFromOptions(A);
113: MatLoad(A,viewer);
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],viewer);
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: }
129: if (constantnullspace) {
130: MatNullSpace constant;
131: MatNullSpaceCreate(PETSC_COMM_WORLD,PETSC_TRUE,0,NULL,&constant);
132: MatSetNullSpace(A,constant);
133: MatNullSpaceDestroy(&constant);
134: }
135: flg = PETSC_FALSE;
136: PetscOptionsGetString(NULL,NULL,"-rhs",file[2],PETSC_MAX_PATH_LEN,&flg);
137: VecCreate(PETSC_COMM_WORLD,&b);
138: if (flg) { /* rhs is stored in a separate file */
139: if (file[2][0] == '0' || file[2][0] == 0) {
140: PetscInt m;
141: PetscScalar one = 1.0;
142: PetscInfo(0,"Using vector of ones for RHS\n");
143: MatGetLocalSize(A,&m,NULL);
144: VecSetSizes(b,m,PETSC_DECIDE);
145: VecSetFromOptions(b);
146: VecSet(b,one);
147: } else {
148: PetscViewerDestroy(&viewer);
149: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[2],FILE_MODE_READ,&viewer);
150: VecSetFromOptions(b);
151: VecLoad(b,viewer);
152: }
153: } else { /* rhs is stored in the same file as matrix */
154: VecSetFromOptions(b);
155: VecLoad(b,viewer);
156: }
157: PetscViewerDestroy(&viewer);
159: /* Make A singular for testing zero-pivot of ilu factorization */
160: /* Example: ./ex72 -f0 <datafile> -test_zeropivot -pc_factor_shift_type <shift_type> */
161: flg = PETSC_FALSE;
162: PetscOptionsGetBool(NULL,NULL, "-test_zeropivot", &flg,NULL);
163: if (flg) { /* set a row as zeros */
164: PetscInt row=0;
165: MatSetOption(A,MAT_KEEP_NONZERO_PATTERN,PETSC_TRUE);
166: MatZeroRows(A,1,&row,0.0,NULL,NULL);
167: }
169: /* Check whether A is symmetric, then set A->symmetric option */
170: flg = PETSC_FALSE;
171: PetscOptionsGetBool(NULL,NULL, "-check_symmetry", &flg,NULL);
172: if (flg) {
173: MatIsSymmetric(A,0.0,&isSymmetric);
174: if (!isSymmetric) {
175: PetscPrintf(PETSC_COMM_WORLD,"Warning: A is non-symmetric \n");
176: }
177: }
179: /*
180: If the loaded matrix is larger than the vector (due to being padded
181: to match the block size of the system), then create a new padded vector.
182: */
184: MatGetLocalSize(A,&m,&n);
185: /* if (m != n) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_SIZ, "This example is not intended for rectangular matrices (%d, %d)", m, n);*/
186: MatGetSize(A,&M,NULL);
187: VecGetSize(b,&m);
188: if (M != m) { /* Create a new vector b by padding the old one */
189: PetscInt j,mvec,start,end,indx;
190: Vec tmp;
191: PetscScalar *bold;
193: VecCreate(PETSC_COMM_WORLD,&tmp);
194: VecSetSizes(tmp,n,PETSC_DECIDE);
195: VecSetFromOptions(tmp);
196: VecGetOwnershipRange(b,&start,&end);
197: VecGetLocalSize(b,&mvec);
198: VecGetArray(b,&bold);
199: for (j=0; j<mvec; j++) {
200: indx = start+j;
201: VecSetValues(tmp,1,&indx,bold+j,INSERT_VALUES);
202: }
203: VecRestoreArray(b,&bold);
204: VecDestroy(&b);
205: VecAssemblyBegin(tmp);
206: VecAssemblyEnd(tmp);
207: b = tmp;
208: }
210: MatCreateVecs(A,&x,NULL);
211: VecDuplicate(b,&u);
212: if (initialguessfile) {
213: PetscViewerBinaryOpen(PETSC_COMM_WORLD,initialguessfilename,FILE_MODE_READ,&viewer);
214: VecLoad(x,viewer);
215: PetscViewerDestroy(&viewer);
216: initialguess = PETSC_TRUE;
217: } else if (initialguess) {
218: VecSet(x,1.0);
219: } else {
220: VecSet(x,0.0);
221: }
224: /* Check scaling in A */
225: flg = PETSC_FALSE;
226: PetscOptionsGetBool(NULL,NULL, "-check_scaling", &flg,NULL);
227: if (flg) {
228: Vec max, min;
229: PetscInt idx;
230: PetscReal val;
232: VecDuplicate(x, &max);
233: VecDuplicate(x, &min);
234: MatGetRowMaxAbs(A, max, NULL);
235: MatGetRowMinAbs(A, min, NULL);
236: {
237: PetscViewerASCIIOpen(PETSC_COMM_WORLD, "max.data", &viewer);
238: VecView(max, viewer);
239: PetscViewerDestroy(&viewer);
240: PetscViewerASCIIOpen(PETSC_COMM_WORLD, "min.data", &viewer);
241: VecView(min, viewer);
242: PetscViewerDestroy(&viewer);
243: }
244: VecView(max, PETSC_VIEWER_DRAW_WORLD);
245: VecMax(max, &idx, &val);
246: PetscPrintf(PETSC_COMM_WORLD, "Largest max row element %g at row %D\n", (double)val, idx);
247: VecView(min, PETSC_VIEWER_DRAW_WORLD);
248: VecMin(min, &idx, &val);
249: PetscPrintf(PETSC_COMM_WORLD, "Smallest min row element %g at row %D\n", (double)val, idx);
250: VecMin(max, &idx, &val);
251: PetscPrintf(PETSC_COMM_WORLD, "Smallest max row element %g at row %D\n", (double)val, idx);
252: VecPointwiseDivide(max, max, min);
253: VecMax(max, &idx, &val);
254: PetscPrintf(PETSC_COMM_WORLD, "Largest row ratio %g at row %D\n", (double)val, idx);
255: VecView(max, PETSC_VIEWER_DRAW_WORLD);
256: VecDestroy(&max);
257: VecDestroy(&min);
258: }
260: /* MatView(A,PETSC_VIEWER_STDOUT_WORLD); */
261: /* - - - - - - - - - - - New Stage - - - - - - - - - - - - -
262: Setup solve for system
263: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
264: /*
265: Conclude profiling last stage; begin profiling next stage.
266: */
267: PetscPreLoadStage("KSPSetUpSolve");
269: /*
270: Create linear solver; set operators; set runtime options.
271: */
272: KSPCreate(PETSC_COMM_WORLD,&ksp);
273: KSPSetInitialGuessNonzero(ksp,initialguess);
274: num_numfac = 1;
275: PetscOptionsGetInt(NULL,NULL,"-num_numfac",&num_numfac,NULL);
276: while (num_numfac--) {
277: PetscBool lsqr;
278: char str[32];
279: PetscOptionsGetString(NULL,NULL,"-ksp_type",str,32,&lsqr);
280: if (lsqr) {
281: PetscStrcmp("lsqr",str,&lsqr);
282: }
283: if (lsqr) {
284: Mat BtB;
285: MatTransposeMatMult(A,A,MAT_INITIAL_MATRIX,4,&BtB);
286: KSPSetOperators(ksp,A,BtB);
287: MatDestroy(&BtB);
288: } else {
289: KSPSetOperators(ksp,A,A);
290: }
291: KSPSetFromOptions(ksp);
293: /*
294: Here we explicitly call KSPSetUp() and KSPSetUpOnBlocks() to
295: enable more precise profiling of setting up the preconditioner.
296: These calls are optional, since both will be called within
297: KSPSolve() if they haven't been called already.
298: */
299: KSPSetUp(ksp);
300: KSPSetUpOnBlocks(ksp);
302: /* - - - - - - - - - - - New Stage - - - - - - - - - - - - -
303: Solve system
304: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
306: /*
307: Solve linear system;
308: */
309: if (trans) {
310: KSPSolveTranspose(ksp,b,x);
311: KSPGetIterationNumber(ksp,&its);
312: } else {
313: PetscInt num_rhs=1;
314: PetscOptionsGetInt(NULL,NULL,"-num_rhs",&num_rhs,NULL);
315: cknorm = PETSC_FALSE;
316: PetscOptionsGetBool(NULL,NULL,"-cknorm",&cknorm,NULL);
317: while (num_rhs--) {
318: if (num_rhs == 1) VecSet(x,0.0);
319: KSPSolve(ksp,b,x);
320: }
321: KSPGetIterationNumber(ksp,&its);
322: if (cknorm) { /* Check error for each rhs */
323: if (trans) {
324: MatMultTranspose(A,x,u);
325: } else {
326: MatMult(A,x,u);
327: }
328: VecAXPY(u,-1.0,b);
329: VecNorm(u,NORM_2,&norm);
330: PetscPrintf(PETSC_COMM_WORLD," Number of iterations = %3D\n",its);
331: if (!PetscIsNanScalar(norm)) {
332: if (norm < 1.e-12) {
333: PetscPrintf(PETSC_COMM_WORLD," Residual norm < 1.e-12\n");
334: } else {
335: PetscPrintf(PETSC_COMM_WORLD," Residual norm %g\n",(double)norm);
336: }
337: }
338: }
339: } /* while (num_rhs--) */
341: /* - - - - - - - - - - - New Stage - - - - - - - - - - - - -
342: Check error, print output, free data structures.
343: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
345: /*
346: Check error
347: */
348: if (trans) {
349: MatMultTranspose(A,x,u);
350: } else {
351: MatMult(A,x,u);
352: }
353: VecAXPY(u,-1.0,b);
354: VecNorm(u,NORM_2,&norm);
355: /*
356: Write output (optinally using table for solver details).
357: - PetscPrintf() handles output for multiprocessor jobs
358: by printing from only one processor in the communicator.
359: - KSPView() prints information about the linear solver.
360: */
361: if (table) {
362: char *matrixname,kspinfo[120];
364: /*
365: Open a string viewer; then write info to it.
366: */
367: PetscViewerStringOpen(PETSC_COMM_WORLD,kspinfo,120,&viewer);
368: KSPView(ksp,viewer);
369: PetscStrrchr(file[PetscPreLoadIt],'/',&matrixname);
370: PetscPrintf(PETSC_COMM_WORLD,"%-8.8s %3D %2.0e %s \n",matrixname,its,norm,kspinfo);
372: /*
373: Destroy the viewer
374: */
375: PetscViewerDestroy(&viewer);
376: } else {
377: PetscPrintf(PETSC_COMM_WORLD,"Number of iterations = %3D\n",its);
378: if (!PetscIsNanScalar(norm)) {
379: if (norm < 1.e-12 && !PetscIsNanScalar((PetscScalar)norm)) {
380: PetscPrintf(PETSC_COMM_WORLD," Residual norm < 1.e-12\n");
381: } else {
382: PetscPrintf(PETSC_COMM_WORLD,"Residual norm %g\n",(double)norm);
383: }
384: }
385: }
386: PetscOptionsGetString(NULL,NULL,"-solution",file[3],PETSC_MAX_PATH_LEN,&flg);
387: if (flg) {
388: Vec xstar;
389: PetscReal norm;
391: PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[3],FILE_MODE_READ,&viewer);
392: VecCreate(PETSC_COMM_WORLD,&xstar);
393: VecLoad(xstar,viewer);
394: VecAXPY(xstar, -1.0, x);
395: VecNorm(xstar, NORM_2, &norm);
396: PetscPrintf(PETSC_COMM_WORLD, "Error norm %g\n", (double)norm);
397: VecDestroy(&xstar);
398: PetscViewerDestroy(&viewer);
399: }
400: if (outputSoln) {
401: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"solution.petsc",FILE_MODE_WRITE,&viewer);
402: VecView(x, viewer);
403: PetscViewerDestroy(&viewer);
404: }
406: flg = PETSC_FALSE;
407: PetscOptionsGetBool(NULL,NULL, "-ksp_reason", &flg,NULL);
408: if (flg) {
409: KSPConvergedReason reason;
410: KSPGetConvergedReason(ksp,&reason);
411: PetscPrintf(PETSC_COMM_WORLD,"KSPConvergedReason: %D\n", reason);
412: }
414: } /* while (num_numfac--) */
416: /*
417: Free work space. All PETSc objects should be destroyed when they
418: are no longer needed.
419: */
420: MatDestroy(&A); VecDestroy(&b);
421: VecDestroy(&u); VecDestroy(&x);
422: KSPDestroy(&ksp);
423: PetscPreLoadEnd();
424: /* -----------------------------------------------------------
425: End of linear solver loop
426: ----------------------------------------------------------- */
428: PetscFinalize();
429: return ierr;
430: }
433: /*TEST
435: build:
436: requires: !complex
438: testset:
439: suffix: 1
440: nsize: 2
441: args: -f0 ${wPETSC_DIR}/share/petsc/datafiles/matrices/spd-real-int@PETSC_INDEX_SIZE@-float@PETSC_SCALAR_SIZE@
442: requires: !__float128
444: testset:
445: suffix: 1a
446: args: -f0 ${wPETSC_DIR}/share/petsc/datafiles/matrices/spd-real-int@PETSC_INDEX_SIZE@-float@PETSC_SCALAR_SIZE@
447: requires: !__float128
449: testset:
450: nsize: 2
451: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
452: args: -f0 ${DATAFILESPATH}/matrices/medium
453: args: -ksp_type bicg
454: test:
455: suffix: 2
457: testset:
458: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
459: args: -f0 ${DATAFILESPATH}/matrices/medium
460: args: -ksp_type bicg
461: test:
462: suffix: 4
463: args: -pc_type lu
464: test:
465: suffix: 5
467: testset:
468: suffix: 6
469: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
470: args: -f0 ${DATAFILESPATH}/matrices/fem1
471: args: -pc_factor_levels 2 -pc_factor_fill 1.73 -ksp_gmres_cgs_refinement_type refine_always
473: testset:
474: suffix: 7
475: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
476: args: -f0 ${DATAFILESPATH}/matrices/medium
477: args: -viewer_binary_skip_info -mat_type seqbaij
478: args: -matload_block_size {{2 3 4 5 6 7 8}separate output}
479: args: -ksp_max_it 100 -ksp_gmres_cgs_refinement_type refine_always
480: args: -ksp_rtol 1.0e-15 -ksp_monitor_short
481: test:
482: suffix: a
483: test:
484: suffix: b
485: args: -pc_factor_mat_ordering_type nd
486: test:
487: suffix: c
488: args: -pc_factor_levels 1
490: testset:
491: suffix: 7_d
492: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
493: args: -f0 ${DATAFILESPATH}/matrices/medium
494: args: -viewer_binary_skip_info -mat_type seqbaij
495: args: -matload_block_size {{2 3 4 5 6 7 8}shared output}
496: args: -ksp_type preonly -pc_type lu
498: testset:
499: suffix: 8
500: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
501: args: -f0 ${DATAFILESPATH}/matrices/medium
502: args: -ksp_diagonal_scale -pc_type eisenstat -ksp_monitor_short -ksp_diagonal_scale_fix -ksp_gmres_cgs_refinement_type refine_always -mat_no_inode
504: testset:
505: suffix: 9
506: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
507: args: -f0 ${DATAFILESPATH}/matrices/medium
508: 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
509: test:
510: suffix: a
511: args: -mat_type seqbaij
512: test:
513: suffix: b
514: args: -mat_type seqbaij -trans
515: test:
516: suffix: c
517: nsize: 2
518: args: -mat_type mpibaij
519: test:
520: suffix: d
521: nsize: 2
522: args: -mat_type mpibaij -trans
523: test:
524: suffix: e
525: nsize: 3
526: args: -mat_type mpibaij
527: test:
528: suffix: f
529: nsize: 3
530: args: -mat_type mpibaij -trans
533: testset:
534: suffix: 10
535: nsize: 2
536: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
537: args: -ksp_type fgmres -pc_type ksp -f0 ${DATAFILESPATH}/matrices/medium -ksp_fgmres_modifypcksp -ksp_monitor_short
539: testset:
540: TODO: Need to determine goal of this test
541: suffix: 11
542: nsize: 2
543: args: -f0 http://ftp.mcs.anl.gov/pub/petsc/Datafiles/matrices/testmatrix.gz
545: testset:
546: suffix: 12
547: requires: matlab
548: args: -pc_type lu -pc_factor_mat_solver_type matlab -f0 ${DATAFILESPATH}/matrices/arco1
550: testset:
551: suffix: 13
552: requires: lusol
553: args: -f0 ${DATAFILESPATH}/matrices/arco1
554: args: -mat_type lusol -pc_type lu
556: testset:
557: nsize: 3
558: args: -f0 ${DATAFILESPATH}/matrices/medium
559: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
560: test:
561: suffix: 14
562: requires: spai
563: args: -pc_type spai
564: test:
565: suffix: 15
566: requires: hypre
567: args: -pc_type hypre -pc_hypre_type pilut
568: test:
569: suffix: 16
570: requires: hypre
571: args: -pc_type hypre -pc_hypre_type parasails
572: test:
573: suffix: 17
574: requires: hypre
575: args: -pc_type hypre -pc_hypre_type boomeramg
577: testset:
578: suffix: 19
579: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
580: args: -f0 ${DATAFILESPATH}/matrices/poisson1
581: args: -ksp_type cg -pc_type icc
582: args: -pc_factor_levels {{0 2 4}separate output}
583: test:
584: test:
585: args: -mat_type seqsbaij
588: testset:
589: suffix: ILU
590: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
591: args: -f0 ${DATAFILESPATH}/matrices/small
592: args: -pc_factor_levels 1
593: test:
594: test:
595: # This is tested against regular ILU (used to be denoted ILUBAIJ)
596: args: -mat_type baij
598: testset:
599: suffix: aijcusparse
600: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) cusparse
601: args: -f0 ${DATAFILESPATH}/matrices/medium -ksp_monitor_short -ksp_view -mat_view ascii::ascii_info -mat_type aijcusparse -pc_factor_mat_solver_type cusparse -pc_type ilu
603: testset:
604: TODO: No output file. Need to determine if deprecated
605: suffix: asm_viennacl
606: nsize: 2
607: requires: viennacl
608: 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}
610: testset:
611: nsize: 2
612: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) hypre
613: args: -f0 ${DATAFILESPATH}/matrices/poisson2.gz -ksp_monitor_short -ksp_rtol 1.E-9 -pc_type hypre -pc_hypre_type boomeramg
614: test:
615: suffix: boomeramg_euclid
616: 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
617: TODO: Need to determine if deprecated
618: test:
619: suffix: boomeramg_euclid_bj
620: 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
621: TODO: Need to determine if deprecated
622: test:
623: suffix: boomeramg_parasails
624: args: -pc_hypre_boomeramg_smooth_type ParaSails -pc_hypre_boomeramg_smooth_num_levels 2
625: test:
626: suffix: boomeramg_pilut
627: args: -pc_hypre_boomeramg_smooth_type Pilut -pc_hypre_boomeramg_smooth_num_levels 2
628: test:
629: suffix: boomeramg_schwarz
630: args: -pc_hypre_boomeramg_smooth_type Schwarz-smoothers
632: testset:
633: suffix: cg_singlereduction
634: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
635: args: -f0 ${DATAFILESPATH}/matrices/small
636: args: -mat_type mpisbaij -ksp_type cg -pc_type eisenstat -ksp_monitor_short -ksp_converged_reason
637: test:
638: test:
639: args: -ksp_cg_single_reduction
641: testset:
642: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
643: args: -f0 ${DATAFILESPATH}/matrices/poisson2.gz
644: args: -ksp_monitor_short -pc_type icc
645: test:
646: suffix: cr
647: args: -ksp_type cr
648: test:
649: suffix: lcd
650: args: -ksp_type lcd
652: testset:
653: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
654: args: -f0 ${DATAFILESPATH}/matrices/small
655: args: -ksp_monitor_short -ksp_view -mat_view ascii::ascii_info
656: test:
657: suffix: seqaijcrl
658: args: -mat_type seqaijcrl
659: test:
660: suffix: seqaijperm
661: args: -mat_type seqaijperm
663: testset:
664: nsize: 2
665: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
666: args: -f0 ${DATAFILESPATH}/matrices/small
667: args: -ksp_monitor_short -ksp_view
668: # Different output files
669: test:
670: suffix: mpiaijcrl
671: args: -mat_type mpiaijcrl
672: test:
673: suffix: mpiaijperm
674: args: -mat_type mpiaijperm
676: testset:
677: nsize: 8
678: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES)
679: args: -ksp_monitor_short -ksp_view
680: test:
681: suffix: xxt
682: args: -f0 ${DATAFILESPATH}/matrices/poisson1 -check_symmetry -ksp_type cg -pc_type tfs
683: test:
684: suffix: xyt
685: args: -f0 ${DATAFILESPATH}/matrices/arco1 -ksp_type gmres -pc_type tfs
687: testset:
688: # The output file here is the same as mumps
689: suffix: mumps_cholesky
690: output_file: output/ex72_mumps.out
691: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) mumps
692: args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type cholesky -pc_factor_mat_solver_type mumps -num_numfac 2 -num_rhs 2
693: nsize: {{1 2}}
694: test:
695: args: -mat_type sbaij -mat_ignore_lower_triangular
696: test:
697: args: -mat_type aij
698: test:
699: args: -mat_type aij -matload_spd
701: testset:
702: # The output file here is the same as mumps
703: suffix: mumps_lu
704: output_file: output/ex72_mumps.out
705: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) mumps
706: args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -pc_factor_mat_solver_type mumps -num_numfac 2 -num_rhs 2
707: test:
708: args: -mat_type seqaij
709: test:
710: nsize: 2
711: args: -mat_type mpiaij
712: test:
713: args: -mat_type seqbaij -matload_block_size 2
714: test:
715: nsize: 2
716: args: -mat_type mpibaij -matload_block_size 2
717: test:
718: args: -mat_type aij -mat_mumps_icntl_7 5
719: TODO: Need to determine if deprecated
720: test:
721: nsize: 2
722: args: -mat_type mpiaij -mat_mumps_icntl_28 2 -mat_mumps_icntl_29 2
723: TODO: Need to determine if deprecated
725: testset:
726: # The output file here is the same as mumps
727: suffix: mumps_redundant
728: output_file: output/ex72_mumps_redundant.out
729: nsize: 8
730: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) mumps
731: 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_type mumps -num_numfac 2 -num_rhs 2
733: testset:
734: suffix: pastix_cholesky
735: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) pastix
736: output_file: output/ex72_mumps.out
737: nsize: {{1 2}}
738: args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_factor_mat_solver_type pastix -num_numfac 2 -num_rhs 2 -pc_type cholesky -mat_type sbaij -mat_ignore_lower_triangular
740: testset:
741: suffix: pastix_lu
742: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) pastix
743: args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -pc_factor_mat_solver_type pastix -num_numfac 2 -num_rhs 2
744: output_file: output/ex72_mumps.out
745: test:
746: args: -mat_type seqaij
747: test:
748: nsize: 2
749: args: -mat_type mpiaij
751: testset:
752: suffix: pastix_redundant
753: output_file: output/ex72_mumps_redundant.out
754: nsize: 8
755: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) pastix
756: 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_type pastix -num_numfac 2 -num_rhs 2
759: testset:
760: suffix: superlu_dist_lu
761: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) superlu_dist
762: output_file: output/ex72_mumps.out
763: args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -pc_factor_mat_solver_type superlu_dist -num_numfac 2 -num_rhs 2
764: nsize: {{1 2}}
766: testset:
767: suffix: superlu_dist_redundant
768: nsize: 8
769: output_file: output/ex72_mumps_redundant.out
770: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) superlu_dist
771: 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_type superlu_dist -num_numfac 2 -num_rhs 2
773: testset:
774: suffix: superlu_lu
775: output_file: output/ex72_mumps.out
776: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) superlu
777: args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -pc_factor_mat_solver_type superlu -num_numfac 2 -num_rhs 2
779: testset:
780: suffix: umfpack
781: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) suitesparse
782: args: -f0 ${DATAFILESPATH}/matrices/small -ksp_type preonly -pc_type lu -mat_type seqaij -pc_factor_mat_solver_type umfpack -num_numfac 2 -num_rhs 2
785: testset:
786: suffix: zeropivot
787: requires: datafilespath double !define(PETSC_USE_64BIT_INDICES) mumps
788: args: -f0 ${DATAFILESPATH}/matrices/small -test_zeropivot -ksp_converged_reason -ksp_type fgmres -pc_type ksp
789: test:
790: nsize: 3
791: args: -ksp_pc_type bjacobi
792: test:
793: nsize: 2
794: args: -ksp_ksp_type cg -ksp_pc_type bjacobi -ksp_pc_bjacobi_blocks 1
795: #test:
796: #nsize: 3
797: #args: -ksp_ksp_converged_reason -ksp_pc_type bjacobi -ksp_sub_ksp_converged_reason
798: #TODO: Need to determine if deprecated
799: TEST*/