Actual source code: svddefault.c

slepc-3.16.2 2022-02-01
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2021, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    Simple default routines for common SVD operations
 12: */

 14: #include <slepc/private/svdimpl.h>

 16: /*
 17:   SVDConvergedAbsolute - Checks convergence absolutely.
 18: */
 19: PetscErrorCode SVDConvergedAbsolute(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
 20: {
 22:   *errest = res;
 23:   return(0);
 24: }

 26: /*
 27:   SVDConvergedRelative - Checks convergence relative to the singular value.
 28: */
 29: PetscErrorCode SVDConvergedRelative(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
 30: {
 32:   *errest = res/sigma;
 33:   return(0);
 34: }

 36: /*
 37:   SVDConvergedNorm - Checks convergence relative to the matrix norms.
 38: */
 39: PetscErrorCode SVDConvergedNorm(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
 40: {
 42:   *errest = res/SlepcAbs(svd->nrma,svd->nrmb);
 43:   return(0);
 44: }

 46: /*
 47:   SVDConvergedMaxIt - Always returns Inf to force reaching the maximum number of iterations.
 48: */
 49: PetscErrorCode SVDConvergedMaxIt(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
 50: {
 52:   *errest = PETSC_MAX_REAL;
 53:   return(0);
 54: }

 56: /*@C
 57:    SVDStoppingBasic - Default routine to determine whether the outer singular value
 58:    solver iteration must be stopped.

 60:    Collective on svd

 62:    Input Parameters:
 63: +  svd    - singular value solver context obtained from SVDCreate()
 64: .  its    - current number of iterations
 65: .  max_it - maximum number of iterations
 66: .  nconv  - number of currently converged singular triplets
 67: .  nsv    - number of requested singular triplets
 68: -  ctx    - context (not used here)

 70:    Output Parameter:
 71: .  reason - result of the stopping test

 73:    Notes:
 74:    A positive value of reason indicates that the iteration has finished successfully
 75:    (converged), and a negative value indicates an error condition (diverged). If
 76:    the iteration needs to be continued, reason must be set to SVD_CONVERGED_ITERATING
 77:    (zero).

 79:    SVDStoppingBasic() will stop if all requested singular values are converged, or if
 80:    the maximum number of iterations has been reached.

 82:    Use SVDSetStoppingTest() to provide your own test instead of using this one.

 84:    Level: advanced

 86: .seealso: SVDSetStoppingTest(), SVDConvergedReason, SVDGetConvergedReason()
 87: @*/
 88: PetscErrorCode SVDStoppingBasic(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)
 89: {

 93:   *reason = SVD_CONVERGED_ITERATING;
 94:   if (nconv >= nsv) {
 95:     PetscInfo2(svd,"Singular value solver finished successfully: %D singular triplets converged at iteration %D\n",nconv,its);
 96:     *reason = SVD_CONVERGED_TOL;
 97:   } else if (its >= max_it) {
 98:     if (svd->conv == SVD_CONV_MAXIT) *reason = SVD_CONVERGED_MAXIT;
 99:     else {
100:       *reason = SVD_DIVERGED_ITS;
101:       PetscInfo1(svd,"Singular value solver iteration reached maximum number of iterations (%D)\n",its);
102:     }
103:   }
104:   return(0);
105: }

107: /*@
108:    SVDSetWorkVecs - Sets a number of work vectors into an SVD object.

110:    Collective on svd

112:    Input Parameters:
113: +  svd    - singular value solver context
114: .  nleft  - number of work vectors of dimension equal to left singular vector
115: -  nright - number of work vectors of dimension equal to right singular vector

117:    Developers Note:
118:    This is SLEPC_EXTERN because it may be required by user plugin SVD
119:    implementations.

121:    Level: developer
122: @*/
123: PetscErrorCode SVDSetWorkVecs(SVD svd,PetscInt nleft,PetscInt nright)
124: {
126:   Vec            t;

132:   if (nleft <= 0) SETERRQ1(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nleft must be > 0: nleft = %D",nleft);
133:   if (nright <= 0) SETERRQ1(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nright must be > 0: nright = %D",nright);
134:   if (svd->nworkl < nleft) {
135:     VecDestroyVecs(svd->nworkl,&svd->workl);
136:     svd->nworkl = nleft;
137:     if (svd->isgeneralized) { SVDCreateLeftTemplate(svd,&t); }
138:     else { MatCreateVecsEmpty(svd->OP,NULL,&t); }
139:     VecDuplicateVecs(t,nleft,&svd->workl);
140:     VecDestroy(&t);
141:     PetscLogObjectParents(svd,nleft,svd->workl);
142:   }
143:   if (svd->nworkr < nright) {
144:     VecDestroyVecs(svd->nworkr,&svd->workr);
145:     svd->nworkr = nright;
146:     MatCreateVecsEmpty(svd->OP,&t,NULL);
147:     VecDuplicateVecs(t,nright,&svd->workr);
148:     VecDestroy(&t);
149:     PetscLogObjectParents(svd,nright,svd->workr);
150:   }
151:   return(0);
152: }