Actual source code: itcl.c
petsc-3.7.7 2017-09-25
2: /*
3: Code for setting KSP options from the options database.
4: */
6: #include <petsc/private/kspimpl.h> /*I "petscksp.h" I*/
7: #include <petscdraw.h>
9: extern PetscBool KSPRegisterAllCalled;
13: /*@C
14: KSPSetOptionsPrefix - Sets the prefix used for searching for all
15: KSP options in the database.
17: Logically Collective on KSP
19: Input Parameters:
20: + ksp - the Krylov context
21: - prefix - the prefix string to prepend to all KSP option requests
23: Notes:
24: A hyphen (-) must NOT be given at the beginning of the prefix name.
25: The first character of all runtime options is AUTOMATICALLY the
26: hyphen.
28: For example, to distinguish between the runtime options for two
29: different KSP contexts, one could call
30: .vb
31: KSPSetOptionsPrefix(ksp1,"sys1_")
32: KSPSetOptionsPrefix(ksp2,"sys2_")
33: .ve
35: This would enable use of different options for each system, such as
36: .vb
37: -sys1_ksp_type gmres -sys1_ksp_rtol 1.e-3
38: -sys2_ksp_type bcgs -sys2_ksp_rtol 1.e-4
39: .ve
41: Level: advanced
43: .keywords: KSP, set, options, prefix, database
45: .seealso: KSPAppendOptionsPrefix(), KSPGetOptionsPrefix()
46: @*/
47: PetscErrorCode KSPSetOptionsPrefix(KSP ksp,const char prefix[])
48: {
53: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
54: PCSetOptionsPrefix(ksp->pc,prefix);
55: PetscObjectSetOptionsPrefix((PetscObject)ksp,prefix);
56: return(0);
57: }
61: /*@C
62: KSPAppendOptionsPrefix - Appends to the prefix used for searching for all
63: KSP options in the database.
65: Logically Collective on KSP
67: Input Parameters:
68: + ksp - the Krylov context
69: - prefix - the prefix string to prepend to all KSP option requests
71: Notes:
72: A hyphen (-) must NOT be given at the beginning of the prefix name.
73: The first character of all runtime options is AUTOMATICALLY the hyphen.
75: Level: advanced
77: .keywords: KSP, append, options, prefix, database
79: .seealso: KSPSetOptionsPrefix(), KSPGetOptionsPrefix()
80: @*/
81: PetscErrorCode KSPAppendOptionsPrefix(KSP ksp,const char prefix[])
82: {
87: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
88: PCAppendOptionsPrefix(ksp->pc,prefix);
89: PetscObjectAppendOptionsPrefix((PetscObject)ksp,prefix);
90: return(0);
91: }
95: /*@
96: KSPGetTabLevel - Gets the number of tabs that ASCII output used by ksp.
98: Not Collective
100: Input Parameter:
101: . ksp - a KSP object.
103: Output Parameter:
104: . tab - the number of tabs
106: Level: developer
108: Notes: this is used in conjunction with KSPSetTabLevel() to manage the output from the KSP and its PC coherently.
111: .seealso: KSPSetTabLevel()
113: @*/
114: PetscErrorCode KSPGetTabLevel(KSP ksp,PetscInt *tab)
115: {
120: PetscObjectGetTabLevel((PetscObject)ksp, tab);
121: return(0);
122: }
126: /*@
127: KSPSetTabLevel - Sets the number of tabs that ASCII output for the ksp andn its pc will use.
129: Not Collective
131: Input Parameters:
132: + ksp - a KSP object
133: - tab - the number of tabs
135: Level: developer
137: Notes: this is used to manage the output from KSP and PC objects that are imbedded in other objects,
138: for example, the KSP object inside a SNES object. By indenting each lower level further the heirarchy
139: of objects is very clear. By setting the KSP object's tab level with KSPSetTabLevel() its PC object
140: automatically receives the same tab level, so that whatever objects the pc might create are tabbed
141: appropriately, too.
143: .seealso: KSPGetTabLevel()
144: @*/
145: PetscErrorCode KSPSetTabLevel(KSP ksp, PetscInt tab)
146: {
151: PetscObjectSetTabLevel((PetscObject)ksp, tab);
152: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
153: /* Do we need a PCSetTabLevel()? */
154: PetscObjectSetTabLevel((PetscObject)ksp->pc, tab);
155: return(0);
156: }
160: /*@C
161: KSPSetUseFischerGuess - Use the Paul Fischer algorithm, see KSPFischerGuessCreate()
163: Logically Collective on KSP
165: Input Parameters:
166: + ksp - the Krylov context
167: . model - use model 1, model 2 or 0 to turn it off
168: - size - size of subspace used to generate initial guess
170: Options Database:
171: . -ksp_fischer_guess <model,size> - uses the Fischer initial guess generator for repeated linear solves
173: Level: advanced
175: .keywords: KSP, set, options, prefix, database
177: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix(), KSPSetUseFischerGuess(), KSPSetFischerGuess(), KSPGetFischerInitialGuess()
178: @*/
179: PetscErrorCode KSPSetUseFischerGuess(KSP ksp,PetscInt model,PetscInt size)
180: {
187: KSPFischerGuessDestroy(&ksp->guess);
188: if (model == 1 || model == 2) {
189: KSPFischerGuessCreate(ksp,model,size,&ksp->guess);
190: KSPFischerGuessSetFromOptions(ksp->guess);
191: } else if (model != 0) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE,"Model must be 1 or 2 (or 0 to turn off guess generation)");
192: return(0);
193: }
197: /*@C
198: KSPSetFischerGuess - Use the Paul Fischer algorithm created by KSPFischerGuessCreate()
200: Logically Collective on KSP
202: Input Parameters:
203: + ksp - the Krylov context
204: - guess - the object created with KSPFischerGuessCreate()
206: Level: advanced
208: Notes: this allows a single KSP to be used with several different initial guess generators (likely for different linear
209: solvers, see KSPSetPC()).
211: This increases the reference count of the guess object, you must destroy the object with KSPFischerGuessDestroy()
212: before the end of the program.
214: .keywords: KSP, set, options, prefix, database
216: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix(), KSPSetUseFischerGuess(), KSPSetFischerGuess(), KSPGetFischerGuess()
217: @*/
218: PetscErrorCode KSPSetFischerGuess(KSP ksp,KSPFischerGuess guess)
219: {
224: KSPFischerGuessDestroy(&ksp->guess);
225: ksp->guess = guess;
226: if (guess) guess->refcnt++;
227: return(0);
228: }
232: /*@C
233: KSPGetFischerGuess - Gets the initial guess generator set with either KSPSetFischerGuess() or KSPCreateFischerGuess()/KSPSetFischerGuess()
235: Not Collective
237: Input Parameters:
238: . ksp - the Krylov context
240: Output Parameters:
241: . guess - the object
243: Level: developer
245: .keywords: KSP, set, options, prefix, database
247: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix(), KSPSetUseFischerGuess(), KSPSetFischerGuess()
248: @*/
249: PetscErrorCode KSPGetFischerGuess(KSP ksp,KSPFischerGuess *guess)
250: {
252: *guess = ksp->guess;
253: return(0);
254: }
258: /*@C
259: KSPGetOptionsPrefix - Gets the prefix used for searching for all
260: KSP options in the database.
262: Not Collective
264: Input Parameters:
265: . ksp - the Krylov context
267: Output Parameters:
268: . prefix - pointer to the prefix string used is returned
270: Notes: On the fortran side, the user should pass in a string 'prifix' of
271: sufficient length to hold the prefix.
273: Level: advanced
275: .keywords: KSP, set, options, prefix, database
277: .seealso: KSPSetOptionsPrefix(), KSPAppendOptionsPrefix()
278: @*/
279: PetscErrorCode KSPGetOptionsPrefix(KSP ksp,const char *prefix[])
280: {
285: PetscObjectGetOptionsPrefix((PetscObject)ksp,prefix);
286: return(0);
287: }
291: /*@C
292: KSPMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type indicated by the user
294: Collective on KSP
296: Input Parameters:
297: + ksp - KSP object you wish to monitor
298: . name - the monitor type one is seeking
299: . help - message indicating what monitoring is done
300: . manual - manual page for the monitor
301: - monitor - the monitor function, the context for this object is a PetscViewerAndFormat
303: Level: developer
305: .seealso: PetscOptionsGetViewer(), PetscOptionsGetReal(), PetscOptionsHasName(), PetscOptionsGetString(),
306: PetscOptionsGetIntArray(), PetscOptionsGetRealArray(), PetscOptionsBool()
307: PetscOptionsInt(), PetscOptionsString(), PetscOptionsReal(), PetscOptionsBool(),
308: PetscOptionsName(), PetscOptionsBegin(), PetscOptionsEnd(), PetscOptionsHead(),
309: PetscOptionsStringArray(),PetscOptionsRealArray(), PetscOptionsScalar(),
310: PetscOptionsBoolGroupBegin(), PetscOptionsBoolGroup(), PetscOptionsBoolGroupEnd(),
311: PetscOptionsFList(), PetscOptionsEList()
312: @*/
313: PetscErrorCode KSPMonitorSetFromOptions(KSP ksp,const char name[],const char help[], const char manual[],PetscErrorCode (*monitor)(KSP,PetscInt,PetscReal,PetscViewerAndFormat*))
314: {
315: PetscErrorCode ierr;
316: PetscBool flg;
317: PetscViewer viewer;
318: PetscViewerFormat format;
321: PetscOptionsGetViewer(PetscObjectComm((PetscObject)ksp),((PetscObject)ksp)->prefix,name,&viewer,&format,&flg);
322: if (flg) {
323: PetscViewerAndFormat *vf;
324: PetscViewerAndFormatCreate(viewer,format,&vf);
325: PetscObjectDereference((PetscObject)viewer);
326: KSPMonitorSet(ksp,(PetscErrorCode (*)(KSP,PetscInt,PetscReal,void*))monitor,vf,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
327: }
328: return(0);
329: }
333: /*@
334: KSPSetFromOptions - Sets KSP options from the options database.
335: This routine must be called before KSPSetUp() if the user is to be
336: allowed to set the Krylov type.
338: Collective on KSP
340: Input Parameters:
341: . ksp - the Krylov space context
343: Options Database Keys:
344: + -ksp_max_it - maximum number of linear iterations
345: . -ksp_rtol rtol - relative tolerance used in default determination of convergence, i.e.
346: if residual norm decreases by this factor than convergence is declared
347: . -ksp_atol abstol - absolute tolerance used in default convergence test, i.e. if residual
348: norm is less than this then convergence is declared
349: . -ksp_divtol tol - if residual norm increases by this factor than divergence is declared
350: . -ksp_converged_use_initial_residual_norm - see KSPConvergedDefaultSetUIRNorm()
351: . -ksp_converged_use_min_initial_residual_norm - see KSPConvergedDefaultSetUMIRNorm()
352: . -ksp_norm_type - none - skip norms used in convergence tests (useful only when not using
353: convergence test (say you always want to run with 5 iterations) to
354: save on communication overhead
355: preconditioned - default for left preconditioning
356: unpreconditioned - see KSPSetNormType()
357: natural - see KSPSetNormType()
358: . -ksp_check_norm_iteration it - do not compute residual norm until iteration number it (does compute at 0th iteration)
359: works only for PCBCGS, PCIBCGS and and PCCG
360: . -ksp_lag_norm - compute the norm of the residual for the ith iteration on the i+1 iteration; this means that one can use
361: the norm of the residual for convergence test WITHOUT an extra MPI_Allreduce() limiting global synchronizations.
362: This will require 1 more iteration of the solver than usual.
363: . -ksp_fischer_guess <model,size> - uses the Fischer initial guess generator for repeated linear solves
364: . -ksp_constant_null_space - assume the operator (matrix) has the constant vector in its null space
365: . -ksp_test_null_space - tests the null space set with MatSetNullSpace() to see if it truly is a null space
366: . -ksp_knoll - compute initial guess by applying the preconditioner to the right hand side
367: . -ksp_monitor_cancel - cancel all previous convergene monitor routines set
368: . -ksp_monitor <optional filename> - print residual norm at each iteration
369: . -ksp_monitor_lg_residualnorm - plot residual norm at each iteration
370: . -ksp_monitor_solution [ascii binary or draw][:filename][:format option] - plot solution at each iteration
371: - -ksp_monitor_singular_value - monitor extreme singular values at each iteration
373: Notes:
374: To see all options, run your program with the -help option
375: or consult Users-Manual: Chapter 4 KSP: Linear Equations Solvers
377: Level: beginner
379: .keywords: KSP, set, from, options, database
381: .seealso: KSPSetUseFischerGuess()
383: @*/
384: PetscErrorCode KSPSetFromOptions(KSP ksp)
385: {
387: PetscInt indx;
388: const char *convtests[] = {"default","skip"};
389: char type[256], monfilename[PETSC_MAX_PATH_LEN];
390: PetscBool flg,flag,reuse,set;
391: PetscInt model[2]={0,0},nmax;
392: KSPNormType normtype;
393: PCSide pcside;
394: void *ctx;
398: if (!ksp->skippcsetfromoptions) {
399: if (!ksp->pc) {KSPGetPC(ksp,&ksp->pc);}
400: PCSetFromOptions(ksp->pc);
401: }
403: KSPRegisterAll();
404: PetscObjectOptionsBegin((PetscObject)ksp);
405: PetscOptionsFList("-ksp_type","Krylov method","KSPSetType",KSPList,(char*)(((PetscObject)ksp)->type_name ? ((PetscObject)ksp)->type_name : KSPGMRES),type,256,&flg);
406: if (flg) {
407: KSPSetType(ksp,type);
408: }
409: /*
410: Set the type if it was never set.
411: */
412: if (!((PetscObject)ksp)->type_name) {
413: KSPSetType(ksp,KSPGMRES);
414: }
416: PetscObjectTypeCompare((PetscObject)ksp,KSPPREONLY,&flg);
417: if (flg) goto skipoptions;
419: PetscOptionsInt("-ksp_max_it","Maximum number of iterations","KSPSetTolerances",ksp->max_it,&ksp->max_it,NULL);
420: PetscOptionsReal("-ksp_rtol","Relative decrease in residual norm","KSPSetTolerances",ksp->rtol,&ksp->rtol,NULL);
421: PetscOptionsReal("-ksp_atol","Absolute value of residual norm","KSPSetTolerances",ksp->abstol,&ksp->abstol,NULL);
422: PetscOptionsReal("-ksp_divtol","Residual norm increase cause divergence","KSPSetTolerances",ksp->divtol,&ksp->divtol,NULL);
424: PetscOptionsBool("-ksp_converged_use_initial_residual_norm","Use initial residual residual norm for computing relative convergence","KSPConvergedDefaultSetUIRNorm",PETSC_FALSE,&flag,&set);
425: if (set && flag) {KSPConvergedDefaultSetUIRNorm(ksp);}
426: PetscOptionsBool("-ksp_converged_use_min_initial_residual_norm","Use minimum of initial residual norm and b for computing relative convergence","KSPConvergedDefaultSetUMIRNorm",PETSC_FALSE,&flag,&set);
427: if (set && flag) {KSPConvergedDefaultSetUMIRNorm(ksp);}
428: PetscOptionsBool("-ksp_initial_guess_nonzero","Use the contents of the solution vector for initial guess","KSPSetInitialNonzero",ksp->guess_zero ? PETSC_FALSE : PETSC_TRUE,&flag,&flg);
429: if (flg) {
430: KSPSetInitialGuessNonzero(ksp,flag);
431: }
432: PCGetReusePreconditioner(ksp->pc,&reuse);
433: PetscOptionsBool("-ksp_reuse_preconditioner","Use initial preconditioner and don't ever compute a new one ","KSPReusePreconditioner",reuse,&reuse,NULL);
434: KSPSetReusePreconditioner(ksp,reuse);
436: PetscOptionsBool("-ksp_knoll","Use preconditioner applied to b for initial guess","KSPSetInitialGuessKnoll",ksp->guess_knoll,&ksp->guess_knoll,NULL);
437: PetscOptionsBool("-ksp_error_if_not_converged","Generate error if solver does not converge","KSPSetErrorIfNotConverged",ksp->errorifnotconverged,&ksp->errorifnotconverged,NULL);
438: nmax = 2;
439: PetscOptionsIntArray("-ksp_fischer_guess","Use Paul Fischer's algorithm for initial guess","KSPSetUseFischerGuess",model,&nmax,&flag);
440: if (flag) {
441: if (nmax != 2) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE,"Must pass in model,size as arguments");
442: KSPSetUseFischerGuess(ksp,model[0],model[1]);
443: }
445: PetscOptionsEList("-ksp_convergence_test","Convergence test","KSPSetConvergenceTest",convtests,2,"default",&indx,&flg);
446: if (flg) {
447: switch (indx) {
448: case 0:
449: KSPConvergedDefaultCreate(&ctx);
450: KSPSetConvergenceTest(ksp,KSPConvergedDefault,ctx,KSPConvergedDefaultDestroy);
451: break;
452: case 1: KSPSetConvergenceTest(ksp,KSPConvergedSkip,NULL,NULL); break;
453: }
454: }
456: KSPSetUpNorms_Private(ksp,&normtype,&pcside);
457: PetscOptionsEnum("-ksp_norm_type","KSP Norm type","KSPSetNormType",KSPNormTypes,(PetscEnum)normtype,(PetscEnum*)&normtype,&flg);
458: if (flg) { KSPSetNormType(ksp,normtype); }
460: PetscOptionsInt("-ksp_check_norm_iteration","First iteration to compute residual norm","KSPSetCheckNormIteration",ksp->chknorm,&ksp->chknorm,NULL);
462: PetscOptionsBool("-ksp_lag_norm","Lag the calculation of the residual norm","KSPSetLagNorm",ksp->lagnorm,&flag,&flg);
463: if (flg) {
464: KSPSetLagNorm(ksp,flag);
465: }
467: KSPGetDiagonalScale(ksp,&flag);
468: PetscOptionsBool("-ksp_diagonal_scale","Diagonal scale matrix before building preconditioner","KSPSetDiagonalScale",flag,&flag,&flg);
469: if (flg) {
470: KSPSetDiagonalScale(ksp,flag);
471: }
472: KSPGetDiagonalScaleFix(ksp,&flag);
473: PetscOptionsBool("-ksp_diagonal_scale_fix","Fix diagonally scaled matrix after solve","KSPSetDiagonalScaleFix",flag,&flag,&flg);
474: if (flg) {
475: KSPSetDiagonalScaleFix(ksp,flag);
476: }
478: PetscOptionsBool("-ksp_constant_null_space","Add constant null space to Krylov solver matrix","MatSetNullSpace",PETSC_FALSE,&flg,&set);
479: if (set && flg) {
480: MatNullSpace nsp;
481: Mat Amat;
483: MatNullSpaceCreate(PetscObjectComm((PetscObject)ksp),PETSC_TRUE,0,NULL,&nsp);
484: PCGetOperators(ksp->pc,&Amat,NULL);
485: if (Amat) {
486: MatSetNullSpace(Amat,nsp);
487: MatNullSpaceDestroy(&nsp);
488: } else SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_WRONGSTATE,"Cannot set nullspace, matrix has not yet been provided");
489: }
491: PetscOptionsBool("-ksp_monitor_cancel","Remove any hardwired monitor routines","KSPMonitorCancel",PETSC_FALSE,&flg,&set);
492: /* -----------------------------------------------------------------------*/
493: /*
494: Cancels all monitors hardwired into code before call to KSPSetFromOptions()
495: */
496: if (set && flg) {
497: KSPMonitorCancel(ksp);
498: }
499: KSPMonitorSetFromOptions(ksp,"-ksp_monitor","Monitor the (preconditioned) residual norm","KSPMonitorDefault",KSPMonitorDefault);
500: KSPMonitorSetFromOptions(ksp,"-ksp_monitor_range","Monitor the percentage of large entries in the residual","KSPMonitorRange",KSPMonitorRange);
501: KSPMonitorSetFromOptions(ksp,"-ksp_monitor_true_residual","Monitor the unprecondiitoned residual norm","KSPMOnitorTrueResidual",KSPMonitorTrueResidualNorm);
502: KSPMonitorSetFromOptions(ksp,"-ksp_monitor_max","Monitor the maximum norm of the residual","KSPMonitorTrueResidualMaxNorm",KSPMonitorTrueResidualMaxNorm);
503: KSPMonitorSetFromOptions(ksp,"-ksp_monitor_short","Monitor preconditioned residual norm with fewer digits","KSPMonitorDefaultShort",KSPMonitorDefaultShort);
504: KSPMonitorSetFromOptions(ksp,"-ksp_monitor_solution","Monitor the solution","KSPMonitorSolution",KSPMonitorSolution);
505: KSPMonitorSetFromOptions(ksp,"-ksp_monitor_singular_value","Monitor singular values","KSPMonitorSingularValue",KSPMonitorSingularValue);
506: PetscOptionsHasName(NULL,((PetscObject)ksp)->prefix,"-ksp_monitor_singular_value",&flg);
507: if (flg) {
508: KSPSetComputeSingularValues(ksp,PETSC_TRUE);
509: }
510: PetscObjectTypeCompare((PetscObject)ksp->pc,PCKSP,&flg);
511: PetscObjectTypeCompare((PetscObject)ksp->pc,PCBJACOBI,&flag);
513: if (flg || flag) {
514: /* A hack for using dynamic tolerance in preconditioner */
515: PetscOptionsString("-sub_ksp_dynamic_tolerance","Use dynamic tolerance for PC if PC is a KSP","KSPMonitorDynamicTolerance","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
516: if (flg) {
517: KSPDynTolCtx *scale;
518: PetscMalloc1(1,&scale);
519: scale->bnrm = -1.0;
520: scale->coef = 1.0;
521: PetscOptionsReal("-sub_ksp_dynamic_tolerance_param","Parameter of dynamic tolerance for inner PCKSP","KSPMonitorDynamicToleranceParam",scale->coef,&scale->coef,&flg);
522: KSPMonitorSet(ksp,KSPMonitorDynamicTolerance,scale,KSPMonitorDynamicToleranceDestroy);
523: }
524: }
525:
527: /*
528: Calls Python function
529: */
530: PetscOptionsString("-ksp_monitor_python","Use Python function","KSPMonitorSet",0,monfilename,PETSC_MAX_PATH_LEN,&flg);
531: if (flg) {PetscPythonMonitorSet((PetscObject)ksp,monfilename);}
532: /*
533: Graphically plots preconditioned residual norm
534: */
535: PetscOptionsBool("-ksp_monitor_lg_residualnorm","Monitor graphically preconditioned residual norm","KSPMonitorSet",PETSC_FALSE,&flg,&set);
536: if (set && flg) {
537: PetscDrawLG ctx;
539: KSPMonitorLGResidualNormCreate(PetscObjectComm((PetscObject)ksp),NULL,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&ctx);
540: KSPMonitorSet(ksp,KSPMonitorLGResidualNorm,ctx,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
541: }
542: /*
543: Graphically plots preconditioned and true residual norm
544: */
545: PetscOptionsBool("-ksp_monitor_lg_true_residualnorm","Monitor graphically true residual norm","KSPMonitorSet",PETSC_FALSE,&flg,&set);
546: if (set && flg) {
547: PetscDrawLG ctx;
549: KSPMonitorLGTrueResidualNormCreate(PetscObjectComm((PetscObject)ksp),NULL,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&ctx);
550: KSPMonitorSet(ksp,KSPMonitorLGTrueResidualNorm,ctx,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
551: }
552: /*
553: Graphically plots preconditioned residual norm and range of residual element values
554: */
555: PetscOptionsBool("-ksp_monitor_lg_range","Monitor graphically range of preconditioned residual norm","KSPMonitorSet",PETSC_FALSE,&flg,&set);
556: if (set && flg) {
557: PetscViewer ctx;
559: PetscViewerDrawOpen(PetscObjectComm((PetscObject)ksp),NULL,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&ctx);
560: KSPMonitorSet(ksp,KSPMonitorLGRange,ctx,(PetscErrorCode (*)(void**))PetscViewerDestroy);
561: }
563: #if defined(PETSC_HAVE_SAWS)
564: /*
565: Publish convergence information using AMS
566: */
567: PetscOptionsBool("-ksp_monitor_saws","Publish KSP progress using SAWs","KSPMonitorSet",PETSC_FALSE,&flg,&set);
568: if (set && flg) {
569: void *ctx;
570: KSPMonitorSAWsCreate(ksp,&ctx);
571: KSPMonitorSet(ksp,KSPMonitorSAWs,ctx,KSPMonitorSAWsDestroy);
572: KSPSetComputeSingularValues(ksp,PETSC_TRUE);
573: }
574: #endif
576: /* -----------------------------------------------------------------------*/
577: KSPSetUpNorms_Private(ksp,&normtype,&pcside);
578: PetscOptionsEnum("-ksp_pc_side","KSP preconditioner side","KSPSetPCSide",PCSides,(PetscEnum)pcside,(PetscEnum*)&pcside,&flg);
579: if (flg) {KSPSetPCSide(ksp,pcside);}
581: PetscOptionsBool("-ksp_compute_singularvalues","Compute singular values of preconditioned operator","KSPSetComputeSingularValues",ksp->calc_sings,&flg,&set);
582: if (set) { KSPSetComputeSingularValues(ksp,flg); }
583: PetscOptionsBool("-ksp_compute_eigenvalues","Compute eigenvalues of preconditioned operator","KSPSetComputeSingularValues",ksp->calc_sings,&flg,&set);
584: if (set) { KSPSetComputeSingularValues(ksp,flg); }
585: PetscOptionsBool("-ksp_plot_eigenvalues","Scatter plot extreme eigenvalues","KSPSetComputeSingularValues",PETSC_FALSE,&flg,&set);
586: if (set) { KSPSetComputeSingularValues(ksp,flg); }
588: #if defined(PETSC_HAVE_SAWS)
589: {
590: PetscBool set;
591: flg = PETSC_FALSE;
592: PetscOptionsBool("-ksp_saws_block","Block for SAWs at end of KSPSolve","PetscObjectSAWsBlock",((PetscObject)ksp)->amspublishblock,&flg,&set);
593: if (set) {
594: PetscObjectSAWsSetBlock((PetscObject)ksp,flg);
595: }
596: }
597: #endif
599: if (ksp->ops->setfromoptions) {
600: (*ksp->ops->setfromoptions)(PetscOptionsObject,ksp);
601: }
602: skipoptions:
603: /* process any options handlers added with PetscObjectAddOptionsHandler() */
604: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)ksp);
605: PetscOptionsEnd();
606: return(0);
607: }