Actual source code: iterativ.c

  1: #define PETSCKSP_DLL

  3: /*
  4:    This file contains some simple default routines.  
  5:    These routines should be SHORT, since they will be included in every
  6:    executable image that uses the iterative routines (note that, through
  7:    the registry system, we provide a way to load only the truely necessary
  8:    files) 
  9:  */
 10:  #include include/private/kspimpl.h

 14: /*
 15:   KSPDefaultFreeWork - Free work vectors

 17:   Input Parameters:
 18: . ksp  - iterative context
 19:  */
 20: PetscErrorCode KSPDefaultFreeWork(KSP ksp)
 21: {
 25:   if (ksp->work)  {
 26:     VecDestroyVecs(ksp->work,ksp->nwork);
 27:     ksp->work = PETSC_NULL;
 28:   }
 29:   return(0);
 30: }

 34: /*@
 35:    KSPGetResidualNorm - Gets the last (approximate preconditioned)
 36:    residual norm that has been computed.
 37:  
 38:    Not Collective

 40:    Input Parameters:
 41: .  ksp - the iterative context

 43:    Output Parameters:
 44: .  rnorm - residual norm

 46:    Level: intermediate

 48: .keywords: KSP, get, residual norm

 50: .seealso: KSPBuildResidual()
 51: @*/
 52: PetscErrorCode  KSPGetResidualNorm(KSP ksp,PetscReal *rnorm)
 53: {
 57:   *rnorm = ksp->rnorm;
 58:   return(0);
 59: }

 63: /*@
 64:    KSPGetIterationNumber - Gets the current iteration number; if the 
 65:          KSPSolve() is complete, returns the number of iterations
 66:          used.
 67:  
 68:    Not Collective

 70:    Input Parameters:
 71: .  ksp - the iterative context

 73:    Output Parameters:
 74: .  its - number of iterations

 76:    Level: intermediate

 78:    Notes:
 79:       During the ith iteration this returns i-1
 80: .keywords: KSP, get, residual norm

 82: .seealso: KSPBuildResidual(), KSPGetResidualNorm()
 83: @*/
 84: PetscErrorCode  KSPGetIterationNumber(KSP ksp,PetscInt *its)
 85: {
 89:   *its = ksp->its;
 90:   return(0);
 91: }

 95: /*@C
 96:     KSPMonitorSingularValue - Prints the two norm of the true residual and
 97:     estimation of the extreme singular values of the preconditioned problem
 98:     at each iteration.
 99:  
100:     Collective on KSP

102:     Input Parameters:
103: +   ksp - the iterative context
104: .   n  - the iteration
105: -   rnorm - the two norm of the residual

107:     Options Database Key:
108: .   -ksp_monitor_singular_value - Activates KSPMonitorSingularValue()

110:     Notes:
111:     The CG solver uses the Lanczos technique for eigenvalue computation, 
112:     while GMRES uses the Arnoldi technique; other iterative methods do
113:     not currently compute singular values.

115:     Level: intermediate

117: .keywords: KSP, CG, default, monitor, extreme, singular values, Lanczos, Arnoldi

119: .seealso: KSPComputeExtremeSingularValues()
120: @*/
121: PetscErrorCode  KSPMonitorSingularValue(KSP ksp,PetscInt n,PetscReal rnorm,void *dummy)
122: {
123:   PetscReal               emin,emax,c;
124:   PetscErrorCode          ierr;
125:   PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;

129:   if (!dummy) {PetscViewerASCIIMonitorCreate(ksp->comm,"stdout",0,&viewer);}
130:   if (!ksp->calc_sings) {
131:     PetscViewerASCIIMonitorPrintf(viewer,"%3D KSP Residual norm %14.12e \n",n,rnorm);
132:   } else {
133:     KSPComputeExtremeSingularValues(ksp,&emax,&emin);
134:     c = emax/emin;
135:     PetscViewerASCIIMonitorPrintf(viewer,"%3D KSP Residual norm %14.12e %% max %G min %G max/min %G\n",n,rnorm,emax,emin,c);
136:   }
137:   if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
138:   return(0);
139: }

143: /*@C
144:    KSPMonitorSolution - Monitors progress of the KSP solvers by calling 
145:    VecView() for the approximate solution at each iteration.

147:    Collective on KSP

149:    Input Parameters:
150: +  ksp - the KSP context
151: .  its - iteration number
152: .  fgnorm - 2-norm of residual (or gradient)
153: -  dummy - either a viewer or PETSC_NULL

155:    Level: intermediate

157:    Notes:
158:     For some Krylov methods such as GMRES constructing the solution at
159:   each iteration is expensive, hence using this will slow the code.

161: .keywords: KSP, nonlinear, vector, monitor, view

163: .seealso: KSPMonitorSet(), KSPMonitorDefault(), VecView()
164: @*/
165: PetscErrorCode  KSPMonitorSolution(KSP ksp,PetscInt its,PetscReal fgnorm,void *dummy)
166: {
168:   Vec            x;
169:   PetscViewer    viewer = (PetscViewer) dummy;

172:   KSPBuildSolution(ksp,PETSC_NULL,&x);
173:   if (!viewer) {
174:     MPI_Comm comm;
175:     PetscObjectGetComm((PetscObject)ksp,&comm);
176:     viewer = PETSC_VIEWER_DRAW_(comm);
177:   }
178:   VecView(x,viewer);

180:   return(0);
181: }

185: /*@C
186:    KSPMonitorDefault - Print the residual norm at each iteration of an
187:    iterative solver.

189:    Collective on KSP

191:    Input Parameters:
192: +  ksp   - iterative context
193: .  n     - iteration number
194: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).  
195: -  dummy - unused monitor context 

197:    Level: intermediate

199: .keywords: KSP, default, monitor, residual

201: .seealso: KSPMonitorSet(), KSPMonitorTrueResidualNorm(), KSPMonitorLGCreate()
202: @*/
203: PetscErrorCode  KSPMonitorDefault(KSP ksp,PetscInt n,PetscReal rnorm,void *dummy)
204: {
205:   PetscErrorCode          ierr;
206:   PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;

209:   if (!dummy) {PetscViewerASCIIMonitorCreate(ksp->comm,"stdout",0,&viewer);}
210:   PetscViewerASCIIMonitorPrintf(viewer,"%3D KSP Residual norm %14.12e \n",n,rnorm);
211:   if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
212:   return(0);
213: }

217: /*@C
218:    KSPMonitorTrueResidualNorm - Prints the true residual norm as well as the preconditioned
219:    residual norm at each iteration of an iterative solver.

221:    Collective on KSP

223:    Input Parameters:
224: +  ksp   - iterative context
225: .  n     - iteration number
226: .  rnorm - 2-norm (preconditioned) residual value (may be estimated).  
227: -  dummy - unused monitor context 

229:    Options Database Key:
230: .  -ksp_monitor_true_residual - Activates KSPMonitorTrueResidualNorm()

232:    Notes:
233:    When using right preconditioning, these values are equivalent.

235:    When using either ICC or ILU preconditioners in BlockSolve95 
236:    (via MATMPIROWBS matrix format), then use this monitor will
237:    print both the residual norm associated with the original
238:    (unscaled) matrix.

240:    Level: intermediate

242: .keywords: KSP, default, monitor, residual

244: .seealso: KSPMonitorSet(), KSPMonitorDefault(), KSPMonitorLGCreate()
245: @*/
246: PetscErrorCode  KSPMonitorTrueResidualNorm(KSP ksp,PetscInt n,PetscReal rnorm,void *dummy)
247: {
248:   PetscErrorCode          ierr;
249:   Vec                     resid,work;
250:   PetscReal               scnorm,bnorm;
251:   PC                      pc;
252:   Mat                     A,B;
253:   PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;
254: 
256:   if (!dummy) {PetscViewerASCIIMonitorCreate(ksp->comm,"stdout",0,&viewer);}
257:   VecDuplicate(ksp->vec_rhs,&work);
258:   KSPBuildResidual(ksp,0,work,&resid);

260:   /*
261:      Unscale the residual if the matrix is, for example, a BlockSolve matrix
262:     but only if both matrices are the same matrix, since only then would 
263:     they be scaled.
264:   */
265:   VecCopy(resid,work);
266:   KSPGetPC(ksp,&pc);
267:   PCGetOperators(pc,&A,&B,PETSC_NULL);
268:   if (A == B) {
269:     MatUnScaleSystem(A,work,PETSC_NULL);
270:   }
271:   VecNorm(work,NORM_2,&scnorm);
272:   VecDestroy(work);
273:   VecNorm(ksp->vec_rhs,NORM_2,&bnorm);
274:   PetscViewerASCIIMonitorPrintf(viewer,"%3D KSP preconditioned resid norm %14.12e true resid norm %14.12e ||Ae||/||Ax|| %14.12e\n",n,rnorm,scnorm,scnorm/bnorm);
275:   if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
276:   return(0);
277: }

281: /*
282:   Default (short) KSP Monitor, same as KSPMonitorDefault() except
283:   it prints fewer digits of the residual as the residual gets smaller.
284:   This is because the later digits are meaningless and are often 
285:   different on different machines; by using this routine different 
286:   machines will usually generate the same output.
287: */
288: PetscErrorCode  KSPMonitorDefaultShort(KSP ksp,PetscInt its,PetscReal fnorm,void *dummy)
289: {
290:   PetscErrorCode          ierr;
291:   PetscViewerASCIIMonitor viewer = (PetscViewerASCIIMonitor) dummy;

294:   if (!dummy) {PetscViewerASCIIMonitorCreate(ksp->comm,"stdout",0,&viewer);}

296:   if (fnorm > 1.e-9) {
297:     PetscViewerASCIIMonitorPrintf(viewer,"%3D KSP Residual norm %G \n",its,fnorm);
298:   } else if (fnorm > 1.e-11){
299:     PetscViewerASCIIMonitorPrintf(viewer,"%3D KSP Residual norm %5.3e \n",its,fnorm);
300:   } else {
301:     PetscViewerASCIIMonitorPrintf(viewer,"%3D KSP Residual norm < 1.e-11\n",its);
302:   }
303:   if (!dummy) {PetscViewerASCIIMonitorDestroy(viewer);}
304:   return(0);
305: }

309: /*@C
310:    KSPSkipConverged - Convergence test that do not return as converged
311:    until the maximum number of iterations is reached.

313:    Collective on KSP

315:    Input Parameters:
316: +  ksp   - iterative context
317: .  n     - iteration number
318: .  rnorm - 2-norm residual value (may be estimated)
319: -  dummy - unused convergence context 

321:    Returns:
322: .  reason - KSP_CONVERGED_ITERATING, KSP_CONVERGED_ITS

324:    Notes: 
325:    This should be used as the convergence test with the option
326:    KSPSetNormType(ksp,KSP_NORM_NO), since norms of the residual are
327:    not computed. Convergence is then declared after the maximum number
328:    of iterations have been reached. Useful when one is using CG or
329:    BiCGStab as a smoother.
330:                     
331:    Level: advanced

333: .keywords: KSP, default, convergence, residual

335: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPSetNormType()
336: @*/
337: PetscErrorCode  KSPSkipConverged(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *dummy)
338: {
342:   *reason = KSP_CONVERGED_ITERATING;
343:   if (n >= ksp->max_it) *reason = KSP_CONVERGED_ITS;
344:   return(0);
345: }

349: /*@C
350:    KSPDefaultConvergedSetUIRNorm - makes the default convergence test use || B*(b - A*(initial guess))||
351:       instead of || B*b ||. In the case of right preconditioner or if KSPSetNormType(ksp,KSP_NORM_UNPRECONDIITONED)
352:       is used there is no B in the above formula. UIRNorm is short for Use Initial Residual Norm.

354:    Collective on KSP

356:    Input Parameters:
357: .  ksp   - iterative context

359:    Options Database:
360: .   -ksp_converged_use_initial_residual_norm

362:    Use KSPSetTolerances() to alter the defaults for rtol, abstol, dtol.

364:    The precise values of reason are macros such as KSP_CONVERGED_RTOL, which
365:    are defined in petscksp.h.

367:    Level: intermediate

369: .keywords: KSP, default, convergence, residual

371: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPSkipConverged(), KSPConvergedReason, KSPGetConvergedReason(), KSPDefaultConvergedSetUMIRNorm()
372: @*/
373: PetscErrorCode  KSPDefaultConvergedSetUIRNorm(KSP ksp)
374: {
377:   if (ksp->defaultconvergedmininitialrtol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot use KSPDefaultConvergedSetUIRNorm() and KSPDefaultConvergedSetUMIRNorm() together");
378:   ksp->defaultconvergedinitialrtol = PETSC_TRUE;
379:   return(0);
380: }

384: /*@C
385:    KSPDefaultConvergedSetUMIRNorm - makes the default convergence test use min(|| B*(b - A*(initial guess))||,|| B*b ||)
386:       In the case of right preconditioner or if KSPSetNormType(ksp,KSP_NORM_UNPRECONDIITONED)
387:       is used there is no B in the above formula. UMIRNorm is short for Use Minimum Initial Residual Norm.

389:    Collective on KSP

391:    Input Parameters:
392: .  ksp   - iterative context

394:    Options Database:
395: .   -ksp_converged_use_min_initial_residual_norm

397:    Use KSPSetTolerances() to alter the defaults for rtol, abstol, dtol.

399:    The precise values of reason are macros such as KSP_CONVERGED_RTOL, which
400:    are defined in petscksp.h.

402:    Level: intermediate

404: .keywords: KSP, default, convergence, residual

406: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPSkipConverged(), KSPConvergedReason, KSPGetConvergedReason(), KSPDefaultConvergedSetUIRNorm()
407: @*/
408: PetscErrorCode  KSPDefaultConvergedSetUMIRNorm(KSP ksp)
409: {
412:   if (ksp->defaultconvergedinitialrtol) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"Cannot use KSPDefaultConvergedSetUIRNorm() and KSPDefaultConvergedSetUMIRNorm() together");
413:   ksp->defaultconvergedmininitialrtol = PETSC_TRUE;
414:   return(0);
415: }

419: /*@C
420:    KSPDefaultConverged - Determines convergence of
421:    the iterative solvers (default code).

423:    Collective on KSP

425:    Input Parameters:
426: +  ksp   - iterative context
427: .  n     - iteration number
428: .  rnorm - 2-norm residual value (may be estimated)
429: -  dummy - unused convergence context 

431:    Returns:
432: +   positive - if the iteration has converged;
433: .   negative - if residual norm exceeds divergence threshold;
434: -   0 - otherwise.

436:    Notes:
437:    KSPDefaultConverged() reaches convergence when
438: $      rnorm < MAX (rtol * rnorm_0, abstol);
439:    Divergence is detected if
440: $      rnorm > dtol * rnorm_0,

442:    where 
443: +     rtol = relative tolerance,
444: .     abstol = absolute tolerance.
445: .     dtol = divergence tolerance,
446: -     rnorm_0 is the two norm of the right hand side. When initial guess is non-zero you 
447:           can call KSPDefaultConvergedSetUIRNorm() to use the norm of (b - A*(initial guess))
448:           as the starting point for relative norm convergence testing.

450:    Use KSPSetTolerances() to alter the defaults for rtol, abstol, dtol.

452:    The precise values of reason are macros such as KSP_CONVERGED_RTOL, which
453:    are defined in petscksp.h.

455:    Level: intermediate

457: .keywords: KSP, default, convergence, residual

459: .seealso: KSPSetConvergenceTest(), KSPSetTolerances(), KSPSkipConverged(), KSPConvergedReason, KSPGetConvergedReason(),
460:           KSPDefaultConvergedSetUIRNorm(), KSPDefaultConvergedSetUMIRNorm().
461: @*/
462: PetscErrorCode  KSPDefaultConverged(KSP ksp,PetscInt n,PetscReal rnorm,KSPConvergedReason *reason,void *dummy)
463: {

469:   *reason = KSP_CONVERGED_ITERATING;

471:   if (!n) {
472:     /* if user gives initial guess need to compute norm of b */
473:     if (!ksp->guess_zero && !ksp->defaultconvergedinitialrtol) {
474:       PetscReal      snorm;
475:       if (ksp->normtype == KSP_NORM_UNPRECONDITIONED || ksp->pc_side == PC_RIGHT) {
476:         PetscInfo(ksp,"user has provided nonzero initial guess, computing 2-norm of RHS\n");
477:         VecNorm(ksp->vec_rhs,NORM_2,&snorm);        /*     <- b'*b */
478:       } else {
479:         Vec z;
480:         VecDuplicate(ksp->vec_rhs,&z);
481:         KSP_PCApply(ksp,ksp->vec_rhs,z);
482:         if (ksp->normtype == KSP_NORM_PRECONDITIONED) {
483:           PetscInfo(ksp,"user has provided nonzero initial guess, computing 2-norm of preconditioned RHS\n");
484:           VecNorm(z,NORM_2,&snorm);                 /*    dp <- b'*B'*B*b */
485:         } else if (ksp->normtype == KSP_NORM_NATURAL) {
486:           PetscScalar norm;
487:           PetscInfo(ksp,"user has provided nonzero initial guess, computing natural norm of RHS\n");
488:           VecDot(ksp->vec_rhs,z,&norm);
489:           snorm = sqrt(PetscAbsScalar(norm));                            /*    dp <- b'*B*b */
490:         }
491:         VecDestroy(z);
492:       }
493:       /* handle special case of zero RHS and nonzero guess */
494:       if (!snorm) {
495:         PetscInfo(ksp,"Special case, user has provided nonzero initial guess and zero RHS\n");
496:         snorm = rnorm;
497:       }
498:       if (ksp->defaultconvergedmininitialrtol) {
499:         ksp->rnorm0 = PetscMin(snorm,rnorm);
500:       } else {
501:         ksp->rnorm0 = snorm;
502:       }
503:     } else {
504:       ksp->rnorm0 = rnorm;
505:     }
506:     ksp->ttol   = PetscMax(ksp->rtol*ksp->rnorm0,ksp->abstol);
507:   }

509:   if (rnorm != rnorm) {
510:     PetscInfo(ksp,"Linear solver has created a not a number (NaN) as the residual norm, declaring divergence \n");
511:     *reason = KSP_DIVERGED_NAN;
512:   } else if (rnorm <= ksp->ttol) {
513:     if (rnorm < ksp->abstol) {
514:       PetscInfo3(ksp,"Linear solver has converged. Residual norm %G is less than absolute tolerance %G at iteration %D\n",rnorm,ksp->abstol,n);
515:       *reason = KSP_CONVERGED_ATOL;
516:     } else {
517:       if (ksp->defaultconvergedinitialrtol) {
518:         PetscInfo4(ksp,"Linear solver has converged. Residual norm %G is less than relative tolerance %G times initial residual norm %G at iteration %D\n",rnorm,ksp->rtol,ksp->rnorm0,n);
519:       } else {
520:         PetscInfo4(ksp,"Linear solver has converged. Residual norm %G is less than relative tolerance %G times initial right hand side norm %G at iteration %D\n",rnorm,ksp->rtol,ksp->rnorm0,n);
521:       }
522:       *reason = KSP_CONVERGED_RTOL;
523:     }
524:   } else if (rnorm >= ksp->divtol*ksp->rnorm0) {
525:     PetscInfo3(ksp,"Linear solver is diverging. Initial right hand size norm %G, current residual norm %G at iteration %D\n",ksp->rnorm0,rnorm,n);
526:     *reason = KSP_DIVERGED_DTOL;
527:   }
528:   return(0);
529: }

533: /*
534:    KSPDefaultBuildSolution - Default code to create/move the solution.

536:    Input Parameters:
537: +  ksp - iterative context
538: -  v   - pointer to the user's vector  

540:    Output Parameter:
541: .  V - pointer to a vector containing the solution

543:    Level: advanced

545: .keywords:  KSP, build, solution, default

547: .seealso: KSPGetSolution(), KSPDefaultBuildResidual()
548: */
549: PetscErrorCode KSPDefaultBuildSolution(KSP ksp,Vec v,Vec *V)
550: {
553:   if (ksp->pc_side == PC_RIGHT) {
554:     if (ksp->pc) {
555:       if (v) {KSP_PCApply(ksp,ksp->vec_sol,v); *V = v;}
556:       else {SETERRQ(PETSC_ERR_SUP,"Not working with right preconditioner");}
557:     } else {
558:       if (v) {VecCopy(ksp->vec_sol,v); *V = v;}
559:       else { *V = ksp->vec_sol;}
560:     }
561:   } else if (ksp->pc_side == PC_SYMMETRIC) {
562:     if (ksp->pc) {
563:       if (ksp->transpose_solve) SETERRQ(PETSC_ERR_SUP,"Not working with symmetric preconditioner and transpose solve");
564:       if (v) {PCApplySymmetricRight(ksp->pc,ksp->vec_sol,v); *V = v;}
565:       else {SETERRQ(PETSC_ERR_SUP,"Not working with symmetric preconditioner");}
566:     } else  {
567:       if (v) {VecCopy(ksp->vec_sol,v); *V = v;}
568:       else { *V = ksp->vec_sol;}
569:     }
570:   } else {
571:     if (v) {VecCopy(ksp->vec_sol,v); *V = v;}
572:     else { *V = ksp->vec_sol; }
573:   }
574:   return(0);
575: }

579: /*
580:    KSPDefaultBuildResidual - Default code to compute the residual.

582:    Input Parameters:
583: .  ksp - iterative context
584: .  t   - pointer to temporary vector
585: .  v   - pointer to user vector  

587:    Output Parameter:
588: .  V - pointer to a vector containing the residual

590:    Level: advanced

592: .keywords:  KSP, build, residual, default

594: .seealso: KSPDefaultBuildSolution()
595: */
596: PetscErrorCode KSPDefaultBuildResidual(KSP ksp,Vec t,Vec v,Vec *V)
597: {
599:   MatStructure   pflag;
600:   Mat            Amat,Pmat;

603:   PCGetOperators(ksp->pc,&Amat,&Pmat,&pflag);
604:   KSPBuildSolution(ksp,t,PETSC_NULL);
605:   KSP_MatMult(ksp,Amat,t,v);
606:   VecAYPX(v,-1.0,ksp->vec_rhs);
607:   *V = v;
608:   return(0);
609: }

613: /*@C
614:   KSPGetVecs - Gets a number of work vectors.

616:   Input Parameters:
617: + ksp  - iterative context
618: . rightn  - number of right work vectors
619: - leftn   - number of left work vectors to allocate

621:   Output Parameter:
622: +  right - the array of vectors created
623: -  left - the array of left vectors

625:    Note: The right vector has as many elements as the matrix has columns. The left
626:      vector has as many elements as the matrix has rows.

628:    Level: advanced

630: .seealso:   MatGetVecs()

632: @*/
633: PetscErrorCode KSPGetVecs(KSP ksp,PetscInt rightn, Vec **right,PetscInt leftn,Vec **left)
634: {
636:   Vec            vecr,vecl;

639:   if (rightn) {
640:     if (!right) SETERRQ(PETSC_ERR_ARG_INCOMP,"You asked for right vectors but did not pass a pointer to hold them");
641:     if (ksp->vec_sol) vecr = ksp->vec_sol;
642:     else {
643:       Mat pmat;
644:       PCGetOperators(ksp->pc,PETSC_NULL,&pmat,PETSC_NULL);
645:       MatGetVecs(pmat,&vecr,PETSC_NULL);
646:     }
647:     VecDuplicateVecs(vecr,rightn,right);
648:     if (!ksp->vec_sol) {
649:       VecDestroy(vecr);
650:     }
651:   }
652:   if (leftn) {
653:     if (!left) SETERRQ(PETSC_ERR_ARG_INCOMP,"You asked for left vectors but did not pass a pointer to hold them");
654:     if (ksp->vec_rhs) vecl = ksp->vec_rhs;
655:     else {
656:       Mat pmat;
657:       PCGetOperators(ksp->pc,PETSC_NULL,&pmat,PETSC_NULL);
658:       MatGetVecs(pmat,PETSC_NULL,&vecl);
659:     }
660:     VecDuplicateVecs(vecl,leftn,left);
661:     if (!ksp->vec_rhs) {
662:       VecDestroy(vecl);
663:     }
664:   }
665:   return(0);
666: }

670: /*
671:   KSPDefaultGetWork - Gets a number of work vectors.

673:   Input Parameters:
674: . ksp  - iterative context
675: . nw   - number of work vectors to allocate

677:   Notes:
678:   Call this only if no work vectors have been allocated 
679:  */
680: PetscErrorCode KSPDefaultGetWork(KSP ksp,PetscInt nw)
681: {

685:   if (ksp->work) {KSPDefaultFreeWork(ksp);}
686:   ksp->nwork = nw;
687:   KSPGetVecs(ksp,nw,&ksp->work,0,PETSC_NULL);
688:   PetscLogObjectParents(ksp,nw,ksp->work);
689:   return(0);
690: }

694: /*
695:   KSPDefaultDestroy - Destroys a iterative context variable for methods with
696:   no separate context.  Preferred calling sequence KSPDestroy().

698:   Input Parameter: 
699: . ksp - the iterative context
700: */
701: PetscErrorCode KSPDefaultDestroy(KSP ksp)
702: {

707:   /* free work vectors */
708:   KSPDefaultFreeWork(ksp);
709:   /* free private data structure */
710:   PetscFree(ksp->data);
711:   return(0);
712: }

716: /*@
717:    KSPGetConvergedReason - Gets the reason the KSP iteration was stopped.

719:    Not Collective

721:    Input Parameter:
722: .  ksp - the KSP context

724:    Output Parameter:
725: .  reason - negative value indicates diverged, positive value converged, see KSPConvergedReason

727:    Possible values for reason:
728: +  KSP_CONVERGED_RTOL (residual 2-norm decreased by a factor of rtol, from 2-norm of right hand side)
729: .  KSP_CONVERGED_ATOL (residual 2-norm less than abstol)
730: .  KSP_CONVERGED_ITS (used by the preonly preconditioner that always uses ONE iteration) 
731: .  KSP_CONVERGED_CG_NEG_CURVE
732: .  KSP_CONVERGED_CG_CONSTRAINED
733: .  KSP_CONVERGED_STEP_LENGTH
734: .  KSP_DIVERGED_ITS  (required more than its to reach convergence)
735: .  KSP_DIVERGED_DTOL (residual norm increased by a factor of divtol)
736: .  KSP_DIVERGED_NAN (residual norm became Not-a-number likely do to 0/0)
737: .  KSP_DIVERGED_BREAKDOWN (generic breakdown in method)
738: -  KSP_DIVERGED_BREAKDOWN_BICG (Initial residual is orthogonal to preconditioned initial
739:                                 residual. Try a different preconditioner, or a different initial guess.)
740:  

742:    Level: beginner

744:    Notes: Can only be called after the call the KSPSolve() is complete.

746: .keywords: KSP, nonlinear, set, convergence, test

748: .seealso: KSPSetConvergenceTest(), KSPDefaultConverged(), KSPSetTolerances(), KSPConvergedReason
749: @*/
750: PetscErrorCode  KSPGetConvergedReason(KSP ksp,KSPConvergedReason *reason)
751: {
755:   *reason = ksp->reason;
756:   return(0);
757: }