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: }