1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, 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: {
 21:   PetscFunctionBegin;
 22:   *errest = res;
 23:   PetscFunctionReturn(PETSC_SUCCESS);
 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: {
 31:   PetscFunctionBegin;
 32:   *errest = (sigma!=0.0)? res/sigma: PETSC_MAX_REAL;
 33:   PetscFunctionReturn(PETSC_SUCCESS);
 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: {
 41:   PetscFunctionBegin;
 42:   *errest = res/PetscMax(svd->nrma,svd->nrmb);
 43:   PetscFunctionReturn(PETSC_SUCCESS);
 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: {
 51:   PetscFunctionBegin;
 52:   *errest = PETSC_MAX_REAL;
 53:   PetscFunctionReturn(PETSC_SUCCESS);
 54: }

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

 60:    Collective

 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: {
 90:   PetscFunctionBegin;
 91:   *reason = SVD_CONVERGED_ITERATING;
 92:   if (nconv >= nsv) {
 93:     PetscCall(PetscInfo(svd,"Singular value solver finished successfully: %" PetscInt_FMT " singular triplets converged at iteration %" PetscInt_FMT "\n",nconv,its));
 94:     *reason = SVD_CONVERGED_TOL;
 95:   } else if (its >= max_it) {
 96:     if (svd->conv == SVD_CONV_MAXIT) *reason = SVD_CONVERGED_MAXIT;
 97:     else {
 98:       *reason = SVD_DIVERGED_ITS;
 99:       PetscCall(PetscInfo(svd,"Singular value solver iteration reached maximum number of iterations (%" PetscInt_FMT ")\n",its));
100:     }
101:   }
102:   PetscFunctionReturn(PETSC_SUCCESS);
103: }

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

108:    Collective

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

115:    Developer Notes:
116:    This is SLEPC_EXTERN because it may be required by user plugin SVD
117:    implementations.

119:    Level: developer

121: .seealso: SVDSetUp()
122: @*/
123: PetscErrorCode SVDSetWorkVecs(SVD svd,PetscInt nleft,PetscInt nright)
124: {
125:   Vec            t;

127:   PetscFunctionBegin;
131:   PetscCheck(nleft>=0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nleft must be >= 0: nleft = %" PetscInt_FMT,nleft);
132:   PetscCheck(nright>=0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nright must be >= 0: nright = %" PetscInt_FMT,nright);
133:   PetscCheck(nleft>0 || nright>0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nleft and nright cannot be both zero");
134:   if (svd->nworkl < nleft) {
135:     PetscCall(VecDestroyVecs(svd->nworkl,&svd->workl));
136:     svd->nworkl = nleft;
137:     if (svd->isgeneralized) PetscCall(SVDCreateLeftTemplate(svd,&t));
138:     else PetscCall(MatCreateVecsEmpty(svd->OP,NULL,&t));
139:     PetscCall(VecDuplicateVecs(t,nleft,&svd->workl));
140:     PetscCall(VecDestroy(&t));
141:   }
142:   if (svd->nworkr < nright) {
143:     PetscCall(VecDestroyVecs(svd->nworkr,&svd->workr));
144:     svd->nworkr = nright;
145:     PetscCall(MatCreateVecsEmpty(svd->OP,&t,NULL));
146:     PetscCall(VecDuplicateVecs(t,nright,&svd->workr));
147:     PetscCall(VecDestroy(&t));
148:   }
149:   PetscFunctionReturn(PETSC_SUCCESS);
150: }