Actual source code: parms.c
petsc-3.3-p7 2013-05-11
1: #define PETSCKSP_DLL
3: /*
4: Provides an interface to pARMS.
5: Requires pARMS 3.2 or later.
6: */
8: #include "petsc-private/pcimpl.h" /*I "petscpc.h" I*/
10: #ifdef PETSC_USE_COMPLEX
11: #define DBL_CMPLX
12: #else
13: #define DBL
14: #endif
15: #define USE_MPI
16: #define REAL double
17: #define HAS_BLAS
18: #define FORTRAN_UNDERSCORE
19: #include "parms_sys.h"
20: #undef FLOAT
21: #define FLOAT PetscScalar
22: #include "parms.h"
24: /*
25: Private context (data structure) for the preconditioner.
26: */
27: typedef struct {
28: parms_Map map;
29: parms_Mat A;
30: parms_PC pc;
31: PCPARMSGlobalType global;
32: PCPARMSLocalType local;
33: PetscInt levels, blocksize, maxdim, maxits, lfil[7];
34: PetscBool nonsymperm, meth[8];
35: PetscReal solvetol, indtol, droptol[7];
36: PetscScalar *lvec0, *lvec1;
37: } PC_PARMS;
42: static PetscErrorCode PCSetUp_PARMS(PC pc)
43: {
44: Mat pmat;
45: PC_PARMS *parms = (PC_PARMS*)pc->data;
46: const PetscInt *mapptr0;
47: PetscInt n, lsize, low, high, i, pos, ncols, length;
48: int *maptmp, *mapptr, *ia, *ja, *ja1, *im;
49: PetscScalar *aa, *aa1;
50: const PetscInt *cols;
51: PetscInt meth[8];
52: const PetscScalar *values;
53: PetscErrorCode ierr;
54: MatInfo matinfo;
55: PetscMPIInt rank, npro;
59: /* Get preconditioner matrix from PETSc and setup pARMS structs */
60: PCGetOperators(pc,PETSC_NULL,&pmat,PETSC_NULL);
61: MPI_Comm_size(((PetscObject)pmat)->comm,&npro);
62: MPI_Comm_rank(((PetscObject)pmat)->comm,&rank);
64: MatGetSize(pmat,&n,PETSC_NULL);
65: PetscMalloc((npro+1)*sizeof(int),&mapptr);
66: PetscMalloc(n*sizeof(int),&maptmp);
67: MatGetOwnershipRanges(pmat,&mapptr0);
68: low = mapptr0[rank];
69: high = mapptr0[rank+1];
70: lsize = high - low;
72: for (i=0; i<npro+1; i++)
73: mapptr[i] = mapptr0[i]+1;
74: for (i = 0; i<n; i++)
75: maptmp[i] = i+1;
77: /* if created, destroy the previous map */
78: if(parms->map) {
79: parms_MapFree(&parms->map);
80: parms->map = PETSC_NULL;
81: }
82:
83: /* create pARMS map object */
84: parms_MapCreateFromPtr(&parms->map,(int)n,maptmp,mapptr,((PetscObject)pmat)->comm,1,NONINTERLACED);
86: /* if created, destroy the previous pARMS matrix */
87: if(parms->A) {
88: parms_MatFree(&parms->A);
89: parms->A = PETSC_NULL;
90: }
92: /* create pARMS mat object */
93: parms_MatCreate(&parms->A,parms->map);
95: /* setup and copy csr data structure for pARMS */
96: PetscMalloc((lsize+1)*sizeof(int),&ia);
97: ia[0] = 1;
98: MatGetInfo(pmat,MAT_LOCAL,&matinfo);
99: length = matinfo.nz_used;
100: PetscMalloc(length*sizeof(int),&ja);
101: PetscMalloc(length*sizeof(PetscScalar),&aa);
103: for (i = low; i<high; i++) {
104: pos = ia[i-low]-1;
105: MatGetRow(pmat,i,&ncols,&cols,&values);
106: ia[i-low+1] = ia[i-low] + ncols;
108: if (ia[i-low+1] >= length) {
109: length += ncols;
110: PetscMalloc(length*sizeof(int),&ja1);
111: PetscMemcpy(ja1,ja,(ia[i-low]-1)*sizeof(int));
112: PetscFree(ja);
113: ja = ja1;
114: PetscMalloc(length*sizeof(PetscScalar),&aa1);
115: PetscMemcpy(aa1,aa,(ia[i-low]-1)*sizeof(PetscScalar));
116: PetscFree(aa);
117: aa = aa1;
118: }
119: PetscMemcpy(&ja[pos],cols,ncols*sizeof(int));
120: PetscMemcpy(&aa[pos],values,ncols*sizeof(PetscScalar));
121: MatRestoreRow(pmat,i,&ncols,&cols,&values);
122: }
124: /* csr info is for local matrix so initialize im[] locally */
125: PetscMalloc(lsize*sizeof(int),&im);
126: PetscMemcpy(im,&maptmp[mapptr[rank]-1],lsize*sizeof(int));
128: /* 1-based indexing */
129: for(i=0; i<ia[lsize]-1; i++)
130: ja[i] = ja[i]+1;
132: /* Now copy csr matrix to parms_mat object */
133: parms_MatSetValues(parms->A,(int)lsize,im,ia,ja,aa,INSERT);
135: /* free memory */
136: PetscFree(maptmp);
137: PetscFree(mapptr);
138: PetscFree(aa);
139: PetscFree(ja);
140: PetscFree(ia);
141: PetscFree(im);
143: /* setup parms matrix */
144: parms_MatSetup(parms->A);
146: /* if created, destroy the previous pARMS pc */
147: if(parms->pc) {
148: parms_PCFree(&parms->pc);
149: parms->pc = PETSC_NULL;
150: }
151:
152: /* Now create pARMS preconditioner object based on A */
153: parms_PCCreate(&parms->pc,parms->A);
154:
155: /* Transfer options from PC to pARMS */
156: switch(parms->global) {
157: case 0: parms_PCSetType(parms->pc, PCRAS); break;
158: case 1: parms_PCSetType(parms->pc, PCSCHUR); break;
159: case 2: parms_PCSetType(parms->pc, PCBJ); break;
160: }
161: switch(parms->local) {
162: case 0: parms_PCSetILUType(parms->pc, PCILU0); break;
163: case 1: parms_PCSetILUType(parms->pc, PCILUK); break;
164: case 2: parms_PCSetILUType(parms->pc, PCILUT); break;
165: case 3: parms_PCSetILUType(parms->pc, PCARMS); break;
166: }
167: parms_PCSetInnerEps(parms->pc, parms->solvetol);
168: parms_PCSetNlevels(parms->pc, parms->levels);
169: parms_PCSetPermType(parms->pc, parms->nonsymperm?1:0);
170: parms_PCSetBsize(parms->pc, parms->blocksize);
171: parms_PCSetTolInd(parms->pc, parms->indtol);
172: parms_PCSetInnerKSize(parms->pc, parms->maxdim);
173: parms_PCSetInnerMaxits(parms->pc, parms->maxits);
174: for(i=0; i<8; i++) meth[i] = parms->meth[i]?1:0;
175: parms_PCSetPermScalOptions(parms->pc, &meth[0], 1);
176: parms_PCSetPermScalOptions(parms->pc, &meth[4], 0);
177: parms_PCSetFill(parms->pc, parms->lfil);
178: parms_PCSetTol(parms->pc, parms->droptol);
180: parms_PCSetup(parms->pc);
182: /* Allocate two auxiliary vector of length lsize */
183: if (parms->lvec0) { PetscFree(parms->lvec0); }
184: PetscMalloc(lsize*sizeof(PetscScalar), &parms->lvec0);
185: if (parms->lvec1) { PetscFree(parms->lvec1); }
186: PetscMalloc(lsize*sizeof(PetscScalar), &parms->lvec1);
187: return(0);
188: }
192: static PetscErrorCode PCView_PARMS(PC pc,PetscViewer viewer)
193: {
194: PetscErrorCode ierr;
195: PetscBool iascii;
196: PC_PARMS *parms = (PC_PARMS*)pc->data;
197: char *str;
198: double fill_fact;
199:
201: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
202: if (iascii) {
203: parms_PCGetName(parms->pc,&str);
204: PetscViewerASCIIPrintf(viewer," global preconditioner: %s\n",str);
205: parms_PCILUGetName(parms->pc,&str);
206: PetscViewerASCIIPrintf(viewer," local preconditioner: %s\n",str);
207: parms_PCGetRatio(parms->pc,&fill_fact);
208: PetscViewerASCIIPrintf(viewer," non-zero elements/original non-zero entries: %-4.2f\n",fill_fact);
209: PetscViewerASCIIPrintf(viewer," Tolerance for local solve: %g\n",parms->solvetol);
210: PetscViewerASCIIPrintf(viewer," Number of levels: %d\n",parms->levels);
211: if (parms->nonsymperm) {
212: PetscViewerASCIIPrintf(viewer," Using nonsymmetric permutation\n");
213: }
214: PetscViewerASCIIPrintf(viewer," Block size: %d\n",parms->blocksize);
215: PetscViewerASCIIPrintf(viewer," Tolerance for independent sets: %g\n",parms->indtol);
216: PetscViewerASCIIPrintf(viewer," Inner Krylov dimension: %d\n",parms->maxdim);
217: PetscViewerASCIIPrintf(viewer," Maximum number of inner iterations: %d\n",parms->maxits);
218: if (parms->meth[0]) {
219: PetscViewerASCIIPrintf(viewer," Using nonsymmetric permutation for interlevel blocks\n");
220: }
221: if (parms->meth[1]) {
222: PetscViewerASCIIPrintf(viewer," Using column permutation for interlevel blocks\n");
223: }
224: if (parms->meth[2]) {
225: PetscViewerASCIIPrintf(viewer," Using row scaling for interlevel blocks\n");
226: }
227: if (parms->meth[3]) {
228: PetscViewerASCIIPrintf(viewer," Using column scaling for interlevel blocks\n");
229: }
230: if (parms->meth[4]) {
231: PetscViewerASCIIPrintf(viewer," Using nonsymmetric permutation for last level blocks\n");
232: }
233: if (parms->meth[5]) {
234: PetscViewerASCIIPrintf(viewer," Using column permutation for last level blocks\n");
235: }
236: if (parms->meth[6]) {
237: PetscViewerASCIIPrintf(viewer," Using row scaling for last level blocks\n");
238: }
239: if (parms->meth[7]) {
240: PetscViewerASCIIPrintf(viewer," Using column scaling for last level blocks\n");
241: }
242: PetscViewerASCIIPrintf(viewer," amount of fill-in for ilut, iluk and arms: %d\n",parms->lfil[0]);
243: PetscViewerASCIIPrintf(viewer," amount of fill-in for schur: %d\n",parms->lfil[4]);
244: PetscViewerASCIIPrintf(viewer," amount of fill-in for ILUT L and U: %d\n",parms->lfil[5]);
245: PetscViewerASCIIPrintf(viewer," drop tolerance for L, U, L^{-1}F and EU^{-1}: %g\n",parms->droptol[0]);
246: PetscViewerASCIIPrintf(viewer," drop tolerance for schur complement at each level: %g\n",parms->droptol[4]);
247: PetscViewerASCIIPrintf(viewer," drop tolerance for ILUT in last level schur complement: %g\n",parms->droptol[5]);
248: }
249:
250: return(0);
251: }
255: static PetscErrorCode PCDestroy_PARMS(PC pc)
256: {
257: PC_PARMS *parms = (PC_PARMS*)pc->data;
261: if(parms->map) {
262: parms_MapFree(&parms->map);
263: }
264: if(parms->A) {
265: parms_MatFree(&parms->A);
266: }
267: if(parms->pc) {
268: parms_PCFree(&parms->pc);
269: }
270: if(parms->lvec0) {
271: PetscFree(parms->lvec0);
272: }
273: if(parms->lvec1) {
274: PetscFree(parms->lvec1);
275: }
276: PetscFree(pc->data);
278: PetscObjectChangeTypeName((PetscObject)pc,0);
279: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetGlobal_C","",PETSC_NULL);
280: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetLocal_C","",PETSC_NULL);
281: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetSolveTolerances_C","",PETSC_NULL);
282: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetSolveRestart_C","",PETSC_NULL);
283: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetNonsymPerm_C","",PETSC_NULL);
284: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetFill_C","",PETSC_NULL);
285: return(0);
286: }
290: static PetscErrorCode PCSetFromOptions_PARMS(PC pc)
291: {
292: PC_PARMS *parms = (PC_PARMS*)pc->data;
293: PetscBool flag;
294: PCPARMSGlobalType global;
295: PCPARMSLocalType local;
296: PetscErrorCode ierr;
299: PetscOptionsHead("PARMS Options");
300: PetscOptionsEnum("-pc_parms_global","Global preconditioner","PCPARMSSetGlobal",PCPARMSGlobalTypes,(PetscEnum)parms->global,(PetscEnum*)&global,&flag);
301: if (flag) { PCPARMSSetGlobal(pc,global); }
302: PetscOptionsEnum("-pc_parms_local","Local preconditioner","PCPARMSSetLocal",PCPARMSLocalTypes,(PetscEnum)parms->local,(PetscEnum*)&local,&flag);
303: if (flag) { PCPARMSSetLocal(pc,local); }
304: PetscOptionsReal("-pc_parms_solve_tol","Tolerance for local solve","PCPARMSSetSolveTolerances",parms->solvetol,&parms->solvetol,&flag);
305: PetscOptionsInt("-pc_parms_levels","Number of levels","None",parms->levels,&parms->levels,&flag);
306: PetscOptionsBool("-pc_parms_nonsymmetric_perm","Use nonsymmetric permutation","PCPARMSSetNonsymPerm",parms->nonsymperm,&parms->nonsymperm,&flag);
307: PetscOptionsInt("-pc_parms_blocksize","Block size","None",parms->blocksize,&parms->blocksize,&flag);
308: PetscOptionsReal("-pc_parms_ind_tol","Tolerance for independent sets","None",parms->indtol,&parms->indtol,&flag);
309: PetscOptionsInt("-pc_parms_max_dim","Inner Krylov dimension","PCPARMSSetSolveRestart",parms->maxdim,&parms->maxdim,&flag);
310: PetscOptionsInt("-pc_parms_max_it","Maximum number of inner iterations","PCPARMSSetSolveTolerances",parms->maxits,&parms->maxits,&flag);
311: PetscOptionsBool("-pc_parms_inter_nonsymmetric_perm","nonsymmetric permutation for interlevel blocks","None",parms->meth[0],&parms->meth[0],&flag);
312: PetscOptionsBool("-pc_parms_inter_column_perm","column permutation for interlevel blocks","None",parms->meth[1],&parms->meth[1],&flag);
313: PetscOptionsBool("-pc_parms_inter_row_scaling","row scaling for interlevel blocks","None",parms->meth[2],&parms->meth[2],&flag);
314: PetscOptionsBool("-pc_parms_inter_column_scaling","column scaling for interlevel blocks","None",parms->meth[3],&parms->meth[3],&flag);
315: PetscOptionsBool("-pc_parms_last_nonsymmetric_perm","nonsymmetric permutation for last level blocks","None",parms->meth[4],&parms->meth[4],&flag);
316: PetscOptionsBool("-pc_parms_last_column_perm","column permutation for last level blocks","None",parms->meth[5],&parms->meth[5],&flag);
317: PetscOptionsBool("-pc_parms_last_row_scaling","row scaling for last level blocks","None",parms->meth[6],&parms->meth[6],&flag);
318: PetscOptionsBool("-pc_parms_last_column_scaling","column scaling for last level blocks","None",parms->meth[7],&parms->meth[7],&flag);
319: PetscOptionsInt("-pc_parms_lfil_ilu_arms","amount of fill-in for ilut, iluk and arms","PCPARMSSetFill",parms->lfil[0],&parms->lfil[0],&flag);
320: if(flag) parms->lfil[1] = parms->lfil[2] = parms->lfil[3] = parms->lfil[0];
321: PetscOptionsInt("-pc_parms_lfil_schur","amount of fill-in for schur","PCPARMSSetFill",parms->lfil[4],&parms->lfil[4],&flag);
322: PetscOptionsInt("-pc_parms_lfil_ilut_L_U","amount of fill-in for ILUT L and U","PCPARMSSetFill",parms->lfil[5],&parms->lfil[5],&flag);
323: if(flag) parms->lfil[6] = parms->lfil[5];
324: PetscOptionsReal("-pc_parms_droptol_factors","drop tolerance for L, U, L^{-1}F and EU^{-1}","None",parms->droptol[0],&parms->droptol[0],&flag);
325: PetscOptionsReal("-pc_parms_droptol_schur_compl","drop tolerance for schur complement at each level","None",parms->droptol[4],&parms->droptol[4],&flag);
326: if(flag) parms->droptol[1] = parms->droptol[2] = parms->droptol[3] = parms->droptol[0];
327: PetscOptionsReal("-pc_parms_droptol_last_schur","drop tolerance for ILUT in last level schur complement","None",parms->droptol[5],&parms->droptol[5],&flag);
328: if(flag) parms->droptol[6] = parms->droptol[5];
329: PetscOptionsTail();
330: return(0);
331: }
335: static PetscErrorCode PCApply_PARMS(PC pc,Vec b,Vec x)
336: {
337: PetscErrorCode ierr;
338: PC_PARMS *parms = (PC_PARMS*)pc->data;
339: const PetscScalar *b1;
340: PetscScalar *x1;
343: VecGetArrayRead(b,&b1);
344: VecGetArray(x,&x1);
345: parms_VecPermAux((PetscScalar*)b1,parms->lvec0,parms->map);
346: parms_PCApply(parms->pc,parms->lvec0,parms->lvec1);
347: parms_VecInvPermAux(parms->lvec1,x1,parms->map);
348: VecRestoreArrayRead(b,&b1);
349: VecRestoreArray(x,&x1);
350: return(0);
351: }
353: EXTERN_C_BEGIN
356: PetscErrorCode PCPARMSSetGlobal_PARMS(PC pc,PCPARMSGlobalType type)
357: {
358: PC_PARMS *parms = (PC_PARMS*)pc->data;
361: if (type != parms->global) {
362: parms->global = type;
363: pc->setupcalled = 0;
364: }
365: return(0);
366: }
367: EXTERN_C_END
371: /*@
372: PCPARMSSetGlobal - Sets the global preconditioner to be used in PARMS.
374: Collective on PC
376: Input Parameters:
377: + pc - the preconditioner context
378: - type - the global preconditioner type, one of
379: .vb
380: PC_PARMS_GLOBAL_RAS - Restricted additive Schwarz
381: PC_PARMS_GLOBAL_SCHUR - Schur complement
382: PC_PARMS_GLOBAL_BJ - Block Jacobi
383: .ve
385: Options Database Keys:
386: -pc_parms_global [ras,schur,bj] - Sets global preconditioner
387:
388: Level: intermediate
390: Notes:
391: See the pARMS function parms_PCSetType for more information.
393: .seealso: PCPARMS, PCPARMSSetLocal()
394: @*/
395: PetscErrorCode PCPARMSSetGlobal(PC pc,PCPARMSGlobalType type)
396: {
402: PetscTryMethod(pc,"PCPARMSSetGlobal_C",(PC,PCPARMSGlobalType),(pc,type));
403: return(0);
404: }
406: EXTERN_C_BEGIN
409: PetscErrorCode PCPARMSSetLocal_PARMS(PC pc,PCPARMSLocalType type)
410: {
411: PC_PARMS *parms = (PC_PARMS*)pc->data;
414: if (type != parms->local) {
415: parms->local = type;
416: pc->setupcalled = 0;
417: }
418: return(0);
419: }
420: EXTERN_C_END
424: /*@
425: PCPARMSSetLocal - Sets the local preconditioner to be used in PARMS.
427: Collective on PC
429: Input Parameters:
430: + pc - the preconditioner context
431: - type - the local preconditioner type, one of
432: .vb
433: PC_PARMS_LOCAL_ILU0 - ILU0 preconditioner
434: PC_PARMS_LOCAL_ILUK - ILU(k) preconditioner
435: PC_PARMS_LOCAL_ILUT - ILUT preconditioner
436: PC_PARMS_LOCAL_ARMS - ARMS preconditioner
437: .ve
439: Options Database Keys:
440: -pc_parms_local [ilu0,iluk,ilut,arms] - Sets local preconditioner
441:
442: Level: intermediate
444: Notes:
445: For the ARMS preconditioner, one can use either the symmetric ARMS or the non-symmetric
446: variant (ARMS-ddPQ) by setting the permutation type with PCPARMSSetNonsymPerm().
448: See the pARMS function parms_PCILUSetType for more information.
450: .seealso: PCPARMS, PCPARMSSetGlobal(), PCPARMSSetNonsymPerm()
452: @*/
453: PetscErrorCode PCPARMSSetLocal(PC pc,PCPARMSLocalType type)
454: {
460: PetscTryMethod(pc,"PCPARMSSetLocal_C",(PC,PCPARMSLocalType),(pc,type));
461: return(0);
462: }
464: EXTERN_C_BEGIN
467: PetscErrorCode PCPARMSSetSolveTolerances_PARMS(PC pc,PetscReal tol,PetscInt maxits)
468: {
469: PC_PARMS *parms = (PC_PARMS*)pc->data;
473: if (tol != parms->solvetol) {
474: parms->solvetol = tol;
475: pc->setupcalled = 0;
476: }
477: if (maxits != parms->maxits) {
478: parms->maxits = maxits;
479: pc->setupcalled = 0;
480: }
481:
482: return(0);
483: }
484: EXTERN_C_END
488: /*@
489: PCPARMSSetSolveTolerances - Sets the convergence tolerance and the maximum iterations for the
490: inner GMRES solver, when the Schur global preconditioner is used.
492: Collective on PC
494: Input Parameters:
495: + pc - the preconditioner context
496: . tol - the convergence tolerance
497: - maxits - the maximum number of iterations to use
499: Options Database Keys:
500: + -pc_parms_solve_tol - set the tolerance for local solve
501: - -pc_parms_max_it - set the maximum number of inner iterations
502:
503: Level: intermediate
505: Notes:
506: See the pARMS functions parms_PCSetInnerEps and parms_PCSetInnerMaxits for more information.
508: .seealso: PCPARMS, PCPARMSSetSolveRestart()
509: @*/
510: PetscErrorCode PCPARMSSetSolveTolerances(PC pc,PetscReal tol,PetscInt maxits)
511: {
516: PetscTryMethod(pc,"PCPARMSSetSolveTolerances_C",(PC,PetscReal,PetscInt),(pc,tol,maxits));
517: return(0);
518: }
520: EXTERN_C_BEGIN
523: PetscErrorCode PCPARMSSetSolveRestart_PARMS(PC pc,PetscInt restart)
524: {
525: PC_PARMS *parms = (PC_PARMS*)pc->data;
529: if (restart != parms->maxdim) {
530: parms->maxdim = restart;
531: pc->setupcalled = 0;
532: }
533:
534: return(0);
535: }
536: EXTERN_C_END
540: /*@
541: PCPARMSSetSolveRestart - Sets the number of iterations at which the
542: inner GMRES solver restarts.
544: Collective on PC
546: Input Parameters:
547: + pc - the preconditioner context
548: - restart - maximum dimension of the Krylov subspace
550: Options Database Keys:
551: . -pc_parms_max_dim - sets the inner Krylov dimension
552:
553: Level: intermediate
555: Notes:
556: See the pARMS function parms_PCSetInnerKSize for more information.
558: .seealso: PCPARMS, PCPARMSSetSolveTolerances()
559: @*/
560: PetscErrorCode PCPARMSSetSolveRestart(PC pc,PetscInt restart)
561: {
566: PetscTryMethod(pc,"PCPARMSSetSolveRestart_C",(PC,PetscInt),(pc,restart));
567: return(0);
568: }
570: EXTERN_C_BEGIN
573: PetscErrorCode PCPARMSSetNonsymPerm_PARMS(PC pc,PetscBool nonsym)
574: {
575: PC_PARMS *parms = (PC_PARMS*)pc->data;
578: if ((nonsym && !parms->nonsymperm) || (!nonsym && parms->nonsymperm)) {
579: parms->nonsymperm = nonsym;
580: pc->setupcalled = 0;
581: }
582: return(0);
583: }
584: EXTERN_C_END
588: /*@
589: PCPARMSSetNonsymPerm - Sets the type of permutation for the ARMS preconditioner: the standard
590: symmetric ARMS or the non-symmetric ARMS (ARMS-ddPQ).
592: Collective on PC
594: Input Parameters:
595: + pc - the preconditioner context
596: - nonsym - PETSC_TRUE indicates the non-symmetric ARMS is used;
597: PETSC_FALSE indicates the symmetric ARMS is used
599: Options Database Keys:
600: . -pc_parms_nonsymmetric_perm - sets the use of nonsymmetric permutation
601:
602: Level: intermediate
604: Notes:
605: See the pARMS function parms_PCSetPermType for more information.
607: .seealso: PCPARMS
608: @*/
609: PetscErrorCode PCPARMSSetNonsymPerm(PC pc,PetscBool nonsym)
610: {
615: PetscTryMethod(pc,"PCPARMSSetNonsymPerm_C",(PC,PetscBool),(pc,nonsym));
616: return(0);
617: }
619: EXTERN_C_BEGIN
622: PetscErrorCode PCPARMSSetFill_PARMS(PC pc,PetscInt lfil0,PetscInt lfil1,PetscInt lfil2)
623: {
624: PC_PARMS *parms = (PC_PARMS*)pc->data;
627: if (lfil0 != parms->lfil[0] || lfil0 != parms->lfil[1] || lfil0 != parms->lfil[2] || lfil0 != parms->lfil[3]) {
628: parms->lfil[1] = parms->lfil[2] = parms->lfil[3] = parms->lfil[0] = lfil0;
629: pc->setupcalled = 0;
630: }
631: if (lfil1 != parms->lfil[4]) {
632: parms->lfil[4] = lfil1;
633: pc->setupcalled = 0;
634: }
635: if (lfil2 != parms->lfil[5] || lfil2 != parms->lfil[6]) {
636: parms->lfil[5] = parms->lfil[6] = lfil2;
637: pc->setupcalled = 0;
638: }
639: return(0);
640: }
641: EXTERN_C_END
645: /*@
646: PCPARMSSetFill - Sets the fill-in parameters for ILUT, ILUK and ARMS preconditioners.
647: Consider the original matrix A = [B F; E C] and the approximate version
648: M = [LB 0; E/UB I]*[UB LB\F; 0 S].
650: Collective on PC
652: Input Parameters:
653: + pc - the preconditioner context
654: . fil0 - the level of fill-in kept in LB, UB, E/UB and LB\F
655: . fil1 - the level of fill-in kept in S
656: - fil2 - the level of fill-in kept in the L and U parts of the LU factorization of S
658: Options Database Keys:
659: + -pc_parms_lfil_ilu_arms - set the amount of fill-in for ilut, iluk and arms
660: . -pc_parms_lfil_schur - set the amount of fill-in for schur
661: - -pc_parms_lfil_ilut_L_U - set the amount of fill-in for ILUT L and U
662:
663: Level: intermediate
665: Notes:
666: See the pARMS function parms_PCSetFill for more information.
668: .seealso: PCPARMS
669: @*/
670: PetscErrorCode PCPARMSSetFill(PC pc,PetscInt lfil0,PetscInt lfil1,PetscInt lfil2)
671: {
676: PetscTryMethod(pc,"PCPARMSSetFill_C",(PC,PetscInt,PetscInt,PetscInt),(pc,lfil0,lfil1,lfil2));
677: return(0);
678: }
680: /*MC
681: PCPARMS - Allows the use of the parallel Algebraic Recursive Multilevel Solvers
682: available in the package pARMS
684: Options Database Keys:
685: + -pc_parms_global - one of ras, schur, bj
686: . -pc_parms_local - one of ilu0, iluk, ilut, arms
687: . -pc_parms_solve_tol - set the tolerance for local solve
688: . -pc_parms_levels - set the number of levels
689: . -pc_parms_nonsymmetric_perm - set the use of nonsymmetric permutation
690: . -pc_parms_blocksize - set the block size
691: . -pc_parms_ind_tol - set the tolerance for independent sets
692: . -pc_parms_max_dim - set the inner krylov dimension
693: . -pc_parms_max_it - set the maximum number of inner iterations
694: . -pc_parms_inter_nonsymmetric_perm - set the use of nonsymmetric permutation for interlevel blocks
695: . -pc_parms_inter_column_perm - set the use of column permutation for interlevel blocks
696: . -pc_parms_inter_row_scaling - set the use of row scaling for interlevel blocks
697: . -pc_parms_inter_column_scaling - set the use of column scaling for interlevel blocks
698: . -pc_parms_last_nonsymmetric_perm - set the use of nonsymmetric permutation for last level blocks
699: . -pc_parms_last_column_perm - set the use of column permutation for last level blocks
700: . -pc_parms_last_row_scaling - set the use of row scaling for last level blocks
701: . -pc_parms_last_column_scaling - set the use of column scaling for last level blocks
702: . -pc_parms_lfil_ilu_arms - set the amount of fill-in for ilut, iluk and arms
703: . -pc_parms_lfil_schur - set the amount of fill-in for schur
704: . -pc_parms_lfil_ilut_L_U - set the amount of fill-in for ILUT L and U
705: . -pc_parms_droptol_factors - set the drop tolerance for L, U, L^{-1}F and EU^{-1}
706: . -pc_parms_droptol_schur_compl - set the drop tolerance for schur complement at each level
707: - -pc_parms_droptol_last_schur - set the drop tolerance for ILUT in last level schur complement
709: IMPORTANT:
710: Unless configured appropriately, this preconditioner performs an inexact solve
711: as part of the preconditioner application. Therefore, it must be used in combination
712: with flexible variants of iterative solvers, such as KSPFGMRES or KSPCGR.
714: Level: intermediate
716: .seealso: PCCreate(), PCSetType(), PCType (for list of available types), PC
717: M*/
719: EXTERN_C_BEGIN
722: PetscErrorCode PCCreate_PARMS(PC pc)
723: {
724: PC_PARMS *parms;
725: PetscErrorCode ierr;
728: PetscNewLog(pc,PC_PARMS,&parms);
729: parms->map = 0;
730: parms->A = 0;
731: parms->pc = 0;
732: parms->global = PC_PARMS_GLOBAL_RAS;
733: parms->local = PC_PARMS_LOCAL_ARMS;
734: parms->levels = 10;
735: parms->nonsymperm = PETSC_TRUE;
736: parms->blocksize = 250;
737: parms->maxdim = 0;
738: parms->maxits = 0;
739: parms->meth[0] = PETSC_FALSE;
740: parms->meth[1] = PETSC_FALSE;
741: parms->meth[2] = PETSC_FALSE;
742: parms->meth[3] = PETSC_FALSE;
743: parms->meth[4] = PETSC_FALSE;
744: parms->meth[5] = PETSC_FALSE;
745: parms->meth[6] = PETSC_FALSE;
746: parms->meth[7] = PETSC_FALSE;
747: parms->solvetol = 0.01;
748: parms->indtol = 0.4;
749: parms->lfil[0] = parms->lfil[1] = parms->lfil[2] = parms->lfil[3] = 20;
750: parms->lfil[4] = parms->lfil[5] = parms->lfil[6] = 20;
751: parms->droptol[0] = parms->droptol[1] = parms->droptol[2] = parms->droptol[3] = 0.00001;
752: parms->droptol[4] = 0.001;
753: parms->droptol[5] = parms->droptol[6] = 0.001;
754: pc->data = parms;
755: pc->ops->destroy = PCDestroy_PARMS;
756: pc->ops->setfromoptions = PCSetFromOptions_PARMS;
757: pc->ops->setup = PCSetUp_PARMS;
758: pc->ops->apply = PCApply_PARMS;
759: pc->ops->view = PCView_PARMS;
760: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetGlobal_C","PCPARMSSetGlobal_PARMS",PCPARMSSetGlobal_PARMS);
761: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetLocal_C","PCPARMSSetLocal_PARMS",PCPARMSSetLocal_PARMS);
762: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetSolveTolerances_C","PCPARMSSetSolveTolerances_PARMS",PCPARMSSetSolveTolerances_PARMS);
763: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetSolveRestart_C","PCPARMSSetSolveRestart_PARMS",PCPARMSSetSolveRestart_PARMS);
764: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetNonsymPerm_C","PCPARMSSetNonsymPerm_PARMS",PCPARMSSetNonsymPerm_PARMS);
765: PetscObjectComposeFunctionDynamic((PetscObject)pc,"PCPARMSSetFill_C","PCPARMSSetFill_PARMS",PCPARMSSetFill_PARMS);
766:
767: return(0);
768: }
769: EXTERN_C_END