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