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:    Level: advanced

141: .keywords: KSP, set, options, prefix, database

143: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix(), KSPSetUseFischerGuess(), KSPSetFischerGuess(), KSPGetFischerInitialGuess()
144: @*/
145: PetscErrorCode  KSPSetUseFischerGuess(KSP ksp,PetscInt model,PetscInt size)
146: {
150:   if (ksp->guess) {
151:     KSPFischerGuessDestroy(ksp->guess);
152:     ksp->guess = PETSC_NULL;
153:   }
154:   if (model == 1 || model == 2) {
155:     KSPFischerGuessCreate(ksp,model,size,&ksp->guess);
156:     KSPFischerGuessSetFromOptions(ksp->guess);
157:   } else if (model != 0) {
158:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Model must be 1 or 2 (or 0 to turn off guess generation)");
159:   }
160:   return(0);
161: }

165: /*@C
166:    KSPSetFischerGuess - Use the Paul Fischer algorithm created by KSPFischerGuessCreate()

168:    Collective on KSP

170:    Input Parameters:
171: +  ksp - the Krylov context
172: -  guess - the object created with KSPFischerGuessCreate()

174:    Level: advanced

176:    Notes: this allows a single KSP to be used with several different initial guess generators (likely for different linear
177:           solvers, see KSPSetPC()).

179:           This increases the reference count of the guess object, you must destroy the object with KSPFischerGuessDestroy()
180:           before the end of the program.

182: .keywords: KSP, set, options, prefix, database

184: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix(), KSPSetUseFischerGuess(), KSPSetFischerGuess(), KSPGetFischerGuess()
185: @*/
186: PetscErrorCode  KSPSetFischerGuess(KSP ksp,KSPFischerGuess guess)
187: {
191:   if (ksp->guess) {
192:     KSPFischerGuessDestroy(ksp->guess);
193:   }
194:   ksp->guess = guess;
195:   if (guess) guess->refcnt++;
196:   return(0);
197: }

201: /*@C
202:    KSPGetFischerGuess - Gets the initial guess generator set with either KSPSetFischerGuess() or KSPCreateFischerGuess()/KSPSetFischerGuess()

204:    Collective on KSP

206:    Input Parameters:
207: .  ksp - the Krylov context

209:    Output Parameters:
210: .   guess - the object

212:    Level: developer

214: .keywords: KSP, set, options, prefix, database

216: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix(), KSPSetUseFischerGuess(), KSPSetFischerGuess()
217: @*/
218: PetscErrorCode  KSPGetFischerGuess(KSP ksp,KSPFischerGuess *guess)
219: {
221:   *guess = ksp->guess;
222:   return(0);
223: }

227: /*@C
228:    KSPGetOptionsPrefix - Gets the prefix used for searching for all 
229:    KSP options in the database.

231:    Not Collective

233:    Input Parameters:
234: .  ksp - the Krylov context

236:    Output Parameters:
237: .  prefix - pointer to the prefix string used is returned

239:    Notes: On the fortran side, the user should pass in a string 'prifix' of
240:    sufficient length to hold the prefix.

242:    Level: advanced

244: .keywords: KSP, set, options, prefix, database

246: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix()
247: @*/
248: PetscErrorCode  KSPGetOptionsPrefix(KSP ksp,const char *prefix[])
249: {
253:   PetscObjectGetOptionsPrefix((PetscObject)ksp,prefix);
254:   return(0);
255: }

259: /*@
260:    KSPSetFromOptions - Sets KSP options from the options database.
261:    This routine must be called before KSPSetUp() if the user is to be 
262:    allowed to set the Krylov type. 

264:    Collective on KSP

266:    Input Parameters:
267: .  ksp - the Krylov space context

269:    Options Database Keys:
270: +   -ksp_max_it - maximum number of linear iterations
271: .   -ksp_rtol rtol - relative tolerance used in default determination of convergence, i.e.
272:                 if residual norm decreases by this factor than convergence is declared
273: .   -ksp_atol abstol - absolute tolerance used in default convergence test, i.e. if residual 
274:                 norm is less than this then convergence is declared
275: .   -ksp_divtol tol - if residual norm increases by this factor than divergence is declared
276: .   -ksp_converged_use_initial_residual_norm - see KSPDefaultConvergedSetUIRNorm()
277: .   -ksp_converged_use_min_initial_residual_norm - see KSPDefaultConvergedSetUMIRNorm()
278: .   -ksp_norm_type - none - skip norms used in convergence tests (useful only when not using 
279: $                       convergence test (say you always want to run with 5 iterations) to 
280: $                       save on communication overhead
281: $                    preconditioned - default for left preconditioning 
282: $                    unpreconditioned - see KSPSetNormType()
283: $                    natural - see KSPSetNormType()
284: .   -ksp_check_norm_iteration it - do not compute residual norm until iteration number it (does compute at 0th iteration)
285: $       works only for PCBCGS, PCIBCGS and and PCCG
286: .   -ksp_fischer_guess <model,size> - uses the Fischer initial guess generator for repeated linear solves
287: .   -ksp_constant_null_space - assume the operator (matrix) has the constant vector in its null space
288: .   -ksp_test_null_space - tests the null space set with KSPSetNullSpace() to see if it truly is a null space
289: .   -ksp_knoll - compute initial guess by applying the preconditioner to the right hand side
290: .   -ksp_monitor_cancel - cancel all previous convergene monitor routines set
291: .   -ksp_monitor <optional filename> - print residual norm at each iteration
292: .   -ksp_monitor_draw - plot residual norm at each iteration
293: .   -ksp_monitor_solution - plot solution at each iteration
294: -   -ksp_monitor_singular_value - monitor extremem singular values at each iteration

296:    Notes:  
297:    To see all options, run your program with the -help option
298:    or consult the users manual.

300:    Level: beginner

302: .keywords: KSP, set, from, options, database

304: .seealso: KSPSetUseFischerInitialGuess()

306: @*/
307: PetscErrorCode  KSPSetFromOptions(KSP ksp)
308: {
309:   PetscErrorCode          ierr;
310:   PetscInt                indx;
311:   const char             *convtests[] = {"default","skip"};
312:   char                    type[256], monfilename[PETSC_MAX_PATH_LEN];
313:   PetscViewerASCIIMonitor monviewer;
314:   PetscTruth              flg,flag;
315:   PetscInt                i,model[2],nmax = 2;
316:   void                    *ctx;

320:   if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
321:   PCSetFromOptions(ksp->pc);

323:   if (!KSPRegisterAllCalled) {KSPRegisterAll(PETSC_NULL);}
324:   PetscOptionsBegin(((PetscObject)ksp)->comm,((PetscObject)ksp)->prefix,"Krylov Method (KSP) Options","KSP");
325:     PetscOptionsList("-ksp_type","Krylov method","KSPSetType",KSPList,(char*)(((PetscObject)ksp)->type_name?((PetscObject)ksp)->type_name:KSPGMRES),type,256,&flg);
326:     if (flg) {
327:       KSPSetType(ksp,type);
328:     }
329:     /*
330:       Set the type if it was never set.
331:     */
332:     if (!((PetscObject)ksp)->type_name) {
333:       KSPSetType(ksp,KSPGMRES);
334:     }

336:     PetscOptionsInt("-ksp_max_it","Maximum number of iterations","KSPSetTolerances",ksp->max_it,&ksp->max_it,PETSC_NULL);
337:     PetscOptionsReal("-ksp_rtol","Relative decrease in residual norm","KSPSetTolerances",ksp->rtol,&ksp->rtol,PETSC_NULL);
338:     PetscOptionsReal("-ksp_atol","Absolute value of residual norm","KSPSetTolerances",ksp->abstol,&ksp->abstol,PETSC_NULL);
339:     PetscOptionsReal("-ksp_divtol","Residual norm increase cause divergence","KSPSetTolerances",ksp->divtol,&ksp->divtol,PETSC_NULL);

341:     PetscOptionsName("-ksp_converged_use_initial_residual_norm","Use initial residual residual norm for computing relative convergence","KSPDefaultConvergedSetUIRNorm",&flag);
342:     if (flag) {KSPDefaultConvergedSetUIRNorm(ksp);}
343:     PetscOptionsName("-ksp_converged_use_min_initial_residual_norm","Use minimum of initial residual norm and b for computing relative convergence","KSPDefaultConvergedSetUMIRNorm",&flag);
344:     if (flag) {KSPDefaultConvergedSetUMIRNorm(ksp);}

346:     PetscOptionsTruth("-ksp_knoll","Use preconditioner applied to b for initial guess","KSPSetInitialGuessKnoll",ksp->guess_knoll,&ksp->guess_knoll,PETSC_NULL);
347:     PetscOptionsIntArray("-ksp_fischer_guess","Use Paul Fischer's algorihtm for initial guess","KSPSetUseFischerGuess",model,&nmax,&flag);
348:     if (flag) {
349:       if (nmax != 2) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Must pass in model,size as arguments");
350:       KSPSetUseFischerGuess(ksp,model[0],model[1]);
351:     }

353:     PetscOptionsEList("-ksp_convergence_test","Convergence test","KSPSetConvergenceTest",convtests,2,"default",&indx,&flg);
354:     if (flg) {
355:       switch (indx) {
356:       case 0:
357:         KSPDefaultConvergedCreate(&ctx);
358:         KSPSetConvergenceTest(ksp,KSPDefaultConverged,ctx,KSPDefaultConvergedDestroy);
359:         break;
360:       case 1: KSPSetConvergenceTest(ksp,KSPSkipConverged,PETSC_NULL,PETSC_NULL);    break;
361:       }
362:     }

364:     PetscOptionsEList("-ksp_norm_type","KSP Norm type","KSPSetNormType",KSPNormTypes,4,"preconditioned",&indx,&flg);
365:     if (flg) { KSPSetNormType(ksp,(KSPNormType)indx); }

367:     PetscOptionsInt("-ksp_check_norm_iteration","First iteration to compute residual norm","KSPSetCheckNormIteration",ksp->chknorm,&ksp->chknorm,PETSC_NULL);

369:     PetscOptionsName("-ksp_lag_norm","Lag the calculation of the residual norm","KSPSetLagNorm",&flg);
370:     if (flg) {
371:       KSPSetLagNorm(ksp,PETSC_TRUE);
372:     }

374:     PetscOptionsName("-ksp_diagonal_scale","Diagonal scale matrix before building preconditioner","KSPSetDiagonalScale",&flg);
375:     if (flg) {
376:       KSPSetDiagonalScale(ksp,PETSC_TRUE);
377:     }
378:     PetscOptionsName("-ksp_diagonal_scale_fix","Fix diagonaled scaled matrix after solve","KSPSetDiagonalScaleFix",&flg);
379:     if (flg) {
380:       KSPSetDiagonalScaleFix(ksp,PETSC_TRUE);
381:     }


384:     PetscOptionsName("-ksp_constant_null_space","Add constant null space to Krylov solver","KSPSetNullSpace",&flg);
385:     if (flg) {
386:       MatNullSpace nsp;

388:       MatNullSpaceCreate(((PetscObject)ksp)->comm,PETSC_TRUE,0,0,&nsp);
389:       KSPSetNullSpace(ksp,nsp);
390:       MatNullSpaceDestroy(nsp);
391:     }

393:     /* option is actually checked in KSPSetUp() */
394:     if (ksp->nullsp) {
395:       PetscOptionsName("-ksp_test_null_space","Is provided null space correct","None",&flg);
396:     }

398:     /*
399:       Prints reason for convergence or divergence of each linear solve
400:     */
401:     PetscOptionsName("-ksp_converged_reason","Print reason for converged or diverged","KSPSolve",&flg);
402:     if (flg) {
403:       ksp->printreason = PETSC_TRUE;
404:     }

406:     PetscOptionsName("-ksp_monitor_cancel","Remove any hardwired monitor routines","KSPMonitorCancel",&flg);
407:     /* -----------------------------------------------------------------------*/
408:     /*
409:       Cancels all monitors hardwired into code before call to KSPSetFromOptions()
410:     */
411:     if (flg) {
412:       KSPMonitorCancel(ksp);
413:     }
414:     /*
415:       Prints preconditioned residual norm at each iteration
416:     */
417:     PetscOptionsString("-ksp_monitor","Monitor preconditioned residual norm","KSPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
418:     if (flg) {
419:       PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
420:       KSPMonitorSet(ksp,KSPMonitorDefault,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
421:     }
422:     /*
423:       Prints preconditioned residual norm at each iteration
424:     */
425:     PetscOptionsString("-ksp_monitor_range","Monitor percent of residual entries more than 10 percent of max","KSPMonitorRange","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
426:     if (flg) {
427:       PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
428:       KSPMonitorSet(ksp,KSPMonitorRange,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
429:     }
430:     /*
431:       Plots the vector solution 
432:     */
433:     PetscOptionsName("-ksp_monitor_solution","Monitor solution graphically","KSPMonitorSet",&flg);
434:     if (flg) {
435:       KSPMonitorSet(ksp,KSPMonitorSolution,PETSC_NULL,PETSC_NULL);
436:     }
437:     /*
438:       Prints preconditioned and true residual norm at each iteration
439:     */
440:     PetscOptionsString("-ksp_monitor_true_residual","Monitor true residual norm","KSPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
441:     if (flg) {
442:       PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
443:       KSPMonitorSet(ksp,KSPMonitorTrueResidualNorm,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
444:     }
445:     /*
446:       Prints extreme eigenvalue estimates at each iteration
447:     */
448:     PetscOptionsString("-ksp_monitor_singular_value","Monitor singular values","KSPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
449:     if (flg) {
450:       KSPSetComputeSingularValues(ksp,PETSC_TRUE);
451:       PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
452:       KSPMonitorSet(ksp,KSPMonitorSingularValue,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
453:     }
454:     /*
455:       Prints preconditioned residual norm with fewer digits
456:     */
457:     PetscOptionsString("-ksp_monitor_short","Monitor preconditioned residual norm with fewer digits","KSPMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
458:     if (flg) {
459:       PetscViewerASCIIMonitorCreate(((PetscObject)ksp)->comm,monfilename,((PetscObject)ksp)->tablevel,&monviewer);
460:       KSPMonitorSet(ksp,KSPMonitorDefaultShort,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
461:     }
462:     /*
463:       Graphically plots preconditioned residual norm
464:     */
465:     PetscOptionsName("-ksp_monitor_draw","Monitor graphically preconditioned residual norm","KSPMonitorSet",&flg);
466:     if (flg) {
467:       KSPMonitorSet(ksp,KSPMonitorLG,PETSC_NULL,PETSC_NULL);
468:     }
469:     /*
470:       Graphically plots preconditioned and true residual norm
471:     */
472:     PetscOptionsName("-ksp_monitor_draw_true_residual","Monitor graphically true residual norm","KSPMonitorSet",&flg);
473:     if (flg){
474:       KSPMonitorSet(ksp,KSPMonitorLGTrueResidualNorm,PETSC_NULL,PETSC_NULL);
475:     }
476:     /*
477:       Graphically plots preconditioned residual norm and range of residual element values
478:     */
479:     PetscOptionsName("-ksp_monitor_range_draw","Monitor graphically preconditioned residual norm","KSPMonitorSet",&flg);
480:     if (flg) {
481:       KSPMonitorSet(ksp,KSPMonitorLGRange,PETSC_NULL,PETSC_NULL);
482:     }

484:     /* -----------------------------------------------------------------------*/

486:     PetscOptionsTruthGroupBegin("-ksp_left_pc","Use left preconditioning","KSPSetPreconditionerSide",&flg);
487:     if (flg) { KSPSetPreconditionerSide(ksp,PC_LEFT); }
488:     PetscOptionsTruthGroup("-ksp_right_pc","Use right preconditioning","KSPSetPreconditionerSide",&flg);
489:     if (flg) { KSPSetPreconditionerSide(ksp,PC_RIGHT);}
490:     PetscOptionsTruthGroupEnd("-ksp_symmetric_pc","Use symmetric (factorized) preconditioning","KSPSetPreconditionerSide",&flg);
491:     if (flg) { KSPSetPreconditionerSide(ksp,PC_SYMMETRIC);}

493:     PetscOptionsName("-ksp_compute_singularvalues","Compute singular values of preconditioned operator","KSPSetComputeSingularValues",&flg);
494:     if (flg) { KSPSetComputeSingularValues(ksp,PETSC_TRUE); }
495:     PetscOptionsName("-ksp_compute_eigenvalues","Compute eigenvalues of preconditioned operator","KSPSetComputeSingularValues",&flg);
496:     if (flg) { KSPSetComputeSingularValues(ksp,PETSC_TRUE); }
497:     PetscOptionsName("-ksp_plot_eigenvalues","Scatter plot extreme eigenvalues","KSPSetComputeSingularValues",&flg);
498:     if (flg) { KSPSetComputeSingularValues(ksp,PETSC_TRUE); }

500:     for(i = 0; i < numberofsetfromoptions; i++) {
501:       (*othersetfromoptions[i])(ksp);
502:     }

504:     if (ksp->ops->setfromoptions) {
505:       (*ksp->ops->setfromoptions)(ksp);
506:     }
507:     PetscOptionsName("-ksp_view","View linear solver parameters","KSPView",&flg);
508:   PetscOptionsEnd();
509:   return(0);
510: }