Actual source code: itfunc.c

  1: /*$Id: itfunc.c,v 1.159 2001/08/07 03:03:45 balay Exp $*/
  2: /*
  3:       Interface KSP routines that the user calls.
  4: */

 6:  #include src/sles/ksp/kspimpl.h

 10: /*@
 11:    KSPComputeExtremeSingularValues - Computes the extreme singular values
 12:    for the preconditioned operator. Called after or during KSPSolve()
 13:    (SLESSolve()).

 15:    Not Collective

 17:    Input Parameter:
 18: .  ksp - iterative context obtained from KSPCreate()

 20:    Output Parameters:
 21: .  emin, emax - extreme singular values

 23:    Notes:
 24:    One must call KSPSetComputeSingularValues() before calling KSPSetUp() 
 25:    (or use the option -ksp_compute_eigenvalues) in order for this routine to work correctly.

 27:    Many users may just want to use the monitoring routine
 28:    KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
 29:    to print the extreme singular values at each iteration of the linear solve.

 31:    Level: advanced

 33: .keywords: KSP, compute, extreme, singular, values

 35: .seealso: KSPSetComputeSingularValues(), KSPSingularValueMonitor(), KSPComputeEigenvalues()
 36: @*/
 37: int KSPComputeExtremeSingularValues(KSP ksp,PetscReal *emax,PetscReal *emin)
 38: {

 45:   if (!ksp->calc_sings) {
 46:     SETERRQ(4,"Singular values not requested before KSPSetUp()");
 47:   }

 49:   if (ksp->ops->computeextremesingularvalues) {
 50:     (*ksp->ops->computeextremesingularvalues)(ksp,emax,emin);
 51:   } else {
 52:     *emin = -1.0;
 53:     *emax = -1.0;
 54:   }
 55:   return(0);
 56: }

 60: /*@
 61:    KSPComputeEigenvalues - Computes the extreme eigenvalues for the
 62:    preconditioned operator. Called after or during KSPSolve() (SLESSolve()).

 64:    Not Collective

 66:    Input Parameter:
 67: +  ksp - iterative context obtained from KSPCreate()
 68: -  n - size of arrays r and c. The number of eigenvalues computed (neig) will, in 
 69:        general, be less than this.

 71:    Output Parameters:
 72: +  r - real part of computed eigenvalues
 73: .  c - complex part of computed eigenvalues
 74: -  neig - number of eigenvalues computed (will be less than or equal to n)

 76:    Options Database Keys:
 77: +  -ksp_compute_eigenvalues - Prints eigenvalues to stdout
 78: -  -ksp_plot_eigenvalues - Plots eigenvalues in an x-window display

 80:    Notes:
 81:    The number of eigenvalues estimated depends on the size of the Krylov space
 82:    generated during the KSPSolve() (that is the SLESSolve); for example, with 
 83:    CG it corresponds to the number of CG iterations, for GMRES it is the number 
 84:    of GMRES iterations SINCE the last restart. Any extra space in r[] and c[]
 85:    will be ignored.

 87:    KSPComputeEigenvalues() does not usually provide accurate estimates; it is
 88:    intended only for assistance in understanding the convergence of iterative 
 89:    methods, not for eigenanalysis. 

 91:    One must call KSPSetComputeEigenvalues() before calling KSPSetUp() 
 92:    in order for this routine to work correctly.

 94:    Many users may just want to use the monitoring routine
 95:    KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
 96:    to print the singular values at each iteration of the linear solve.

 98:    Level: advanced

100: .keywords: KSP, compute, extreme, singular, values

102: .seealso: KSPSetComputeSingularValues(), KSPSingularValueMonitor(), KSPComputeExtremeSingularValues()
103: @*/
104: int KSPComputeEigenvalues(KSP ksp,int n,PetscReal *r,PetscReal *c,int *neig)
105: {

112:   if (!ksp->calc_sings) {
113:     SETERRQ(4,"Eigenvalues not requested before KSPSetUp()");
114:   }

116:   if (ksp->ops->computeeigenvalues) {
117:     (*ksp->ops->computeeigenvalues)(ksp,n,r,c,neig);
118:   } else {
119:     *neig = 0;
120:   }
121:   return(0);
122: }

126: /*@
127:    KSPSetUp - Sets up the internal data structures for the
128:    later use of an iterative solver.

130:    Collective on KSP

132:    Input Parameter:
133: .  ksp   - iterative context obtained from KSPCreate()

135:    Level: developer

137: .keywords: KSP, setup

139: .seealso: KSPCreate(), KSPSolve(), KSPDestroy()
140: @*/
141: int KSPSetUp(KSP ksp)
142: {


148:   /* reset the convergence flag from the previous solves */
149:   ksp->reason = KSP_CONVERGED_ITERATING;

151:   if (!ksp->type_name){
152:     KSPSetType(ksp,KSPGMRES);
153:   }

155:   if (ksp->setupcalled) return(0);
156:   ksp->setupcalled = 1;
157:   (*ksp->ops->setup)(ksp);
158:   return(0);
159: }


164: /*@
165:    KSPSolve - Solves linear system; usually not called directly, rather 
166:    it is called by a call to SLESSolve().

168:    Collective on KSP

170:    Input Parameter:
171: .  ksp - iterative context obtained from KSPCreate()

173:    Output Parameter:
174: .  its - number of iterations required

176:    Options Database Keys:
177: +  -ksp_compute_eigenvalues - compute preconditioned operators eigenvalues
178: .  -ksp_plot_eigenvalues - plot the computed eigenvalues in an X-window
179: .  -ksp_compute_eigenvalues_explicitly - compute the eigenvalues by forming the 
180:       dense operator and useing LAPACK
181: -  -ksp_plot_eigenvalues_explicitly - plot the explicitly computing eigenvalues

183:    Notes:
184:    On return, the parameter "its" contains either the iteration
185:    number at which convergence was successfully reached or failure was detected.

187:    Call KSPGetConvergedReason() to determine if the solver converged or failed and 
188:    why.
189:    
190:    If using a direct method (e.g., via the KSP solver
191:    KSPPREONLY and a preconditioner such as PCLU/PCILU),
192:    then its=1.  See KSPSetTolerances() and KSPDefaultConverged()
193:    for more details.

195:    Understanding Convergence:
196:    The routines KSPSetMonitor(), KSPComputeEigenvalues(), and
197:    KSPComputeEigenvaluesExplicitly() provide information on additional
198:    options to monitor convergence and print eigenvalue information.

200:    Level: developer

202: .keywords: KSP, solve, linear system

204: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPDefaultConverged(),
205:           SLESSolve(), KSPSolveTranspose(), SLESGetKSP()
206: @*/
207: int KSPSolve(KSP ksp,int *its)
208: {
209:   int          ierr,rank,nits;
210:   PetscTruth   flag1,flag2;
211:   PetscScalar  zero = 0.0;


216:   if (!ksp->setupcalled){ KSPSetUp(ksp);}
217:   if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
218:   if (ksp->guess_knoll) {
219:     PCApply(ksp->B,ksp->vec_rhs,ksp->vec_sol,PC_LEFT);
220:     ksp->guess_zero = PETSC_FALSE;
221:   }

223:   /* reset the residual history list if requested */
224:   if (ksp->res_hist_reset) ksp->res_hist_len = 0;

226:   ksp->transpose_solve = PETSC_FALSE;
227:   (*ksp->ops->solve)(ksp,&nits);
228:   if (!ksp->reason) {
229:     SETERRQ(1,"Internal error, solver returned without setting converged reason");
230:   }
231:   if (its) *its = nits;

233:   MPI_Comm_rank(ksp->comm,&rank);

235:   PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues",&flag1);
236:   PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues",&flag2);
237:   if (flag1 || flag2) {
238:     int       n = nits + 2,i,neig;
239:     PetscReal *r,*c;

241:     if (!n) {
242:       PetscPrintf(ksp->comm,"Zero iterations in solver, cannot approximate any eigenvalues\n");
243:     } else {
244:       PetscMalloc(2*n*sizeof(PetscReal),&r);
245:       c = r + n;
246:       KSPComputeEigenvalues(ksp,n,r,c,&neig);
247:       if (flag1) {
248:         PetscPrintf(ksp->comm,"Iteratively computed eigenvalues\n");
249:         for (i=0; i<neig; i++) {
250:           if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gi\n",r[i],c[i]);}
251:           else             {PetscPrintf(ksp->comm,"%g - %gi\n",r[i],-c[i]);}
252:         }
253:       }
254:       if (flag2 && !rank) {
255:         PetscViewer viewer;
256:         PetscDraw   draw;
257:         PetscDrawSP drawsp;

259:         PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Iteratively Computed Eigenvalues",
260:                                PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
261:         PetscViewerDrawGetDraw(viewer,0,&draw);
262:         PetscDrawSPCreate(draw,1,&drawsp);
263:         for (i=0; i<neig; i++) {
264:           PetscDrawSPAddPoint(drawsp,r+i,c+i);
265:         }
266:         PetscDrawSPDraw(drawsp);
267:         PetscDrawSPDestroy(drawsp);
268:         PetscViewerDestroy(viewer);
269:       }
270:       PetscFree(r);
271:     }
272:   }

274:   PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues_explicitly",&flag1);
275:   PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues_explicitly",&flag2);
276:   if (flag1 || flag2) {
277:     int       n,i;
278:     PetscReal *r,*c;
279:     VecGetSize(ksp->vec_sol,&n);
280:     PetscMalloc(2*n*sizeof(PetscReal),&r);
281:     c = r + n;
282:     KSPComputeEigenvaluesExplicitly(ksp,n,r,c);
283:     if (flag1) {
284:       PetscPrintf(ksp->comm,"Explicitly computed eigenvalues\n");
285:       for (i=0; i<n; i++) {
286:         if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gi\n",r[i],c[i]);}
287:         else             {PetscPrintf(ksp->comm,"%g - %gi\n",r[i],-c[i]);}
288:       }
289:     }
290:     if (flag2 && !rank) {
291:       PetscViewer viewer;
292:       PetscDraw   draw;
293:       PetscDrawSP drawsp;

295:       PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Explicitly Computed Eigenvalues",0,320,300,300,&viewer);
296:       PetscViewerDrawGetDraw(viewer,0,&draw);
297:       PetscDrawSPCreate(draw,1,&drawsp);
298:       for (i=0; i<n; i++) {
299:         PetscDrawSPAddPoint(drawsp,r+i,c+i);
300:       }
301:       PetscDrawSPDraw(drawsp);
302:       PetscDrawSPDestroy(drawsp);
303:       PetscViewerDestroy(viewer);
304:     }
305:     PetscFree(r);
306:   }

308:   PetscOptionsHasName(ksp->prefix,"-ksp_view_operator",&flag2);
309:   if (flag2) {
310:     Mat A,B;
311:     PCGetOperators(ksp->B,&A,PETSC_NULL,PETSC_NULL);
312:     MatComputeExplicitOperator(A,&B);
313:     PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(ksp->comm),PETSC_VIEWER_ASCII_MATLAB);
314:     MatView(B,PETSC_VIEWER_STDOUT_(ksp->comm));
315:     PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(ksp->comm));
316:     MatDestroy(B);
317:   }
318:   PetscOptionsHasName(ksp->prefix,"-ksp_view_operator_binary",&flag2);
319:   if (flag2) {
320:     Mat A,B;
321:     PCGetOperators(ksp->B,&A,PETSC_NULL,PETSC_NULL);
322:     MatComputeExplicitOperator(A,&B);
323:     MatView(B,PETSC_VIEWER_BINARY_(ksp->comm));
324:     MatDestroy(B);
325:   }
326:   PetscOptionsHasName(ksp->prefix,"-ksp_view_preconditioned_operator_binary",&flag2);
327:   if (flag2) {
328:     Mat B;
329:     KSPComputeExplicitOperator(ksp,&B);
330:     MatView(B,PETSC_VIEWER_BINARY_(ksp->comm));
331:     MatDestroy(B);
332:   }
333:   return(0);
334: }

338: /*@
339:    KSPSolveTranspose - Solves the transpose of a linear system. Usually
340:    accessed through SLESSolveTranspose().

342:    Collective on KSP

344:    Input Parameter:
345: .  ksp - iterative context obtained from KSPCreate()

347:    Output Parameter:
348: .  its - number of iterations required

350:    Notes:
351:    On return, the parameter "its" contains either the iteration
352:    number at which convergence was successfully reached, or the
353:    negative of the iteration at which divergence or breakdown was detected.

355:    Currently only supported by KSPType of KSPPREONLY. This routine is usally 
356:    only used internally by the BiCG solver on the subblocks in BJacobi and ASM.

358:    Level: developer

360: .keywords: KSP, solve, linear system

362: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPDefaultConverged(),
363:           SLESSolve(), SLESGetKSP()
364: @*/
365: int KSPSolveTranspose(KSP ksp,int *its)
366: {
367:   int           ierr;
368:   PetscScalar   zero = 0.0;


373:   if (!ksp->setupcalled){ KSPSetUp(ksp);}
374:   if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
375:   ksp->transpose_solve = PETSC_TRUE;
376:   (*ksp->ops->solve)(ksp,its);
377:   return(0);
378: }

382: /*@C
383:    KSPDestroy - Destroys KSP context.

385:    Collective on KSP

387:    Input Parameter:
388: .  ksp - iterative context obtained from KSPCreate()

390:    Level: developer

392: .keywords: KSP, destroy

394: .seealso: KSPCreate(), KSPSetUp(), KSPSolve()
395: @*/
396: int KSPDestroy(KSP ksp)
397: {
398:   int i,ierr;

402:   if (--ksp->refct > 0) return(0);

404:   /* if memory was published with AMS then destroy it */
405:   PetscObjectDepublish(ksp);

407:   if (ksp->ops->destroy) {
408:     (*ksp->ops->destroy)(ksp);
409:   }
410:   for (i=0; i<ksp->numbermonitors; i++) {
411:     if (ksp->monitordestroy[i]) {
412:       (*ksp->monitordestroy[i])(ksp->monitorcontext[i]);
413:     }
414:   }
415:   PetscLogObjectDestroy(ksp);
416:   PetscHeaderDestroy(ksp);
417:   return(0);
418: }

422: /*@
423:     KSPSetPreconditionerSide - Sets the preconditioning side.

425:     Collective on KSP

427:     Input Parameter:
428: .   ksp - iterative context obtained from KSPCreate()

430:     Output Parameter:
431: .   side - the preconditioning side, where side is one of
432: .vb
433:       PC_LEFT - left preconditioning (default)
434:       PC_RIGHT - right preconditioning
435:       PC_SYMMETRIC - symmetric preconditioning
436: .ve

438:     Options Database Keys:
439: +   -ksp_left_pc - Sets left preconditioning
440: .   -ksp_right_pc - Sets right preconditioning
441: -   -ksp_symmetric_pc - Sets symmetric preconditioning

443:     Notes:
444:     Left preconditioning is used by default.  Symmetric preconditioning is
445:     currently available only for the KSPQCG method. Note, however, that
446:     symmetric preconditioning can be emulated by using either right or left
447:     preconditioning and a pre or post processing step.

449:     Level: intermediate

451: .keywords: KSP, set, right, left, symmetric, side, preconditioner, flag

453: .seealso: KSPGetPreconditionerSide()
454: @*/
455: int KSPSetPreconditionerSide(KSP ksp,PCSide side)
456: {
459:   ksp->pc_side = side;
460:   return(0);
461: }

465: /*@C
466:     KSPGetPreconditionerSide - Gets the preconditioning side.

468:     Not Collective

470:     Input Parameter:
471: .   ksp - iterative context obtained from KSPCreate()

473:     Output Parameter:
474: .   side - the preconditioning side, where side is one of
475: .vb
476:       PC_LEFT - left preconditioning (default)
477:       PC_RIGHT - right preconditioning
478:       PC_SYMMETRIC - symmetric preconditioning
479: .ve

481:     Level: intermediate

483: .keywords: KSP, get, right, left, symmetric, side, preconditioner, flag

485: .seealso: KSPSetPreconditionerSide()
486: @*/
487: int KSPGetPreconditionerSide(KSP ksp,PCSide *side)
488: {
491:   *side = ksp->pc_side;
492:   return(0);
493: }

497: /*@
498:    KSPGetTolerances - Gets the relative, absolute, divergence, and maximum
499:    iteration tolerances used by the default KSP convergence tests. 

501:    Not Collective

503:    Input Parameter:
504: .  ksp - the Krylov subspace context
505:   
506:    Output Parameters:
507: +  rtol - the relative convergence tolerance
508: .  atol - the absolute convergence tolerance
509: .  dtol - the divergence tolerance
510: -  maxits - maximum number of iterations

512:    Notes:
513:    The user can specify PETSC_NULL for any parameter that is not needed.

515:    Level: intermediate

517: .keywords: KSP, get, tolerance, absolute, relative, divergence, convergence,
518:            maximum, iterations

520: .seealso: KSPSetTolerances()
521: @*/
522: int KSPGetTolerances(KSP ksp,PetscReal *rtol,PetscReal *atol,PetscReal *dtol,int *maxits)
523: {
526:   if (atol)   *atol   = ksp->atol;
527:   if (rtol)   *rtol   = ksp->rtol;
528:   if (dtol)   *dtol   = ksp->divtol;
529:   if (maxits) *maxits = ksp->max_it;
530:   return(0);
531: }

535: /*@
536:    KSPSetTolerances - Sets the relative, absolute, divergence, and maximum
537:    iteration tolerances used by the default KSP convergence testers. 

539:    Collective on KSP

541:    Input Parameters:
542: +  ksp - the Krylov subspace context
543: .  rtol - the relative convergence tolerance
544:    (relative decrease in the residual norm)
545: .  atol - the absolute convergence tolerance 
546:    (absolute size of the residual norm)
547: .  dtol - the divergence tolerance
548:    (amount residual can increase before KSPDefaultConverged()
549:    concludes that the method is diverging)
550: -  maxits - maximum number of iterations to use

552:    Options Database Keys:
553: +  -ksp_atol <atol> - Sets atol
554: .  -ksp_rtol <rtol> - Sets rtol
555: .  -ksp_divtol <dtol> - Sets dtol
556: -  -ksp_max_it <maxits> - Sets maxits

558:    Notes:
559:    Use PETSC_DEFAULT to retain the default value of any of the tolerances.

561:    See KSPDefaultConverged() for details on the use of these parameters
562:    in the default convergence test.  See also KSPSetConvergenceTest() 
563:    for setting user-defined stopping criteria.

565:    Level: intermediate

567: .keywords: KSP, set, tolerance, absolute, relative, divergence, 
568:            convergence, maximum, iterations

570: .seealso: KSPGetTolerances(), KSPDefaultConverged(), KSPSetConvergenceTest()
571: @*/
572: int KSPSetTolerances(KSP ksp,PetscReal rtol,PetscReal atol,PetscReal dtol,int maxits)
573: {
576:   if (atol != PETSC_DEFAULT)   ksp->atol   = atol;
577:   if (rtol != PETSC_DEFAULT)   ksp->rtol   = rtol;
578:   if (dtol != PETSC_DEFAULT)   ksp->divtol = dtol;
579:   if (maxits != PETSC_DEFAULT) ksp->max_it = maxits;
580:   return(0);
581: }

585: /*@
586:    KSPSetInitialGuessNonzero - Tells the iterative solver that the 
587:    initial guess is nonzero; otherwise KSP assumes the initial guess
588:    is to be zero (and thus zeros it out before solving).

590:    Collective on KSP

592:    Input Parameters:
593: +  ksp - iterative context obtained from SLESGetKSP() or KSPCreate()
594: -  flg - PETSC_TRUE indicates the guess is non-zero, PETSC_FALSE indicates the guess is zero

596:    Level: beginner

598:    Notes:
599:     If this is not called the X vector is zeroed in the call to 
600: SLESSolve() (or KSPSolve()).

602: .keywords: KSP, set, initial guess, nonzero

604: .seealso: KSPGetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
605: @*/
606: int KSPSetInitialGuessNonzero(KSP ksp,PetscTruth flg)
607: {
609:   ksp->guess_zero   = (PetscTruth)!(int)flg;
610:   return(0);
611: }

615: /*@
616:    KSPGetInitialGuessNonzero - Determines whether the KSP solver is using
617:    a zero initial guess.

619:    Not Collective

621:    Input Parameter:
622: .  ksp - iterative context obtained from KSPCreate()

624:    Output Parameter:
625: .  flag - PETSC_TRUE if guess is nonzero, else PETSC_FALSE

627:    Level: intermediate

629: .keywords: KSP, set, initial guess, nonzero

631: .seealso: KSPSetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
632: @*/
633: int KSPGetInitialGuessNonzero(KSP ksp,PetscTruth *flag)
634: {
636:   if (ksp->guess_zero) *flag = PETSC_FALSE;
637:   else                 *flag = PETSC_TRUE;
638:   return(0);
639: }

643: /*@
644:    KSPSetInitialGuessKnoll - Tells the iterative solver to use PCApply(pc,b,..) to compute the initial guess (The Knoll trick)

646:    Collective on KSP

648:    Input Parameters:
649: +  ksp - iterative context obtained from SLESGetKSP() or KSPCreate()
650: -  flg - PETSC_TRUE or PETSC_FALSE

652:    Level: advanced


655: .keywords: KSP, set, initial guess, nonzero

657: .seealso: KSPGetInitialGuessKnoll(), KSPSetInitialGuessNonzero(), KSPGetInitialGuessNonzero()
658: @*/
659: int KSPSetInitialGuessKnoll(KSP ksp,PetscTruth flg)
660: {
662:   ksp->guess_knoll   = flg;
663:   return(0);
664: }

668: /*@
669:    KSPGetInitialGuessKnoll - Determines whether the KSP solver is using the Knoll trick (using PCApply(pc,b,...) to compute
670:      the initial guess

672:    Not Collective

674:    Input Parameter:
675: .  ksp - iterative context obtained from KSPCreate()

677:    Output Parameter:
678: .  flag - PETSC_TRUE if using Knoll trick, else PETSC_FALSE

680:    Level: advanced

682: .keywords: KSP, set, initial guess, nonzero

684: .seealso: KSPSetInitialGuessKnoll(), KSPSetInitialGuessNonzero(), KSPGetInitialGuessNonzero()
685: @*/
686: int KSPGetInitialGuessKnoll(KSP ksp,PetscTruth *flag)
687: {
689:   *flag = ksp->guess_knoll;
690:   return(0);
691: }

695: /*@
696:    KSPSetComputeSingularValues - Sets a flag so that the extreme singular 
697:    values will be calculated via a Lanczos or Arnoldi process as the linear 
698:    system is solved.

700:    Collective on KSP

702:    Input Parameters:
703: +  ksp - iterative context obtained from KSPCreate()
704: -  flg - PETSC_TRUE or PETSC_FALSE

706:    Options Database Key:
707: .  -ksp_singmonitor - Activates KSPSetComputeSingularValues()

709:    Notes:
710:    Currently this option is not valid for all iterative methods.

712:    Many users may just want to use the monitoring routine
713:    KSPSingularValueMonitor() (which can be set with option -ksp_singmonitor)
714:    to print the singular values at each iteration of the linear solve.

716:    Level: advanced

718: .keywords: KSP, set, compute, singular values

720: .seealso: KSPComputeExtremeSingularValues(), KSPSingularValueMonitor()
721: @*/
722: int KSPSetComputeSingularValues(KSP ksp,PetscTruth flg)
723: {
726:   ksp->calc_sings  = flg;
727:   return(0);
728: }

732: /*@
733:    KSPSetComputeEigenvalues - Sets a flag so that the extreme eigenvalues
734:    values will be calculated via a Lanczos or Arnoldi process as the linear 
735:    system is solved.

737:    Collective on KSP

739:    Input Parameters:
740: +  ksp - iterative context obtained from KSPCreate()
741: -  flg - PETSC_TRUE or PETSC_FALSE

743:    Notes:
744:    Currently this option is not valid for all iterative methods.

746:    Level: advanced

748: .keywords: KSP, set, compute, eigenvalues

750: .seealso: KSPComputeEigenvalues(), KSPComputeEigenvaluesExplicitly()
751: @*/
752: int KSPSetComputeEigenvalues(KSP ksp,PetscTruth flg)
753: {
756:   ksp->calc_sings  = flg;
757:   return(0);
758: }

762: /*@
763:    KSPSetRhs - Sets the right-hand-side vector for the linear system to
764:    be solved.

766:    Collective on KSP and Vec

768:    Input Parameters:
769: +  ksp - iterative context obtained from KSPCreate()
770: -  b   - right-hand-side vector

772:    Level: developer

774: .keywords: KSP, set, right-hand-side, rhs

776: .seealso: KSPGetRhs(), KSPSetSolution()
777: @*/
778: int KSPSetRhs(KSP ksp,Vec b)
779: {
783:   ksp->vec_rhs    = (b);
784:   return(0);
785: }

789: /*@C
790:    KSPGetRhs - Gets the right-hand-side vector for the linear system to
791:    be solved.

793:    Not Collective

795:    Input Parameter:
796: .  ksp - iterative context obtained from KSPCreate()

798:    Output Parameter:
799: .  r - right-hand-side vector

801:    Level: developer

803: .keywords: KSP, get, right-hand-side, rhs

805: .seealso: KSPSetRhs(), KSPGetSolution()
806: @*/
807: int KSPGetRhs(KSP ksp,Vec *r)
808: {
811:   *r = ksp->vec_rhs;
812:   return(0);
813: }

817: /*@
818:    KSPSetSolution - Sets the location of the solution for the 
819:    linear system to be solved.

821:    Collective on KSP and Vec

823:    Input Parameters:
824: +  ksp - iterative context obtained from KSPCreate()
825: -  x   - solution vector

827:    Level: developer

829: .keywords: KSP, set, solution

831: .seealso: KSPSetRhs(), KSPGetSolution()
832: @*/
833: int KSPSetSolution(KSP ksp,Vec x)
834: {
838:   ksp->vec_sol    = (x);
839:   return(0);
840: }

844: /*@C
845:    KSPGetSolution - Gets the location of the solution for the 
846:    linear system to be solved.  Note that this may not be where the solution
847:    is stored during the iterative process; see KSPBuildSolution().

849:    Not Collective

851:    Input Parameters:
852: .  ksp - iterative context obtained from KSPCreate()

854:    Output Parameters:
855: .  v - solution vector

857:    Level: developer

859: .keywords: KSP, get, solution

861: .seealso: KSPGetRhs(), KSPSetSolution(), KSPBuildSolution()
862: @*/
863: int KSPGetSolution(KSP ksp,Vec *v)
864: {
867:   *v = ksp->vec_sol;
868:   return(0);
869: }

873: /*@
874:    KSPSetPC - Sets the preconditioner to be used to calculate the 
875:    application of the preconditioner on a vector. 

877:    Collective on KSP

879:    Input Parameters:
880: +  ksp - iterative context obtained from KSPCreate()
881: -  B   - the preconditioner object

883:    Notes:
884:    Use KSPGetPC() to retrieve the preconditioner context (for example,
885:    to free it at the end of the computations).

887:    Level: developer

889: .keywords: KSP, set, precondition, Binv

891: .seealso: KSPGetPC()
892: @*/
893: int KSPSetPC(KSP ksp,PC B)
894: {
899:   ksp->B = B;
900:   return(0);
901: }

905: /*@C
906:    KSPGetPC - Returns a pointer to the preconditioner context
907:    set with KSPSetPC().

909:    Not Collective

911:    Input Parameters:
912: .  ksp - iterative context obtained from KSPCreate()

914:    Output Parameter:
915: .  B - preconditioner context

917:    Level: developer

919: .keywords: KSP, get, preconditioner, Binv

921: .seealso: KSPSetPC()
922: @*/
923: int KSPGetPC(KSP ksp,PC *B)
924: {
927:   *B = ksp->B;
928:   return(0);
929: }

933: /*@C
934:    KSPSetMonitor - Sets an ADDITIONAL function to be called at every iteration to monitor 
935:    the residual/error etc.
936:       
937:    Collective on KSP

939:    Input Parameters:
940: +  ksp - iterative context obtained from KSPCreate()
941: .  monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring
942: .  mctx    - [optional] context for private data for the
943:              monitor routine (use PETSC_NULL if no context is desired)
944: -  monitordestroy - [optional] routine that frees monitor context
945:           (may be PETSC_NULL)

947:    Calling Sequence of monitor:
948: $     monitor (KSP ksp, int it, PetscReal rnorm, void *mctx)

950: +  ksp - iterative context obtained from KSPCreate()
951: .  it - iteration number
952: .  rnorm - (estimated) 2-norm of (preconditioned) residual
953: -  mctx  - optional monitoring context, as set by KSPSetMonitor()

955:    Options Database Keys:
956: +    -ksp_monitor        - sets KSPDefaultMonitor()
957: .    -ksp_truemonitor    - sets KSPTrueMonitor()
958: .    -ksp_xmonitor       - sets line graph monitor,
959:                            uses KSPLGMonitorCreate()
960: .    -ksp_xtruemonitor   - sets line graph monitor,
961:                            uses KSPLGMonitorCreate()
962: .    -ksp_singmonitor    - sets KSPSingularValueMonitor()
963: -    -ksp_cancelmonitors - cancels all monitors that have
964:                           been hardwired into a code by 
965:                           calls to KSPSetMonitor(), but
966:                           does not cancel those set via
967:                           the options database.

969:    Notes:  
970:    The default is to do nothing.  To print the residual, or preconditioned 
971:    residual if KSPSetNormType(ksp,KSP_PRECONDITIONED_NORM) was called, use 
972:    KSPDefaultMonitor() as the monitoring routine, with a null monitoring 
973:    context. 

975:    Several different monitoring routines may be set by calling
976:    KSPSetMonitor() multiple times; all will be called in the 
977:    order in which they were set.

979:    Level: beginner

981: .keywords: KSP, set, monitor

983: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPClearMonitor()
984: @*/
985: int KSPSetMonitor(KSP ksp,int (*monitor)(KSP,int,PetscReal,void*),void *mctx,int (*monitordestroy)(void*))
986: {
989:   if (ksp->numbermonitors >= MAXKSPMONITORS) {
990:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many KSP monitors set");
991:   }
992:   ksp->monitor[ksp->numbermonitors]           = monitor;
993:   ksp->monitordestroy[ksp->numbermonitors]    = monitordestroy;
994:   ksp->monitorcontext[ksp->numbermonitors++]  = (void*)mctx;
995:   return(0);
996: }

1000: /*@
1001:    KSPClearMonitor - Clears all monitors for a KSP object.

1003:    Collective on KSP

1005:    Input Parameters:
1006: .  ksp - iterative context obtained from KSPCreate()

1008:    Options Database Key:
1009: .  -ksp_cancelmonitors - Cancels all monitors that have
1010:     been hardwired into a code by calls to KSPSetMonitor(), 
1011:     but does not cancel those set via the options database.

1013:    Level: intermediate

1015: .keywords: KSP, set, monitor

1017: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPSetMonitor()
1018: @*/
1019: int KSPClearMonitor(KSP ksp)
1020: {
1023:   ksp->numbermonitors = 0;
1024:   return(0);
1025: }

1029: /*@C
1030:    KSPGetMonitorContext - Gets the monitoring context, as set by 
1031:    KSPSetMonitor() for the FIRST monitor only.

1033:    Not Collective

1035:    Input Parameter:
1036: .  ksp - iterative context obtained from KSPCreate()

1038:    Output Parameter:
1039: .  ctx - monitoring context

1041:    Level: intermediate

1043: .keywords: KSP, get, monitor, context

1045: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate()
1046: @*/
1047: int KSPGetMonitorContext(KSP ksp,void **ctx)
1048: {
1051:   *ctx =      (ksp->monitorcontext[0]);
1052:   return(0);
1053: }

1057: /*@
1058:    KSPSetResidualHistory - Sets the array used to hold the residual history.
1059:    If set, this array will contain the residual norms computed at each
1060:    iteration of the solver.

1062:    Not Collective

1064:    Input Parameters:
1065: +  ksp - iterative context obtained from KSPCreate()
1066: .  a   - array to hold history
1067: .  na  - size of a
1068: -  reset - PETSC_TRUE indicates the history counter is reset to zero
1069:            for each new linear solve

1071:    Level: advanced

1073:    Notes: The array is NOT freed by PETSc so the user needs to keep track of 
1074:            it and destroy once the KSP object is destroyed.

1076: .keywords: KSP, set, residual, history, norm

1078: .seealso: KSPGetResidualHistory()

1080: @*/
1081: int KSPSetResidualHistory(KSP ksp,PetscReal a[],int na,PetscTruth reset)
1082: {

1087:   if (na != PETSC_DECIDE && a != PETSC_NULL) {
1088:     ksp->res_hist        = a;
1089:     ksp->res_hist_max    = na;
1090:   } else {
1091:     ksp->res_hist_max    = 1000;
1092:     PetscMalloc(ksp->res_hist_max*sizeof(PetscReal),&ksp->res_hist);
1093:   }
1094:   ksp->res_hist_len    = 0;
1095:   ksp->res_hist_reset  = reset;


1098:   return(0);
1099: }

1103: /*@C
1104:    KSPGetResidualHistory - Gets the array used to hold the residual history
1105:    and the number of residuals it contains.

1107:    Not Collective

1109:    Input Parameter:
1110: .  ksp - iterative context obtained from KSPCreate()

1112:    Output Parameters:
1113: +  a   - pointer to array to hold history (or PETSC_NULL)
1114: -  na  - number of used entries in a (or PETSC_NULL)

1116:    Level: advanced

1118:    Notes:
1119:      Can only be called after a KSPSetResidualHistory() otherwise a and na are set to zero

1121:      The Fortran version of this routine has a calling sequence
1122: $   call KSPGetResidualHistory(KSP ksp, integer na, integer ierr)

1124: .keywords: KSP, get, residual, history, norm

1126: .seealso: KSPGetResidualHistory()

1128: @*/
1129: int KSPGetResidualHistory(KSP ksp,PetscReal *a[],int *na)
1130: {
1133:   if (a)  *a = ksp->res_hist;
1134:   if (na) *na = ksp->res_hist_len;
1135:   return(0);
1136: }

1140: /*@C
1141:    KSPSetConvergenceTest - Sets the function to be used to determine
1142:    convergence.  

1144:    Collective on KSP

1146:    Input Parameters:
1147: +  ksp - iterative context obtained from KSPCreate()
1148: .  converge - pointer to int function
1149: -  cctx    - context for private data for the convergence routine (may be null)

1151:    Calling sequence of converge:
1152: $     converge (KSP ksp, int it, PetscReal rnorm, KSPConvergedReason *reason,void *mctx)

1154: +  ksp - iterative context obtained from KSPCreate()
1155: .  it - iteration number
1156: .  rnorm - (estimated) 2-norm of (preconditioned) residual
1157: .  reason - the reason why it has converged or diverged
1158: -  cctx  - optional convergence context, as set by KSPSetConvergenceTest()


1161:    Notes:
1162:    Must be called after the KSP type has been set so put this after
1163:    a call to KSPSetType(), or KSPSetFromOptions().

1165:    The default convergence test, KSPDefaultConverged(), aborts if the 
1166:    residual grows to more than 10000 times the initial residual.

1168:    The default is a combination of relative and absolute tolerances.  
1169:    The residual value that is tested may be an approximation; routines 
1170:    that need exact values should compute them.

1172:    Level: advanced

1174: .keywords: KSP, set, convergence, test, context

1176: .seealso: KSPDefaultConverged(), KSPGetConvergenceContext()
1177: @*/
1178: int KSPSetConvergenceTest(KSP ksp,int (*converge)(KSP,int,PetscReal,KSPConvergedReason*,void*),void *cctx)
1179: {
1182:   ksp->converged = converge;
1183:   ksp->cnvP      = (void*)cctx;
1184:   return(0);
1185: }

1189: /*@C
1190:    KSPGetConvergenceContext - Gets the convergence context set with 
1191:    KSPSetConvergenceTest().  

1193:    Not Collective

1195:    Input Parameter:
1196: .  ksp - iterative context obtained from KSPCreate()

1198:    Output Parameter:
1199: .  ctx - monitoring context

1201:    Level: advanced

1203: .keywords: KSP, get, convergence, test, context

1205: .seealso: KSPDefaultConverged(), KSPSetConvergenceTest()
1206: @*/
1207: int KSPGetConvergenceContext(KSP ksp,void **ctx)
1208: {
1211:   *ctx = ksp->cnvP;
1212:   return(0);
1213: }

1217: /*@C
1218:    KSPBuildSolution - Builds the approximate solution in a vector provided.
1219:    This routine is NOT commonly needed (see SLESSolve()).

1221:    Collective on KSP

1223:    Input Parameter:
1224: .  ctx - iterative context obtained from KSPCreate()

1226:    Output Parameter: 
1227:    Provide exactly one of
1228: +  v - location to stash solution.   
1229: -  V - the solution is returned in this location. This vector is created 
1230:        internally. This vector should NOT be destroyed by the user with
1231:        VecDestroy().

1233:    Notes:
1234:    This routine can be used in one of two ways
1235: .vb
1236:       KSPBuildSolution(ksp,PETSC_NULL,&V);
1237:    or
1238:       KSPBuildSolution(ksp,v,PETSC_NULL); 
1239: .ve
1240:    In the first case an internal vector is allocated to store the solution
1241:    (the user cannot destroy this vector). In the second case the solution
1242:    is generated in the vector that the user provides. Note that for certain 
1243:    methods, such as KSPCG, the second case requires a copy of the solution,
1244:    while in the first case the call is essentially free since it simply 
1245:    returns the vector where the solution already is stored.

1247:    Level: advanced

1249: .keywords: KSP, build, solution

1251: .seealso: KSPGetSolution(), KSPBuildResidual()
1252: @*/
1253: int KSPBuildSolution(KSP ksp,Vec v,Vec *V)
1254: {

1259:   if (!V && !v) SETERRQ(PETSC_ERR_ARG_WRONG,"Must provide either v or V");
1260:   if (!V) V = &v;
1261:   (*ksp->ops->buildsolution)(ksp,v,V);
1262:   return(0);
1263: }

1267: /*@C
1268:    KSPBuildResidual - Builds the residual in a vector provided.

1270:    Collective on KSP

1272:    Input Parameter:
1273: .  ksp - iterative context obtained from KSPCreate()

1275:    Output Parameters:
1276: +  v - optional location to stash residual.  If v is not provided,
1277:        then a location is generated.
1278: .  t - work vector.  If not provided then one is generated.
1279: -  V - the residual

1281:    Notes:
1282:    Regardless of whether or not v is provided, the residual is 
1283:    returned in V.

1285:    Level: advanced

1287: .keywords: KSP, build, residual

1289: .seealso: KSPBuildSolution()
1290: @*/
1291: int KSPBuildResidual(KSP ksp,Vec t,Vec v,Vec *V)
1292: {
1293:   int flag = 0,ierr;
1294:   Vec w = v,tt = t;

1298:   if (!w) {
1299:     VecDuplicate(ksp->vec_rhs,&w);
1300:     PetscLogObjectParent((PetscObject)ksp,w);
1301:   }
1302:   if (!tt) {
1303:     VecDuplicate(ksp->vec_rhs,&tt); flag = 1;
1304:     PetscLogObjectParent((PetscObject)ksp,tt);
1305:   }
1306:   (*ksp->ops->buildresidual)(ksp,tt,w,V);
1307:   if (flag) {VecDestroy(tt);}
1308:   return(0);
1309: }