Actual source code: itfunc.c

  1: /*
  2:       Interface KSP routines that the user calls.
  3: */

 5:  #include src/ksp/ksp/kspimpl.h

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

 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: PetscErrorCode 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:   } else {
 50:     *emin = -1.0;
 51:     *emax = -1.0;
 52:   }
 53:   return(0);
 54: }

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

 62:    Not Collective

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

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

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

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

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

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

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

 96:    Level: advanced

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

100: .seealso: KSPSetComputeSingularValues(), KSPSingularValueMonitor(), KSPComputeExtremeSingularValues()
101: @*/
102: PetscErrorCode KSPComputeEigenvalues(KSP ksp,PetscInt n,PetscReal *r,PetscReal *c,PetscInt *neig)
103: {

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

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

125: /*@
126:    KSPSetUpOnBlocks - Sets up the preconditioner for each block in
127:    the block Jacobi, block Gauss-Seidel, and overlapping Schwarz 
128:    methods.

130:    Collective on KSP

132:    Input Parameter:
133: .  ksp - the KSP context

135:    Notes:
136:    KSPSetUpOnBlocks() is a routine that the user can optinally call for
137:    more precise profiling (via -log_summary) of the setup phase for these
138:    block preconditioners.  If the user does not call KSPSetUpOnBlocks(),
139:    it will automatically be called from within KSPSolve().
140:    
141:    Calling KSPSetUpOnBlocks() is the same as calling PCSetUpOnBlocks()
142:    on the PC context within the KSP context.

144:    Level: advanced

146: .keywords: KSP, setup, blocks

148: .seealso: PCSetUpOnBlocks(), KSPSetUp(), PCSetUp()
149: @*/
150: PetscErrorCode KSPSetUpOnBlocks(KSP ksp)
151: {

156:   PCSetUpOnBlocks(ksp->pc);
157:   return(0);
158: }

162: /*@
163:    KSPSetUp - Sets up the internal data structures for the
164:    later use of an iterative solver.

166:    Collective on KSP

168:    Input Parameter:
169: .  ksp   - iterative context obtained from KSPCreate()

171:    Level: developer

173: .keywords: KSP, setup

175: .seealso: KSPCreate(), KSPSolve(), KSPDestroy()
176: @*/
177: PetscErrorCode KSPSetUp(KSP ksp)
178: {


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

187:   if (!ksp->type_name){
188:     KSPSetType(ksp,KSPGMRES);
189:   }

191:   if (ksp->setupcalled == 2) return(0);

193:   PetscLogEventBegin(KSP_SetUp,ksp,ksp->vec_rhs,ksp->vec_sol,0);

195:   if (!ksp->setupcalled) {
196:     (*ksp->ops->setup)(ksp);
197:   }

199:   /* scale the matrix if requested */
200:   if (ksp->dscale) {
201:     Mat mat,pmat;
202:     PCGetOperators(ksp->pc,&mat,&pmat,PETSC_NULL);
203:     if (mat == pmat) {
204:       PetscScalar  *xx;
205:       PetscInt          i,n;
206:       PetscTruth   zeroflag = PETSC_FALSE;

208:       if (!ksp->diagonal) { /* allocate vector to hold diagonal */
209:         MatGetVecs(pmat,&ksp->diagonal,0);
210:       }
211:       MatGetDiagonal(mat,ksp->diagonal);
212:       VecGetLocalSize(ksp->diagonal,&n);
213:       VecGetArray(ksp->diagonal,&xx);
214:       for (i=0; i<n; i++) {
215:         if (xx[i] != 0.0) xx[i] = 1.0/sqrt(PetscAbsScalar(xx[i]));
216:         else {
217:           xx[i]     = 1.0;
218:           zeroflag  = PETSC_TRUE;
219:         }
220:       }
221:       VecRestoreArray(ksp->diagonal,&xx);
222:       if (zeroflag) {
223:         PetscLogInfo(ksp,"KSPSetUp:Zero detected in diagonal of matrix, using 1 at those locations\n");
224:       }
225:       MatDiagonalScale(mat,ksp->diagonal,ksp->diagonal);
226:       ksp->dscalefix2 = PETSC_FALSE;
227:     } else {
228:       SETERRQ(PETSC_ERR_SUP,"No support for diagonal scaling of linear system if preconditioner matrix not actual matrix");
229:     }
230:   }
231:   PetscLogEventEnd(KSP_SetUp,ksp,ksp->vec_rhs,ksp->vec_sol,0);
232:   PCSetUp(ksp->pc);
233:   if (ksp->nullsp) {
234:     PetscTruth test;
235:     PetscOptionsHasName(ksp->prefix,"-ksp_test_null_space",&test);
236:     if (test) {
237:       Mat mat;
238:       PCGetOperators(ksp->pc,&mat,PETSC_NULL,PETSC_NULL);
239:       MatNullSpaceTest(ksp->nullsp,mat);
240:     }
241:   }
242:   ksp->setupcalled = 2;
243:   return(0);
244: }

246: /*
247:     The numerical values for these is set in include/petscksp.h changes
248:    there must be made here.
249: */
250: static const char *convergedreasons[] = {"preconditioner is indefinite",                         "matrix or preconditioner is nonsymmetric",
251:                                          "breakdown in BICG",                                    "breakdown",
252:                                          "residual norm increased by dtol",                      "reached maximum number of iterations",
253:                                          "null operator or preconditioner",                      "not used",
254:                                          "never reached",                                        "not used",
255:                                          "residual 2-norm < relative tolerance * 2-norm of RHS", "residual 2-norm below absolute tolerance",
256:                                          "only one iteration requested",                         "negative curvature obtained in QCG",
257:                                          "constrained in QCG",                                   "small step length reached"};

261: /*@
262:    KSPSolve - Solves linear system.

264:    Collective on KSP

266:    Parameter:
267: +  ksp - iterative context obtained from KSPCreate()
268: .  b - the right hand side vector
269: -  x - the solution 

271:    Options Database Keys:
272: +  -ksp_compute_eigenvalues - compute preconditioned operators eigenvalues
273: .  -ksp_plot_eigenvalues - plot the computed eigenvalues in an X-window
274: .  -ksp_compute_eigenvalues_explicitly - compute the eigenvalues by forming the 
275:       dense operator and useing LAPACK
276: .  -ksp_plot_eigenvalues_explicitly - plot the explicitly computing eigenvalues
277: .  -ksp_view_binary - save matrix and right hand side that define linear system to the 
278:                       default binary viewer (can be
279:        read later with src/ksp/examples/tutorials/ex10.c for testing solvers)
280: -  -ksp_view - print the ksp data structure at the end of the system solution

282:    Notes:

284:    The operator is specified with PCSetOperators().

286:    Call KSPGetConvergedReason() to determine if the solver converged or failed and 
287:    why. The number of iterations can be obtained from KSPGetIterationNumber().
288:    
289:    If using a direct method (e.g., via the KSP solver
290:    KSPPREONLY and a preconditioner such as PCLU/PCILU),
291:    then its=1.  See KSPSetTolerances() and KSPDefaultConverged()
292:    for more details.

294:    Understanding Convergence:
295:    The routines KSPSetMonitor(), KSPComputeEigenvalues(), and
296:    KSPComputeEigenvaluesExplicitly() provide information on additional
297:    options to monitor convergence and print eigenvalue information.

299:    Level: beginner

301: .keywords: KSP, solve, linear system

303: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPDefaultConverged(),
304:           KSPSolveTranspose(), KSPGetIterationNumber()
305: @*/
306: PetscErrorCode KSPSolve(KSP ksp,Vec b,Vec x)
307: {
309:   PetscMPIInt    rank;
310:   PetscTruth     flag1,flag2,viewed=PETSC_FALSE,flg;
311:   char           view[10];
312:   PetscScalar    zero = 0.0;

318: 
319:   ksp->vec_rhs = b;
320:   ksp->vec_sol = x;
321:   PetscOptionsHasName(ksp->prefix,"-ksp_view_binary",&flg);
322:   if (flg) {
323:     Mat mat;
324:     PCGetOperators(ksp->pc,&mat,PETSC_NULL,PETSC_NULL);
325:     MatView(mat,PETSC_VIEWER_BINARY_(ksp->comm));
326:     VecView(ksp->vec_rhs,PETSC_VIEWER_BINARY_(ksp->comm));
327:   }

329:   PetscLogEventBegin(KSP_Solve,ksp,ksp->vec_rhs,ksp->vec_sol,0);

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

334:   PetscOptionsGetString(ksp->prefix,"-ksp_view",view,10,&flg);
335:   if (flg) {
336:     PetscStrcmp(view,"before",&viewed);
337:     if (viewed){
338:       KSPView(ksp,PETSC_VIEWER_STDOUT_(ksp->comm));
339:     }
340:   }

342:   /* KSPSetUp() scales the matrix if needed */
343:   KSPSetUp(ksp);
344:   KSPSetUpOnBlocks(ksp);

346:   ksp->transpose_solve = PETSC_FALSE;

348:   /* diagonal scale RHS if called for */
349:   if (ksp->dscale) {
350:     VecPointwiseMult(ksp->diagonal,ksp->vec_rhs,ksp->vec_rhs);
351:     /* second time in, but matrix was scaled back to original */
352:     if (ksp->dscalefix && ksp->dscalefix2) {
353:       Mat mat;

355:       PCGetOperators(ksp->pc,&mat,PETSC_NULL,PETSC_NULL);
356:       MatDiagonalScale(mat,ksp->diagonal,ksp->diagonal);
357:     }
358:   }
359:   PCPreSolve(ksp->pc,ksp);

361:   if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
362:   if (ksp->guess_knoll) {
363:     PCApply(ksp->pc,ksp->vec_rhs,ksp->vec_sol);
364:     KSP_RemoveNullSpace(ksp,ksp->vec_sol);
365:     ksp->guess_zero = PETSC_FALSE;
366:   }
367:   (*ksp->ops->solve)(ksp);
368:   if (!ksp->reason) {
369:     SETERRQ(PETSC_ERR_PLIB,"Internal error, solver returned without setting converged reason");
370:   }
371:   if (ksp->printreason) {
372:     if (ksp->reason > 0) {
373:       PetscPrintf(ksp->comm,"Linear solve converged due to %s\n",convergedreasons[ksp->reason+8]);
374:     } else {
375:       PetscPrintf(ksp->comm,"Linear solve did not converge due to %s\n",convergedreasons[ksp->reason+8]);
376:     }
377:   }

379:   /* diagonal scale solution if called for */
380:   PCPostSolve(ksp->pc,ksp);
381:   if (ksp->dscale) {
382:     VecPointwiseMult(ksp->diagonal,ksp->vec_sol,ksp->vec_sol);
383:     /* unscale right hand side and matrix */
384:     if (ksp->dscalefix) {
385:       Mat mat;

387:       VecReciprocal(ksp->diagonal);
388:       VecPointwiseMult(ksp->diagonal,ksp->vec_rhs,ksp->vec_rhs);
389:       PCGetOperators(ksp->pc,&mat,PETSC_NULL,PETSC_NULL);
390:       MatDiagonalScale(mat,ksp->diagonal,ksp->diagonal);
391:       VecReciprocal(ksp->diagonal);
392:       ksp->dscalefix2 = PETSC_TRUE;
393:     }
394:   }
395:   PetscLogEventEnd(KSP_Solve,ksp,ksp->vec_rhs,ksp->vec_sol,0);

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

399:   PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues",&flag1);
400:   PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues",&flag2);
401:   if (flag1 || flag2) {
402:     PetscInt   nits,n,i,neig;
403:     PetscReal *r,*c;
404: 
405:     KSPGetIterationNumber(ksp,&nits);
406:     n = nits+2;

408:     if (!n) {
409:       PetscPrintf(ksp->comm,"Zero iterations in solver, cannot approximate any eigenvalues\n");
410:     } else {
411:       PetscMalloc(2*n*sizeof(PetscReal),&r);
412:       c = r + n;
413:       KSPComputeEigenvalues(ksp,n,r,c,&neig);
414:       if (flag1) {
415:         PetscPrintf(ksp->comm,"Iteratively computed eigenvalues\n");
416:         for (i=0; i<neig; i++) {
417:           if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gi\n",r[i],c[i]);}
418:           else             {PetscPrintf(ksp->comm,"%g - %gi\n",r[i],-c[i]);}
419:         }
420:       }
421:       if (flag2 && !rank) {
422:         PetscViewer viewer;
423:         PetscDraw   draw;
424:         PetscDrawSP drawsp;

426:         PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Iteratively Computed Eigenvalues",
427:                                PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
428:         PetscViewerDrawGetDraw(viewer,0,&draw);
429:         PetscDrawSPCreate(draw,1,&drawsp);
430:         for (i=0; i<neig; i++) {
431:           PetscDrawSPAddPoint(drawsp,r+i,c+i);
432:         }
433:         PetscDrawSPDraw(drawsp);
434:         PetscDrawSPDestroy(drawsp);
435:         PetscViewerDestroy(viewer);
436:       }
437:       PetscFree(r);
438:     }
439:   }

441:   PetscOptionsHasName(ksp->prefix,"-ksp_compute_eigenvalues_explicitly",&flag1);
442:   PetscOptionsHasName(ksp->prefix,"-ksp_plot_eigenvalues_explicitly",&flag2);
443:   if (flag1 || flag2) {
444:     PetscInt       n,i;
445:     PetscReal *r,*c;
446:     VecGetSize(ksp->vec_sol,&n);
447:     PetscMalloc(2*n*sizeof(PetscReal),&r);
448:     c = r + n;
449:     KSPComputeEigenvaluesExplicitly(ksp,n,r,c);
450:     if (flag1) {
451:       PetscPrintf(ksp->comm,"Explicitly computed eigenvalues\n");
452:       for (i=0; i<n; i++) {
453:         if (c[i] >= 0.0) {PetscPrintf(ksp->comm,"%g + %gi\n",r[i],c[i]);}
454:         else             {PetscPrintf(ksp->comm,"%g - %gi\n",r[i],-c[i]);}
455:       }
456:     }
457:     if (flag2 && !rank) {
458:       PetscViewer viewer;
459:       PetscDraw   draw;
460:       PetscDrawSP drawsp;

462:       PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Explicitly Computed Eigenvalues",0,320,300,300,&viewer);
463:       PetscViewerDrawGetDraw(viewer,0,&draw);
464:       PetscDrawSPCreate(draw,1,&drawsp);
465:       for (i=0; i<n; i++) {
466:         PetscDrawSPAddPoint(drawsp,r+i,c+i);
467:       }
468:       PetscDrawSPDraw(drawsp);
469:       PetscDrawSPDestroy(drawsp);
470:       PetscViewerDestroy(viewer);
471:     }
472:     PetscFree(r);
473:   }

475:   PetscOptionsHasName(ksp->prefix,"-ksp_view_operator",&flag2);
476:   if (flag2) {
477:     Mat A,B;
478:     PCGetOperators(ksp->pc,&A,PETSC_NULL,PETSC_NULL);
479:     MatComputeExplicitOperator(A,&B);
480:     PetscViewerPushFormat(PETSC_VIEWER_STDOUT_(ksp->comm),PETSC_VIEWER_ASCII_MATLAB);
481:     MatView(B,PETSC_VIEWER_STDOUT_(ksp->comm));
482:     PetscViewerPopFormat(PETSC_VIEWER_STDOUT_(ksp->comm));
483:     MatDestroy(B);
484:   }
485:   PetscOptionsHasName(ksp->prefix,"-ksp_view_operator_binary",&flag2);
486:   if (flag2) {
487:     Mat A,B;
488:     PCGetOperators(ksp->pc,&A,PETSC_NULL,PETSC_NULL);
489:     MatComputeExplicitOperator(A,&B);
490:     MatView(B,PETSC_VIEWER_BINARY_(ksp->comm));
491:     MatDestroy(B);
492:   }
493:   PetscOptionsHasName(ksp->prefix,"-ksp_view_preconditioned_operator_binary",&flag2);
494:   if (flag2) {
495:     Mat B;
496:     KSPComputeExplicitOperator(ksp,&B);
497:     MatView(B,PETSC_VIEWER_BINARY_(ksp->comm));
498:     MatDestroy(B);
499:   }
500:   if (!viewed) {
501:     PetscOptionsHasName(ksp->prefix,"-ksp_view",&flg);
502:     if (flg) {
503:       KSPView(ksp,PETSC_VIEWER_STDOUT_(ksp->comm));
504:     }
505:   }
506:   return(0);
507: }

511: /*@
512:    KSPSolveTranspose - Solves the transpose of a linear system. Usually
513:    accessed through KSPSolveTranspose().

515:    Collective on KSP

517:    Input Parameter:
518: +  ksp - iterative context obtained from KSPCreate()
519: .  b - right hand side vector
520: -  x - solution vector

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

526:    Level: developer

528: .keywords: KSP, solve, linear system

530: .seealso: KSPCreate(), KSPSetUp(), KSPDestroy(), KSPSetTolerances(), KSPDefaultConverged(),
531:           KSPSolve()
532: @*/
533: PetscErrorCode KSPSolveTranspose(KSP ksp,Vec b,Vec x)
534: {
536:   PetscScalar    zero = 0.0;


543:   ksp->vec_rhs = b;
544:   ksp->vec_sol = x;
545:   KSPSetUp(ksp);
546:   if (ksp->guess_zero) { VecSet(&zero,ksp->vec_sol);}
547:   ksp->transpose_solve = PETSC_TRUE;
548:   (*ksp->ops->solve)(ksp);
549:   return(0);
550: }

554: /*@C
555:    KSPDestroy - Destroys KSP context.

557:    Collective on KSP

559:    Input Parameter:
560: .  ksp - iterative context obtained from KSPCreate()

562:    Level: beginner

564: .keywords: KSP, destroy

566: .seealso: KSPCreate(), KSPSetUp(), KSPSolve()
567: @*/
568: PetscErrorCode KSPDestroy(KSP ksp)
569: {
571:   PetscInt       i;

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

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

580:   if (ksp->ops->destroy) {
581:     (*ksp->ops->destroy)(ksp);
582:   }
583:   for (i=0; i<ksp->numbermonitors; i++) {
584:     if (ksp->monitordestroy[i]) {
585:       (*ksp->monitordestroy[i])(ksp->monitorcontext[i]);
586:     }
587:   }
588:   PCDestroy(ksp->pc);
589:   if (ksp->diagonal) {VecDestroy(ksp->diagonal);}
590:   if (ksp->nullsp) {MatNullSpaceDestroy(ksp->nullsp);}
591:   PetscLogObjectDestroy(ksp);
592:   PetscHeaderDestroy(ksp);
593:   return(0);
594: }

598: /*@
599:     KSPSetPreconditionerSide - Sets the preconditioning side.

601:     Collective on KSP

603:     Input Parameter:
604: .   ksp - iterative context obtained from KSPCreate()

606:     Output Parameter:
607: .   side - the preconditioning side, where side is one of
608: .vb
609:       PC_LEFT - left preconditioning (default)
610:       PC_RIGHT - right preconditioning
611:       PC_SYMMETRIC - symmetric preconditioning
612: .ve

614:     Options Database Keys:
615: +   -ksp_left_pc - Sets left preconditioning
616: .   -ksp_right_pc - Sets right preconditioning
617: -   -ksp_symmetric_pc - Sets symmetric preconditioning

619:     Notes:
620:     Left preconditioning is used by default.  Symmetric preconditioning is
621:     currently available only for the KSPQCG method. Note, however, that
622:     symmetric preconditioning can be emulated by using either right or left
623:     preconditioning and a pre or post processing step.

625:     Level: intermediate

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

629: .seealso: KSPGetPreconditionerSide()
630: @*/
631: PetscErrorCode KSPSetPreconditionerSide(KSP ksp,PCSide side)
632: {
635:   ksp->pc_side = side;
636:   return(0);
637: }

641: /*@C
642:     KSPGetPreconditionerSide - Gets the preconditioning side.

644:     Not Collective

646:     Input Parameter:
647: .   ksp - iterative context obtained from KSPCreate()

649:     Output Parameter:
650: .   side - the preconditioning side, where side is one of
651: .vb
652:       PC_LEFT - left preconditioning (default)
653:       PC_RIGHT - right preconditioning
654:       PC_SYMMETRIC - symmetric preconditioning
655: .ve

657:     Level: intermediate

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

661: .seealso: KSPSetPreconditionerSide()
662: @*/
663: PetscErrorCode KSPGetPreconditionerSide(KSP ksp,PCSide *side)
664: {
668:   *side = ksp->pc_side;
669:   return(0);
670: }

674: /*@
675:    KSPGetTolerances - Gets the relative, absolute, divergence, and maximum
676:    iteration tolerances used by the default KSP convergence tests. 

678:    Not Collective

680:    Input Parameter:
681: .  ksp - the Krylov subspace context
682:   
683:    Output Parameters:
684: +  rtol - the relative convergence tolerance
685: .  abstol - the absolute convergence tolerance
686: .  dtol - the divergence tolerance
687: -  maxits - maximum number of iterations

689:    Notes:
690:    The user can specify PETSC_NULL for any parameter that is not needed.

692:    Level: intermediate

694: .keywords: KSP, get, tolerance, absolute, relative, divergence, convergence,
695:            maximum, iterations

697: .seealso: KSPSetTolerances()
698: @*/
699: PetscErrorCode KSPGetTolerances(KSP ksp,PetscReal *rtol,PetscReal *abstol,PetscReal *dtol,PetscInt *maxits)
700: {
703:   if (abstol)   *abstol   = ksp->abstol;
704:   if (rtol)   *rtol   = ksp->rtol;
705:   if (dtol)   *dtol   = ksp->divtol;
706:   if (maxits) *maxits = ksp->max_it;
707:   return(0);
708: }

712: /*@
713:    KSPSetTolerances - Sets the relative, absolute, divergence, and maximum
714:    iteration tolerances used by the default KSP convergence testers. 

716:    Collective on KSP

718:    Input Parameters:
719: +  ksp - the Krylov subspace context
720: .  rtol - the relative convergence tolerance
721:    (relative decrease in the residual norm)
722: .  abstol - the absolute convergence tolerance 
723:    (absolute size of the residual norm)
724: .  dtol - the divergence tolerance
725:    (amount residual can increase before KSPDefaultConverged()
726:    concludes that the method is diverging)
727: -  maxits - maximum number of iterations to use

729:    Options Database Keys:
730: +  -ksp_atol <abstol> - Sets abstol
731: .  -ksp_rtol <rtol> - Sets rtol
732: .  -ksp_divtol <dtol> - Sets dtol
733: -  -ksp_max_it <maxits> - Sets maxits

735:    Notes:
736:    Use PETSC_DEFAULT to retain the default value of any of the tolerances.

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

742:    Level: intermediate

744: .keywords: KSP, set, tolerance, absolute, relative, divergence, 
745:            convergence, maximum, iterations

747: .seealso: KSPGetTolerances(), KSPDefaultConverged(), KSPSetConvergenceTest()
748: @*/
749: PetscErrorCode KSPSetTolerances(KSP ksp,PetscReal rtol,PetscReal abstol,PetscReal dtol,PetscInt maxits)
750: {
753:   if (abstol != PETSC_DEFAULT)   ksp->abstol   = abstol;
754:   if (rtol != PETSC_DEFAULT)   ksp->rtol   = rtol;
755:   if (dtol != PETSC_DEFAULT)   ksp->divtol = dtol;
756:   if (maxits != PETSC_DEFAULT) ksp->max_it = maxits;
757:   return(0);
758: }

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

767:    Collective on KSP

769:    Input Parameters:
770: +  ksp - iterative context obtained from KSPCreate()
771: -  flg - PETSC_TRUE indicates the guess is non-zero, PETSC_FALSE indicates the guess is zero

773:    Level: beginner

775:    Notes:
776:     If this is not called the X vector is zeroed in the call to KSPSolve().

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

780: .seealso: KSPGetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
781: @*/
782: PetscErrorCode KSPSetInitialGuessNonzero(KSP ksp,PetscTruth flg)
783: {
785:   ksp->guess_zero   = (PetscTruth)!(int)flg;
786:   return(0);
787: }

791: /*@
792:    KSPGetInitialGuessNonzero - Determines whether the KSP solver is using
793:    a zero initial guess.

795:    Not Collective

797:    Input Parameter:
798: .  ksp - iterative context obtained from KSPCreate()

800:    Output Parameter:
801: .  flag - PETSC_TRUE if guess is nonzero, else PETSC_FALSE

803:    Level: intermediate

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

807: .seealso: KSPSetInitialGuessNonzero(), KSPSetInitialGuessKnoll(), KSPGetInitialGuessKnoll()
808: @*/
809: PetscErrorCode KSPGetInitialGuessNonzero(KSP ksp,PetscTruth *flag)
810: {
812:   if (ksp->guess_zero) *flag = PETSC_FALSE;
813:   else                 *flag = PETSC_TRUE;
814:   return(0);
815: }

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

822:    Collective on KSP

824:    Input Parameters:
825: +  ksp - iterative context obtained from KSPCreate()
826: -  flg - PETSC_TRUE or PETSC_FALSE

828:    Level: advanced


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

833: .seealso: KSPGetInitialGuessKnoll(), KSPSetInitialGuessNonzero(), KSPGetInitialGuessNonzero()
834: @*/
835: PetscErrorCode KSPSetInitialGuessKnoll(KSP ksp,PetscTruth flg)
836: {
838:   ksp->guess_knoll   = flg;
839:   return(0);
840: }

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

848:    Not Collective

850:    Input Parameter:
851: .  ksp - iterative context obtained from KSPCreate()

853:    Output Parameter:
854: .  flag - PETSC_TRUE if using Knoll trick, else PETSC_FALSE

856:    Level: advanced

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

860: .seealso: KSPSetInitialGuessKnoll(), KSPSetInitialGuessNonzero(), KSPGetInitialGuessNonzero()
861: @*/
862: PetscErrorCode KSPGetInitialGuessKnoll(KSP ksp,PetscTruth *flag)
863: {
865:   *flag = ksp->guess_knoll;
866:   return(0);
867: }

871: /*@
872:    KSPSetComputeSingularValues - Sets a flag so that the extreme singular 
873:    values will be calculated via a Lanczos or Arnoldi process as the linear 
874:    system is solved.

876:    Collective on KSP

878:    Input Parameters:
879: +  ksp - iterative context obtained from KSPCreate()
880: -  flg - PETSC_TRUE or PETSC_FALSE

882:    Options Database Key:
883: .  -ksp_singmonitor - Activates KSPSetComputeSingularValues()

885:    Notes:
886:    Currently this option is not valid for all iterative methods.

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

892:    Level: advanced

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

896: .seealso: KSPComputeExtremeSingularValues(), KSPSingularValueMonitor()
897: @*/
898: PetscErrorCode KSPSetComputeSingularValues(KSP ksp,PetscTruth flg)
899: {
902:   ksp->calc_sings  = flg;
903:   return(0);
904: }

908: /*@
909:    KSPSetComputeEigenvalues - Sets a flag so that the extreme eigenvalues
910:    values will be calculated via a Lanczos or Arnoldi process as the linear 
911:    system is solved.

913:    Collective on KSP

915:    Input Parameters:
916: +  ksp - iterative context obtained from KSPCreate()
917: -  flg - PETSC_TRUE or PETSC_FALSE

919:    Notes:
920:    Currently this option is not valid for all iterative methods.

922:    Level: advanced

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

926: .seealso: KSPComputeEigenvalues(), KSPComputeEigenvaluesExplicitly()
927: @*/
928: PetscErrorCode KSPSetComputeEigenvalues(KSP ksp,PetscTruth flg)
929: {
932:   ksp->calc_sings  = flg;
933:   return(0);
934: }

938: /*@C
939:    KSPGetRhs - Gets the right-hand-side vector for the linear system to
940:    be solved.

942:    Not Collective

944:    Input Parameter:
945: .  ksp - iterative context obtained from KSPCreate()

947:    Output Parameter:
948: .  r - right-hand-side vector

950:    Level: developer

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

954: .seealso: KSPGetSolution(), KSPSolve()
955: @*/
956: PetscErrorCode KSPGetRhs(KSP ksp,Vec *r)
957: {
960:   *r = ksp->vec_rhs;
961:   return(0);
962: }

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

971:    Not Collective

973:    Input Parameters:
974: .  ksp - iterative context obtained from KSPCreate()

976:    Output Parameters:
977: .  v - solution vector

979:    Level: developer

981: .keywords: KSP, get, solution

983: .seealso: KSPGetRhs(),  KSPBuildSolution(), KSPSolve()
984: @*/
985: PetscErrorCode KSPGetSolution(KSP ksp,Vec *v)
986: {
990:   *v = ksp->vec_sol;
991:   return(0);
992: }

996: /*@
997:    KSPSetPC - Sets the preconditioner to be used to calculate the 
998:    application of the preconditioner on a vector. 

1000:    Collective on KSP

1002:    Input Parameters:
1003: +  ksp - iterative context obtained from KSPCreate()
1004: -  pc   - the preconditioner object

1006:    Notes:
1007:    Use KSPGetPC() to retrieve the preconditioner context (for example,
1008:    to free it at the end of the computations).

1010:    Level: developer

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

1014: .seealso: KSPGetPC()
1015: @*/
1016: PetscErrorCode KSPSetPC(KSP ksp,PC pc)
1017: {

1024:   if (ksp->pc) {PCDestroy(ksp->pc);}
1025:   ksp->pc = pc;
1026:   PetscObjectReference((PetscObject)ksp->pc);
1027:   return(0);
1028: }

1032: /*@C
1033:    KSPGetPC - Returns a pointer to the preconditioner context
1034:    set with KSPSetPC().

1036:    Not Collective

1038:    Input Parameters:
1039: .  ksp - iterative context obtained from KSPCreate()

1041:    Output Parameter:
1042: .  pc - preconditioner context

1044:    Level: developer

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

1048: .seealso: KSPSetPC()
1049: @*/
1050: PetscErrorCode KSPGetPC(KSP ksp,PC *pc)
1051: {
1055:   *pc = ksp->pc;
1056:   return(0);
1057: }

1061: /*@C
1062:    KSPSetMonitor - Sets an ADDITIONAL function to be called at every iteration to monitor 
1063:    the residual/error etc.
1064:       
1065:    Collective on KSP

1067:    Input Parameters:
1068: +  ksp - iterative context obtained from KSPCreate()
1069: .  monitor - pointer to function (if this is PETSC_NULL, it turns off monitoring
1070: .  mctx    - [optional] context for private data for the
1071:              monitor routine (use PETSC_NULL if no context is desired)
1072: -  monitordestroy - [optional] routine that frees monitor context
1073:           (may be PETSC_NULL)

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

1078: +  ksp - iterative context obtained from KSPCreate()
1079: .  it - iteration number
1080: .  rnorm - (estimated) 2-norm of (preconditioned) residual
1081: -  mctx  - optional monitoring context, as set by KSPSetMonitor()

1083:    Options Database Keys:
1084: +    -ksp_monitor        - sets KSPDefaultMonitor()
1085: .    -ksp_truemonitor    - sets KSPTrueMonitor()
1086: .    -ksp_xmonitor       - sets line graph monitor,
1087:                            uses KSPLGMonitorCreate()
1088: .    -ksp_xtruemonitor   - sets line graph monitor,
1089:                            uses KSPLGMonitorCreate()
1090: .    -ksp_singmonitor    - sets KSPSingularValueMonitor()
1091: -    -ksp_cancelmonitors - cancels all monitors that have
1092:                           been hardwired into a code by 
1093:                           calls to KSPSetMonitor(), but
1094:                           does not cancel those set via
1095:                           the options database.

1097:    Notes:  
1098:    The default is to do nothing.  To print the residual, or preconditioned 
1099:    residual if KSPSetNormType(ksp,KSP_PRECONDITIONED_NORM) was called, use 
1100:    KSPDefaultMonitor() as the monitoring routine, with a null monitoring 
1101:    context. 

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

1107:    Level: beginner

1109: .keywords: KSP, set, monitor

1111: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPClearMonitor()
1112: @*/
1113: PetscErrorCode KSPSetMonitor(KSP ksp,PetscErrorCode (*monitor)(KSP,PetscInt,PetscReal,void*),void *mctx,PetscErrorCode (*monitordestroy)(void*))
1114: {
1117:   if (ksp->numbermonitors >= MAXKSPMONITORS) {
1118:     SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Too many KSP monitors set");
1119:   }
1120:   ksp->monitor[ksp->numbermonitors]           = monitor;
1121:   ksp->monitordestroy[ksp->numbermonitors]    = monitordestroy;
1122:   ksp->monitorcontext[ksp->numbermonitors++]  = (void*)mctx;
1123:   return(0);
1124: }

1128: /*@
1129:    KSPClearMonitor - Clears all monitors for a KSP object.

1131:    Collective on KSP

1133:    Input Parameters:
1134: .  ksp - iterative context obtained from KSPCreate()

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

1141:    Level: intermediate

1143: .keywords: KSP, set, monitor

1145: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate(), KSPSetMonitor()
1146: @*/
1147: PetscErrorCode KSPClearMonitor(KSP ksp)
1148: {
1151:   ksp->numbermonitors = 0;
1152:   return(0);
1153: }

1157: /*@C
1158:    KSPGetMonitorContext - Gets the monitoring context, as set by 
1159:    KSPSetMonitor() for the FIRST monitor only.

1161:    Not Collective

1163:    Input Parameter:
1164: .  ksp - iterative context obtained from KSPCreate()

1166:    Output Parameter:
1167: .  ctx - monitoring context

1169:    Level: intermediate

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

1173: .seealso: KSPDefaultMonitor(), KSPLGMonitorCreate()
1174: @*/
1175: PetscErrorCode KSPGetMonitorContext(KSP ksp,void **ctx)
1176: {
1179:   *ctx =      (ksp->monitorcontext[0]);
1180:   return(0);
1181: }

1185: /*@
1186:    KSPSetResidualHistory - Sets the array used to hold the residual history.
1187:    If set, this array will contain the residual norms computed at each
1188:    iteration of the solver.

1190:    Not Collective

1192:    Input Parameters:
1193: +  ksp - iterative context obtained from KSPCreate()
1194: .  a   - array to hold history
1195: .  na  - size of a
1196: -  reset - PETSC_TRUE indicates the history counter is reset to zero
1197:            for each new linear solve

1199:    Level: advanced

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

1204:    If 'na' is PETSC_DECIDE or 'a' is PETSC_NULL, then a default array of
1205:    length 1000 is allocated.

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

1209: .seealso: KSPGetResidualHistory()

1211: @*/
1212: PetscErrorCode KSPSetResidualHistory(KSP ksp,PetscReal a[],PetscInt na,PetscTruth reset)
1213: {

1218:   if (na != PETSC_DECIDE && a != PETSC_NULL) {
1219:     ksp->res_hist        = a;
1220:     ksp->res_hist_max    = na;
1221:   } else {
1222:     ksp->res_hist_max    = 1000;
1223:     PetscMalloc(ksp->res_hist_max*sizeof(PetscReal),&ksp->res_hist);
1224:   }
1225:   ksp->res_hist_len    = 0;
1226:   ksp->res_hist_reset  = reset;


1229:   return(0);
1230: }

1234: /*@C
1235:    KSPGetResidualHistory - Gets the array used to hold the residual history
1236:    and the number of residuals it contains.

1238:    Not Collective

1240:    Input Parameter:
1241: .  ksp - iterative context obtained from KSPCreate()

1243:    Output Parameters:
1244: +  a   - pointer to array to hold history (or PETSC_NULL)
1245: -  na  - number of used entries in a (or PETSC_NULL)

1247:    Level: advanced

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

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

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

1257: .seealso: KSPGetResidualHistory()

1259: @*/
1260: PetscErrorCode KSPGetResidualHistory(KSP ksp,PetscReal *a[],PetscInt *na)
1261: {
1264:   if (a)  *a = ksp->res_hist;
1265:   if (na) *na = ksp->res_hist_len;
1266:   return(0);
1267: }

1271: /*@C
1272:    KSPSetConvergenceTest - Sets the function to be used to determine
1273:    convergence.  

1275:    Collective on KSP

1277:    Input Parameters:
1278: +  ksp - iterative context obtained from KSPCreate()
1279: .  converge - pointer to int function
1280: -  cctx    - context for private data for the convergence routine (may be null)

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

1285: +  ksp - iterative context obtained from KSPCreate()
1286: .  it - iteration number
1287: .  rnorm - (estimated) 2-norm of (preconditioned) residual
1288: .  reason - the reason why it has converged or diverged
1289: -  cctx  - optional convergence context, as set by KSPSetConvergenceTest()


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

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

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

1303:    In the default PETSc convergence test, the precise values of reason
1304:    are macros such as KSP_CONVERGED_RTOL, which are defined in petscksp.h.

1306:    Level: advanced

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

1310: .seealso: KSPDefaultConverged(), KSPGetConvergenceContext()
1311: @*/
1312: PetscErrorCode KSPSetConvergenceTest(KSP ksp,PetscErrorCode (*converge)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*),void *cctx)
1313: {
1316:   ksp->converged = converge;
1317:   ksp->cnvP      = (void*)cctx;
1318:   return(0);
1319: }

1323: /*@C
1324:    KSPGetConvergenceContext - Gets the convergence context set with 
1325:    KSPSetConvergenceTest().  

1327:    Not Collective

1329:    Input Parameter:
1330: .  ksp - iterative context obtained from KSPCreate()

1332:    Output Parameter:
1333: .  ctx - monitoring context

1335:    Level: advanced

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

1339: .seealso: KSPDefaultConverged(), KSPSetConvergenceTest()
1340: @*/
1341: PetscErrorCode KSPGetConvergenceContext(KSP ksp,void **ctx)
1342: {
1345:   *ctx = ksp->cnvP;
1346:   return(0);
1347: }

1351: /*@C
1352:    KSPBuildSolution - Builds the approximate solution in a vector provided.
1353:    This routine is NOT commonly needed (see KSPSolve()).

1355:    Collective on KSP

1357:    Input Parameter:
1358: .  ctx - iterative context obtained from KSPCreate()

1360:    Output Parameter: 
1361:    Provide exactly one of
1362: +  v - location to stash solution.   
1363: -  V - the solution is returned in this location. This vector is created 
1364:        internally. This vector should NOT be destroyed by the user with
1365:        VecDestroy().

1367:    Notes:
1368:    This routine can be used in one of two ways
1369: .vb
1370:       KSPBuildSolution(ksp,PETSC_NULL,&V);
1371:    or
1372:       KSPBuildSolution(ksp,v,PETSC_NULL); 
1373: .ve
1374:    In the first case an internal vector is allocated to store the solution
1375:    (the user cannot destroy this vector). In the second case the solution
1376:    is generated in the vector that the user provides. Note that for certain 
1377:    methods, such as KSPCG, the second case requires a copy of the solution,
1378:    while in the first case the call is essentially free since it simply 
1379:    returns the vector where the solution already is stored.

1381:    Level: advanced

1383: .keywords: KSP, build, solution

1385: .seealso: KSPGetSolution(), KSPBuildResidual()
1386: @*/
1387: PetscErrorCode KSPBuildSolution(KSP ksp,Vec v,Vec *V)
1388: {

1393:   if (!V && !v) SETERRQ(PETSC_ERR_ARG_WRONG,"Must provide either v or V");
1394:   if (!V) V = &v;
1395:   (*ksp->ops->buildsolution)(ksp,v,V);
1396:   return(0);
1397: }

1401: /*@C
1402:    KSPBuildResidual - Builds the residual in a vector provided.

1404:    Collective on KSP

1406:    Input Parameter:
1407: .  ksp - iterative context obtained from KSPCreate()

1409:    Output Parameters:
1410: +  v - optional location to stash residual.  If v is not provided,
1411:        then a location is generated.
1412: .  t - work vector.  If not provided then one is generated.
1413: -  V - the residual

1415:    Notes:
1416:    Regardless of whether or not v is provided, the residual is 
1417:    returned in V.

1419:    Level: advanced

1421: .keywords: KSP, build, residual

1423: .seealso: KSPBuildSolution()
1424: @*/
1425: PetscErrorCode KSPBuildResidual(KSP ksp,Vec t,Vec v,Vec *V)
1426: {
1428:   PetscTruth     flag = PETSC_FALSE;
1429:   Vec            w = v,tt = t;

1433:   if (!w) {
1434:     VecDuplicate(ksp->vec_rhs,&w);
1435:     PetscLogObjectParent((PetscObject)ksp,w);
1436:   }
1437:   if (!tt) {
1438:     VecDuplicate(ksp->vec_rhs,&tt); flag = PETSC_TRUE;
1439:     PetscLogObjectParent((PetscObject)ksp,tt);
1440:   }
1441:   (*ksp->ops->buildresidual)(ksp,tt,w,V);
1442:   if (flag) {VecDestroy(tt);}
1443:   return(0);
1444: }

1448: /*@
1449:    KSPSetDiagonalScale - Tells KSP to symmetrically diagonally scale the system
1450:      before solving. This actually CHANGES the matrix (and right hand side).

1452:    Collective on KSP

1454:    Input Parameter:
1455: +  ksp - the KSP context
1456: -  scale - PETSC_TRUE or PETSC_FALSE

1458:    Options Database Key:
1459: +   -ksp_diagonal_scale - 
1460: -   -ksp_diagonal_scale_fix - scale the matrix back AFTER the solve 


1463:     BE CAREFUL with this routine: it actually scales the matrix and right 
1464:     hand side that define the system. After the system is solved the matrix
1465:     and right hand side remain scaled.

1467:     This routine is only used if the matrix and preconditioner matrix are
1468:     the same thing.
1469:  
1470:     If you use this with the PCType Eisenstat preconditioner than you can 
1471:     use the PCEisenstatNoDiagonalScaling() option, or -pc_eisenstat_no_diagonal_scaling
1472:     to save some unneeded, redundant flops.

1474:    Level: intermediate

1476: .keywords: KSP, set, options, prefix, database

1478: .seealso: KSPGetDiagonalScale(), KSPSetDiagonalScaleFix()
1479: @*/
1480: PetscErrorCode KSPSetDiagonalScale(KSP ksp,PetscTruth scale)
1481: {
1484:   ksp->dscale = scale;
1485:   return(0);
1486: }

1490: /*@C
1491:    KSPGetDiagonalScale - Checks if KSP solver scales the matrix and
1492:                           right hand side

1494:    Not Collective

1496:    Input Parameter:
1497: .  ksp - the KSP context

1499:    Output Parameter:
1500: .  scale - PETSC_TRUE or PETSC_FALSE

1502:    Notes:
1503:     BE CAREFUL with this routine: it actually scales the matrix and right 
1504:     hand side that define the system. After the system is solved the matrix
1505:     and right hand side remain scaled.

1507:     This routine is only used if the matrix and preconditioner matrix are
1508:     the same thing.

1510:    Level: intermediate

1512: .keywords: KSP, set, options, prefix, database

1514: .seealso: KSPSetDiagonalScale(), KSPSetDiagonalScaleFix()
1515: @*/
1516: PetscErrorCode KSPGetDiagonalScale(KSP ksp,PetscTruth *scale)
1517: {
1521:   *scale = ksp->dscale;
1522:   return(0);
1523: }

1527: /*@
1528:    KSPSetDiagonalScaleFix - Tells KSP to diagonally scale the system
1529:      back after solving.

1531:    Collective on KSP

1533:    Input Parameter:
1534: +  ksp - the KSP context
1535: -  fix - PETSC_TRUE to scale back after the system solve, PETSC_FALSE to not 
1536:          rescale (default)

1538:    Notes:
1539:      Must be called after KSPSetDiagonalScale()

1541:      Using this will slow things down, because it rescales the matrix before and
1542:      after each linear solve. This is intended mainly for testing to allow one
1543:      to easily get back the original system to make sure the solution computed is
1544:      accurate enough.

1546:     This routine is only used if the matrix and preconditioner matrix are
1547:     the same thing.

1549:    Level: intermediate

1551: .keywords: KSP, set, options, prefix, database

1553: .seealso: KSPGetDiagonalScale(), KSPSetDiagonalScale(), KSPGetDiagonalScaleFix()
1554: @*/
1555: PetscErrorCode KSPSetDiagonalScaleFix(KSP ksp,PetscTruth fix)
1556: {
1559:   ksp->dscalefix = fix;
1560:   return(0);
1561: }

1565: /*@
1566:    KSPGetDiagonalScaleFix - Determines if KSP diagonally scales the system
1567:      back after solving.

1569:    Collective on KSP

1571:    Input Parameter:
1572: .  ksp - the KSP context

1574:    Output Parameter:
1575: .  fix - PETSC_TRUE to scale back after the system solve, PETSC_FALSE to not 
1576:          rescale (default)

1578:    Notes:
1579:      Must be called after KSPSetDiagonalScale()

1581:      If PETSC_TRUE will slow things down, because it rescales the matrix before and
1582:      after each linear solve. This is intended mainly for testing to allow one
1583:      to easily get back the original system to make sure the solution computed is
1584:      accurate enough.

1586:     This routine is only used if the matrix and preconditioner matrix are
1587:     the same thing.

1589:    Level: intermediate

1591: .keywords: KSP, set, options, prefix, database

1593: .seealso: KSPGetDiagonalScale(), KSPSetDiagonalScale(), KSPSetDiagonalScaleFix()
1594: @*/
1595: PetscErrorCode KSPGetDiagonalScaleFix(KSP ksp,PetscTruth *fix)
1596: {
1600:   *fix = ksp->dscalefix;
1601:   return(0);
1602: }