Actual source code: svddefault.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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> /*I "slepcsvd.h" I*/
16: /*
17: SVDConvergedRelative - Checks convergence relative to the eigenvalue.
18: */
19: PetscErrorCode SVDConvergedRelative(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
20: {
22: *errest = res/sigma;
23: return(0);
24: }
26: /*
27: SVDConvergedAbsolute - Checks convergence absolutely.
28: */
29: PetscErrorCode SVDConvergedAbsolute(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
30: {
32: *errest = res;
33: return(0);
34: }
36: /*@C
37: SVDStoppingBasic - Default routine to determine whether the outer singular value
38: solver iteration must be stopped.
40: Collective on SVD
42: Input Parameters:
43: + svd - singular value solver context obtained from SVDCreate()
44: . its - current number of iterations
45: . max_it - maximum number of iterations
46: . nconv - number of currently converged singular triplets
47: . nsv - number of requested singular triplets
48: - ctx - context (not used here)
50: Output Parameter:
51: . reason - result of the stopping test
53: Notes:
54: A positive value of reason indicates that the iteration has finished successfully
55: (converged), and a negative value indicates an error condition (diverged). If
56: the iteration needs to be continued, reason must be set to SVD_CONVERGED_ITERATING
57: (zero).
59: SVDStoppingBasic() will stop if all requested singular values are converged, or if
60: the maximum number of iterations has been reached.
62: Use SVDSetStoppingTest() to provide your own test instead of using this one.
64: Level: advanced
66: .seealso: SVDSetStoppingTest(), SVDConvergedReason, SVDGetConvergedReason()
67: @*/
68: PetscErrorCode SVDStoppingBasic(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)
69: {
73: *reason = SVD_CONVERGED_ITERATING;
74: if (nconv >= nsv) {
75: PetscInfo2(svd,"Singular value solver finished successfully: %D singular triplets converged at iteration %D\n",nconv,its);
76: *reason = SVD_CONVERGED_TOL;
77: } else if (its >= max_it) {
78: *reason = SVD_DIVERGED_ITS;
79: PetscInfo1(svd,"Singular value solver iteration reached maximum number of iterations (%D)\n",its);
80: }
81: return(0);
82: }