Actual source code: itcreate.c
1: /*$Id: itcreate.c,v 1.204 2001/04/10 19:36:24 bsmith Exp bsmith $*/
2: /*
3: The basic KSP routines, Create, View etc. are here.
4: */
5: #include src/sles/ksp/kspimpl.h
6: #include petscsys.h
8: PetscTruth KSPRegisterAllCalled = PETSC_FALSE;
10: /*@C
11: KSPView - Prints the KSP data structure.
13: Collective on KSP
15: Input Parameters:
16: + ksp - the Krylov space context
17: - viewer - visualization context
19: Note:
20: The available visualization contexts include
21: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
22: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
23: output where only the first processor opens
24: the file. All other processors send their
25: data to the first processor to print.
27: The user can open an alternative visualization context with
28: PetscViewerASCIIOpen() - output to a specified file.
30: Level: developer
32: .keywords: KSP, view
34: .seealso: PCView(), PetscViewerASCIIOpen()
35: @*/
36: int KSPView(KSP ksp,PetscViewer viewer)
37: {
38: char *type;
39: int ierr;
40: PetscTruth isascii;
44: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(ksp->comm);
48: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
49: if (isascii) {
50: KSPGetType(ksp,&type);
51: PetscViewerASCIIPrintf(viewer,"KSP Object:n");
52: if (type) {
53: PetscViewerASCIIPrintf(viewer," type: %sn",type);
54: } else {
55: PetscViewerASCIIPrintf(viewer," type: not yet setn");
56: }
57: if (ksp->ops->view) {
58: PetscViewerASCIIPushTab(viewer);
59: (*ksp->ops->view)(ksp,viewer);
60: PetscViewerASCIIPopTab(viewer);
61: }
62: if (ksp->guess_zero) {PetscViewerASCIIPrintf(viewer," maximum iterations=%d, initial guess is zeron",ksp->max_it);}
63: else {PetscViewerASCIIPrintf(viewer," maximum iterations=%dn", ksp->max_it);}
64: PetscViewerASCIIPrintf(viewer," tolerances: relative=%g, absolute=%g, divergence=%gn",ksp->rtol,ksp->atol,ksp->divtol);
65: if (ksp->pc_side == PC_RIGHT) {PetscViewerASCIIPrintf(viewer," right preconditioningn");}
66: else if (ksp->pc_side == PC_SYMMETRIC) {PetscViewerASCIIPrintf(viewer," symmetric preconditioningn");}
67: else {PetscViewerASCIIPrintf(viewer," left preconditioningn");}
68: } else {
69: if (ksp->ops->view) {
70: (*ksp->ops->view)(ksp,viewer);
71: }
72: }
73: return(0);
74: }
76: /*
77: Contains the list of registered KSP routines
78: */
79: PetscFList KSPList = 0;
81: /*@C
82: KSPSetAvoidNorms - Sets the KSP solver to avoid computing the residual norm
83: when possible. This, for example, reduces the number of collective operations
84: when using the Krylov method as a smoother.
86: Collective on KSP
88: Input Parameter:
89: . ksp - Krylov solver context
91: Notes:
92: One cannot use the default convergence test routines when this option is
93: set, since these are based on decreases in the residual norms. Thus, this
94: option automatically switches to activate the KSPSkipConverged() test function.
96: Currently only works with the CG, Richardson, Bi-CG-stab, CR, and CGS methods.
98: Level: advanced
100: .keywords: KSP, create, context, norms
102: .seealso: KSPSetUp(), KSPSolve(), KSPDestroy(), KSPSkipConverged()
103: @*/
104: int KSPSetAvoidNorms(KSP ksp)
105: {
110: ksp->avoidnorms = PETSC_TRUE;
111: KSPSetConvergenceTest(ksp,KSPSkipConverged,PETSC_NULL);
112: return(0);
113: }
115: static int KSPPublish_Petsc(PetscObject obj)
116: {
117: #if defined(PETSC_HAVE_AMS)
118: KSP v = (KSP) obj;
119: int ierr;
120: #endif
124: #if defined(PETSC_HAVE_AMS)
125: /* if it is already published then return */
126: if (v->amem >=0) return(0);
128: PetscObjectPublishBaseBegin(obj);
129: AMS_Memory_add_field((AMS_Memory)v->amem,"Iteration",&v->its,1,AMS_INT,AMS_READ,
130: AMS_COMMON,AMS_REDUCT_UNDEF);
131: AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->rnorm,1,AMS_DOUBLE,AMS_READ,
132: AMS_COMMON,AMS_REDUCT_UNDEF);
134: if (v->res_hist_max > 0) {
135: AMS_Memory_add_field((AMS_Memory)v->amem,"ResidualNormsCount",&v->res_hist_len,1,AMS_INT,AMS_READ,
136: AMS_COMMON,AMS_REDUCT_UNDEF);
137: AMS_Memory_add_field((AMS_Memory)v->amem,"ResidualNormsCountMax",&v->res_hist_max,1,AMS_INT,AMS_READ,
138: AMS_COMMON,AMS_REDUCT_UNDEF);
139: AMS_Memory_add_field((AMS_Memory)v->amem,"ResidualNorms",v->res_hist,v->res_hist_max,AMS_DOUBLE,AMS_READ,
140: AMS_COMMON,AMS_REDUCT_UNDEF);
141: }
143: PetscObjectPublishBaseEnd(obj);
144: #endif
146: return(0);
147: }
150: /*@C
151: KSPCreate - Creates the default KSP context.
153: Collective on MPI_Comm
155: Input Parameter:
156: . comm - MPI communicator
158: Output Parameter:
159: . ksp - location to put the KSP context
161: Notes:
162: The default KSP type is GMRES with a restart of 30, using modified Gram-Schmidt
163: orthogonalization.
165: Level: developer
167: .keywords: KSP, create, context
169: .seealso: KSPSetUp(), KSPSolve(), KSPDestroy(), KSP
170: @*/
171: int KSPCreate(MPI_Comm comm,KSP *inksp)
172: {
173: KSP ksp;
177: *inksp = 0;
178: PetscHeaderCreate(ksp,_p_KSP,struct _KSPOps,KSP_COOKIE,-1,"KSP",comm,KSPDestroy,KSPView);
179: PetscLogObjectCreate(ksp);
180: *inksp = ksp;
181: ksp->bops->publish = KSPPublish_Petsc;
183: ksp->type = -1;
184: ksp->max_it = 10000;
185: ksp->pc_side = PC_LEFT;
186: ksp->use_pres = PETSC_FALSE;
187: ksp->rtol = 1.e-5;
188: ksp->atol = 1.e-50;
189: ksp->divtol = 1.e4;
190: ksp->avoidnorms = PETSC_FALSE;
192: ksp->rnorm = 0.0;
193: ksp->its = 0;
194: ksp->guess_zero = PETSC_TRUE;
195: ksp->calc_sings = PETSC_FALSE;
196: ksp->calc_res = PETSC_FALSE;
197: ksp->res_hist = PETSC_NULL;
198: ksp->res_hist_len = 0;
199: ksp->res_hist_max = 0;
200: ksp->res_hist_reset = PETSC_TRUE;
201: ksp->numbermonitors = 0;
202: ksp->converged = KSPDefaultConverged;
203: ksp->ops->buildsolution = KSPDefaultBuildSolution;
204: ksp->ops->buildresidual = KSPDefaultBuildResidual;
206: ksp->ops->setfromoptions = 0;
208: ksp->vec_sol = 0;
209: ksp->vec_rhs = 0;
210: ksp->B = 0;
212: ksp->ops->solve = 0;
213: ksp->ops->setup = 0;
214: ksp->ops->destroy = 0;
216: ksp->data = 0;
217: ksp->nwork = 0;
218: ksp->work = 0;
220: ksp->cnvP = 0;
222: ksp->reason = KSP_CONVERGED_ITERATING;
224: ksp->setupcalled = 0;
225: PetscPublishAll(ksp);
226: return(0);
227: }
228:
229: /*@C
230: KSPSetType - Builds KSP for a particular solver.
232: Collective on KSP
234: Input Parameters:
235: + ksp - the Krylov space context
236: - type - a known method
238: Options Database Key:
239: . -ksp_type <method> - Sets the method; use -help for a list
240: of available methods (for instance, cg or gmres)
242: Notes:
243: See "petsc/include/petscksp.h" for available methods (for instance,
244: KSPCG or KSPGMRES).
246: Normally, it is best to use the SLESSetFromOptions() command and
247: then set the KSP type from the options database rather than by using
248: this routine. Using the options database provides the user with
249: maximum flexibility in evaluating the many different Krylov methods.
250: The KSPSetType() routine is provided for those situations where it
251: is necessary to set the iterative solver independently of the command
252: line or options database. This might be the case, for example, when
253: the choice of iterative solver changes during the execution of the
254: program, and the user's application is taking responsibility for
255: choosing the appropriate method. In other words, this routine is
256: not for beginners.
258: Level: intermediate
260: .keywords: KSP, set, method
262: .seealso: PCSetType(), KSPType
264: @*/
265: int KSPSetType(KSP ksp,KSPType type)
266: {
267: int ierr,(*r)(KSP);
268: PetscTruth match;
274: PetscTypeCompare((PetscObject)ksp,type,&match);
275: if (match) return(0);
277: if (ksp->data) {
278: /* destroy the old private KSP context */
279: (*ksp->ops->destroy)(ksp);
280: ksp->data = 0;
281: }
282: /* Get the function pointers for the iterative method requested */
283: if (!KSPRegisterAllCalled) {KSPRegisterAll(PETSC_NULL);}
285: PetscFListFind(ksp->comm,KSPList,type,(void (**)(void)) &r);
287: if (!r) SETERRQ1(1,"Unknown KSP type given: %s",type);
289: ksp->setupcalled = 0;
290: (*r)(ksp);
292: PetscObjectChangeTypeName((PetscObject)ksp,type);
293: return(0);
294: }
296: /*@C
297: KSPRegisterDestroy - Frees the list of KSP methods that were
298: registered by KSPRegisterDynamic().
300: Not Collective
302: Level: advanced
304: .keywords: KSP, register, destroy
306: .seealso: KSPRegisterDynamic(), KSPRegisterAll()
307: @*/
308: int KSPRegisterDestroy(void)
309: {
313: if (KSPList) {
314: PetscFListDestroy(&KSPList);
315: KSPList = 0;
316: }
317: KSPRegisterAllCalled = PETSC_FALSE;
318: return(0);
319: }
321: /*@C
322: KSPGetType - Gets the KSP type as a string from the KSP object.
324: Not Collective
326: Input Parameter:
327: . ksp - Krylov context
329: Output Parameter:
330: . name - name of KSP method
332: Level: intermediate
334: .keywords: KSP, get, method, name
336: .seealso: KSPSetType()
337: @*/
338: int KSPGetType(KSP ksp,KSPType *type)
339: {
342: *type = ksp->type_name;
343: return(0);
344: }
346: /*@
347: KSPSetFromOptions - Sets KSP options from the options database.
348: This routine must be called before KSPSetUp() if the user is to be
349: allowed to set the Krylov type.
351: Collective on KSP
353: Input Parameters:
354: . ksp - the Krylov space context
356: Options Database Keys:
357: + -ksp_max_it - maximum number of linear iterations
358: . -ksp_rtol rtol - relative tolerance used in default determination of convergence, i.e.
359: if residual norm decreases by this factor than convergence is declared
360: . -ksp_atol atol - absolute tolerance used in default convergence test, i.e. if residual
361: norm is less than this then convergence is declared
362: . -ksp_divtol tol - if residual norm increases by this factor than divergence is declared
363: . -ksp_avoid_norms - skip norms used in convergence tests (useful only when not using
364: convergence test (say you always want to run with 5 iterations) to
365: save on communication overhead
366: . -ksp_cancelmonitors - cancel all previous convergene monitor routines set
367: . -ksp_monitor - print residual norm at each iteration
368: . -ksp_xmonitor - plot residual norm at each iteration
369: . -ksp_vecmonitor - plot solution at each iteration
370: - -ksp_singmonitor - monitor extremem singular values at each iteration
372: Notes:
373: To see all options, run your program with the -help option
374: or consult the users manual.
376: Level: developer
378: .keywords: KSP, set, from, options, database
380: .seealso:
381: @*/
382: int KSPSetFromOptions(KSP ksp)
383: {
384: int ierr;
385: char type[256];
386: PetscTruth flg;
390: if (!KSPRegisterAllCalled) {KSPRegisterAll(PETSC_NULL);}
391: PetscOptionsBegin(ksp->comm,ksp->prefix,"Krylov Method (KSP) Options","KSP");
392: PetscOptionsList("-ksp_type","Krylov method","KSPSetType",KSPList,(char*)(ksp->type_name?ksp->type_name:KSPGMRES),type,256,&flg);
393: if (flg) {
394: KSPSetType(ksp,type);
395: }
396: /*
397: Set the type if it was never set.
398: */
399: if (!ksp->type_name) {
400: KSPSetType(ksp,KSPGMRES);
401: }
403: PetscOptionsInt("-ksp_max_it","Maximum number of iterations","KSPSetTolerances",ksp->max_it,&ksp->max_it,PETSC_NULL);
404: PetscOptionsDouble("-ksp_rtol","Relative decrease in residual norm","KSPSetTolerances",ksp->rtol,&ksp->rtol,PETSC_NULL);
405: PetscOptionsDouble("-ksp_atol","Absolute value of residual norm","KSPSetTolerances",ksp->atol,&ksp->atol,PETSC_NULL);
406: PetscOptionsDouble("-ksp_divtol","Residual norm increase cause divergence","KSPSetTolerances",ksp->divtol,&ksp->divtol,PETSC_NULL);
408: PetscOptionsName("-ksp_avoid_norms","Do not compute norms for convergence tests","KSPSetAvoidNorms",&flg);
409: if (flg) {
410: KSPSetAvoidNorms(ksp);
411: }
413: PetscOptionsName("-ksp_cancelmonitors","Remove any hardwired monitor routines","KSPClearMonitor",&flg);
414: /* -----------------------------------------------------------------------*/
415: /*
416: Cancels all monitors hardwired into code before call to KSPSetFromOptions()
417: */
418: if (flg) {
419: KSPClearMonitor(ksp);
420: }
421: /*
422: Prints preconditioned residual norm at each iteration
423: */
424: PetscOptionsName("-ksp_monitor","Monitor preconditioned residual norm","KSPSetMonitor",&flg);
425: if (flg) {
426: KSPSetMonitor(ksp,KSPDefaultMonitor,PETSC_NULL,PETSC_NULL);
427: }
428: /*
429: Plots the vector solution
430: */
431: PetscOptionsName("-ksp_vecmonitor","Monitor solution graphically","KSPSetMonitor",&flg);
432: if (flg) {
433: KSPSetMonitor(ksp,KSPVecViewMonitor,PETSC_NULL,PETSC_NULL);
434: }
435: /*
436: Prints preconditioned and true residual norm at each iteration
437: */
438: PetscOptionsName("-ksp_truemonitor","Monitor true (unpreconditioned) residual norm","KSPSetMonitor",&flg);
439: if (flg) {
440: KSPSetMonitor(ksp,KSPTrueMonitor,PETSC_NULL,PETSC_NULL);
441: }
442: /*
443: Prints extreme eigenvalue estimates at each iteration
444: */
445: PetscOptionsName("-ksp_singmonitor","Monitor singular values","KSPSetMonitor",&flg);
446: if (flg) {
447: KSPSetComputeSingularValues(ksp);
448: KSPSetMonitor(ksp,KSPSingularValueMonitor,PETSC_NULL,PETSC_NULL);
449: }
450: /*
451: Prints preconditioned residual norm with fewer digits
452: */
453: PetscOptionsName("-ksp_smonitor","Monitor preconditioned residual norm with fewer digitis","KSPSetMonitor",&flg);
454: if (flg) {
455: KSPSetMonitor(ksp,KSPDefaultSMonitor,PETSC_NULL,PETSC_NULL);
456: }
457: /*
458: Graphically plots preconditioned residual norm
459: */
460: PetscOptionsName("-ksp_xmonitor","Monitor graphically preconditioned residual norm","KSPSetMonitor",&flg);
461: if (flg) {
462: KSPSetMonitor(ksp,KSPLGMonitor,PETSC_NULL,PETSC_NULL);
463: }
464: /*
465: Graphically plots preconditioned and true residual norm
466: */
467: PetscOptionsName("-ksp_xtruemonitor","Monitor graphically true residual norm","KSPSetMonitor",&flg);
468: if (flg){
469: KSPSetMonitor(ksp,KSPLGTrueMonitor,PETSC_NULL,PETSC_NULL);
470: }
472: /* -----------------------------------------------------------------------*/
473: PetscOptionsName("-ksp_preres","Use preconditioner residual norm in convergence tests","KSPSetUsePreconditionedResidual",&flg);
474: if (flg) { KSPSetUsePreconditionedResidual(ksp);}
476: PetscOptionsLogicalGroupBegin("-ksp_left_pc","Use left preconditioning","KSPSetPreconditionerSide",&flg);
477: if (flg) { KSPSetPreconditionerSide(ksp,PC_LEFT); }
478: PetscOptionsLogicalGroup("-ksp_right_pc","Use right preconditioning","KSPSetPreconditionerSide",&flg);
479: if (flg) { KSPSetPreconditionerSide(ksp,PC_RIGHT);}
480: PetscOptionsLogicalGroupEnd("-ksp_symmetric_pc","Use symmetric (factorized) preconditioning","KSPSetPreconditionerSide",&flg);
481: if (flg) { KSPSetPreconditionerSide(ksp,PC_SYMMETRIC);}
483: PetscOptionsName("-ksp_compute_singularvalues","Compute singular values of preconditioned operator","KSPSetComputeSingularValues",&flg);
484: if (flg) { KSPSetComputeSingularValues(ksp); }
485: PetscOptionsName("-ksp_compute_eigenvalues","Compute eigenvalues of preconditioned operator","KSPSetComputeSingularValues",&flg);
486: if (flg) { KSPSetComputeSingularValues(ksp); }
487: PetscOptionsName("-ksp_plot_eigenvalues","Scatter plot extreme eigenvalues","KSPSetComputeSingularValues",&flg);
488: if (flg) { KSPSetComputeSingularValues(ksp); }
490: if (ksp->ops->setfromoptions) {
491: (*ksp->ops->setfromoptions)(ksp);
492: }
493: PetscOptionsEnd();
496: return(0);
497: }
499: /*MC
500: KSPRegisterDynamic - Adds a method to the Krylov subspace solver package.
502: Synopsis:
503: int KSPRegisterDynamic(char *name_solver,char *path,char *name_create,int (*routine_create)(KSP))
505: Not Collective
507: Input Parameters:
508: + name_solver - name of a new user-defined solver
509: . path - path (either absolute or relative) the library containing this solver
510: . name_create - name of routine to create method context
511: - routine_create - routine to create method context
513: Notes:
514: KSPRegisterDynamic() may be called multiple times to add several user-defined solvers.
516: If dynamic libraries are used, then the fourth input argument (routine_create)
517: is ignored.
519: Sample usage:
520: .vb
521: KSPRegisterDynamic("my_solver",/home/username/my_lib/lib/libO/solaris/mylib.a,
522: "MySolverCreate",MySolverCreate);
523: .ve
525: Then, your solver can be chosen with the procedural interface via
526: $ KSPSetType(ksp,"my_solver")
527: or at runtime via the option
528: $ -ksp_type my_solver
530: Level: advanced
532: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, ${BOPT},
533: and others of the form ${any_environmental_variable} occuring in pathname will be
534: replaced with appropriate values.
536: .keywords: KSP, register
538: .seealso: KSPRegisterAll(), KSPRegisterDestroy()
540: M*/
542: int KSPRegister(char *sname,char *path,char *name,int (*function)(KSP))
543: {
544: int ierr;
545: char fullname[256];
548: PetscFListConcat(path,name,fullname);
549: PetscFListAdd(&KSPList,sname,fullname,(void (*)())function);
550: return(0);
551: }