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: }