Actual source code: itfunc.c
1: /*$Id: itfunc.c,v 1.159 2001/08/07 03:03:45 balay Exp $*/
2: /*
3: Interface KSP routines that the user calls.
4: */
6: #include src/sles/ksp/kspimpl.h
10: /*@
11: KSPComputeExtremeSingularValues - Computes the extreme singular values
12: for the preconditioned operator. Called after or during KSPSolve()
13: (SLESSolve()).
15: Not Collective
17: Input Parameter:
18: . ksp - iterative context obtained from KSPCreate()
20: Output Parameters:
21: . emin, emax - extreme singular values
23: Notes:
24: One must call KSPSetComputeSingularValues() before calling KSPSetUp()
25: (or use the option -ksp_compute_eigenvalues) in order for this routine to work correctly.
27: Many users may just want to use the monitoring routine
28: KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
29: to print the extreme singular values at each iteration of the linear solve.
31: Level: advanced
33: .keywords: KSP, compute, extreme, singular, values
35: .seealso: KSPSetComputeSingularValues(), KSPSingularValueMonitor(), KSPComputeEigenvalues()
36: @*/
37: int KSPComputeExtremeSingularValues(KSP ksp,PetscReal *emax,PetscReal *emin)
38: {
45: if (!ksp->calc_sings) {
46: SETERRQ(4,"Singular values not requested before KSPSetUp()");
47: }
49: if (ksp->ops->computeextremesingularvalues) {
50: (*ksp->ops->computeextremesingularvalues)(ksp,emax,emin);
51: } else {
52: *emin = -1.0;
53: *emax = -1.0;
54: }
55: return(0);
56: }
60: /*@
61: KSPComputeEigenvalues - Computes the extreme eigenvalues for the
62: preconditioned operator. Called after or during KSPSolve() (SLESSolve()).
64: Not Collective
66: Input Parameter:
67: + ksp - iterative context obtained from KSPCreate()
68: - n - size of arrays r and c. The number of eigenvalues computed (neig) will, in
69: general, be less than this.
71: Output Parameters:
72: + r - real part of computed eigenvalues
73: . c - complex part of computed eigenvalues
74: - neig - number of eigenvalues computed (will be less than or equal to n)
76: Options Database Keys:
77: + -ksp_compute_eigenvalues - Prints eigenvalues to stdout
78: - -ksp_plot_eigenvalues - Plots eigenvalues in an x-window display
80: Notes:
81: The number of eigenvalues estimated depends on the size of the Krylov space
82: generated during the KSPSolve() (that is the SLESSolve); for example, with
83: CG it corresponds to the number of CG iterations, for GMRES it is the number
84: of GMRES iterations SINCE the last restart. Any extra space in r[] and c[]
85: will be ignored.
87: KSPComputeEigenvalues() does not usually provide accurate estimates; it is
88: intended only for assistance in understanding the convergence of iterative
89: methods, not for eigenanalysis.
91: One must call KSPSetComputeEigenvalues() before calling KSPSetUp()
92: in order for this routine to work correctly.
94: Many users may just want to use the monitoring routine
95: KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
96: to print the singular values at each iteration of the linear solve.
98: Level: advanced
100: .keywords: KSP, compute, extreme, singular, values
102: .seealso: KSPSetComputeSingularValues(), KSPSingularValueMonitor(), KSPComputeExtremeSingularValues()
103: @*/
104: int KSPComputeEigenvalues(KSP ksp,int n,PetscReal *r,PetscReal *c,int *neig)
105: {
112: if (!ksp->calc_sings) {
113: SETERRQ(4,"Eigenvalues not requested before KSPSetUp()");
114: }
116: if (ksp->ops->computeeigenvalues) {
117: (*ksp->ops->computeeigenvalues)(ksp,n,r,c,neig);
118: } else {
119: *neig = 0;
120: }
121: return(0);
122: }
126: /*@
127: KSPSetUp - Sets up the internal data structures for the
128: later use of an iterative solver.
130: Collective on KSP
132: Input Parameter:
133: . ksp - iterative context obtained from KSPCreate()
135: Level: developer
137: .keywords: KSP, setup
139: .seealso: KSPCreate(), KSPSolve(), KSPDestroy()
140: @*/
141: int KSPSetUp(KSP ksp)
142: {
148: /* reset the convergence flag from the previous solves */
149: ksp->reason = KSP_CONVERGED_ITERATING;
151: if (!ksp->type_name){
152: KSPSetType(ksp,KSPGMRES);
153: }
155: if (ksp->setupcalled) return(0);
156: ksp->setupcalled = 1;
157: (*ksp->ops->setup)(ksp);
158: return(0);
159: }
164: /*@
165: KSPSolve - Solves linear system; usually not called directly, rather
166: it is called by a call to SLESSolve().
168: Collective on KSP
170: Input Parameter:
171: . ksp - iterative context obtained from KSPCreate()
173: Output Parameter:
174: . its - number of iterations required
176: Options Database Keys:
177: + -ksp_compute_eigenvalues - compute preconditioned operators eigenvalues
178: . -ksp_plot_eigenvalues - plot the computed eigenvalues in an X-window
179: . -ksp_compute_eigenvalues_explicitly - compute the eigenvalues by forming the
180: dense operator and useing LAPACK
181: - -ksp_plot_eigenvalues_explicitly - plot the explicitly computing eigenvalues
183: Notes:
184: On return, the parameter "its" contains either the iteration
185: number at which convergence was successfully reached or failure was detected.
187: Call KSPGetConvergedReason() to determine if the solver converged or failed and
188: why.
189:
190: If using a direct method (e.g., via the KSP solver
191: KSPPREONLY and a preconditioner such as PCLU/PCILU),
192: then its=1. See KSPSetTolerances() and KSPDefaultConverged()
193: for more details.
195: Understanding Convergence:
196: The routines KSPSetMonitor(), KSPComputeEigenvalues(), and
197: KSPComputeEigenvaluesExplicitly() provide information on additional
198: options to monitor convergence and print eigenvalue information.
200: Level: developer
202: .keywords: KSP, solve, linear system
204: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPDefaultConverged(),
205: SLESSolve(), KSPSolveTranspose(), SLESGetKSP()
206: @*/
207: int KSPSolve(KSP ksp,int *its)
208: {
209: int ierr,rank,nits;
210: PetscTruth flag1,flag2;
211: PetscScalar zero = 0.0;
216: if (!ksp->setupcalled){ KSPSetUp(ksp);}
217: if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
218: if (ksp->guess_knoll) {
219: PCApply(ksp->B,ksp->vec_rhs,ksp->vec_sol,PC_LEFT);
220: ksp->guess_zero = PETSC_FALSE;
221: }
223: /* reset the residual history list if requested */
224: if (ksp->res_hist_reset) ksp->res_hist_len = 0;
226: ksp->transpose_solve = PETSC_FALSE;
227: (*ksp->ops->solve)(ksp,&nits);
228: if (!ksp->reason) {
229: SETERRQ(1,"Internal error, solver returned without setting converged reason");
230: }
231: if (its) *its = nits;
233: MPI_Comm_rank(ksp->comm,&rank);
235: PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues",&flag1);
236: PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues",&flag2);
237: if (flag1 || flag2) {
238: int n = nits + 2,i,neig;
239: PetscReal *r,*c;
241: if (!n) {
242: PetscPrintf(ksp->comm,"Zero iterations in solver, cannot approximate any eigenvalues\n");
243: } else {
244: PetscMalloc(2*n*sizeof(PetscReal),&r);
245: c = r + n;
246: KSPComputeEigenvalues(ksp,n,r,c,&neig);
247: if (flag1) {
248: PetscPrintf(ksp->comm,"Iteratively computed eigenvalues\n");
249: for (i=0; i<neig; i++) {
250: if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gi\n",r[i],c[i]);}
251: else {PetscPrintf(ksp->comm,"%g - %gi\n",r[i],-c[i]);}
252: }
253: }
254: if (flag2 && !rank) {
255: PetscViewer viewer;
256: PetscDraw draw;
257: PetscDrawSP drawsp;
259: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Iteratively Computed Eigenvalues",
260: PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
261: PetscViewerDrawGetDraw(viewer,0,&draw);
262: PetscDrawSPCreate(draw,1,&drawsp);
263: for (i=0; i<neig; i++) {
264: PetscDrawSPAddPoint(drawsp,r+i,c+i);
265: }
266: PetscDrawSPDraw(drawsp);
267: PetscDrawSPDestroy(drawsp);
268: PetscViewerDestroy(viewer);
269: }
270: PetscFree(r);
271: }
272: }
274: PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues_explicitly",&flag1);
275: PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues_explicitly",&flag2);
276: if (flag1 || flag2) {
277: int n,i;
278: PetscReal *r,*c;
279: VecGetSize(ksp->vec_sol,&n);
280: PetscMalloc(2*n*sizeof(PetscReal),&r);
281: c = r + n;
282: KSPComputeEigenvaluesExplicitly(ksp,n,r,c);
283: if (flag1) {
284: PetscPrintf(ksp->comm,"Explicitly computed eigenvalues\n");
285: for (i=0; i<n; i++) {
286: if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gi\n",r[i],c[i]);}
287: else {PetscPrintf(ksp->comm,"%g - %gi\n",r[i],-c[i]);}
288: }
289: }
290: if (flag2 && !rank) {
291: PetscViewer viewer;
292: PetscDraw draw;
293: PetscDrawSP drawsp;
295: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Explicitly Computed Eigenvalues",0,320,300,300,&viewer);
296: PetscViewerDrawGetDraw(viewer,0,&draw);
297: PetscDrawSPCreate(draw,1,&drawsp);
298: for (i=0; i<n; i++) {
299: PetscDrawSPAddPoint(drawsp,r+i,c+i);
300: }
301: PetscDrawSPDraw(drawsp);
302: PetscDrawSPDestroy(drawsp);
303: PetscViewerDestroy(viewer);
304: }
305: PetscFree(r);
306: }
308: PetscOptionsHasName(ksp->prefix,"-ksp_view_operator",&flag2);
309: if (flag2) {
310: Mat A,B;
311: PCGetOperators(ksp->B,&A,PETSC_NULL,PETSC_NULL);
312: MatComputeExplicitOperator(A,&B);
313: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(ksp->comm),PETSC_VIEWER_ASCII_MATLAB);
314: MatView(B,PETSC_VIEWER_STDOUT_(ksp->comm));
315: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(ksp->comm));
316: MatDestroy(B);
317: }
318: PetscOptionsHasName(ksp->prefix,"-ksp_view_operator_binary",&flag2);
319: if (flag2) {
320: Mat A,B;
321: PCGetOperators(ksp->B,&A,PETSC_NULL,PETSC_NULL);
322: MatComputeExplicitOperator(A,&B);
323: MatView(B,PETSC_VIEWER_BINARY_(ksp->comm));
324: MatDestroy(B);
325: }
326: PetscOptionsHasName(ksp->prefix,"-ksp_view_preconditioned_operator_binary",&flag2);
327: if (flag2) {
328: Mat B;
329: KSPComputeExplicitOperator(ksp,&B);
330: MatView(B,PETSC_VIEWER_BINARY_(ksp->comm));
331: MatDestroy(B);
332: }
333: return(0);
334: }
338: /*@
339: KSPSolveTranspose - Solves the transpose of a linear system. Usually
340: accessed through SLESSolveTranspose().
342: Collective on KSP
344: Input Parameter:
345: . ksp - iterative context obtained from KSPCreate()
347: Output Parameter:
348: . its - number of iterations required
350: Notes:
351: On return, the parameter "its" contains either the iteration
352: number at which convergence was successfully reached, or the
353: negative of the iteration at which divergence or breakdown was detected.
355: Currently only supported by KSPType of KSPPREONLY. This routine is usally
356: only used internally by the BiCG solver on the subblocks in BJacobi and ASM.
358: Level: developer
360: .keywords: KSP, solve, linear system
362: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPDefaultConverged(),
363: SLESSolve(), SLESGetKSP()
364: @*/
365: int KSPSolveTranspose(KSP ksp,int *its)
366: {
367: int ierr;
368: PetscScalar zero = 0.0;
373: if (!ksp->setupcalled){ KSPSetUp(ksp);}
374: if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
375: ksp->transpose_solve = PETSC_TRUE;
376: (*ksp->ops->solve)(ksp,its);
377: return(0);
378: }
382: /*@C
383: KSPDestroy - Destroys KSP context.
385: Collective on KSP
387: Input Parameter:
388: . ksp - iterative context obtained from KSPCreate()
390: Level: developer
392: .keywords: KSP, destroy
394: .seealso: KSPCreate(), KSPSetUp(), KSPSolve()
395: @*/
396: int KSPDestroy(KSP ksp)
397: {
398: int i,ierr;
402: if (--ksp->refct > 0) return(0);
404: /* if memory was published with AMS then destroy it */
405: PetscObjectDepublish(ksp);
407: if (ksp->ops->destroy) {
408: (*ksp->ops->destroy)(ksp);
409: }
410: for (i=0; i<ksp->numbermonitors; i++) {
411: if (ksp->monitordestroy[i]) {
412: (*ksp->monitordestroy[i])(ksp->monitorcontext[i]);
413: }
414: }
415: PetscLogObjectDestroy(ksp);
416: PetscHeaderDestroy(ksp);
417: return(0);
418: }
422: /*@
423: KSPSetPreconditionerSide - Sets the preconditioning side.
425: Collective on KSP
427: Input Parameter:
428: . ksp - iterative context obtained from KSPCreate()
430: Output Parameter:
431: . side - the preconditioning side, where side is one of
432: .vb
433: PC_LEFT - left preconditioning (default)
434: PC_RIGHT - right preconditioning
435: PC_SYMMETRIC - symmetric preconditioning
436: .ve
438: Options Database Keys:
439: + -ksp_left_pc - Sets left preconditioning
440: . -ksp_right_pc - Sets right preconditioning
441: - -ksp_symmetric_pc - Sets symmetric preconditioning
443: Notes:
444: Left preconditioning is used by default. Symmetric preconditioning is
445: currently available only for the KSPQCG method. Note, however, that
446: symmetric preconditioning can be emulated by using either right or left
447: preconditioning and a pre or post processing step.
449: Level: intermediate
451: .keywords: KSP, set, right, left, symmetric, side, preconditioner, flag
453: .seealso: KSPGetPreconditionerSide()
454: @*/
455: int KSPSetPreconditionerSide(KSP ksp,PCSide side)
456: {
459: ksp->pc_side = side;
460: return(0);
461: }
465: /*@C
466: KSPGetPreconditionerSide - Gets the preconditioning side.
468: Not Collective
470: Input Parameter:
471: . ksp - iterative context obtained from KSPCreate()
473: Output Parameter:
474: . side - the preconditioning side, where side is one of
475: .vb
476: PC_LEFT - left preconditioning (default)
477: PC_RIGHT - right preconditioning
478: PC_SYMMETRIC - symmetric preconditioning
479: .ve
481: Level: intermediate
483: .keywords: KSP, get, right, left, symmetric, side, preconditioner, flag
485: .seealso: KSPSetPreconditionerSide()
486: @*/
487: int KSPGetPreconditionerSide(KSP ksp,PCSide *side)
488: {
491: *side = ksp->pc_side;
492: return(0);
493: }
497: /*@
498: KSPGetTolerances - Gets the relative, absolute, divergence, and maximum
499: iteration tolerances used by the default KSP convergence tests.
501: Not Collective
503: Input Parameter:
504: . ksp - the Krylov subspace context
505:
506: Output Parameters:
507: + rtol - the relative convergence tolerance
508: . atol - the absolute convergence tolerance
509: . dtol - the divergence tolerance
510: - maxits - maximum number of iterations
512: Notes:
513: The user can specify PETSC_NULL for any parameter that is not needed.
515: Level: intermediate
517: .keywords: KSP, get, tolerance, absolute, relative, divergence, convergence,
518: maximum, iterations
520: .seealso: KSPSetTolerances()
521: @*/
522: int KSPGetTolerances(KSP ksp,PetscReal *rtol,PetscReal *atol,PetscReal *dtol,int *maxits)
523: {
526: if (atol) *atol = ksp->atol;
527: if (rtol) *rtol = ksp->rtol;
528: if (dtol) *dtol = ksp->divtol;
529: if (maxits) *maxits = ksp->max_it;
530: return(0);
531: }
535: /*@
536: KSPSetTolerances - Sets the relative, absolute, divergence, and maximum
537: iteration tolerances used by the default KSP convergence testers.
539: Collective on KSP
541: Input Parameters:
542: + ksp - the Krylov subspace context
543: . rtol - the relative convergence tolerance
544: (relative decrease in the residual norm)
545: . atol - the absolute convergence tolerance
546: (absolute size of the residual norm)
547: . dtol - the divergence tolerance
548: (amount residual can increase before KSPDefaultConverged()
549: concludes that the method is diverging)
550: - maxits - maximum number of iterations to use
552: Options Database Keys:
553: + -ksp_atol <atol> - Sets atol
554: . -ksp_rtol <rtol> - Sets rtol
555: . -ksp_divtol <dtol> - Sets dtol
556: - -ksp_max_it <maxits> - Sets maxits
558: Notes:
559: Use PETSC_DEFAULT to retain the default value of any of the tolerances.
561: See KSPDefaultConverged() for details on the use of these parameters
562: in the default convergence test. See also KSPSetConvergenceTest()
563: for setting user-defined stopping criteria.
565: Level: intermediate
567: .keywords: KSP, set, tolerance, absolute, relative, divergence,
568: convergence, maximum, iterations
570: .seealso: KSPGetTolerances(), KSPDefaultConverged(), KSPSetConvergenceTest()
571: @*/
572: int KSPSetTolerances(KSP ksp,PetscReal rtol,PetscReal atol,PetscReal dtol,int maxits)
573: {
576: if (atol != PETSC_DEFAULT) ksp->atol = atol;
577: if (rtol != PETSC_DEFAULT) ksp->rtol = rtol;
578: if (dtol != PETSC_DEFAULT) ksp->divtol = dtol;
579: if (maxits != PETSC_DEFAULT) ksp->max_it = maxits;
580: return(0);
581: }
585: /*@
586: KSPSetInitialGuessNonzero - Tells the iterative solver that the
587: initial guess is nonzero; otherwise KSP assumes the initial guess
588: is to be zero (and thus zeros it out before solving).
590: Collective on KSP
592: Input Parameters:
593: + ksp - iterative context obtained from SLESGetKSP() or KSPCreate()
594: - flg - PETSC_TRUE indicates the guess is non-zero, PETSC_FALSE indicates the guess is zero
596: Level: beginner
598: Notes:
599: If this is not called the X vector is zeroed in the call to
600: SLESSolve() (or KSPSolve()).
602: .keywords: KSP, set, initial guess, nonzero
604: .seealso: KSPGetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
605: @*/
606: int KSPSetInitialGuessNonzero(KSP ksp,PetscTruth flg)
607: {
609: ksp->guess_zero = (PetscTruth)!(int)flg;
610: return(0);
611: }
615: /*@
616: KSPGetInitialGuessNonzero - Determines whether the KSP solver is using
617: a zero initial guess.
619: Not Collective
621: Input Parameter:
622: . ksp - iterative context obtained from KSPCreate()
624: Output Parameter:
625: . flag - PETSC_TRUE if guess is nonzero, else PETSC_FALSE
627: Level: intermediate
629: .keywords: KSP, set, initial guess, nonzero
631: .seealso: KSPSetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
632: @*/
633: int KSPGetInitialGuessNonzero(KSP ksp,PetscTruth *flag)
634: {
636: if (ksp->guess_zero) *flag = PETSC_FALSE;
637: else *flag = PETSC_TRUE;
638: return(0);
639: }
643: /*@
644: KSPSetInitialGuessKnoll - Tells the iterative solver to use PCApply(pc,b,..) to compute the initial guess (The Knoll trick)
646: Collective on KSP
648: Input Parameters:
649: + ksp - iterative context obtained from SLESGetKSP() or KSPCreate()
650: - flg - PETSC_TRUE or PETSC_FALSE
652: Level: advanced
655: .keywords: KSP, set, initial guess, nonzero
657: .seealso: KSPGetInitialGuessKnoll(), KSPSetInitialGuessNonzero(), KSPGetInitialGuessNonzero()
658: @*/
659: int KSPSetInitialGuessKnoll(KSP ksp,PetscTruth flg)
660: {
662: ksp->guess_knoll = flg;
663: return(0);
664: }
668: /*@
669: KSPGetInitialGuessKnoll - Determines whether the KSP solver is using the Knoll trick (using PCApply(pc,b,...) to compute
670: the initial guess
672: Not Collective
674: Input Parameter:
675: . ksp - iterative context obtained from KSPCreate()
677: Output Parameter:
678: . flag - PETSC_TRUE if using Knoll trick, else PETSC_FALSE
680: Level: advanced
682: .keywords: KSP, set, initial guess, nonzero
684: .seealso: KSPSetInitialGuessKnoll(), KSPSetInitialGuessNonzero(), KSPGetInitialGuessNonzero()
685: @*/
686: int KSPGetInitialGuessKnoll(KSP ksp,PetscTruth *flag)
687: {
689: *flag = ksp->guess_knoll;
690: return(0);
691: }
695: /*@
696: KSPSetComputeSingularValues - Sets a flag so that the extreme singular
697: values will be calculated via a Lanczos or Arnoldi process as the linear
698: system is solved.
700: Collective on KSP
702: Input Parameters:
703: + ksp - iterative context obtained from KSPCreate()
704: - flg - PETSC_TRUE or PETSC_FALSE
706: Options Database Key:
707: . -ksp_singmonitor - Activates KSPSetComputeSingularValues()
709: Notes:
710: Currently this option is not valid for all iterative methods.
712: Many users may just want to use the monitoring routine
713: KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
714: to print the singular values at each iteration of the linear solve.
716: Level: advanced
718: .keywords: KSP, set, compute, singular values
720: .seealso: KSPComputeExtremeSingularValues(), KSPSingularValueMonitor()
721: @*/
722: int KSPSetComputeSingularValues(KSP ksp,PetscTruth flg)
723: {
726: ksp->calc_sings = flg;
727: return(0);
728: }
732: /*@
733: KSPSetComputeEigenvalues - Sets a flag so that the extreme eigenvalues
734: values will be calculated via a Lanczos or Arnoldi process as the linear
735: system is solved.
737: Collective on KSP
739: Input Parameters:
740: + ksp - iterative context obtained from KSPCreate()
741: - flg - PETSC_TRUE or PETSC_FALSE
743: Notes:
744: Currently this option is not valid for all iterative methods.
746: Level: advanced
748: .keywords: KSP, set, compute, eigenvalues
750: .seealso: KSPComputeEigenvalues(), KSPComputeEigenvaluesExplicitly()
751: @*/
752: int KSPSetComputeEigenvalues(KSP ksp,PetscTruth flg)
753: {
756: ksp->calc_sings = flg;
757: return(0);
758: }
762: /*@
763: KSPSetRhs - Sets the right-hand-side vector for the linear system to
764: be solved.
766: Collective on KSP and Vec
768: Input Parameters:
769: + ksp - iterative context obtained from KSPCreate()
770: - b - right-hand-side vector
772: Level: developer
774: .keywords: KSP, set, right-hand-side, rhs
776: .seealso: KSPGetRhs(), KSPSetSolution()
777: @*/
778: int KSPSetRhs(KSP ksp,Vec b)
779: {
783: ksp->vec_rhs = (b);
784: return(0);
785: }
789: /*@C
790: KSPGetRhs - Gets the right-hand-side vector for the linear system to
791: be solved.
793: Not Collective
795: Input Parameter:
796: . ksp - iterative context obtained from KSPCreate()
798: Output Parameter:
799: . r - right-hand-side vector
801: Level: developer
803: .keywords: KSP, get, right-hand-side, rhs
805: .seealso: KSPSetRhs(), KSPGetSolution()
806: @*/
807: int KSPGetRhs(KSP ksp,Vec *r)
808: {
811: *r = ksp->vec_rhs;
812: return(0);
813: }
817: /*@
818: KSPSetSolution - Sets the location of the solution for the
819: linear system to be solved.
821: Collective on KSP and Vec
823: Input Parameters:
824: + ksp - iterative context obtained from KSPCreate()
825: - x - solution vector
827: Level: developer
829: .keywords: KSP, set, solution
831: .seealso: KSPSetRhs(), KSPGetSolution()
832: @*/
833: int KSPSetSolution(KSP ksp,Vec x)
834: {
838: ksp->vec_sol = (x);
839: return(0);
840: }
844: /*@C
845: KSPGetSolution - Gets the location of the solution for the
846: linear system to be solved. Note that this may not be where the solution
847: is stored during the iterative process; see KSPBuildSolution().
849: Not Collective
851: Input Parameters:
852: . ksp - iterative context obtained from KSPCreate()
854: Output Parameters:
855: . v - solution vector
857: Level: developer
859: .keywords: KSP, get, solution
861: .seealso: KSPGetRhs(), KSPSetSolution(), KSPBuildSolution()
862: @*/
863: int KSPGetSolution(KSP ksp,Vec *v)
864: {
867: *v = ksp->vec_sol;
868: return(0);
869: }
873: /*@
874: KSPSetPC - Sets the preconditioner to be used to calculate the
875: application of the preconditioner on a vector.
877: Collective on KSP
879: Input Parameters:
880: + ksp - iterative context obtained from KSPCreate()
881: - B - the preconditioner object
883: Notes:
884: Use KSPGetPC() to retrieve the preconditioner context (for example,
885: to free it at the end of the computations).
887: Level: developer
889: .keywords: KSP, set, precondition, Binv
891: .seealso: KSPGetPC()
892: @*/
893: int KSPSetPC(KSP ksp,PC B)
894: {
899: ksp->B = B;
900: return(0);
901: }
905: /*@C
906: KSPGetPC - Returns a pointer to the preconditioner context
907: set with KSPSetPC().
909: Not Collective
911: Input Parameters:
912: . ksp - iterative context obtained from KSPCreate()
914: Output Parameter:
915: . B - preconditioner context
917: Level: developer
919: .keywords: KSP, get, preconditioner, Binv
921: .seealso: KSPSetPC()
922: @*/
923: int KSPGetPC(KSP ksp,PC *B)
924: {
927: *B = ksp->B;
928: return(0);
929: }
933: /*@C
934: KSPSetMonitor - Sets an ADDITIONAL function to be called at every iteration to monitor
935: the residual/error etc.
936:
937: Collective on KSP
939: Input Parameters:
940: + ksp - iterative context obtained from KSPCreate()
941: . monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring
942: . mctx - [optional] context for private data for the
943: monitor routine (use PETSC_NULL if no context is desired)
944: - monitordestroy - [optional] routine that frees monitor context
945: (may be PETSC_NULL)
947: Calling Sequence of monitor:
948: $ monitor (KSP ksp, int it, PetscReal rnorm, void *mctx)
950: + ksp - iterative context obtained from KSPCreate()
951: . it - iteration number
952: . rnorm - (estimated) 2-norm of (preconditioned) residual
953: - mctx - optional monitoring context, as set by KSPSetMonitor()
955: Options Database Keys:
956: + -ksp_monitor - sets KSPDefaultMonitor()
957: . -ksp_truemonitor - sets KSPTrueMonitor()
958: . -ksp_xmonitor - sets line graph monitor,
959: uses KSPLGMonitorCreate()
960: . -ksp_xtruemonitor - sets line graph monitor,
961: uses KSPLGMonitorCreate()
962: . -ksp_singmonitor - sets KSPSingularValueMonitor()
963: - -ksp_cancelmonitors - cancels all monitors that have
964: been hardwired into a code by
965: calls to KSPSetMonitor(), but
966: does not cancel those set via
967: the options database.
969: Notes:
970: The default is to do nothing. To print the residual, or preconditioned
971: residual if KSPSetNormType(ksp,KSP_PRECONDITIONED_NORM) was called, use
972: KSPDefaultMonitor() as the monitoring routine, with a null monitoring
973: context.
975: Several different monitoring routines may be set by calling
976: KSPSetMonitor() multiple times; all will be called in the
977: order in which they were set.
979: Level: beginner
981: .keywords: KSP, set, monitor
983: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPClearMonitor()
984: @*/
985: int KSPSetMonitor(KSP ksp,int (*monitor)(KSP,int,PetscReal,void*),void *mctx,int (*monitordestroy)(void*))
986: {
989: if (ksp->numbermonitors >= MAXKSPMONITORS) {
990: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many KSP monitors set");
991: }
992: ksp->monitor[ksp->numbermonitors] = monitor;
993: ksp->monitordestroy[ksp->numbermonitors] = monitordestroy;
994: ksp->monitorcontext[ksp->numbermonitors++] = (void*)mctx;
995: return(0);
996: }
1000: /*@
1001: KSPClearMonitor - Clears all monitors for a KSP object.
1003: Collective on KSP
1005: Input Parameters:
1006: . ksp - iterative context obtained from KSPCreate()
1008: Options Database Key:
1009: . -ksp_cancelmonitors - Cancels all monitors that have
1010: been hardwired into a code by calls to KSPSetMonitor(),
1011: but does not cancel those set via the options database.
1013: Level: intermediate
1015: .keywords: KSP, set, monitor
1017: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPSetMonitor()
1018: @*/
1019: int KSPClearMonitor(KSP ksp)
1020: {
1023: ksp->numbermonitors = 0;
1024: return(0);
1025: }
1029: /*@C
1030: KSPGetMonitorContext - Gets the monitoring context, as set by
1031: KSPSetMonitor() for the FIRST monitor only.
1033: Not Collective
1035: Input Parameter:
1036: . ksp - iterative context obtained from KSPCreate()
1038: Output Parameter:
1039: . ctx - monitoring context
1041: Level: intermediate
1043: .keywords: KSP, get, monitor, context
1045: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate()
1046: @*/
1047: int KSPGetMonitorContext(KSP ksp,void **ctx)
1048: {
1051: *ctx = (ksp->monitorcontext[0]);
1052: return(0);
1053: }
1057: /*@
1058: KSPSetResidualHistory - Sets the array used to hold the residual history.
1059: If set, this array will contain the residual norms computed at each
1060: iteration of the solver.
1062: Not Collective
1064: Input Parameters:
1065: + ksp - iterative context obtained from KSPCreate()
1066: . a - array to hold history
1067: . na - size of a
1068: - reset - PETSC_TRUE indicates the history counter is reset to zero
1069: for each new linear solve
1071: Level: advanced
1073: Notes: The array is NOT freed by PETSc so the user needs to keep track of
1074: it and destroy once the KSP object is destroyed.
1076: .keywords: KSP, set, residual, history, norm
1078: .seealso: KSPGetResidualHistory()
1080: @*/
1081: int KSPSetResidualHistory(KSP ksp,PetscReal a[],int na,PetscTruth reset)
1082: {
1087: if (na != PETSC_DECIDE && a != PETSC_NULL) {
1088: ksp->res_hist = a;
1089: ksp->res_hist_max = na;
1090: } else {
1091: ksp->res_hist_max = 1000;
1092: PetscMalloc(ksp->res_hist_max*sizeof(PetscReal),&ksp->res_hist);
1093: }
1094: ksp->res_hist_len = 0;
1095: ksp->res_hist_reset = reset;
1098: return(0);
1099: }
1103: /*@C
1104: KSPGetResidualHistory - Gets the array used to hold the residual history
1105: and the number of residuals it contains.
1107: Not Collective
1109: Input Parameter:
1110: . ksp - iterative context obtained from KSPCreate()
1112: Output Parameters:
1113: + a - pointer to array to hold history (or PETSC_NULL)
1114: - na - number of used entries in a (or PETSC_NULL)
1116: Level: advanced
1118: Notes:
1119: Can only be called after a KSPSetResidualHistory() otherwise a and na are set to zero
1121: The Fortran version of this routine has a calling sequence
1122: $ call KSPGetResidualHistory(KSP ksp, integer na, integer ierr)
1124: .keywords: KSP, get, residual, history, norm
1126: .seealso: KSPGetResidualHistory()
1128: @*/
1129: int KSPGetResidualHistory(KSP ksp,PetscReal *a[],int *na)
1130: {
1133: if (a) *a = ksp->res_hist;
1134: if (na) *na = ksp->res_hist_len;
1135: return(0);
1136: }
1140: /*@C
1141: KSPSetConvergenceTest - Sets the function to be used to determine
1142: convergence.
1144: Collective on KSP
1146: Input Parameters:
1147: + ksp - iterative context obtained from KSPCreate()
1148: . converge - pointer to int function
1149: - cctx - context for private data for the convergence routine (may be null)
1151: Calling sequence of converge:
1152: $ converge (KSP ksp, int it, PetscReal rnorm, KSPConvergedReason *reason,void *mctx)
1154: + ksp - iterative context obtained from KSPCreate()
1155: . it - iteration number
1156: . rnorm - (estimated) 2-norm of (preconditioned) residual
1157: . reason - the reason why it has converged or diverged
1158: - cctx - optional convergence context, as set by KSPSetConvergenceTest()
1161: Notes:
1162: Must be called after the KSP type has been set so put this after
1163: a call to KSPSetType(), or KSPSetFromOptions().
1165: The default convergence test, KSPDefaultConverged(), aborts if the
1166: residual grows to more than 10000 times the initial residual.
1168: The default is a combination of relative and absolute tolerances.
1169: The residual value that is tested may be an approximation; routines
1170: that need exact values should compute them.
1172: Level: advanced
1174: .keywords: KSP, set, convergence, test, context
1176: .seealso: KSPDefaultConverged(), KSPGetConvergenceContext()
1177: @*/
1178: int KSPSetConvergenceTest(KSP ksp,int (*converge)(KSP,int,PetscReal,KSPConvergedReason*,void*),void *cctx)
1179: {
1182: ksp->converged = converge;
1183: ksp->cnvP = (void*)cctx;
1184: return(0);
1185: }
1189: /*@C
1190: KSPGetConvergenceContext - Gets the convergence context set with
1191: KSPSetConvergenceTest().
1193: Not Collective
1195: Input Parameter:
1196: . ksp - iterative context obtained from KSPCreate()
1198: Output Parameter:
1199: . ctx - monitoring context
1201: Level: advanced
1203: .keywords: KSP, get, convergence, test, context
1205: .seealso: KSPDefaultConverged(), KSPSetConvergenceTest()
1206: @*/
1207: int KSPGetConvergenceContext(KSP ksp,void **ctx)
1208: {
1211: *ctx = ksp->cnvP;
1212: return(0);
1213: }
1217: /*@C
1218: KSPBuildSolution - Builds the approximate solution in a vector provided.
1219: This routine is NOT commonly needed (see SLESSolve()).
1221: Collective on KSP
1223: Input Parameter:
1224: . ctx - iterative context obtained from KSPCreate()
1226: Output Parameter:
1227: Provide exactly one of
1228: + v - location to stash solution.
1229: - V - the solution is returned in this location. This vector is created
1230: internally. This vector should NOT be destroyed by the user with
1231: VecDestroy().
1233: Notes:
1234: This routine can be used in one of two ways
1235: .vb
1236: KSPBuildSolution(ksp,PETSC_NULL,&V);
1237: or
1238: KSPBuildSolution(ksp,v,PETSC_NULL);
1239: .ve
1240: In the first case an internal vector is allocated to store the solution
1241: (the user cannot destroy this vector). In the second case the solution
1242: is generated in the vector that the user provides. Note that for certain
1243: methods, such as KSPCG, the second case requires a copy of the solution,
1244: while in the first case the call is essentially free since it simply
1245: returns the vector where the solution already is stored.
1247: Level: advanced
1249: .keywords: KSP, build, solution
1251: .seealso: KSPGetSolution(), KSPBuildResidual()
1252: @*/
1253: int KSPBuildSolution(KSP ksp,Vec v,Vec *V)
1254: {
1259: if (!V && !v) SETERRQ(PETSC_ERR_ARG_WRONG,"Must provide either v or V");
1260: if (!V) V = &v;
1261: (*ksp->ops->buildsolution)(ksp,v,V);
1262: return(0);
1263: }
1267: /*@C
1268: KSPBuildResidual - Builds the residual in a vector provided.
1270: Collective on KSP
1272: Input Parameter:
1273: . ksp - iterative context obtained from KSPCreate()
1275: Output Parameters:
1276: + v - optional location to stash residual. If v is not provided,
1277: then a location is generated.
1278: . t - work vector. If not provided then one is generated.
1279: - V - the residual
1281: Notes:
1282: Regardless of whether or not v is provided, the residual is
1283: returned in V.
1285: Level: advanced
1287: .keywords: KSP, build, residual
1289: .seealso: KSPBuildSolution()
1290: @*/
1291: int KSPBuildResidual(KSP ksp,Vec t,Vec v,Vec *V)
1292: {
1293: int flag = 0,ierr;
1294: Vec w = v,tt = t;
1298: if (!w) {
1299: VecDuplicate(ksp->vec_rhs,&w);
1300: PetscLogObjectParent((PetscObject)ksp,w);
1301: }
1302: if (!tt) {
1303: VecDuplicate(ksp->vec_rhs,&tt); flag = 1;
1304: PetscLogObjectParent((PetscObject)ksp,tt);
1305: }
1306: (*ksp->ops->buildresidual)(ksp,tt,w,V);
1307: if (flag) {VecDestroy(tt);}
1308: return(0);
1309: }