Actual source code: svddefault.c
slepc-3.15.0 2021-03-31
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: SVDConvergedMaxIt - Always returns Inf to force reaching the maximum number of iterations.
38: */
39: PetscErrorCode SVDConvergedMaxIt(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
40: {
42: *errest = PETSC_MAX_REAL;
43: return(0);
44: }
46: /*@C
47: SVDStoppingBasic - Default routine to determine whether the outer singular value
48: solver iteration must be stopped.
50: Collective on svd
52: Input Parameters:
53: + svd - singular value solver context obtained from SVDCreate()
54: . its - current number of iterations
55: . max_it - maximum number of iterations
56: . nconv - number of currently converged singular triplets
57: . nsv - number of requested singular triplets
58: - ctx - context (not used here)
60: Output Parameter:
61: . reason - result of the stopping test
63: Notes:
64: A positive value of reason indicates that the iteration has finished successfully
65: (converged), and a negative value indicates an error condition (diverged). If
66: the iteration needs to be continued, reason must be set to SVD_CONVERGED_ITERATING
67: (zero).
69: SVDStoppingBasic() will stop if all requested singular values are converged, or if
70: the maximum number of iterations has been reached.
72: Use SVDSetStoppingTest() to provide your own test instead of using this one.
74: Level: advanced
76: .seealso: SVDSetStoppingTest(), SVDConvergedReason, SVDGetConvergedReason()
77: @*/
78: PetscErrorCode SVDStoppingBasic(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)
79: {
83: *reason = SVD_CONVERGED_ITERATING;
84: if (nconv >= nsv) {
85: PetscInfo2(svd,"Singular value solver finished successfully: %D singular triplets converged at iteration %D\n",nconv,its);
86: *reason = SVD_CONVERGED_TOL;
87: } else if (its >= max_it) {
88: if (svd->conv == SVD_CONV_MAXIT) *reason = SVD_CONVERGED_MAXIT;
89: else {
90: *reason = SVD_DIVERGED_ITS;
91: PetscInfo1(svd,"Singular value solver iteration reached maximum number of iterations (%D)\n",its);
92: }
93: }
94: return(0);
95: }
97: /*@
98: SVDSetWorkVecs - Sets a number of work vectors into an SVD object.
100: Collective on svd
102: Input Parameters:
103: + svd - singular value solver context
104: . nleft - number of work vectors of dimension equal to left singular vector
105: - nright - number of work vectors of dimension equal to right singular vector
107: Developers Note:
108: This is SLEPC_EXTERN because it may be required by user plugin SVD
109: implementations.
111: Level: developer
112: @*/
113: PetscErrorCode SVDSetWorkVecs(SVD svd,PetscInt nleft,PetscInt nright)
114: {
116: Vec t;
122: if (nleft <= 0) SETERRQ1(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nleft must be > 0: nleft = %D",nleft);
123: if (nright <= 0) SETERRQ1(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nright must be > 0: nright = %D",nright);
124: if (svd->nworkl < nleft) {
125: VecDestroyVecs(svd->nworkl,&svd->workl);
126: svd->nworkl = nleft;
127: if (svd->isgeneralized) { SVDCreateLeftTemplate(svd,&t); }
128: else { MatCreateVecsEmpty(svd->OP,NULL,&t); }
129: VecDuplicateVecs(t,nleft,&svd->workl);
130: VecDestroy(&t);
131: PetscLogObjectParents(svd,nleft,svd->workl);
132: }
133: if (svd->nworkr < nright) {
134: VecDestroyVecs(svd->nworkr,&svd->workr);
135: svd->nworkr = nright;
136: MatCreateVecsEmpty(svd->OP,&t,NULL);
137: VecDuplicateVecs(t,nright,&svd->workr);
138: VecDestroy(&t);
139: PetscLogObjectParents(svd,nright,svd->workr);
140: }
141: return(0);
142: }