Actual source code: svddefault.c
slepc-3.16.1 2021-11-17
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: }