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
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: PetscScalar zero = 0.0;
203: if (!ksp->setupcalled){ KSPSetUp(ksp);}
204: if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
205: if (ksp->guess_knoll) {
206: ierr = PCApply(ksp->B,ksp->vec_rhs,ksp->vec_sol);
207: ksp->guess_zero = PETSC_FALSE;
208: }
210: /* reset the residual history list if requested */
211: if (ksp->res_hist_reset) ksp->res_hist_len = 0;
213: ksp->transpose_solve = PETSC_FALSE;
214: (*ksp->ops->solve)(ksp,&nits);
215: if (!ksp->reason) {
216: SETERRQ(1,"Internal error, solver returned without setting converged reason");
217: }
218: if (its) *its = nits;
220: MPI_Comm_rank(ksp->comm,&rank);
222: PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues",&flag1);
223: PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues",&flag2);
224: if (flag1 || flag2) {
225: int n = nits + 2,i,neig;
226: PetscReal *r,*c;
228: if (!n) {
229: PetscPrintf(ksp->comm,"Zero iterations in solver, cannot approximate any eigenvaluesn");
230: } else {
231: PetscMalloc(2*n*sizeof(PetscReal),&r);
232: c = r + n;
233: KSPComputeEigenvalues(ksp,n,r,c,&neig);
234: if (flag1) {
235: PetscPrintf(ksp->comm,"Iteratively computed eigenvaluesn");
236: for (i=0; i<neig; i++) {
237: if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gin",r[i],c[i]);}
238: else {PetscPrintf(ksp->comm,"%g - %gin",r[i],-c[i]);}
239: }
240: }
241: if (flag2 && !rank) {
242: PetscViewer viewer;
243: PetscDraw draw;
244: PetscDrawSP drawsp;
246: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Iteratively Computed Eigenvalues",
247: PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
248: PetscViewerDrawGetDraw(viewer,0,&draw);
249: PetscDrawSPCreate(draw,1,&drawsp);
250: for (i=0; i<neig; i++) {
251: PetscDrawSPAddPoint(drawsp,r+i,c+i);
252: }
253: PetscDrawSPDraw(drawsp);
254: PetscDrawSPDestroy(drawsp);
255: PetscViewerDestroy(viewer);
256: }
257: PetscFree(r);
258: }
259: }
261: PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues_explicitly",&flag1);
262: PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues_explicitly",&flag2);
263: if (flag1 || flag2) {
264: int n,i;
265: PetscReal *r,*c;
266: VecGetSize(ksp->vec_sol,&n);
267: PetscMalloc(2*n*sizeof(PetscReal),&r);
268: c = r + n;
269: KSPComputeEigenvaluesExplicitly(ksp,n,r,c);
270: if (flag1) {
271: PetscPrintf(ksp->comm,"Explicitly computed eigenvaluesn");
272: for (i=0; i<n; i++) {
273: if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gin",r[i],c[i]);}
274: else {PetscPrintf(ksp->comm,"%g - %gin",r[i],-c[i]);}
275: }
276: }
277: if (flag2 && !rank) {
278: PetscViewer viewer;
279: PetscDraw draw;
280: PetscDrawSP drawsp;
282: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Explicitly Computed Eigenvalues",0,320,300,300,&viewer);
283: PetscViewerDrawGetDraw(viewer,0,&draw);
284: PetscDrawSPCreate(draw,1,&drawsp);
285: for (i=0; i<n; i++) {
286: PetscDrawSPAddPoint(drawsp,r+i,c+i);
287: }
288: PetscDrawSPDraw(drawsp);
289: PetscDrawSPDestroy(drawsp);
290: PetscViewerDestroy(viewer);
291: }
292: PetscFree(r);
293: }
295: PetscOptionsHasName(ksp->prefix,"-ksp_view_operator",&flag2);
296: if (flag2) {
297: Mat A,B;
298: PCGetOperators(ksp->B,&A,PETSC_NULL,PETSC_NULL);
299: MatComputeExplicitOperator(A,&B);
300: PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(ksp->comm),PETSC_VIEWER_ASCII_MATLAB);
301: MatView(B,PETSC_VIEWER_STDOUT_(ksp->comm));
302: PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(ksp->comm));
303: MatDestroy(B);
304: }
305: return(0);
306: }
308: /*@
309: KSPSolveTranspose - Solves the transpose of a linear system. Usually
310: accessed through SLESSolveTranspose().
312: Collective on KSP
314: Input Parameter:
315: . ksp - iterative context obtained from KSPCreate()
317: Output Parameter:
318: . its - number of iterations required
320: Notes:
321: On return, the parameter "its" contains either the iteration
322: number at which convergence was successfully reached, or the
323: negative of the iteration at which divergence or breakdown was detected.
325: Currently only supported by KSPType of KSPPREONLY. This routine is usally
326: only used internally by the BiCG solver on the subblocks in BJacobi and ASM.
328: Level: developer
330: .keywords: KSP, solve, linear system
332: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPDefaultConverged(),
333: SLESSolve(), SLESGetKSP()
334: @*/
335: int KSPSolveTranspose(KSP ksp,int *its)
336: {
337: int ierr;
338: PetscScalar zero = 0.0;
343: if (!ksp->setupcalled){ KSPSetUp(ksp);}
344: if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
345: ksp->transpose_solve = PETSC_TRUE;
346: (*ksp->ops->solve)(ksp,its);
347: return(0);
348: }
350: /*@C
351: KSPDestroy - Destroys KSP context.
353: Collective on KSP
355: Input Parameter:
356: . ksp - iterative context obtained from KSPCreate()
358: Level: developer
360: .keywords: KSP, destroy
362: .seealso: KSPCreate(), KSPSetUp(), KSPSolve()
363: @*/
364: int KSPDestroy(KSP ksp)
365: {
366: int i,ierr;
370: if (--ksp->refct > 0) return(0);
372: /* if memory was published with AMS then destroy it */
373: PetscObjectDepublish(ksp);
375: if (ksp->ops->destroy) {
376: (*ksp->ops->destroy)(ksp);
377: }
378: for (i=0; i<ksp->numbermonitors; i++) {
379: if (ksp->monitordestroy[i]) {
380: (*ksp->monitordestroy[i])(ksp->monitorcontext[i]);
381: }
382: }
383: PetscLogObjectDestroy(ksp);
384: PetscHeaderDestroy(ksp);
385: return(0);
386: }
388: /*@
389: KSPSetPreconditionerSide - Sets the preconditioning side.
391: Collective on KSP
393: Input Parameter:
394: . ksp - iterative context obtained from KSPCreate()
396: Output Parameter:
397: . side - the preconditioning side, where side is one of
398: .vb
399: PC_LEFT - left preconditioning (default)
400: PC_RIGHT - right preconditioning
401: PC_SYMMETRIC - symmetric preconditioning
402: .ve
404: Options Database Keys:
405: + -ksp_left_pc - Sets left preconditioning
406: . -ksp_right_pc - Sets right preconditioning
407: - -ksp_symmetric_pc - Sets symmetric preconditioning
409: Notes:
410: Left preconditioning is used by default. Symmetric preconditioning is
411: currently available only for the KSPQCG method. Note, however, that
412: symmetric preconditioning can be emulated by using either right or left
413: preconditioning and a pre or post processing step.
415: Level: intermediate
417: .keywords: KSP, set, right, left, symmetric, side, preconditioner, flag
419: .seealso: KSPGetPreconditionerSide()
420: @*/
421: int KSPSetPreconditionerSide(KSP ksp,PCSide side)
422: {
425: ksp->pc_side = side;
426: return(0);
427: }
429: /*@C
430: KSPGetPreconditionerSide - Gets the preconditioning side.
432: Not Collective
434: Input Parameter:
435: . ksp - iterative context obtained from KSPCreate()
437: Output Parameter:
438: . side - the preconditioning side, where side is one of
439: .vb
440: PC_LEFT - left preconditioning (default)
441: PC_RIGHT - right preconditioning
442: PC_SYMMETRIC - symmetric preconditioning
443: .ve
445: Level: intermediate
447: .keywords: KSP, get, right, left, symmetric, side, preconditioner, flag
449: .seealso: KSPSetPreconditionerSide()
450: @*/
451: int KSPGetPreconditionerSide(KSP ksp,PCSide *side)
452: {
455: *side = ksp->pc_side;
456: return(0);
457: }
459: /*@
460: KSPGetTolerances - Gets the relative, absolute, divergence, and maximum
461: iteration tolerances used by the default KSP convergence tests.
463: Not Collective
465: Input Parameter:
466: . ksp - the Krylov subspace context
467:
468: Output Parameters:
469: + rtol - the relative convergence tolerance
470: . atol - the absolute convergence tolerance
471: . dtol - the divergence tolerance
472: - maxits - maximum number of iterations
474: Notes:
475: The user can specify PETSC_NULL for any parameter that is not needed.
477: Level: intermediate
479: .keywords: KSP, get, tolerance, absolute, relative, divergence, convergence,
480: maximum, iterations
482: .seealso: KSPSetTolerances()
483: @*/
484: int KSPGetTolerances(KSP ksp,PetscReal *rtol,PetscReal *atol,PetscReal *dtol,int *maxits)
485: {
488: if (atol) *atol = ksp->atol;
489: if (rtol) *rtol = ksp->rtol;
490: if (dtol) *dtol = ksp->divtol;
491: if (maxits) *maxits = ksp->max_it;
492: return(0);
493: }
495: /*@
496: KSPSetTolerances - Sets the relative, absolute, divergence, and maximum
497: iteration tolerances used by the default KSP convergence testers.
499: Collective on KSP
501: Input Parameters:
502: + ksp - the Krylov subspace context
503: . rtol - the relative convergence tolerance
504: (relative decrease in the residual norm)
505: . atol - the absolute convergence tolerance
506: (absolute size of the residual norm)
507: . dtol - the divergence tolerance
508: (amount residual can increase before KSPDefaultConverged()
509: concludes that the method is diverging)
510: - maxits - maximum number of iterations to use
512: Options Database Keys:
513: + -ksp_atol <atol> - Sets atol
514: . -ksp_rtol <rtol> - Sets rtol
515: . -ksp_divtol <dtol> - Sets dtol
516: - -ksp_max_it <maxits> - Sets maxits
518: Notes:
519: Use PETSC_DEFAULT to retain the default value of any of the tolerances.
521: See KSPDefaultConverged() for details on the use of these parameters
522: in the default convergence test. See also KSPSetConvergenceTest()
523: for setting user-defined stopping criteria.
525: Level: intermediate
527: .keywords: KSP, set, tolerance, absolute, relative, divergence,
528: convergence, maximum, iterations
530: .seealso: KSPGetTolerances(), KSPDefaultConverged(), KSPSetConvergenceTest()
531: @*/
532: int KSPSetTolerances(KSP ksp,PetscReal rtol,PetscReal atol,PetscReal dtol,int maxits)
533: {
536: if (atol != PETSC_DEFAULT) ksp->atol = atol;
537: if (rtol != PETSC_DEFAULT) ksp->rtol = rtol;
538: if (dtol != PETSC_DEFAULT) ksp->divtol = dtol;
539: if (maxits != PETSC_DEFAULT) ksp->max_it = maxits;
540: return(0);
541: }
543: /*@
544: KSPSetInitialGuessNonzero - Tells the iterative solver that the
545: initial guess is nonzero; otherwise KSP assumes the initial guess
546: is to be zero (and thus zeros it out before solving).
548: Collective on KSP
550: Input Parameters:
551: + ksp - iterative context obtained from SLESGetKSP() or KSPCreate()
552: - flg - PETSC_TRUE indicates the guess is non-zero, PETSC_FALSE indicates the guess is zero
554: Level: beginner
556: Notes:
557: If this is not called the X vector is zeroed in the call to
558: SLESSolve() (or KSPSolve()).
560: .keywords: KSP, set, initial guess, nonzero
562: .seealso: KSPGetIntialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
563: @*/
564: int KSPSetInitialGuessNonzero(KSP ksp,PetscTruth flg)
565: {
567: ksp->guess_zero = (PetscTruth)!(int)flg;
568: return(0);
569: }
571: /*@
572: KSPGetInitialGuessNonzero - Determines whether the KSP solver is using
573: a zero initial guess.
575: Not Collective
577: Input Parameter:
578: . ksp - iterative context obtained from KSPCreate()
580: Output Parameter:
581: . flag - PETSC_TRUE if guess is nonzero, else PETSC_FALSE
583: Level: intermediate
585: .keywords: KSP, set, initial guess, nonzero
587: .seealso: KSPSetIntialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
588: @*/
589: int KSPGetInitialGuessNonzero(KSP ksp,PetscTruth *flag)
590: {
592: if (ksp->guess_zero) *flag = PETSC_FALSE;
593: else *flag = PETSC_TRUE;
594: return(0);
595: }
597: /*@
598: KSPSetInitialGuessKnoll - Tells the iterative solver to use PCApply(pc,b,..) to compute the initial guess (The Knoll trick)
600: Collective on KSP
602: Input Parameters:
603: + ksp - iterative context obtained from SLESGetKSP() or KSPCreate()
604: - flg - PETSC_TRUE or PETSC_FALSE
606: Level: advanced
609: .keywords: KSP, set, initial guess, nonzero
611: .seealso: KSPGetIntialGuessKnoll(), KSPSetInitialGuessNonzero(), KSPGetInitialGuessNonzero()
612: @*/
613: int KSPSetInitialGuessKnoll(KSP ksp,PetscTruth flg)
614: {
616: ksp->guess_knoll = flg;
617: return(0);
618: }
620: /*@
621: KSPGetInitialGuessKnoll - Determines whether the KSP solver is using the Knoll trick (using PCApply(pc,b,...) to compute
622: the initial guess
624: Not Collective
626: Input Parameter:
627: . ksp - iterative context obtained from KSPCreate()
629: Output Parameter:
630: . flag - PETSC_TRUE if using Knoll trick, else PETSC_FALSE
632: Level: advanced
634: .keywords: KSP, set, initial guess, nonzero
636: .seealso: KSPSetIntialGuessKnoll(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
637: @*/
638: int KSPGetInitialGuessKnoll(KSP ksp,PetscTruth *flag)
639: {
641: *flag = ksp->guess_knoll;
642: return(0);
643: }
645: /*@
646: KSPSetComputeSingularValues - Sets a flag so that the extreme singular
647: values will be calculated via a Lanczos or Arnoldi process as the linear
648: system is solved.
650: Collective on KSP
652: Input Parameters:
653: + ksp - iterative context obtained from KSPCreate()
654: - flg - PETSC_TRUE or PETSC_FALSE
656: Options Database Key:
657: . -ksp_singmonitor - Activates KSPSetComputeSingularValues()
659: Notes:
660: Currently this option is not valid for all iterative methods.
662: Many users may just want to use the monitoring routine
663: KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
664: to print the singular values at each iteration of the linear solve.
666: Level: advanced
668: .keywords: KSP, set, compute, singular values
670: .seealso: KSPComputeExtremeSingularValues(), KSPSingularValueMonitor()
671: @*/
672: int KSPSetComputeSingularValues(KSP ksp,PetscTruth flg)
673: {
676: ksp->calc_sings = flg;
677: return(0);
678: }
680: /*@
681: KSPSetComputeEigenvalues - Sets a flag so that the extreme eigenvalues
682: values will be calculated via a Lanczos or Arnoldi process as the linear
683: system is solved.
685: Collective on KSP
687: Input Parameters:
688: + ksp - iterative context obtained from KSPCreate()
689: - flg - PETSC_TRUE or PETSC_FALSE
691: Notes:
692: Currently this option is not valid for all iterative methods.
694: Level: advanced
696: .keywords: KSP, set, compute, eigenvalues
698: .seealso: KSPComputeEigenvalues(), KSPComputeEigenvaluesExplicitly()
699: @*/
700: int KSPSetComputeEigenvalues(KSP ksp,PetscTruth flg)
701: {
704: ksp->calc_sings = flg;
705: return(0);
706: }
708: /*@
709: KSPSetRhs - Sets the right-hand-side vector for the linear system to
710: be solved.
712: Collective on KSP and Vec
714: Input Parameters:
715: + ksp - iterative context obtained from KSPCreate()
716: - b - right-hand-side vector
718: Level: developer
720: .keywords: KSP, set, right-hand-side, rhs
722: .seealso: KSPGetRhs(), KSPSetSolution()
723: @*/
724: int KSPSetRhs(KSP ksp,Vec b)
725: {
729: ksp->vec_rhs = (b);
730: return(0);
731: }
733: /*@C
734: KSPGetRhs - Gets the right-hand-side vector for the linear system to
735: be solved.
737: Not Collective
739: Input Parameter:
740: . ksp - iterative context obtained from KSPCreate()
742: Output Parameter:
743: . r - right-hand-side vector
745: Level: developer
747: .keywords: KSP, get, right-hand-side, rhs
749: .seealso: KSPSetRhs(), KSPGetSolution()
750: @*/
751: int KSPGetRhs(KSP ksp,Vec *r)
752: {
755: *r = ksp->vec_rhs;
756: return(0);
757: }
759: /*@
760: KSPSetSolution - Sets the location of the solution for the
761: linear system to be solved.
763: Collective on KSP and Vec
765: Input Parameters:
766: + ksp - iterative context obtained from KSPCreate()
767: - x - solution vector
769: Level: developer
771: .keywords: KSP, set, solution
773: .seealso: KSPSetRhs(), KSPGetSolution()
774: @*/
775: int KSPSetSolution(KSP ksp,Vec x)
776: {
780: ksp->vec_sol = (x);
781: return(0);
782: }
784: /*@C
785: KSPGetSolution - Gets the location of the solution for the
786: linear system to be solved. Note that this may not be where the solution
787: is stored during the iterative process; see KSPBuildSolution().
789: Not Collective
791: Input Parameters:
792: . ksp - iterative context obtained from KSPCreate()
794: Output Parameters:
795: . v - solution vector
797: Level: developer
799: .keywords: KSP, get, solution
801: .seealso: KSPGetRhs(), KSPSetSolution(), KSPBuildSolution()
802: @*/
803: int KSPGetSolution(KSP ksp,Vec *v)
804: {
807: *v = ksp->vec_sol;
808: return(0);
809: }
811: /*@
812: KSPSetPC - Sets the preconditioner to be used to calculate the
813: application of the preconditioner on a vector.
815: Collective on KSP
817: Input Parameters:
818: + ksp - iterative context obtained from KSPCreate()
819: - B - the preconditioner object
821: Notes:
822: Use KSPGetPC() to retrieve the preconditioner context (for example,
823: to free it at the end of the computations).
825: Level: developer
827: .keywords: KSP, set, precondition, Binv
829: .seealso: KSPGetPC()
830: @*/
831: int KSPSetPC(KSP ksp,PC B)
832: {
837: ksp->B = B;
838: return(0);
839: }
841: /*@C
842: KSPGetPC - Returns a pointer to the preconditioner context
843: set with KSPSetPC().
845: Not Collective
847: Input Parameters:
848: . ksp - iterative context obtained from KSPCreate()
850: Output Parameter:
851: . B - preconditioner context
853: Level: developer
855: .keywords: KSP, get, preconditioner, Binv
857: .seealso: KSPSetPC()
858: @*/
859: int KSPGetPC(KSP ksp,PC *B)
860: {
863: *B = ksp->B;
864: return(0);
865: }
867: /*@C
868: KSPSetMonitor - Sets an ADDITIONAL function to be called at every iteration to monitor
869: the residual/error etc.
870:
871: Collective on KSP
873: Input Parameters:
874: + ksp - iterative context obtained from KSPCreate()
875: . monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring
876: . mctx - [optional] context for private data for the
877: monitor routine (use PETSC_NULL if no context is desired)
878: - monitordestroy - [optional] routine that frees monitor context
879: (may be PETSC_NULL)
881: Calling Sequence of monitor:
882: $ monitor (KSP ksp, int it, PetscReal rnorm, void *mctx)
884: + ksp - iterative context obtained from KSPCreate()
885: . it - iteration number
886: . rnorm - (estimated) 2-norm of (preconditioned) residual
887: - mctx - optional monitoring context, as set by KSPSetMonitor()
889: Options Database Keys:
890: + -ksp_monitor - sets KSPDefaultMonitor()
891: . -ksp_truemonitor - sets KSPTrueMonitor()
892: . -ksp_xmonitor - sets line graph monitor,
893: uses KSPLGMonitorCreate()
894: . -ksp_xtruemonitor - sets line graph monitor,
895: uses KSPLGMonitorCreate()
896: . -ksp_singmonitor - sets KSPSingularValueMonitor()
897: - -ksp_cancelmonitors - cancels all monitors that have
898: been hardwired into a code by
899: calls to KSPSetMonitor(), but
900: does not cancel those set via
901: the options database.
903: Notes:
904: The default is to do nothing. To print the residual, or preconditioned
905: residual if KSPSetNormType(ksp,KSP_PRECONDITIONED_NORM) was called, use
906: KSPDefaultMonitor() as the monitoring routine, with a null monitoring
907: context.
909: Several different monitoring routines may be set by calling
910: KSPSetMonitor() multiple times; all will be called in the
911: order in which they were set.
913: Level: beginner
915: .keywords: KSP, set, monitor
917: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPClearMonitor()
918: @*/
919: int KSPSetMonitor(KSP ksp,int (*monitor)(KSP,int,PetscReal,void*),void *mctx,int (*monitordestroy)(void*))
920: {
923: if (ksp->numbermonitors >= MAXKSPMONITORS) {
924: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many KSP monitors set");
925: }
926: ksp->monitor[ksp->numbermonitors] = monitor;
927: ksp->monitordestroy[ksp->numbermonitors] = monitordestroy;
928: ksp->monitorcontext[ksp->numbermonitors++] = (void*)mctx;
929: return(0);
930: }
932: /*@
933: KSPClearMonitor - Clears all monitors for a KSP object.
935: Collective on KSP
937: Input Parameters:
938: . ksp - iterative context obtained from KSPCreate()
940: Options Database Key:
941: . -ksp_cancelmonitors - Cancels all monitors that have
942: been hardwired into a code by calls to KSPSetMonitor(),
943: but does not cancel those set via the options database.
945: Level: intermediate
947: .keywords: KSP, set, monitor
949: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPSetMonitor()
950: @*/
951: int KSPClearMonitor(KSP ksp)
952: {
955: ksp->numbermonitors = 0;
956: return(0);
957: }
959: /*@C
960: KSPGetMonitorContext - Gets the monitoring context, as set by
961: KSPSetMonitor() for the FIRST monitor only.
963: Not Collective
965: Input Parameter:
966: . ksp - iterative context obtained from KSPCreate()
968: Output Parameter:
969: . ctx - monitoring context
971: Level: intermediate
973: .keywords: KSP, get, monitor, context
975: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate()
976: @*/
977: int KSPGetMonitorContext(KSP ksp,void **ctx)
978: {
981: *ctx = (ksp->monitorcontext[0]);
982: return(0);
983: }
985: /*@
986: KSPSetResidualHistory - Sets the array used to hold the residual history.
987: If set, this array will contain the residual norms computed at each
988: iteration of the solver.
990: Not Collective
992: Input Parameters:
993: + ksp - iterative context obtained from KSPCreate()
994: . a - array to hold history
995: . na - size of a
996: - reset - PETSC_TRUE indicates the history counter is reset to zero
997: for each new linear solve
999: Level: advanced
1001: .keywords: KSP, set, residual, history, norm
1003: .seealso: KSPGetResidualHistory()
1005: @*/
1006: int KSPSetResidualHistory(KSP ksp,PetscReal *a,int na,PetscTruth reset)
1007: {
1012: if (na != PETSC_DECIDE && a != PETSC_NULL) {
1013: ksp->res_hist = a;
1014: ksp->res_hist_max = na;
1015: } else {
1016: ksp->res_hist_max = 1000;
1017: PetscMalloc(ksp->res_hist_max*sizeof(PetscReal),&ksp->res_hist);
1018: }
1019: ksp->res_hist_len = 0;
1020: ksp->res_hist_reset = reset;
1023: return(0);
1024: }
1026: /*@C
1027: KSPGetResidualHistory - Gets the array used to hold the residual history
1028: and the number of residuals it contains.
1030: Not Collective
1032: Input Parameter:
1033: . ksp - iterative context obtained from KSPCreate()
1035: Output Parameters:
1036: + a - pointer to array to hold history (or PETSC_NULL)
1037: - na - number of used entries in a (or PETSC_NULL)
1039: Level: advanced
1041: Notes:
1042: Can only call after a KSPSetResidualHistory() otherwise returns 0.
1044: The Fortran version of this routine has a calling sequence
1045: $ call KSPGetResidualHistory(KSP ksp, integer na, integer ierr)
1047: .keywords: KSP, get, residual, history, norm
1049: .seealso: KSPGetResidualHistory()
1051: @*/
1052: int KSPGetResidualHistory(KSP ksp,PetscReal **a,int *na)
1053: {
1056: if (a) *a = ksp->res_hist;
1057: if (na) *na = ksp->res_hist_len;
1058: return(0);
1059: }
1061: /*@C
1062: KSPSetConvergenceTest - Sets the function to be used to determine
1063: convergence.
1065: Collective on KSP
1067: Input Parameters:
1068: + ksp - iterative context obtained from KSPCreate()
1069: . converge - pointer to int function
1070: - cctx - context for private data for the convergence routine (may be null)
1072: Calling sequence of converge:
1073: $ converge (KSP ksp, int it, PetscReal rnorm, KSPConvergedReason *reason,void *mctx)
1075: + ksp - iterative context obtained from KSPCreate()
1076: . it - iteration number
1077: . rnorm - (estimated) 2-norm of (preconditioned) residual
1078: . reason - the reason why it has converged or diverged
1079: - cctx - optional convergence context, as set by KSPSetConvergenceTest()
1082: Notes:
1083: Must be called after the KSP type has been set so put this after
1084: a call to KSPSetType(), or KSPSetFromOptions().
1086: The default convergence test, KSPDefaultConverged(), aborts if the
1087: residual grows to more than 10000 times the initial residual.
1089: The default is a combination of relative and absolute tolerances.
1090: The residual value that is tested may be an approximation; routines
1091: that need exact values should compute them.
1093: Level: advanced
1095: .keywords: KSP, set, convergence, test, context
1097: .seealso: KSPDefaultConverged(), KSPGetConvergenceContext()
1098: @*/
1099: int KSPSetConvergenceTest(KSP ksp,int (*converge)(KSP,int,PetscReal,KSPConvergedReason*,void*),void *cctx)
1100: {
1103: ksp->converged = converge;
1104: ksp->cnvP = (void*)cctx;
1105: return(0);
1106: }
1108: /*@C
1109: KSPGetConvergenceContext - Gets the convergence context set with
1110: KSPSetConvergenceTest().
1112: Not Collective
1114: Input Parameter:
1115: . ksp - iterative context obtained from KSPCreate()
1117: Output Parameter:
1118: . ctx - monitoring context
1120: Level: advanced
1122: .keywords: KSP, get, convergence, test, context
1124: .seealso: KSPDefaultConverged(), KSPSetConvergenceTest()
1125: @*/
1126: int KSPGetConvergenceContext(KSP ksp,void **ctx)
1127: {
1130: *ctx = ksp->cnvP;
1131: return(0);
1132: }
1134: /*@C
1135: KSPBuildSolution - Builds the approximate solution in a vector provided.
1136: This routine is NOT commonly needed (see SLESSolve()).
1138: Collective on KSP
1140: Input Parameter:
1141: . ctx - iterative context obtained from KSPCreate()
1143: Output Parameter:
1144: Provide exactly one of
1145: + v - location to stash solution.
1146: - V - the solution is returned in this location. This vector is created
1147: internally. This vector should NOT be destroyed by the user with
1148: VecDestroy().
1150: Notes:
1151: This routine can be used in one of two ways
1152: .vb
1153: KSPBuildSolution(ksp,PETSC_NULL,&V);
1154: or
1155: KSPBuildSolution(ksp,v,PETSC_NULL);
1156: .ve
1157: In the first case an internal vector is allocated to store the solution
1158: (the user cannot destroy this vector). In the second case the solution
1159: is generated in the vector that the user provides. Note that for certain
1160: methods, such as KSPCG, the second case requires a copy of the solution,
1161: while in the first case the call is essentially free since it simply
1162: returns the vector where the solution already is stored.
1164: Level: advanced
1166: .keywords: KSP, build, solution
1168: .seealso: KSPGetSolution(), KSPBuildResidual()
1169: @*/
1170: int KSPBuildSolution(KSP ksp,Vec v,Vec *V)
1171: {
1176: if (!V && !v) SETERRQ(PETSC_ERR_ARG_WRONG,"Must provide either v or V");
1177: if (!V) V = &v;
1178: (*ksp->ops->buildsolution)(ksp,v,V);
1179: return(0);
1180: }
1182: /*@C
1183: KSPBuildResidual - Builds the residual in a vector provided.
1185: Collective on KSP
1187: Input Parameter:
1188: . ksp - iterative context obtained from KSPCreate()
1190: Output Parameters:
1191: + v - optional location to stash residual. If v is not provided,
1192: then a location is generated.
1193: . t - work vector. If not provided then one is generated.
1194: - V - the residual
1196: Notes:
1197: Regardless of whether or not v is provided, the residual is
1198: returned in V.
1200: Level: advanced
1202: .keywords: KSP, build, residual
1204: .seealso: KSPBuildSolution()
1205: @*/
1206: int KSPBuildResidual(KSP ksp,Vec t,Vec v,Vec *V)
1207: {
1208: int flag = 0,ierr;
1209: Vec w = v,tt = t;
1213: if (!w) {
1214: VecDuplicate(ksp->vec_rhs,&w);
1215: PetscLogObjectParent((PetscObject)ksp,w);
1216: }
1217: if (!tt) {
1218: VecDuplicate(ksp->vec_rhs,&tt); flag = 1;
1219: PetscLogObjectParent((PetscObject)ksp,tt);
1220: }
1221: (*ksp->ops->buildresidual)(ksp,tt,w,V);
1222: if (flag) {VecDestroy(tt);}
1223: return(0);
1224: }