Actual source code: itcl.c
1: #define PETSCKSP_DLL
3: /*
4: Code for setting KSP options from the options database.
5: */
7: #include private/kspimpl.h
8: #include petscsys.h
10: /*
11: We retain a list of functions that also take KSP command
12: line options. These are called at the end KSPSetFromOptions()
13: */
14: #define MAXSETFROMOPTIONS 5
15: PetscInt numberofsetfromoptions = 0;
16: PetscErrorCode (*othersetfromoptions[MAXSETFROMOPTIONS])(KSP) = {0};
22: /*@C
23: KSPAddOptionsChecker - Adds an additional function to check for KSP options.
25: Not Collective
27: Input Parameter:
28: . kspcheck - function that checks for options
30: Level: developer
32: .keywords: KSP, add, options, checker
34: .seealso: KSPSetFromOptions()
35: @*/
36: PetscErrorCode KSPAddOptionsChecker(PetscErrorCode (*kspcheck)(KSP))
37: {
39: if (numberofsetfromoptions >= MAXSETFROMOPTIONS) {
40: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many options checkers, only 5 allowed");
41: }
43: othersetfromoptions[numberofsetfromoptions++] = kspcheck;
44: return(0);
45: }
49: /*@C
50: KSPSetOptionsPrefix - Sets the prefix used for searching for all
51: KSP options in the database.
53: Collective on KSP
55: Input Parameters:
56: + ksp - the Krylov context
57: - prefix - the prefix string to prepend to all KSP option requests
59: Notes:
60: A hyphen (-) must NOT be given at the beginning of the prefix name.
61: The first character of all runtime options is AUTOMATICALLY the
62: hyphen.
64: For example, to distinguish between the runtime options for two
65: different KSP contexts, one could call
66: .vb
67: KSPSetOptionsPrefix(ksp1,"sys1_")
68: KSPSetOptionsPrefix(ksp2,"sys2_")
69: .ve
71: This would enable use of different options for each system, such as
72: .vb
73: -sys1_ksp_type gmres -sys1_ksp_rtol 1.e-3
74: -sys2_ksp_type bcgs -sys2_ksp_rtol 1.e-4
75: .ve
77: Level: advanced
79: .keywords: KSP, set, options, prefix, database
81: .seealso: KSPAppendOptionsPrefix(), KSPGetOptionsPrefix()
82: @*/
83: PetscErrorCode KSPSetOptionsPrefix(KSP ksp,const char prefix[])
84: {
88: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
89: PCSetOptionsPrefix(ksp->pc,prefix);
90: PetscObjectSetOptionsPrefix((PetscObject)ksp,prefix);
91: return(0);
92: }
93:
96: /*@C
97: KSPAppendOptionsPrefix - Appends to the prefix used for searching for all
98: KSP options in the database.
100: Collective on KSP
102: Input Parameters:
103: + ksp - the Krylov context
104: - prefix - the prefix string to prepend to all KSP option requests
106: Notes:
107: A hyphen (-) must NOT be given at the beginning of the prefix name.
108: The first character of all runtime options is AUTOMATICALLY the hyphen.
110: Level: advanced
112: .keywords: KSP, append, options, prefix, database
114: .seealso: KSPSetOptionsPrefix(), KSPGetOptionsPrefix()
115: @*/
116: PetscErrorCode KSPAppendOptionsPrefix(KSP ksp,const char prefix[])
117: {
121: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
122: PCAppendOptionsPrefix(ksp->pc,prefix);
123: PetscObjectAppendOptionsPrefix((PetscObject)ksp,prefix);
124: return(0);
125: }
129: /*@C
130: KSPSetUseFischerGuess - Use the Paul Fischer algorithm, see KSPFischerGuessCreate()
132: Collective on KSP
134: Input Parameters:
135: + ksp - the Krylov context
136: . model - use model 1, model 2 or 0 to turn it off
137: - size - size of subspace used to generate initial guess
139: Options Database:
140: . -ksp_fischer_guess <model,size> - uses the Fischer initial guess generator for repeated linear solves
142: Level: advanced
144: .keywords: KSP, set, options, prefix, database
146: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix(), KSPSetUseFischerGuess(), KSPSetFischerGuess(), KSPGetFischerInitialGuess()
147: @*/
148: PetscErrorCode KSPSetUseFischerGuess(KSP ksp,PetscInt model,PetscInt size)
149: {
153: if (ksp->guess) {
154: KSPFischerGuessDestroy(ksp->guess);
155: ksp->guess = PETSC_NULL;
156: }
157: if (model == 1 || model == 2) {
158: KSPFischerGuessCreate(ksp,model,size,&ksp->guess);
159: KSPFischerGuessSetFromOptions(ksp->guess);
160: } else if (model != 0) {
161: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Model must be 1 or 2 (or 0 to turn off guess generation)");
162: }
163: return(0);
164: }
168: /*@C
169: KSPSetFischerGuess - Use the Paul Fischer algorithm created by KSPFischerGuessCreate()
171: Collective on KSP
173: Input Parameters:
174: + ksp - the Krylov context
175: - guess - the object created with KSPFischerGuessCreate()
177: Level: advanced
179: Notes: this allows a single KSP to be used with several different initial guess generators (likely for different linear
180: solvers, see KSPSetPC()).
182: This increases the reference count of the guess object, you must destroy the object with KSPFischerGuessDestroy()
183: before the end of the program.
185: .keywords: KSP, set, options, prefix, database
187: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix(), KSPSetUseFischerGuess(), KSPSetFischerGuess(), KSPGetFischerGuess()
188: @*/
189: PetscErrorCode KSPSetFischerGuess(KSP ksp,KSPFischerGuess guess)
190: {
194: if (ksp->guess) {
195: KSPFischerGuessDestroy(ksp->guess);
196: }
197: ksp->guess = guess;
198: if (guess) guess->refcnt++;
199: return(0);
200: }
204: /*@C
205: KSPGetFischerGuess - Gets the initial guess generator set with either KSPSetFischerGuess() or KSPCreateFischerGuess()/KSPSetFischerGuess()
207: Collective on KSP
209: Input Parameters:
210: . ksp - the Krylov context
212: Output Parameters:
213: . guess - the object
215: Level: developer
217: .keywords: KSP, set, options, prefix, database
219: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix(), KSPSetUseFischerGuess(), KSPSetFischerGuess()
220: @*/
221: PetscErrorCode KSPGetFischerGuess(KSP ksp,KSPFischerGuess *guess)
222: {
224: *guess = ksp->guess;
225: return(0);
226: }
230: /*@C
231: KSPGetOptionsPrefix - Gets the prefix used for searching for all
232: KSP options in the database.
234: Not Collective
236: Input Parameters:
237: . ksp - the Krylov context
239: Output Parameters:
240: . prefix - pointer to the prefix string used is returned
242: Notes: On the fortran side, the user should pass in a string 'prifix' of
243: sufficient length to hold the prefix.
245: Level: advanced
247: .keywords: KSP, set, options, prefix, database
249: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix()
250: @*/
251: PetscErrorCode KSPGetOptionsPrefix(KSP ksp,const char *prefix[])
252: {
256: PetscObjectGetOptionsPrefix((PetscObject)ksp,prefix);
257: return(0);
258: }
262: /*@
263: KSPSetFromOptions - Sets KSP options from the options database.
264: This routine must be called before KSPSetUp() if the user is to be
265: allowed to set the Krylov type.
267: Collective on KSP
269: Input Parameters:
270: . ksp - the Krylov space context
272: Options Database Keys:
273: + -ksp_max_it - maximum number of linear iterations
274: . -ksp_rtol rtol - relative tolerance used in default determination of convergence, i.e.
275: if residual norm decreases by this factor than convergence is declared
276: . -ksp_atol abstol - absolute tolerance used in default convergence test, i.e. if residual
277: norm is less than this then convergence is declared
278: . -ksp_divtol tol - if residual norm increases by this factor than divergence is declared
279: . -ksp_converged_use_initial_residual_norm - see KSPDefaultConvergedSetUIRNorm()
280: . -ksp_converged_use_min_initial_residual_norm - see KSPDefaultConvergedSetUMIRNorm()
281: . -ksp_norm_type - none - skip norms used in convergence tests (useful only when not using
282: $ convergence test (say you always want to run with 5 iterations) to
283: $ save on communication overhead
284: $ preconditioned - default for left preconditioning
285: $ unpreconditioned - see KSPSetNormType()
286: $ natural - see KSPSetNormType()
287: . -ksp_check_norm_iteration it - do not compute residual norm until iteration number it (does compute at 0th iteration)
288: $ works only for PCBCGS, PCIBCGS and and PCCG
289: . -ksp_fischer_guess <model,size> - uses the Fischer initial guess generator for repeated linear solves
290: . -ksp_constant_null_space - assume the operator (matrix) has the constant vector in its null space
291: . -ksp_test_null_space - tests the null space set with KSPSetNullSpace() to see if it truly is a null space
292: . -ksp_knoll - compute initial guess by applying the preconditioner to the right hand side
293: . -ksp_monitor_cancel - cancel all previous convergene monitor routines set
294: . -ksp_monitor <optional filename> - print residual norm at each iteration
295: . -ksp_monitor_draw - plot residual norm at each iteration
296: . -ksp_monitor_solution - plot solution at each iteration
297: - -ksp_monitor_singular_value - monitor extremem singular values at each iteration
299: Notes:
300: To see all options, run your program with the -help option
301: or consult the users manual.
303: Level: beginner
305: .keywords: KSP, set, from, options, database
307: .seealso: KSPSetUseFischerInitialGuess()
309: @*/
310: PetscErrorCode KSPSetFromOptions(KSP ksp)
311: {
312: PetscErrorCode ierr;
313: PetscInt indx;
314: const char *convtests[] = {"default","skip"};
315: char type[256], monfilename[PETSC_MAX_PATH_LEN];
316: PetscViewerASCIIMonitor monviewer;
317: PetscTruth flg,flag;
318: PetscInt i,model[2],nmax = 2;
319: void *ctx;
323: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
324: PCSetFromOptions(ksp->pc);
326: if (!KSPRegisterAllCalled) {KSPRegisterAll(PETSC_NULL);}
327: PetscOptionsBegin(((PetscObject)ksp)->comm,((PetscObject)ksp)->prefix,"Krylov Method (KSP) Options","KSP");
328: PetscOptionsList("-ksp_type","Krylov method","KSPSetType",KSPList,(char*)(((PetscObject)ksp)->type_name?((PetscObject)ksp)->type_name:KSPGMRES),type,256,&flg);
329: if (flg) {
330: KSPSetType(ksp,type);
331: }
332: /*
333: Set the type if it was never set.
334: */
335: if (!((PetscObject)ksp)->type_name) {
336: KSPSetType(ksp,KSPGMRES);
337: }
339: PetscOptionsInt("-ksp_max_it","Maximum number of iterations","KSPSetTolerances",ksp->max_it,&ksp->max_it,PETSC_NULL);
340: PetscOptionsReal("-ksp_rtol","Relative decrease in residual norm","KSPSetTolerances",ksp->rtol,&ksp->rtol,PETSC_NULL);
341: PetscOptionsReal("-ksp_atol","Absolute value of residual norm","KSPSetTolerances",ksp->abstol,&ksp->abstol,PETSC_NULL);
342: PetscOptionsReal("-ksp_divtol","Residual norm increase cause divergence","KSPSetTolerances",ksp->divtol,&ksp->divtol,PETSC_NULL);
344: PetscOptionsName("-ksp_converged_use_initial_residual_norm","Use initial residual residual norm for computing relative convergence","KSPDefaultConvergedSetUIRNorm",&flag);
345: if (flag) {KSPDefaultConvergedSetUIRNorm(ksp);}
346: PetscOptionsName("-ksp_converged_use_min_initial_residual_norm","Use minimum of initial residual norm and b for computing relative convergence","KSPDefaultConvergedSetUMIRNorm",&flag);
347: if (flag) {KSPDefaultConvergedSetUMIRNorm(ksp);}
349: PetscOptionsTruth("-ksp_knoll","Use preconditioner applied to b for initial guess","KSPSetInitialGuessKnoll",ksp->guess_knoll,&ksp->guess_knoll,PETSC_NULL);
350: PetscOptionsIntArray("-ksp_fischer_guess","Use Paul Fischer's algorihtm for initial guess","KSPSetUseFischerGuess",model,&nmax,&flag);
351: if (flag) {
352: if (nmax != 2) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Must pass in model,size as arguments");
353: KSPSetUseFischerGuess(ksp,model[0],model[1]);
354: }
356: PetscOptionsEList("-ksp_convergence_test","Convergence test","KSPSetConvergenceTest",convtests,2,"default",&indx,&flg);
357: if (flg) {
358: switch (indx) {
359: case 0:
360: KSPDefaultConvergedCreate(&ctx);
361: KSPSetConvergenceTest(ksp,KSPDefaultConverged,ctx,KSPDefaultConvergedDestroy);
362: break;
363: case 1: KSPSetConvergenceTest(ksp,KSPSkipConverged,PETSC_NULL,PETSC_NULL); break;
364: }
365: }
367: PetscOptionsEList("-ksp_norm_type","KSP Norm type","KSPSetNormType",KSPNormTypes,4,"preconditioned",&indx,&flg);
368: if (flg) { KSPSetNormType(ksp,(KSPNormType)indx); }
370: PetscOptionsInt("-ksp_check_norm_iteration","First iteration to compute residual norm","KSPSetCheckNormIteration",ksp->chknorm,&ksp->chknorm,PETSC_NULL);
372: PetscOptionsName("-ksp_lag_norm","Lag the calculation of the residual norm","KSPSetLagNorm",&flg);
373: if (flg) {
374: KSPSetLagNorm(ksp,PETSC_TRUE);
375: }
377: PetscOptionsName("-ksp_diagonal_scale","Diagonal scale matrix before building preconditioner","KSPSetDiagonalScale",&flg);
378: if (flg) {
379: KSPSetDiagonalScale(ksp,PETSC_TRUE);
380: }
381: PetscOptionsName("-ksp_diagonal_scale_fix","Fix diagonaled scaled matrix after solve","KSPSetDiagonalScaleFix",&flg);
382: if (flg) {
383: KSPSetDiagonalScaleFix(ksp,PETSC_TRUE);
384: }
387: PetscOptionsName("-ksp_constant_null_space","Add constant null space to Krylov solver","KSPSetNullSpace",&flg);
388: if (flg) {
389: MatNullSpace nsp;
391: MatNullSpaceCreate(((PetscObject)ksp)->comm,PETSC_TRUE,0,0,&nsp);
392: KSPSetNullSpace(ksp,nsp);
393: MatNullSpaceDestroy(nsp);
394: }
396: /* option is actually checked in KSPSetUp() */
397: if (ksp->nullsp) {
398: PetscOptionsName("-ksp_test_null_space","Is provided null space correct","None",&flg);
399: }
401: /*
402: Prints reason for convergence or divergence of each linear solve
403: */
404: PetscOptionsName("-ksp_converged_reason","Print reason for converged or diverged","KSPSolve",&flg);
405: if (flg) {
406: ksp->printreason = PETSC_TRUE;
407: }
409: PetscOptionsName("-ksp_monitor_cancel","Remove any hardwired monitor routines","KSPMonitorCancel",&flg);
410: /* -----------------------------------------------------------------------*/
411: /*
412: Cancels all monitors hardwired into code before call to KSPSetFromOptions()
413: */
414: if (flg) {
415: KSPMonitorCancel(ksp);
416: }
417: /*
418: Prints preconditioned residual norm at each iteration
419: */
420: PetscOptionsString("-ksp_monitor","Monitor preconditioned residual norm","KSPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
421: if (flg) {
422: PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
423: KSPMonitorSet(ksp,KSPMonitorDefault,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
424: }
425: /*
426: Prints preconditioned residual norm at each iteration
427: */
428: PetscOptionsString("-ksp_monitor_range","Monitor percent of residual entries more than 10 percent of max","KSPMonitorRange","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
429: if (flg) {
430: PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
431: KSPMonitorSet(ksp,KSPMonitorRange,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
432: }
433: /*
434: Plots the vector solution
435: */
436: PetscOptionsName("-ksp_monitor_solution","Monitor solution graphically","KSPMonitorSet",&flg);
437: if (flg) {
438: KSPMonitorSet(ksp,KSPMonitorSolution,PETSC_NULL,PETSC_NULL);
439: }
440: /*
441: Prints preconditioned and true residual norm at each iteration
442: */
443: PetscOptionsString("-ksp_monitor_true_residual","Monitor true residual norm","KSPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
444: if (flg) {
445: PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
446: KSPMonitorSet(ksp,KSPMonitorTrueResidualNorm,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
447: }
448: /*
449: Prints extreme eigenvalue estimates at each iteration
450: */
451: PetscOptionsString("-ksp_monitor_singular_value","Monitor singular values","KSPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
452: if (flg) {
453: KSPSetComputeSingularValues(ksp,PETSC_TRUE);
454: PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
455: KSPMonitorSet(ksp,KSPMonitorSingularValue,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
456: }
457: /*
458: Prints preconditioned residual norm with fewer digits
459: */
460: PetscOptionsString("-ksp_monitor_short","Monitor preconditioned residual norm with fewer digits","KSPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
461: if (flg) {
462: PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
463: KSPMonitorSet(ksp,KSPMonitorDefaultShort,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
464: }
465: /*
466: Graphically plots preconditioned residual norm
467: */
468: PetscOptionsName("-ksp_monitor_draw","Monitor graphically preconditioned residual norm","KSPMonitorSet",&flg);
469: if (flg) {
470: KSPMonitorSet(ksp,KSPMonitorLG,PETSC_NULL,PETSC_NULL);
471: }
472: /*
473: Graphically plots preconditioned and true residual norm
474: */
475: PetscOptionsName("-ksp_monitor_draw_true_residual","Monitor graphically true residual norm","KSPMonitorSet",&flg);
476: if (flg){
477: KSPMonitorSet(ksp,KSPMonitorLGTrueResidualNorm,PETSC_NULL,PETSC_NULL);
478: }
479: /*
480: Graphically plots preconditioned residual norm and range of residual element values
481: */
482: PetscOptionsName("-ksp_monitor_range_draw","Monitor graphically preconditioned residual norm","KSPMonitorSet",&flg);
483: if (flg) {
484: KSPMonitorSet(ksp,KSPMonitorLGRange,PETSC_NULL,PETSC_NULL);
485: }
487: /* -----------------------------------------------------------------------*/
489: PetscOptionsTruthGroupBegin("-ksp_left_pc","Use left preconditioning","KSPSetPreconditionerSide",&flg);
490: if (flg) { KSPSetPreconditionerSide(ksp,PC_LEFT); }
491: PetscOptionsTruthGroup("-ksp_right_pc","Use right preconditioning","KSPSetPreconditionerSide",&flg);
492: if (flg) { KSPSetPreconditionerSide(ksp,PC_RIGHT);}
493: PetscOptionsTruthGroupEnd("-ksp_symmetric_pc","Use symmetric (factorized) preconditioning","KSPSetPreconditionerSide",&flg);
494: if (flg) { KSPSetPreconditionerSide(ksp,PC_SYMMETRIC);}
496: PetscOptionsName("-ksp_compute_singularvalues","Compute singular values of preconditioned operator","KSPSetComputeSingularValues",&flg);
497: if (flg) { KSPSetComputeSingularValues(ksp,PETSC_TRUE); }
498: PetscOptionsName("-ksp_compute_eigenvalues","Compute eigenvalues of preconditioned operator","KSPSetComputeSingularValues",&flg);
499: if (flg) { KSPSetComputeSingularValues(ksp,PETSC_TRUE); }
500: PetscOptionsName("-ksp_plot_eigenvalues","Scatter plot extreme eigenvalues","KSPSetComputeSingularValues",&flg);
501: if (flg) { KSPSetComputeSingularValues(ksp,PETSC_TRUE); }
503: for(i = 0; i < numberofsetfromoptions; i++) {
504: (*othersetfromoptions[i])(ksp);
505: }
507: if (ksp->ops->setfromoptions) {
508: (*ksp->ops->setfromoptions)(ksp);
509: }
510: PetscOptionsName("-ksp_view","View linear solver parameters","KSPView",&flg);
511: PetscOptionsEnd();
512: return(0);
513: }