Actual source code: svddefault.c
slepc-3.18.1 2022-11-02
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: *errest = res;
22: return 0;
23: }
25: /*
26: SVDConvergedRelative - Checks convergence relative to the singular value.
27: */
28: PetscErrorCode SVDConvergedRelative(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
29: {
30: *errest = res/sigma;
31: return 0;
32: }
34: /*
35: SVDConvergedNorm - Checks convergence relative to the matrix norms.
36: */
37: PetscErrorCode SVDConvergedNorm(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
38: {
39: *errest = res/PetscMax(svd->nrma,svd->nrmb);
40: return 0;
41: }
43: /*
44: SVDConvergedMaxIt - Always returns Inf to force reaching the maximum number of iterations.
45: */
46: PetscErrorCode SVDConvergedMaxIt(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
47: {
48: *errest = PETSC_MAX_REAL;
49: return 0;
50: }
52: /*@C
53: SVDStoppingBasic - Default routine to determine whether the outer singular value
54: solver iteration must be stopped.
56: Collective on svd
58: Input Parameters:
59: + svd - singular value solver context obtained from SVDCreate()
60: . its - current number of iterations
61: . max_it - maximum number of iterations
62: . nconv - number of currently converged singular triplets
63: . nsv - number of requested singular triplets
64: - ctx - context (not used here)
66: Output Parameter:
67: . reason - result of the stopping test
69: Notes:
70: A positive value of reason indicates that the iteration has finished successfully
71: (converged), and a negative value indicates an error condition (diverged). If
72: the iteration needs to be continued, reason must be set to SVD_CONVERGED_ITERATING
73: (zero).
75: SVDStoppingBasic() will stop if all requested singular values are converged, or if
76: the maximum number of iterations has been reached.
78: Use SVDSetStoppingTest() to provide your own test instead of using this one.
80: Level: advanced
82: .seealso: SVDSetStoppingTest(), SVDConvergedReason, SVDGetConvergedReason()
83: @*/
84: PetscErrorCode SVDStoppingBasic(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)
85: {
86: *reason = SVD_CONVERGED_ITERATING;
87: if (nconv >= nsv) {
88: PetscInfo(svd,"Singular value solver finished successfully: %" PetscInt_FMT " singular triplets converged at iteration %" PetscInt_FMT "\n",nconv,its);
89: *reason = SVD_CONVERGED_TOL;
90: } else if (its >= max_it) {
91: if (svd->conv == SVD_CONV_MAXIT) *reason = SVD_CONVERGED_MAXIT;
92: else {
93: *reason = SVD_DIVERGED_ITS;
94: PetscInfo(svd,"Singular value solver iteration reached maximum number of iterations (%" PetscInt_FMT ")\n",its);
95: }
96: }
97: return 0;
98: }
100: /*@
101: SVDSetWorkVecs - Sets a number of work vectors into an SVD object.
103: Collective on svd
105: Input Parameters:
106: + svd - singular value solver context
107: . nleft - number of work vectors of dimension equal to left singular vector
108: - nright - number of work vectors of dimension equal to right singular vector
110: Developer Notes:
111: This is SLEPC_EXTERN because it may be required by user plugin SVD
112: implementations.
114: Level: developer
116: .seealso: SVDSetUp()
117: @*/
118: PetscErrorCode SVDSetWorkVecs(SVD svd,PetscInt nleft,PetscInt nright)
119: {
120: Vec t;
128: if (svd->nworkl < nleft) {
129: VecDestroyVecs(svd->nworkl,&svd->workl);
130: svd->nworkl = nleft;
131: if (svd->isgeneralized) SVDCreateLeftTemplate(svd,&t);
132: else MatCreateVecsEmpty(svd->OP,NULL,&t);
133: VecDuplicateVecs(t,nleft,&svd->workl);
134: VecDestroy(&t);
135: }
136: if (svd->nworkr < nright) {
137: VecDestroyVecs(svd->nworkr,&svd->workr);
138: svd->nworkr = nright;
139: MatCreateVecsEmpty(svd->OP,&t,NULL);
140: VecDuplicateVecs(t,nright,&svd->workr);
141: VecDestroy(&t);
142: }
143: return 0;
144: }