Actual source code: hyppilut.c
1: /*$Id: bvec2.c,v 1.202 2001/09/12 03:26:24 bsmith Exp $*/
2: /*
4: */
6: #include src/sles/pc/pcimpl.h
7: EXTERN_C_BEGIN
8: #include "HYPRE.h"
9: #include "IJ_mv.h"
10: #include "HYPRE_parcsr_ls.h"
11: EXTERN_C_END
13: extern int MatHYPRE_IJMatrixCreate(Mat,HYPRE_IJMatrix*);
14: extern int MatHYPRE_IJMatrixCopy(Mat,HYPRE_IJMatrix);
15: extern int VecHYPRE_IJVectorCreate(Vec,HYPRE_IJVector*);
17: /*
18: Private context (data structure) for the preconditioner.
19: */
20: typedef struct {
21: HYPRE_Solver hsolver;
22: HYPRE_IJMatrix ij;
23: HYPRE_IJVector b,x;
25: int (*destroy)(HYPRE_Solver);
26: int (*solve)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector);
27: int (*setup)(HYPRE_Solver,HYPRE_ParCSRMatrix,HYPRE_ParVector,HYPRE_ParVector);
29: /* options for pilut and BoomerAMG*/
30: int maxiter;
31: double tol;
33: /* options for pilut */
34: int factorrowsize;
36: /* options for parasails */
37: int nlevels;
38: double threshhold;
39: double filter;
40: int sym;
41: double loadbal;
42: int logging;
43: int ruse;
44: int symt;
46: /* options for euclid */
47: PetscTruth bjilu;
48: int levels;
50: /* options for euclid and BoomerAMG */
51: PetscTruth printstatistics;
53: /* options for BoomerAMG */
54: int maxlevels;
55: double strongthreshold;
56: double maxrowsum;
57: int *gridsweeps;
58: int coarsentype;
59: int measuretype;
60: int *relaxtype;
61: int **gridrelaxpoints;
62: } PC_HYPRE;
65: #undef __FUNCT__
67: static int PCSetUp_HYPRE(PC pc)
68: {
69: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
70: int ierr;
71: HYPRE_ParCSRMatrix hmat;
72: HYPRE_ParVector bv,xv;
75: if (!jac->ij) { /* create the matrix the first time through */
76: MatHYPRE_IJMatrixCreate(pc->pmat,&jac->ij);
77: }
78: if (!jac->b) {
79: VecHYPRE_IJVectorCreate(pc->vec,&jac->b);
80: VecHYPRE_IJVectorCreate(pc->vec,&jac->x);
81: }
82: MatHYPRE_IJMatrixCopy(pc->pmat,jac->ij);
83: HYPRE_IJMatrixGetObject(jac->ij,(void**)&hmat);
84: HYPRE_IJVectorGetObject(jac->b,(void**)&bv);
85: HYPRE_IJVectorGetObject(jac->x,(void**)&xv);
86: (*jac->setup)(jac->hsolver,hmat,bv,xv);
87: if (ierr) SETERRQ1(1,"Error in HYPRE setup, error code %d",ierr);
88: return(0);
89: }
91: /*
92: Replaces the address where the HYPRE vector points to its data with the address of
93: PETSc's data. Saves the old address so it can be reset when we are finished with it.
94: Allows use to get the data into a HYPRE vector without the cost of memcopies
95: */
96: #define HYPREReplacePointer(b,newvalue,savedvalue) {
97: hypre_ParVector *par_vector = (hypre_ParVector *)hypre_IJVectorObject(((hypre_IJVector*)b));
98: hypre_Vector *local_vector = hypre_ParVectorLocalVector(par_vector);
99: savedvalue = local_vector->data;
100: local_vector->data = newvalue;}
102: #undef __FUNCT__
104: static int PCApply_HYPRE(PC pc,Vec b,Vec x)
105: {
106: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
107: int ierr;
108: HYPRE_ParCSRMatrix hmat;
109: PetscScalar *bv,*xv;
110: HYPRE_ParVector jbv,jxv;
111: PetscScalar *sbv,*sxv;
114: VecGetArray(b,&bv);
115: VecGetArray(x,&xv);
116: HYPREReplacePointer(jac->b,bv,sbv);
117: HYPREReplacePointer(jac->x,xv,sxv);
119: HYPRE_IJMatrixGetObject(jac->ij,(void**)&hmat);
120: HYPRE_IJVectorGetObject(jac->b,(void**)&jbv);
121: HYPRE_IJVectorGetObject(jac->x,(void**)&jxv);
122: (*jac->solve)(jac->hsolver,hmat,jbv,jxv);
123: /* error code of 1 in boomerAMG merely means convergence not achieved */
124: if (ierr && (ierr != 1 || jac->solve != HYPRE_BoomerAMGSolve)) SETERRQ1(1,"Error in HYPRE solver, error code %d",ierr);
125:
127: HYPREReplacePointer(jac->b,sbv,bv);
128: HYPREReplacePointer(jac->x,sxv,xv);
129: VecRestoreArray(x,&xv);
130: VecRestoreArray(b,&bv);
132: return(0);
133: }
135: #undef __FUNCT__
137: static int PCDestroy_HYPRE(PC pc)
138: {
139: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
140: int ierr;
143: HYPRE_IJMatrixDestroy(jac->ij);
144: HYPRE_IJVectorDestroy(jac->b);
145: HYPRE_IJVectorDestroy(jac->x);
146: (*jac->destroy)(jac->hsolver);
147: PetscFree(jac);
148: return(0);
149: }
151: /* --------------------------------------------------------------------------------------------*/
152: #undef __FUNCT__
154: static int PCSetFromOptions_HYPRE_Pilut(PC pc)
155: {
156: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
157: int ierr;
158: PetscTruth flag;
161: PetscOptionsHead("HYPRE Pilut Options");
162: PetscOptionsInt("-pc_hypre_pilut_maxiter","Number of iterations","None",jac->maxiter,&jac->maxiter,&flag);
163: if (flag) {
164: HYPRE_ParCSRPilutSetMaxIter(jac->hsolver,jac->maxiter);
165: }
166: PetscOptionsReal("-pc_hypre_pilut_tol","Drop tolerance","None",jac->tol,&jac->tol,&flag);
167: if (flag) {
168: HYPRE_ParCSRPilutSetDropTolerance(jac->hsolver,jac->tol);
169: }
170: PetscOptionsInt("-pc_hypre_pilut_factorrowsize","FactorRowSize","None",jac->factorrowsize,&jac->factorrowsize,&flag);
171: if (flag) {
172: HYPRE_ParCSRPilutSetFactorRowSize(jac->hsolver,jac->factorrowsize);
173: }
174: PetscOptionsTail();
175: return(0);
176: }
178: #undef __FUNCT__
180: static int PCView_HYPRE_Pilut(PC pc,PetscViewer viewer)
181: {
182: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
183: int ierr;
184: PetscTruth isascii;
187: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
188: if (isascii) {
189: PetscViewerASCIIPrintf(viewer," HYPRE Pilut preconditioningn");
190: if (jac->maxiter != PETSC_DEFAULT) {
191: PetscViewerASCIIPrintf(viewer," HYPRE Pilut: maximum number of iterations %dn",jac->maxiter);
192: } else {
193: PetscViewerASCIIPrintf(viewer," HYPRE Pilut: default maximum number of iterations n");
194: }
195: if (jac->tol != PETSC_DEFAULT) {
196: PetscViewerASCIIPrintf(viewer," HYPRE Pilut: drop tolerance %gn",jac->tol);
197: } else {
198: PetscViewerASCIIPrintf(viewer," HYPRE Pilut: default drop tolerance n");
199: }
200: if (jac->factorrowsize != PETSC_DEFAULT) {
201: PetscViewerASCIIPrintf(viewer," HYPRE Pilut: factor row size %dn",jac->factorrowsize);
202: } else {
203: PetscViewerASCIIPrintf(viewer," HYPRE Pilut: default factor row size n");
204: }
205: }
206: return(0);
207: }
209: /* --------------------------------------------------------------------------------------------*/
210: #undef __FUNCT__
212: static int PCSetFromOptions_HYPRE_Euclid(PC pc)
213: {
214: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
215: int ierr;
216: PetscTruth flag;
217: char *args[2];
220: jac->bjilu = PETSC_FALSE;
221: jac->levels = 1;
223: PetscOptionsHead("HYPRE Euclid Options");
224: PetscOptionsInt("-pc_hypre_euclid_levels","Number of levels of fill ILU(k)","None",jac->levels,&jac->levels,&flag);
225: if (flag) {
226: char levels[16];
227: if (jac->levels < 0) SETERRQ1(1,"Number of levels %d must be nonegative",jac->levels);
228: sprintf(levels,"%d",jac->levels);
229: args[0] = "-level"; args[1] = levels;
230: HYPRE_EuclidSetParams(jac->hsolver,2,args);
231: }
232: PetscOptionsLogical("-pc_hypre_euclid_bj","Use block Jacobi ILU(k)","None",jac->bjilu,&jac->bjilu,PETSC_NULL);
233: if (jac->bjilu) {
234: args[0] = "-bj"; args[1] = "1";
235: HYPRE_EuclidSetParams(jac->hsolver,2,args);
236: }
237:
238: PetscOptionsLogical("-pc_hypre_euclid_print_statistics","Print statistics","None",jac->printstatistics,&jac->printstatistics,PETSC_NULL);
239: if (jac->printstatistics) {
240: args[0] = "-eu_stats"; args[1] = "1";
241: HYPRE_EuclidSetParams(jac->hsolver,2,args);
242: args[0] = "-eu_mem"; args[1] = "1";
243: HYPRE_EuclidSetParams(jac->hsolver,2,args);
244: }
245: PetscOptionsTail();
246: return(0);
247: }
249: #undef __FUNCT__
251: static int PCView_HYPRE_Euclid(PC pc,PetscViewer viewer)
252: {
253: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
254: int ierr;
255: PetscTruth isascii;
258: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
259: if (isascii) {
260: PetscViewerASCIIPrintf(viewer," HYPRE Euclid preconditioningn");
261: PetscViewerASCIIPrintf(viewer," HYPRE Euclid: number of levels %dn",jac->levels);
262: if (jac->bjilu) {
263: PetscViewerASCIIPrintf(viewer," HYPRE Euclid: Using block Jacobi ILU instead of parallel ILUn");
264: }
265: }
266: return(0);
267: }
269: /* --------------------------------------------------------------------------------------------*/
271: static char *HYPREBoomerAMGCoarsenType[] = {"CLJP","Ruge-Stueben","","modifiedRuge-Stueben","","","Falgout"};
272: static char *HYPREBoomerAMGMeasureType[] = {"local","global"};
273: static char *HYPREBoomerAMGRelaxType[] = {"Jacobi","sequential-Gauss-Seidel","","Gauss-Seidel/Jacobi","","","symmetric-Gauss-Seidel/Jacobi",
274: "","","Gaussian-elimination"};
275: #undef __FUNCT__
277: static int PCSetFromOptions_HYPRE_BoomerAMG(PC pc)
278: {
279: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
280: int ierr,n = 4,i;
281: PetscTruth flg;
282: char result[32];
285: jac->maxlevels = 25;
286: jac->maxiter = 20;
287: jac->tol = 1.e-7;
288: jac->strongthreshold = .25;
289: jac->maxrowsum = .9;
290: jac->coarsentype = 6;
291: jac->measuretype = 0;
292:
294: /* this is terrible; HYPRE frees this array so we have to malloc it */
295: jac->gridsweeps = (int*)malloc(4*sizeof(int));
296: jac->gridsweeps[0] = jac->gridsweeps[1] = jac->gridsweeps[2] = 2;
297: jac->gridsweeps[3] = 1;
299: jac->relaxtype = (int*)malloc(4*sizeof(int));
300: jac->relaxtype[0] = jac->relaxtype[1] = jac->relaxtype[2] = 3;
301: jac->relaxtype[3] = 9;
303: PetscOptionsHead("HYPRE BoomerAMG Options");
304: PetscOptionsInt("-pc_hypre_boomeramg_max_levels","Number of levels (of grids) allowed","None",jac->maxlevels,&jac->maxlevels,&flg);
305: if (flg) {
306: if (jac->maxlevels < 2) SETERRQ1(1,"Number of levels %d must be at least two",jac->maxlevels);
307: HYPRE_BoomerAMGSetMaxLevels(jac->hsolver,jac->maxlevels);
308: }
309: PetscOptionsInt("-pc_hypre_boomeramg_max_iter","Maximum iterations used","None",jac->maxiter,&jac->maxiter,&flg);
310: if (flg) {
311: if (jac->maxiter < 1) SETERRQ1(1,"Number of iterations %d must be at least one",jac->maxiter);
312: HYPRE_BoomerAMGSetMaxIter(jac->hsolver,jac->maxiter);
313: }
314: PetscOptionsScalar("-pc_hypre_boomeramg_tol","Convergence tolerance","None",jac->tol,&jac->tol,&flg);
315: if (flg) {
316: if (jac->tol < 0.0) SETERRQ1(1,"Tolerance %g must be great than or equal zero",jac->tol);
317: HYPRE_BoomerAMGSetTol(jac->hsolver,jac->tol);
318: }
319: PetscOptionsScalar("-pc_hypre_boomeramg_strong_threshold","Threshold for being strongly connected","None",jac->strongthreshold,&jac->strongthreshold,&flg);
320: if (flg) {
321: if (jac->strongthreshold < 0.0) SETERRQ1(1,"Strong threshold %g must be great than or equal zero",jac->strongthreshold);
322: HYPRE_BoomerAMGSetStrongThreshold(jac->hsolver,jac->strongthreshold);
323: }
324: PetscOptionsScalar("-pc_hypre_boomeramg_max_row_sum","Maximum row sum","None",jac->maxrowsum,&jac->maxrowsum,&flg);
325: if (flg) {
326: if (jac->maxrowsum < 0.0) SETERRQ1(1,"Maximum row sum %g must be greater than zero",jac->maxrowsum);
327: if (jac->maxrowsum > 1.0) SETERRQ1(1,"Maximum row sum %g must be less than or equal one",jac->maxrowsum);
328: HYPRE_BoomerAMGSetMaxRowSum(jac->hsolver,jac->maxrowsum);
329: }
330:
331: n = 4;
332: PetscOptionsIntArray("-pc_hypre_boomeramg_grid_sweeps","Grid sweeps for fine,down,up,coarse","None",jac->gridsweeps,&n,&flg);
333: if (flg) {
334: if (n == 1) {
335: jac->gridsweeps[1] = jac->gridsweeps[2] = jac->gridsweeps[3] = jac->gridsweeps[0];
336: n = 4;
337: }
338: if (n != 4) SETERRQ1(1,"You must provide either 1 or 4 values seperated by commas, you provided %d",n);
339: HYPRE_BoomerAMGSetNumGridSweeps(jac->hsolver,jac->gridsweeps);
340: CHKMEMQ;
341: }
343: /*
344: Suggested by QUANDALLE Philippe <Philippe.QUANDALLE@ifp.fr>
346: gridrelaxpoints[i][j] are for i=0,1,2,3 (fine,down,up,coarse) and j=sweep number
347: 0 indicates smooth all points
348: 1 indicates smooth coarse points
349: -1 indicates smooth fine points
351: Here when j=1 it first smooths all the coarse points, then all the fine points.
352: */
353: jac->gridrelaxpoints = (int**)malloc(4*sizeof(int*));
354: if(jac->gridsweeps[0]>0) jac->gridrelaxpoints[0] = (int*)malloc(jac->gridsweeps[0]*sizeof(int));
355: if(jac->gridsweeps[1]>0) jac->gridrelaxpoints[1] = (int*)malloc(jac->gridsweeps[1]*sizeof(int));
356: if(jac->gridsweeps[2]>0) jac->gridrelaxpoints[2] = (int*)malloc(jac->gridsweeps[2]*sizeof(int));
357: if(jac->gridsweeps[3]>0) jac->gridrelaxpoints[3] = (int*)malloc(jac->gridsweeps[3]*sizeof(int));
359: PetscOptionsLogical("-pc_hypre_boomeramg_sweep_all","Sweep all points","None",PETSC_FALSE,&flg,0);
360: if(jac->gridsweeps[0] == 1) jac->gridrelaxpoints[0][0] = 0;
361: else if(jac->gridsweeps[0] == 2) {
362: if (flg) {
363: jac->gridrelaxpoints[0][0] = 0; jac->gridrelaxpoints[0][1] = 0;
364: } else {
365: jac->gridrelaxpoints[0][0] = 1; jac->gridrelaxpoints[0][1] = -1;
366: }
367: } else if (jac->gridsweeps[0] > 2) {
368: SETERRQ(1,"Grid sweeps can only be 0, 1, or 2");
369: }
370:
371: if(jac->gridsweeps[1] == 1) jac->gridrelaxpoints[1][0] = 0;
372: else if(jac->gridsweeps[1] == 2) {
373: if (flg) {
374: jac->gridrelaxpoints[1][0] = 0; jac->gridrelaxpoints[1][1] = 0;
375: } else {
376: jac->gridrelaxpoints[1][0] = 1; jac->gridrelaxpoints[1][1] = -1;
377: }
378: } else if (jac->gridsweeps[1] > 2) {
379: SETERRQ(1,"Grid sweeps can only be 0, 1, or 2");
380: }
382: if(jac->gridsweeps[2] == 1) jac->gridrelaxpoints[2][0] = 0;
383: else if(jac->gridsweeps[2] == 2) {
384: if (flg) {
385: jac->gridrelaxpoints[2][0] = 0; jac->gridrelaxpoints[2][1] = 0;
386: } else {
387: jac->gridrelaxpoints[2][0] = -1; jac->gridrelaxpoints[2][1] = 1;
388: }
389: } else if (jac->gridsweeps[2] > 2) {
390: SETERRQ(1,"Grid sweeps can only be 0, 1, or 2");
391: }
393: for (i=0; i<jac->gridsweeps[3]; i++) {
394: jac->gridrelaxpoints[3][i] = 0;
395: }
396: HYPRE_BoomerAMGSetGridRelaxPoints(jac->hsolver,jac->gridrelaxpoints);
399: PetscOptionsEList("-pc_hypre_boomeramg_measure_type","Measure type","None",HYPREBoomerAMGMeasureType,2,HYPREBoomerAMGMeasureType[0],result,16,&flg);
400: if (flg) {
401: int type = -1;
402: for (i=0; i<2; i++) {
403: PetscStrcmp(result,HYPREBoomerAMGMeasureType[i],&flg);
404: if (flg) {
405: type = i;
406: break;
407: }
408: }
409: if (type == -1) SETERRQ1(1,"Unknown measure type %s",result);
410: jac->measuretype = type;
411: HYPRE_BoomerAMGSetMeasureType(jac->hsolver,jac->measuretype);
412: }
413: PetscOptionsEList("-pc_hypre_boomeramg_coarsen_type","Coarsen type","None",HYPREBoomerAMGCoarsenType,7,HYPREBoomerAMGCoarsenType[6],result,25,&flg);
414: if (flg) {
415: int type = -1;
416: for (i=0; i<7; i++) {
417: PetscStrcmp(result,HYPREBoomerAMGCoarsenType[i],&flg);
418: if (flg) {
419: type = i;
420: break;
421: }
422: }
423: if (type == -1) SETERRQ1(1,"Unknown coarsen type %s",result);
424: jac->coarsentype = type;
425: HYPRE_BoomerAMGSetCoarsenType(jac->hsolver,jac->coarsentype);
426: }
427: PetscOptionsEList("-pc_hypre_boomeramg_relax_type","Relax type","None",HYPREBoomerAMGRelaxType,10,HYPREBoomerAMGRelaxType[3],result,32,&flg);
428: if (flg) {
429: int type = -1;
430: for (i=0; i<10; i++) {
431: PetscStrcmp(result,HYPREBoomerAMGRelaxType[i],&flg);
432: if (flg) {
433: type = i;
434: break;
435: }
436: }
437: if (type == -1) SETERRQ1(1,"Unknown relax type %s",result);
438: jac->relaxtype[0] = jac->relaxtype[1] = jac->relaxtype[2] = type;
439: }
440: PetscOptionsEList("-pc_hypre_boomeramg_relax_type_coarse","Relax type on coarse grid","None",HYPREBoomerAMGRelaxType,10,HYPREBoomerAMGRelaxType[3],result,32,&flg);
441: if (flg) {
442: int type = -1;
443: for (i=0; i<10; i++) {
444: PetscStrcmp(result,HYPREBoomerAMGRelaxType[i],&flg);
445: if (flg) {
446: type = i;
447: break;
448: }
449: }
450: if (type == -1) SETERRQ1(1,"Unknown relax type %s",result);
451: jac->relaxtype[3] = type;
452: }
453: HYPRE_BoomerAMGSetGridRelaxType(jac->hsolver,jac->relaxtype);
454:
455: PetscOptionsLogical("-pc_hypre_boomeramg_print_statistics","Print statistics","None",jac->printstatistics,&jac->printstatistics,PETSC_NULL);
456: if (jac->printstatistics) {
457: HYPRE_BoomerAMGSetPrintLevel(jac->hsolver,3);
458: }
459: PetscOptionsTail();
460: return(0);
461: }
463: #undef __FUNCT__
465: static int PCView_HYPRE_BoomerAMG(PC pc,PetscViewer viewer)
466: {
467: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
468: int ierr;
469: PetscTruth isascii;
472: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
473: if (isascii) {
474: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG preconditioningn");
475: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Maximum number of levels %dn",jac->maxlevels);
476: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Maximum number of iterations %dn",jac->maxiter);
477: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Convergence tolerance %gn",jac->tol);
478: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Threshold for strong coupling %gn",jac->strongthreshold);
479: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Maximum row sums %gn",jac->maxrowsum);
480: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Sweeps on fine grid %dn",jac->gridsweeps[0]);
481: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Sweeps down %dn",jac->gridsweeps[1]);
482: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Sweeps up %dn",jac->gridsweeps[2]);
483: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Sweeps on coarse %dn",jac->gridsweeps[3]);
485: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Relax on fine grid %sn",HYPREBoomerAMGRelaxType[jac->relaxtype[0]]);
486: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Relax down %sn",HYPREBoomerAMGRelaxType[jac->relaxtype[1]]);
487: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Relax up %sn",HYPREBoomerAMGRelaxType[jac->relaxtype[2]]);
488: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Relax on coarse %sn",HYPREBoomerAMGRelaxType[jac->relaxtype[3]]);
490: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Measure type %sn",HYPREBoomerAMGMeasureType[jac->measuretype]);
491: PetscViewerASCIIPrintf(viewer," HYPRE BoomerAMG: Coarsen type %sn",HYPREBoomerAMGCoarsenType[jac->coarsentype]);
492: }
493: return(0);
494: }
496: /* --------------------------------------------------------------------------------------------*/
497: #undef __FUNCT__
499: static int PCSetFromOptions_HYPRE_ParaSails(PC pc)
500: {
501: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
502: int ierr;
503: PetscTruth flag;
504: char *symtlist[] = {"nonsymmetric","SPD","nonsymmetric,SPD"},buff[32];
507: jac->nlevels = 1;
508: jac->threshhold = .1;
509: jac->filter = .1;
510: jac->loadbal = 0;
511: if (PetscLogPrintInfo) {
512: jac->logging = (int) PETSC_TRUE;
513: } else {
514: jac->logging = (int) PETSC_FALSE;
515: }
516: jac->ruse = (int) PETSC_FALSE;
517: jac->symt = 0;
519: PetscOptionsHead("HYPRE ParaSails Options");
520: PetscOptionsInt("-pc_hypre_parasails_nlevels","Number of number of levels","None",jac->nlevels,&jac->nlevels,0);
521: PetscOptionsReal("-pc_hypre_parasails_thresh","Threshold","None",jac->threshhold,&jac->threshhold,0);
522: HYPRE_ParaSailsSetParams(jac->hsolver,jac->threshhold,jac->nlevels);
524: PetscOptionsReal("-pc_hypre_parasails_filter","filter","None",jac->filter,&jac->filter,0);
525: HYPRE_ParaSailsSetFilter(jac->hsolver,jac->filter);
527: PetscOptionsReal("-pc_hypre_parasails_loadbal","Load balance","None",jac->loadbal,&jac->loadbal,0);
528: HYPRE_ParaSailsSetLoadbal(jac->hsolver,jac->loadbal);
530: PetscOptionsLogical("-pc_hypre_parasails_logging","Print info to screen","None",(PetscTruth)jac->logging,(PetscTruth*)&jac->logging,0);
531: HYPRE_ParaSailsSetLogging(jac->hsolver,jac->logging);
533: PetscOptionsLogical("-pc_hypre_parasails_reuse","Reuse nonzero pattern in preconditioner","None",(PetscTruth)jac->ruse,(PetscTruth*)&jac->ruse,0);
534: HYPRE_ParaSailsSetReuse(jac->hsolver,jac->ruse);
536: PetscOptionsEList("-pc_hypre_parasails_sym","Symmetry of matrix and preconditioner","None",symtlist,3,symtlist[0],buff,32,&flag);
537: while (flag) {
538: PetscStrcmp(symtlist[0],buff,&flag);
539: if (flag) {
540: jac->symt = 0;
541: break;
542: }
543: PetscStrcmp(symtlist[1],buff,&flag);
544: if (flag) {
545: jac->symt = 1;
546: break;
547: }
548: PetscStrcmp(symtlist[2],buff,&flag);
549: if (flag) {
550: jac->symt = 2;
551: break;
552: }
553: SETERRQ1(1,"Unknown HYPRE ParaSails Sym option %s",buff);
554: }
555: HYPRE_ParaSailsSetSym(jac->hsolver,jac->symt);
557: PetscOptionsTail();
558: return(0);
559: }
561: #undef __FUNCT__
563: static int PCView_HYPRE_ParaSails(PC pc,PetscViewer viewer)
564: {
565: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
566: int ierr;
567: PetscTruth isascii;
568: char *symt;
571: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
572: if (isascii) {
573: PetscViewerASCIIPrintf(viewer," HYPRE ParaSails preconditioningn");
574: PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: nlevels %dn",jac->nlevels);
575: PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: threshold %gn",jac->threshhold);
576: PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: filter %gn",jac->filter);
577: PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: load balance %gn",jac->loadbal);
578: PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: reuse nonzero structure %sn",jac->ruse ? "true" : "false");
579: PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: print info to screen %sn",jac->logging ? "true" : "false");
580: if (jac->symt == 0) {
581: symt = "nonsymmetric matrix and preconditioner";
582: } else if (jac->symt == 1) {
583: symt = "SPD matrix and preconditioner";
584: } else if (jac->symt == 2) {
585: symt = "nonsymmetric matrix but SPD preconditioner";
586: } else {
587: SETERRQ1(1,"Unknown HYPRE ParaSails sym option %d",jac->symt);
588: }
589: PetscViewerASCIIPrintf(viewer," HYPRE ParaSails: %sn",symt);
590: }
591: return(0);
592: }
593: /* ---------------------------------------------------------------------------------*/
595: #undef __FUNCT__
597: static int PCHYPRESetType_HYPRE(PC pc,char *name)
598: {
599: PC_HYPRE *jac = (PC_HYPRE*)pc->data;
600: int ierr;
601: PetscTruth flag;
604: if (pc->ops->setup) {
605: SETERRQ(1,"Cannot set the HYPRE preconditioner type once it has been set");
606: }
608: pc->ops->setup = PCSetUp_HYPRE;
609: pc->ops->apply = PCApply_HYPRE;
610: pc->ops->destroy = PCDestroy_HYPRE;
612: jac->maxiter = PETSC_DEFAULT;
613: jac->tol = PETSC_DEFAULT;
614: jac->printstatistics = PetscLogPrintInfo;
616: PetscStrcmp("pilut",name,&flag);
617: if (flag) {
618: ierr = HYPRE_ParCSRPilutCreate(pc->comm,&jac->hsolver);
619: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Pilut;
620: pc->ops->view = PCView_HYPRE_Pilut;
621: jac->destroy = HYPRE_ParCSRPilutDestroy;
622: jac->setup = HYPRE_ParCSRPilutSetup;
623: jac->solve = HYPRE_ParCSRPilutSolve;
624: jac->factorrowsize = PETSC_DEFAULT;
625: return(0);
626: }
627: PetscStrcmp("parasails",name,&flag);
628: if (flag) {
629: ierr = HYPRE_ParaSailsCreate(pc->comm,&jac->hsolver);
630: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_ParaSails;
631: pc->ops->view = PCView_HYPRE_ParaSails;
632: jac->destroy = HYPRE_ParaSailsDestroy;
633: jac->setup = HYPRE_ParaSailsSetup;
634: jac->solve = HYPRE_ParaSailsSolve;
635: return(0);
636: }
637: PetscStrcmp("euclid",name,&flag);
638: if (flag) {
639: ierr = HYPRE_EuclidCreate(pc->comm,&jac->hsolver);
640: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_Euclid;
641: pc->ops->view = PCView_HYPRE_Euclid;
642: jac->destroy = HYPRE_EuclidDestroy;
643: jac->setup = HYPRE_EuclidSetup;
644: jac->solve = HYPRE_EuclidSolve;
645: return(0);
646: }
647: PetscStrcmp("boomeramg",name,&flag);
648: if (flag) {
649: ierr = HYPRE_BoomerAMGCreate(&jac->hsolver);
650: pc->ops->setfromoptions = PCSetFromOptions_HYPRE_BoomerAMG;
651: pc->ops->view = PCView_HYPRE_BoomerAMG;
652: jac->destroy = HYPRE_BoomerAMGDestroy;
653: jac->setup = HYPRE_BoomerAMGSetup;
654: jac->solve = HYPRE_BoomerAMGSolve;
655: return(0);
656: }
657: SETERRQ1(1,"Unknown HYPRE preconditioner %s; Choices are pilut, parasails, euclid, boomeramg",name);
658: return(0);
659: }
661: /*
662: It only gets here if the HYPRE type has not been set before the call to
663: ...SetFromOptions() which actually is most of the time
664: */
665: #undef __FUNCT__
667: static int PCSetFromOptions_HYPRE(PC pc)
668: {
669: int ierr;
670: char buff[32],*type[] = {"pilut","parasails","boomerAMG","euclid"};
671: PetscTruth flg;
674: PetscOptionsHead("HYPRE preconditioner options");
675: PetscOptionsEList("-pc_hypre_type","HYPRE preconditioner type","PCHYPRESetType",type,4,"pilut",buff,32,&flg);
677:
678: if (PetscOptionsPublishCount) { /* force the default if it was not yet set and user did not set with option */
679: if (!flg && !pc->ops->apply) {
680: flg = PETSC_TRUE;
681: PetscStrcpy(buff,"pilut");
682: }
683: }
685: if (flg) {
686: PCHYPRESetType_HYPRE(pc,buff);
687: }
688: if (pc->ops->setfromoptions) {
689: pc->ops->setfromoptions(pc);
690: }
691: PetscOptionsTail();
692: return(0);
693: }
695: #undef __FUNCT__
697: /*@
698: PCHYPRESetType - Sets which hypre preconditioner you wish to use
700: Input Parameters:
701: + pc - the preconditioner context
702: - name - either pilut, parasails, boomerAMG, euclid
704: Options Database Keys:
705: -pc_hypre_type - One of pilut, parasails, boomerAMG, euclid
706:
707: Level: intermediate
709: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
710: PCHYPRE
712: @*/
713: int PCHYPRESetType(PC pc,char *name)
714: {
715: int ierr,(*f)(PC,char*);
719: PetscObjectQueryFunction((PetscObject)pc,"PCHYPRESetType_C",(void (**)(void))&f);
720: if (f) {
721: (*f)(pc,name);
722: }
723: return(0);
724: }
726: /*MC
727: PCHYPRE - Allows you to use the matrix element based preconditioners in the LLNL package hypre
729: Options Database Keys:
730: + -pc_hypre_type - One of pilut, parasails, boomerAMG, euclid
731: - Too many others to list, run with -pc_type hypre -pc_hypre_type XXX to see options for the XXX
732: preconditioner
733:
734: Level: intermediate
736: Notes: The many hypre options can ONLY be set via the options database (e.g. the command line
737: or with PetscOptionsSetValue(), there are no functions to set them)
739: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC,
740: PCHYPRESetType()
742: M*/
744: EXTERN_C_BEGIN
745: #undef __FUNCT__
747: int PCCreate_HYPRE(PC pc)
748: {
749: PC_HYPRE *jac;
750: int ierr;
753: ierr = PetscNew(PC_HYPRE,&jac);
754: ierr = PetscMemzero(jac,sizeof(PC_HYPRE));
755: pc->data = jac;
756: pc->ops->setfromoptions = PCSetFromOptions_HYPRE;
757: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCHYPRESetType_C",
758: "PCHYPRESetType_HYPRE",PCHYPRESetType_HYPRE);
759: return(0);
760: }
761: EXTERN_C_END