1: /*
2: This file contains some simple default routines for common operations.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include <slepc/private/svdimpl.h> /*I "slepcsvd.h" I*/
28: /*
29: SVDConvergedRelative - Checks convergence relative to the eigenvalue.
30: */
31: PetscErrorCode SVDConvergedRelative(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx) 32: {
34: *errest = res/sigma;
35: return(0);
36: }
40: /*
41: SVDConvergedAbsolute - Checks convergence absolutely.
42: */
43: PetscErrorCode SVDConvergedAbsolute(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx) 44: {
46: *errest = res;
47: return(0);
48: }
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: {
89: *reason = SVD_CONVERGED_ITERATING;
90: if (nconv >= nsv) {
91: PetscInfo2(svd,"Singular value solver finished successfully: %D singular triplets converged at iteration %D\n",nconv,its);
92: *reason = SVD_CONVERGED_TOL;
93: } else if (its >= max_it) {
94: *reason = SVD_DIVERGED_ITS;
95: PetscInfo1(svd,"Singular value solver iteration reached maximum number of iterations (%D)\n",its);
96: }
97: return(0);
98: }