Actual source code: itfunc.c
1: /*$Id: itfunc.c,v 1.156 2001/04/10 19:36:24 bsmith Exp $*/
2: /*
3: Interface KSP routines that the user calls.
4: */
6: #include src/sles/ksp/kspimpl.h
8: /*@
9: KSPComputeExtremeSingularValues - Computes the extreme singular values
10: for the preconditioned operator. Called after or during KSPSolve()
11: (SLESSolve()).
13: Not Collective
15: Input Parameter:
16: . ksp - iterative context obtained from KSPCreate()
18: Output Parameters:
19: . emin, emax - extreme singular values
21: Notes:
22: One must call KSPSetComputeSingularValues() before calling KSPSetUp()
23: (or use the option -ksp_compute_eigenvalues) in order for this routine to work correctly.
25: Many users may just want to use the monitoring routine
26: KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
27: to print the extreme singular values at each iteration of the linear solve.
29: Level: advanced
31: .keywords: KSP, compute, extreme, singular, values
33: .seealso: KSPSetComputeSingularValues(), KSPSingularValueMonitor(), KSPComputeEigenvalues()
34: @*/
35: int KSPComputeExtremeSingularValues(KSP ksp,PetscReal *emax,PetscReal *emin)
36: {
43: if (!ksp->calc_sings) {
44: SETERRQ(4,"Singular values not requested before KSPSetUp()");
45: }
47: if (ksp->ops->computeextremesingularvalues) {
48: (*ksp->ops->computeextremesingularvalues)(ksp,emax,emin);
49: }
50: return(0);
51: }
53: /*@
54: KSPComputeEigenvalues - Computes the extreme eigenvalues for the
55: preconditioned operator. Called after or during KSPSolve() (SLESSolve()).
57: Not Collective
59: Input Parameter:
60: + ksp - iterative context obtained from KSPCreate()
61: - n - size of arrays r and c. The number of eigenvalues computed (neig) will, in
62: general, be less than this.
64: Output Parameters:
65: + r - real part of computed eigenvalues
66: . c - complex part of computed eigenvalues
67: - neig - number of eigenvalues computed (will be less than or equal to n)
69: Options Database Keys:
70: + -ksp_compute_eigenvalues - Prints eigenvalues to stdout
71: - -ksp_plot_eigenvalues - Plots eigenvalues in an x-window display
73: Notes:
74: The number of eigenvalues estimated depends on the size of the Krylov space
75: generated during the KSPSolve() (that is the SLESSolve); for example, with
76: CG it corresponds to the number of CG iterations, for GMRES it is the number
77: of GMRES iterations SINCE the last restart. Any extra space in r[] and c[]
78: will be ignored.
80: KSPComputeEigenvalues() does not usually provide accurate estimates; it is
81: intended only for assistance in understanding the convergence of iterative
82: methods, not for eigenanalysis.
84: One must call KSPSetComputeEigenvalues() before calling KSPSetUp()
85: in order for this routine to work correctly.
87: Many users may just want to use the monitoring routine
88: KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
89: to print the singular values at each iteration of the linear solve.
91: Level: advanced
93: .keywords: KSP, compute, extreme, singular, values
95: .seealso: KSPSetComputeSingularValues(), KSPSingularValueMonitor(), KSPComputeExtremeSingularValues()
96: @*/
97: int KSPComputeEigenvalues(KSP ksp,int n,PetscReal *r,PetscReal *c,int *neig)
98: {
105: if (!ksp->calc_sings) {
106: SETERRQ(4,"Eigenvalues not requested before KSPSetUp()");
107: }
109: if (ksp->ops->computeeigenvalues) {
110: (*ksp->ops->computeeigenvalues)(ksp,n,r,c,neig);
111: }
112: return(0);
113: }
115: /*@
116: KSPSetUp - Sets up the internal data structures for the
117: later use of an iterative solver.
119: Collective on KSP
121: Input Parameter:
122: . ksp - iterative context obtained from KSPCreate()
124: Level: developer
126: .keywords: KSP, setup
128: .seealso: KSPCreate(), KSPSolve(), KSPDestroy()
129: @*/
130: int KSPSetUp(KSP ksp)
131: {
137: /* reset the convergence flag from the previous solves */
138: ksp->reason = KSP_CONVERGED_ITERATING;
140: if (!ksp->type_name){
141: KSPSetType(ksp,KSPGMRES);
142: }
144: if (ksp->setupcalled) return(0);
145: ksp->setupcalled = 1;
146: (*ksp->ops->setup)(ksp);
147: return(0);
148: }
151: /*@
152: KSPSolve - Solves linear system; usually not called directly, rather
153: it is called by a call to SLESSolve().
155: Collective on KSP
157: Input Parameter:
158: . ksp - iterative context obtained from KSPCreate()
160: Output Parameter:
161: . its - number of iterations required
163: Options Database:
164: + -ksp_compute_eigenvalues - compute preconditioned operators eigenvalues
165: . -ksp_plot_eigenvalues - plot the computed eigenvalues in an X-window
166: . -ksp_compute_eigenvalues_explicitly - compute the eigenvalues by forming the
167: dense operator and useing LAPACK
168: - -ksp_plot_eigenvalues_explicitly - plot the explicitly computing eigenvalues
170: Notes:
171: On return, the parameter "its" contains either the iteration
172: number at which convergence was successfully reached or failure was detected.
174: Call KSPGetConvergedReason() to determine if the solver converged or failed and
175: why.
176:
177: If using a direct method (e.g., via the KSP solver
178: KSPPREONLY and a preconditioner such as PCLU/PCILU),
179: then its=1. See KSPSetTolerances() and KSPDefaultConverged()
180: for more details.
182: Understanding Convergence:
183: The routines KSPSetMonitor(), KSPComputeEigenvalues(), and
184: KSPComputeEigenvaluesExplicitly() provide information on additional
185: options to monitor convergence and print eigenvalue information.
187: Level: developer
189: .keywords: KSP, solve, linear system
191: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPDefaultConverged(),
192: SLESSolve(), KSPSolveTranspose(), SLESGetKSP()
193: @*/
194: int KSPSolve(KSP ksp,int *its)
195: {
196: int ierr,rank,nits;
197: PetscTruth flag1,flag2;
198: Scalar zero = 0.0;
203: if (!ksp->setupcalled){ KSPSetUp(ksp);}
204: if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
205: /* reset the residual history list if requested */
206: if (ksp->res_hist_reset) ksp->res_hist_len = 0;
208: ksp->transpose_solve = PETSC_FALSE;
209: (*ksp->ops->solve)(ksp,&nits);
210: if (!ksp->reason) {
211: SETERRQ(1,"Internal error, solver returned without setting converged reason");
212: }
213: if (its) *its = nits;
215: MPI_Comm_rank(ksp->comm,&rank);
217: PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues",&flag1);
218: PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues",&flag2);
219: if (flag1 || flag2) {
220: int n = nits + 2,i,neig;
221: PetscReal *r,*c;
223: if (!n) {
224: PetscPrintf(ksp->comm,"Zero iterations in solver, cannot approximate any eigenvaluesn");
225: } else {
226: PetscMalloc(2*n*sizeof(PetscReal),&r);
227: c = r + n;
228: KSPComputeEigenvalues(ksp,n,r,c,&neig);
229: if (flag1) {
230: PetscPrintf(ksp->comm,"Iteratively computed eigenvaluesn");
231: for (i=0; i<neig; i++) {
232: if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gin",r[i],c[i]);}
233: else {PetscPrintf(ksp->comm,"%g - %gin",r[i],-c[i]);}
234: }
235: }
236: if (flag2 && !rank) {
237: PetscViewer viewer;
238: PetscDraw draw;
239: PetscDrawSP drawsp;
241: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Iteratively Computed Eigenvalues",
242: PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
243: PetscViewerDrawGetDraw(viewer,0,&draw);
244: PetscDrawSPCreate(draw,1,&drawsp);
245: for (i=0; i<neig; i++) {
246: PetscDrawSPAddPoint(drawsp,r+i,c+i);
247: }
248: PetscDrawSPDraw(drawsp);
249: PetscDrawSPDestroy(drawsp);
250: PetscViewerDestroy(viewer);
251: }
252: PetscFree(r);
253: }
254: }
256: PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues_explicitly",&flag1);
257: PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues_explicitly",&flag2);
258: if (flag1 || flag2) {
259: int n,i;
260: PetscReal *r,*c;
261: VecGetSize(ksp->vec_sol,&n);
262: PetscMalloc(2*n*sizeof(PetscReal),&r);
263: c = r + n;
264: KSPComputeEigenvaluesExplicitly(ksp,n,r,c);
265: if (flag1) {
266: PetscPrintf(ksp->comm,"Explicitly computed eigenvaluesn");
267: for (i=0; i<n; i++) {
268: if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gin",r[i],c[i]);}
269: else {PetscPrintf(ksp->comm,"%g - %gin",r[i],-c[i]);}
270: }
271: }
272: if (flag2 && !rank) {
273: PetscViewer viewer;
274: PetscDraw draw;
275: PetscDrawSP drawsp;
277: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Explicitly Computed Eigenvalues",0,320,300,300,&viewer);
278: PetscViewerDrawGetDraw(viewer,0,&draw);
279: PetscDrawSPCreate(draw,1,&drawsp);
280: for (i=0; i<n; i++) {
281: PetscDrawSPAddPoint(drawsp,r+i,c+i);
282: }
283: PetscDrawSPDraw(drawsp);
284: PetscDrawSPDestroy(drawsp);
285: PetscViewerDestroy(viewer);
286: }
287: PetscFree(r);
288: }
290: PetscOptionsHasName(ksp->prefix,"-ksp_view_operator",&flag2);
291: if (flag2) {
292: Mat A,B;
293: PCGetOperators(ksp->B,&A,PETSC_NULL,PETSC_NULL);
294: MatComputeExplicitOperator(A,&B);
295: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(ksp->comm),PETSC_VIEWER_ASCII_MATLAB);
296: MatView(B,PETSC_VIEWER_STDOUT_(ksp->comm));
297: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(ksp->comm));
298: MatDestroy(B);
299: }
300: return(0);
301: }
303: /*@
304: KSPSolveTranspose - Solves the transpose of a linear system. Usually
305: accessed through SLESSolveTranspose().
307: Collective on KSP
309: Input Parameter:
310: . ksp - iterative context obtained from KSPCreate()
312: Output Parameter:
313: . its - number of iterations required
315: Notes:
316: On return, the parameter "its" contains either the iteration
317: number at which convergence was successfully reached, or the
318: negative of the iteration at which divergence or breakdown was detected.
320: Currently only supported by KSPType of KSPPREONLY. This routine is usally
321: only used internally by the BiCG solver on the subblocks in BJacobi and ASM.
323: Level: developer
325: .keywords: KSP, solve, linear system
327: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPDefaultConverged(),
328: SLESSolve(), SLESGetKSP()
329: @*/
330: int KSPSolveTranspose(KSP ksp,int *its)
331: {
332: int ierr;
333: Scalar zero = 0.0;
339: if (!ksp->setupcalled){ KSPSetUp(ksp);}
340: if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
341: ksp->transpose_solve = PETSC_TRUE;
342: (*ksp->ops->solve)(ksp,its);
343: return(0);
344: }
346: /*@C
347: KSPDestroy - Destroys KSP context.
349: Collective on KSP
351: Input Parameter:
352: . ksp - iterative context obtained from KSPCreate()
354: Level: developer
356: .keywords: KSP, destroy
358: .seealso: KSPCreate(), KSPSetUp(), KSPSolve()
359: @*/
360: int KSPDestroy(KSP ksp)
361: {
362: int i,ierr;
366: if (--ksp->refct > 0) return(0);
368: /* if memory was published with AMS then destroy it */
369: PetscObjectDepublish(ksp);
371: if (ksp->ops->destroy) {
372: (*ksp->ops->destroy)(ksp);
373: }
374: for (i=0; i<ksp->numbermonitors; i++) {
375: if (ksp->monitordestroy[i]) {
376: (*ksp->monitordestroy[i])(ksp->monitorcontext[i]);
377: }
378: }
379: PetscLogObjectDestroy(ksp);
380: PetscHeaderDestroy(ksp);
381: return(0);
382: }
384: /*@
385: KSPSetPreconditionerSide - Sets the preconditioning side.
387: Collective on KSP
389: Input Parameter:
390: . ksp - iterative context obtained from KSPCreate()
392: Output Parameter:
393: . side - the preconditioning side, where side is one of
394: .vb
395: PC_LEFT - left preconditioning (default)
396: PC_RIGHT - right preconditioning
397: PC_SYMMETRIC - symmetric preconditioning
398: .ve
400: Options Database Keys:
401: + -ksp_left_pc - Sets left preconditioning
402: . -ksp_right_pc - Sets right preconditioning
403: - -ksp_symmetric_pc - Sets symmetric preconditioning
405: Notes:
406: Left preconditioning is used by default. Symmetric preconditioning is
407: currently available only for the KSPQCG method. Note, however, that
408: symmetric preconditioning can be emulated by using either right or left
409: preconditioning and a pre or post processing step.
411: Level: intermediate
413: .keywords: KSP, set, right, left, symmetric, side, preconditioner, flag
415: .seealso: KSPGetPreconditionerSide()
416: @*/
417: int KSPSetPreconditionerSide(KSP ksp,PCSide side)
418: {
421: ksp->pc_side = side;
422: return(0);
423: }
425: /*@C
426: KSPGetPreconditionerSide - Gets the preconditioning side.
428: Not Collective
430: Input Parameter:
431: . ksp - iterative context obtained from KSPCreate()
433: Output Parameter:
434: . side - the preconditioning side, where side is one of
435: .vb
436: PC_LEFT - left preconditioning (default)
437: PC_RIGHT - right preconditioning
438: PC_SYMMETRIC - symmetric preconditioning
439: .ve
441: Level: intermediate
443: .keywords: KSP, get, right, left, symmetric, side, preconditioner, flag
445: .seealso: KSPSetPreconditionerSide()
446: @*/
447: int KSPGetPreconditionerSide(KSP ksp,PCSide *side)
448: {
451: *side = ksp->pc_side;
452: return(0);
453: }
455: /*@
456: KSPGetTolerances - Gets the relative, absolute, divergence, and maximum
457: iteration tolerances used by the default KSP convergence tests.
459: Not Collective
461: Input Parameter:
462: . ksp - the Krylov subspace context
463:
464: Output Parameters:
465: + rtol - the relative convergence tolerance
466: . atol - the absolute convergence tolerance
467: . dtol - the divergence tolerance
468: - maxits - maximum number of iterations
470: Notes:
471: The user can specify PETSC_NULL for any parameter that is not needed.
473: Level: intermediate
475: .keywords: KSP, get, tolerance, absolute, relative, divergence, convergence,
476: maximum, iterations
478: .seealso: KSPSetTolerances()
479: @*/
480: int KSPGetTolerances(KSP ksp,PetscReal *rtol,PetscReal *atol,PetscReal *dtol,int *maxits)
481: {
484: if (atol) *atol = ksp->atol;
485: if (rtol) *rtol = ksp->rtol;
486: if (dtol) *dtol = ksp->divtol;
487: if (maxits) *maxits = ksp->max_it;
488: return(0);
489: }
491: /*@
492: KSPSetTolerances - Sets the relative, absolute, divergence, and maximum
493: iteration tolerances used by the default KSP convergence testers.
495: Collective on KSP
497: Input Parameters:
498: + ksp - the Krylov subspace context
499: . rtol - the relative convergence tolerance
500: (relative decrease in the residual norm)
501: . atol - the absolute convergence tolerance
502: (absolute size of the residual norm)
503: . dtol - the divergence tolerance
504: (amount residual can increase before KSPDefaultConverged()
505: concludes that the method is diverging)
506: - maxits - maximum number of iterations to use
508: Options Database Keys:
509: + -ksp_atol <atol> - Sets atol
510: . -ksp_rtol <rtol> - Sets rtol
511: . -ksp_divtol <dtol> - Sets dtol
512: - -ksp_max_it <maxits> - Sets maxits
514: Notes:
515: Use PETSC_DEFAULT to retain the default value of any of the tolerances.
517: See KSPDefaultConverged() for details on the use of these parameters
518: in the default convergence test. See also KSPSetConvergenceTest()
519: for setting user-defined stopping criteria.
521: Level: intermediate
523: .keywords: KSP, set, tolerance, absolute, relative, divergence,
524: convergence, maximum, iterations
526: .seealso: KSPGetTolerances(), KSPDefaultConverged(), KSPSetConvergenceTest()
527: @*/
528: int KSPSetTolerances(KSP ksp,PetscReal rtol,PetscReal atol,PetscReal dtol,int maxits)
529: {
532: if (atol != PETSC_DEFAULT) ksp->atol = atol;
533: if (rtol != PETSC_DEFAULT) ksp->rtol = rtol;
534: if (dtol != PETSC_DEFAULT) ksp->divtol = dtol;
535: if (maxits != PETSC_DEFAULT) ksp->max_it = maxits;
536: return(0);
537: }
539: /*@
540: KSPSetComputeResidual - Sets a flag to indicate whether the two norm
541: of the residual is calculated at each iteration.
543: Collective on KSP
545: Input Parameters:
546: + ksp - iterative context obtained from KSPCreate()
547: - flag - PETSC_TRUE or PETSC_FALSE
549: Notes:
550: Most Krylov methods do not yet take advantage of flag = PETSC_FALSE.
552: Level: advanced
554: .keywords: KSP, set, residual, norm, calculate, flag
555: @*/
556: int KSPSetComputeResidual(KSP ksp,PetscTruth flag)
557: {
560: ksp->calc_res = flag;
561: return(0);
562: }
564: /*@
565: KSPSetUsePreconditionedResidual - Sets a flag so that the two norm of the
566: preconditioned residual is used rather than the true residual, in the
567: default convergence tests.
569: Collective on KSP
571: Input Parameter:
572: . ksp - iterative context obtained from KSPCreate()
574: Notes:
575: Currently only CG, CHEBYCHEV, and RICHARDSON use this with left
576: preconditioning. All other methods always used the preconditioned
577: residual. With right preconditioning this flag is ignored, since
578: the preconditioned residual and true residual are the same.
580: Options Database Key:
581: . -ksp_preres - Activates KSPSetUsePreconditionedResidual()
583: Level: advanced
585: .keywords: KSP, set, residual, precondition, flag
586: @*/
587: int KSPSetUsePreconditionedResidual(KSP ksp)
588: {
591: ksp->use_pres = PETSC_TRUE;
592: return(0);
593: }
595: /*@
596: KSPSetInitialGuessNonzero - Tells the iterative solver that the
597: initial guess is nonzero; otherwise KSP assumes the initial guess
598: is to be zero (and thus zeros it out before solving).
600: Collective on KSP
602: Input Parameters:
603: . ksp - iterative context obtained from SLESGetKSP() or KSPCreate()
605: Level: beginner
607: Notes:
608: If this is not called the X vector is zeroed in the call to
609: SLESSolve() (or KSPSolve()).
611: .keywords: KSP, set, initial guess, nonzero
613: .seealso: KSPGetIntialGuessNonzero()
614: @*/
615: int KSPSetInitialGuessNonzero(KSP ksp)
616: {
618: ksp->guess_zero = PETSC_FALSE;
619: return(0);
620: }
622: /*@
623: KSPGetInitialGuessNonzero - Determines whether the KSP solver is using
624: a zero initial guess.
626: Not Collective
628: Input Parameter:
629: . ksp - iterative context obtained from KSPCreate()
631: Output Parameter:
632: . flag - PETSC_TRUE if guess is nonzero, else PETSC_FALSE
634: Level: intermediate
636: .keywords: KSP, set, initial guess, nonzero
638: .seealso: KSPSetIntialGuessNonzero()
639: @*/
640: int KSPGetInitialGuessNonzero(KSP ksp,PetscTruth *flag)
641: {
643: if (ksp->guess_zero) *flag = PETSC_FALSE;
644: else *flag = PETSC_TRUE;
645: return(0);
646: }
648: /*@
649: KSPSetComputeSingularValues - Sets a flag so that the extreme singular
650: values will be calculated via a Lanczos or Arnoldi process as the linear
651: system is solved.
653: Collective on KSP
655: Input Parameters:
656: . ksp - iterative context obtained from KSPCreate()
658: Options Database Key:
659: . -ksp_singmonitor - Activates KSPSetComputeSingularValues()
661: Notes:
662: Currently this option is not valid for all iterative methods.
664: Many users may just want to use the monitoring routine
665: KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
666: to print the singular values at each iteration of the linear solve.
668: Level: advanced
670: .keywords: KSP, set, compute, singular values
672: .seealso: KSPComputeExtremeSingularValues(), KSPSingularValueMonitor()
673: @*/
674: int KSPSetComputeSingularValues(KSP ksp)
675: {
678: ksp->calc_sings = PETSC_TRUE;
679: return(0);
680: }
682: /*@
683: KSPSetComputeEigenvalues - Sets a flag so that the extreme eigenvalues
684: values will be calculated via a Lanczos or Arnoldi process as the linear
685: system is solved.
687: Collective on KSP
689: Input Parameters:
690: . ksp - iterative context obtained from KSPCreate()
692: Notes:
693: Currently this option is not valid for all iterative methods.
695: Level: advanced
697: .keywords: KSP, set, compute, eigenvalues
699: .seealso: KSPComputeEigenvalues(), KSPComputeEigenvaluesExplicitly()
700: @*/
701: int KSPSetComputeEigenvalues(KSP ksp)
702: {
705: ksp->calc_sings = PETSC_TRUE;
706: return(0);
707: }
709: /*@
710: KSPSetRhs - Sets the right-hand-side vector for the linear system to
711: be solved.
713: Collective on KSP and Vec
715: Input Parameters:
716: + ksp - iterative context obtained from KSPCreate()
717: - b - right-hand-side vector
719: Level: developer
721: .keywords: KSP, set, right-hand-side, rhs
723: .seealso: KSPGetRhs(), KSPSetSolution()
724: @*/
725: int KSPSetRhs(KSP ksp,Vec b)
726: {
730: ksp->vec_rhs = (b);
731: return(0);
732: }
734: /*@C
735: KSPGetRhs - Gets the right-hand-side vector for the linear system to
736: be solved.
738: Not Collective
740: Input Parameter:
741: . ksp - iterative context obtained from KSPCreate()
743: Output Parameter:
744: . r - right-hand-side vector
746: Level: developer
748: .keywords: KSP, get, right-hand-side, rhs
750: .seealso: KSPSetRhs(), KSPGetSolution()
751: @*/
752: int KSPGetRhs(KSP ksp,Vec *r)
753: {
756: *r = ksp->vec_rhs;
757: return(0);
758: }
760: /*@
761: KSPSetSolution - Sets the location of the solution for the
762: linear system to be solved.
764: Collective on KSP and Vec
766: Input Parameters:
767: + ksp - iterative context obtained from KSPCreate()
768: - x - solution vector
770: Level: developer
772: .keywords: KSP, set, solution
774: .seealso: KSPSetRhs(), KSPGetSolution()
775: @*/
776: int KSPSetSolution(KSP ksp,Vec x)
777: {
781: ksp->vec_sol = (x);
782: return(0);
783: }
785: /*@C
786: KSPGetSolution - Gets the location of the solution for the
787: linear system to be solved. Note that this may not be where the solution
788: is stored during the iterative process; see KSPBuildSolution().
790: Not Collective
792: Input Parameters:
793: . ksp - iterative context obtained from KSPCreate()
795: Output Parameters:
796: . v - solution vector
798: Level: developer
800: .keywords: KSP, get, solution
802: .seealso: KSPGetRhs(), KSPSetSolution(), KSPBuildSolution()
803: @*/
804: int KSPGetSolution(KSP ksp,Vec *v)
805: {
808: *v = ksp->vec_sol;
809: return(0);
810: }
812: /*@
813: KSPSetPC - Sets the preconditioner to be used to calculate the
814: application of the preconditioner on a vector.
816: Collective on KSP
818: Input Parameters:
819: + ksp - iterative context obtained from KSPCreate()
820: - B - the preconditioner object
822: Notes:
823: Use KSPGetPC() to retrieve the preconditioner context (for example,
824: to free it at the end of the computations).
826: Level: developer
828: .keywords: KSP, set, precondition, Binv
830: .seealso: KSPGetPC()
831: @*/
832: int KSPSetPC(KSP ksp,PC B)
833: {
838: ksp->B = B;
839: return(0);
840: }
842: /*@C
843: KSPGetPC - Returns a pointer to the preconditioner context
844: set with KSPSetPC().
846: Not Collective
848: Input Parameters:
849: . ksp - iterative context obtained from KSPCreate()
851: Output Parameter:
852: . B - preconditioner context
854: Level: developer
856: .keywords: KSP, get, preconditioner, Binv
858: .seealso: KSPSetPC()
859: @*/
860: int KSPGetPC(KSP ksp,PC *B)
861: {
864: *B = ksp->B;
865: return(0);
866: }
868: /*@C
869: KSPSetMonitor - Sets an ADDITIONAL function to be called at every iteration to monitor
870: the residual/error etc.
871:
872: Collective on KSP
874: Input Parameters:
875: + ksp - iterative context obtained from KSPCreate()
876: . monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring
877: . mctx - [optional] context for private data for the
878: monitor routine (use PETSC_NULL if no context is desired)
879: - monitordestroy - [optional] routine that frees monitor context
880: (may be PETSC_NULL)
882: Calling Sequence of monitor:
883: $ monitor (KSP ksp, int it, PetscReal rnorm, void *mctx)
885: + ksp - iterative context obtained from KSPCreate()
886: . it - iteration number
887: . rnorm - (estimated) 2-norm of (preconditioned) residual
888: - mctx - optional monitoring context, as set by KSPSetMonitor()
890: Options Database Keys:
891: + -ksp_monitor - sets KSPDefaultMonitor()
892: . -ksp_truemonitor - sets KSPTrueMonitor()
893: . -ksp_xmonitor - sets line graph monitor,
894: uses KSPLGMonitorCreate()
895: . -ksp_xtruemonitor - sets line graph monitor,
896: uses KSPLGMonitorCreate()
897: . -ksp_singmonitor - sets KSPSingularValueMonitor()
898: - -ksp_cancelmonitors - cancels all monitors that have
899: been hardwired into a code by
900: calls to KSPSetMonitor(), but
901: does not cancel those set via
902: the options database.
904: Notes:
905: The default is to do nothing. To print the residual, or preconditioned
906: residual if KSPSetUsePreconditionedResidual() was called, use
907: KSPDefaultMonitor() as the monitoring routine, with a null monitoring
908: context.
910: Several different monitoring routines may be set by calling
911: KSPSetMonitor() multiple times; all will be called in the
912: order in which they were set.
914: Level: beginner
916: .keywords: KSP, set, monitor
918: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPClearMonitor()
919: @*/
920: int KSPSetMonitor(KSP ksp,int (*monitor)(KSP,int,PetscReal,void*),void *mctx,int (*monitordestroy)(void*))
921: {
924: if (ksp->numbermonitors >= MAXKSPMONITORS) {
925: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many KSP monitors set");
926: }
927: ksp->monitor[ksp->numbermonitors] = monitor;
928: ksp->monitordestroy[ksp->numbermonitors] = monitordestroy;
929: ksp->monitorcontext[ksp->numbermonitors++] = (void*)mctx;
930: return(0);
931: }
933: /*@
934: KSPClearMonitor - Clears all monitors for a KSP object.
936: Collective on KSP
938: Input Parameters:
939: . ksp - iterative context obtained from KSPCreate()
941: Options Database Key:
942: . -ksp_cancelmonitors - Cancels all monitors that have
943: been hardwired into a code by calls to KSPSetMonitor(),
944: but does not cancel those set via the options database.
946: Level: intermediate
948: .keywords: KSP, set, monitor
950: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPSetMonitor()
951: @*/
952: int KSPClearMonitor(KSP ksp)
953: {
956: ksp->numbermonitors = 0;
957: return(0);
958: }
960: /*@C
961: KSPGetMonitorContext - Gets the monitoring context, as set by
962: KSPSetMonitor() for the FIRST monitor only.
964: Not Collective
966: Input Parameter:
967: . ksp - iterative context obtained from KSPCreate()
969: Output Parameter:
970: . ctx - monitoring context
972: Level: intermediate
974: .keywords: KSP, get, monitor, context
976: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate()
977: @*/
978: int KSPGetMonitorContext(KSP ksp,void **ctx)
979: {
982: *ctx = (ksp->monitorcontext[0]);
983: return(0);
984: }
986: /*@
987: KSPSetResidualHistory - Sets the array used to hold the residual history.
988: If set, this array will contain the residual norms computed at each
989: iteration of the solver.
991: Not Collective
993: Input Parameters:
994: + ksp - iterative context obtained from KSPCreate()
995: . a - array to hold history
996: . na - size of a
997: - reset - PETSC_TRUE indicates the history counter is reset to zero
998: for each new linear solve
1000: Level: advanced
1002: .keywords: KSP, set, residual, history, norm
1004: .seealso: KSPGetResidualHistory()
1006: @*/
1007: int KSPSetResidualHistory(KSP ksp,PetscReal *a,int na,PetscTruth reset)
1008: {
1013: if (na != PETSC_DECIDE && a != PETSC_NULL) {
1014: ksp->res_hist = a;
1015: ksp->res_hist_max = na;
1016: } else {
1017: ksp->res_hist_max = 1000;
1018: PetscMalloc(ksp->res_hist_max*sizeof(PetscReal),&ksp->res_hist);
1019: }
1020: ksp->res_hist_len = 0;
1021: ksp->res_hist_reset = reset;
1024: return(0);
1025: }
1027: /*@C
1028: KSPGetResidualHistory - Gets the array used to hold the residual history
1029: and the number of residuals it contains.
1031: Not Collective
1033: Input Parameter:
1034: . ksp - iterative context obtained from KSPCreate()
1036: Output Parameters:
1037: + a - pointer to array to hold history (or PETSC_NULL)
1038: - na - number of used entries in a (or PETSC_NULL)
1040: Level: advanced
1042: Notes:
1043: Can only call after a KSPSetResidualHistory() otherwise returns 0.
1045: The Fortran version of this routine has a calling sequence
1046: $ call KSPGetResidualHistory(KSP ksp, integer na, integer ierr)
1048: .keywords: KSP, get, residual, history, norm
1050: .seealso: KSPGetResidualHistory()
1052: @*/
1053: int KSPGetResidualHistory(KSP ksp,PetscReal **a,int *na)
1054: {
1057: if (a) *a = ksp->res_hist;
1058: if (na) *na = ksp->res_hist_len;
1059: return(0);
1060: }
1062: /*@C
1063: KSPSetConvergenceTest - Sets the function to be used to determine
1064: convergence.
1066: Collective on KSP
1068: Input Parameters:
1069: + ksp - iterative context obtained from KSPCreate()
1070: . converge - pointer to int function
1071: - cctx - context for private data for the convergence routine (may be null)
1073: Calling sequence of converge:
1074: $ converge (KSP ksp, int it, PetscReal rnorm, KSPConvergedReason *reason,void *mctx)
1076: + ksp - iterative context obtained from KSPCreate()
1077: . it - iteration number
1078: . rnorm - (estimated) 2-norm of (preconditioned) residual
1079: . reason - the reason why it has converged or diverged
1080: - cctx - optional convergence context, as set by KSPSetConvergenceTest()
1083: Notes:
1084: Must be called after the KSP type has been set so put this after
1085: a call to KSPSetType(), KSPSetTypeFromOptions() or KSPSetFromOptions().
1087: The default convergence test, KSPDefaultConverged(), aborts if the
1088: residual grows to more than 10000 times the initial residual.
1090: The default is a combination of relative and absolute tolerances.
1091: The residual value that is tested may be an approximation; routines
1092: that need exact values should compute them.
1094: Level: advanced
1096: .keywords: KSP, set, convergence, test, context
1098: .seealso: KSPDefaultConverged(), KSPGetConvergenceContext()
1099: @*/
1100: int KSPSetConvergenceTest(KSP ksp,int (*converge)(KSP,int,PetscReal,KSPConvergedReason*,void*),void *cctx)
1101: {
1104: ksp->converged = converge;
1105: ksp->cnvP = (void*)cctx;
1106: return(0);
1107: }
1109: /*@C
1110: KSPGetConvergenceContext - Gets the convergence context set with
1111: KSPSetConvergenceTest().
1113: Not Collective
1115: Input Parameter:
1116: . ksp - iterative context obtained from KSPCreate()
1118: Output Parameter:
1119: . ctx - monitoring context
1121: Level: advanced
1123: .keywords: KSP, get, convergence, test, context
1125: .seealso: KSPDefaultConverged(), KSPSetConvergenceTest()
1126: @*/
1127: int KSPGetConvergenceContext(KSP ksp,void **ctx)
1128: {
1131: *ctx = ksp->cnvP;
1132: return(0);
1133: }
1135: /*@C
1136: KSPBuildSolution - Builds the approximate solution in a vector provided.
1137: This routine is NOT commonly needed (see SLESSolve()).
1139: Collective on KSP
1141: Input Parameter:
1142: . ctx - iterative context obtained from KSPCreate()
1144: Output Parameter:
1145: Provide exactly one of
1146: + v - location to stash solution.
1147: - V - the solution is returned in this location. This vector is created
1148: internally. This vector should NOT be destroyed by the user with
1149: VecDestroy().
1151: Notes:
1152: This routine can be used in one of two ways
1153: .vb
1154: KSPBuildSolution(ksp,PETSC_NULL,&V);
1155: or
1156: KSPBuildSolution(ksp,v,PETSC_NULL);
1157: .ve
1158: In the first case an internal vector is allocated to store the solution
1159: (the user cannot destroy this vector). In the second case the solution
1160: is generated in the vector that the user provides. Note that for certain
1161: methods, such as KSPCG, the second case requires a copy of the solution,
1162: while in the first case the call is essentially free since it simply
1163: returns the vector where the solution already is stored.
1165: Level: advanced
1167: .keywords: KSP, build, solution
1169: .seealso: KSPGetSolution(), KSPBuildResidual()
1170: @*/
1171: int KSPBuildSolution(KSP ksp,Vec v,Vec *V)
1172: {
1177: if (!V && !v) SETERRQ(PETSC_ERR_ARG_WRONG,"Must provide either v or V");
1178: if (!V) V = &v;
1179: (*ksp->ops->buildsolution)(ksp,v,V);
1180: return(0);
1181: }
1183: /*@C
1184: KSPBuildResidual - Builds the residual in a vector provided.
1186: Collective on KSP
1188: Input Parameter:
1189: . ksp - iterative context obtained from KSPCreate()
1191: Output Parameters:
1192: + v - optional location to stash residual. If v is not provided,
1193: then a location is generated.
1194: . t - work vector. If not provided then one is generated.
1195: - V - the residual
1197: Notes:
1198: Regardless of whether or not v is provided, the residual is
1199: returned in V.
1201: Level: advanced
1203: .keywords: KSP, build, residual
1205: .seealso: KSPBuildSolution()
1206: @*/
1207: int KSPBuildResidual(KSP ksp,Vec t,Vec v,Vec *V)
1208: {
1209: int flag = 0,ierr;
1210: Vec w = v,tt = t;
1214: if (!w) {
1215: VecDuplicate(ksp->vec_rhs,&w);
1216: PetscLogObjectParent((PetscObject)ksp,w);
1217: }
1218: if (!tt) {
1219: VecDuplicate(ksp->vec_rhs,&tt); flag = 1;
1220: PetscLogObjectParent((PetscObject)ksp,tt);
1221: }
1222: (*ksp->ops->buildresidual)(ksp,tt,w,V);
1223: if (flag) {VecDestroy(tt);}
1224: return(0);
1225: }