Line data Source code
1 : /*
2 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3 : SLEPc - Scalable Library for Eigenvalue Problem Computations
4 : Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5 :
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 : */
13 :
14 : #include <slepc/private/svdimpl.h> /*I "slepcsvd.h" I*/
15 :
16 : /*
17 : SVDConvergedAbsolute - Checks convergence absolutely.
18 : */
19 1042 : PetscErrorCode SVDConvergedAbsolute(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
20 : {
21 1042 : PetscFunctionBegin;
22 1042 : *errest = res;
23 1042 : PetscFunctionReturn(PETSC_SUCCESS);
24 : }
25 :
26 : /*
27 : SVDConvergedRelative - Checks convergence relative to the singular value.
28 : */
29 2751 : PetscErrorCode SVDConvergedRelative(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
30 : {
31 2751 : PetscFunctionBegin;
32 2751 : *errest = (sigma!=0.0)? res/sigma: PETSC_MAX_REAL;
33 2751 : PetscFunctionReturn(PETSC_SUCCESS);
34 : }
35 :
36 : /*
37 : SVDConvergedNorm - Checks convergence relative to the matrix norms.
38 : */
39 20 : PetscErrorCode SVDConvergedNorm(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
40 : {
41 20 : PetscFunctionBegin;
42 20 : *errest = res/PetscMax(svd->nrma,svd->nrmb);
43 20 : PetscFunctionReturn(PETSC_SUCCESS);
44 : }
45 :
46 : /*
47 : SVDConvergedMaxIt - Always returns Inf to force reaching the maximum number of iterations.
48 : */
49 2 : PetscErrorCode SVDConvergedMaxIt(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)
50 : {
51 2 : PetscFunctionBegin;
52 2 : *errest = PETSC_MAX_REAL;
53 2 : PetscFunctionReturn(PETSC_SUCCESS);
54 : }
55 :
56 : /*@C
57 : SVDStoppingBasic - Default routine to determine whether the outer singular value
58 : solver iteration must be stopped.
59 :
60 : Collective
61 :
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)
69 :
70 : Output Parameter:
71 : . reason - result of the stopping test
72 :
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).
78 :
79 : SVDStoppingBasic() will stop if all requested singular values are converged, or if
80 : the maximum number of iterations has been reached.
81 :
82 : Use SVDSetStoppingTest() to provide your own test instead of using this one.
83 :
84 : Level: advanced
85 :
86 : .seealso: SVDSetStoppingTest(), SVDStoppingThreshold(), SVDConvergedReason, SVDGetConvergedReason()
87 : @*/
88 2080 : PetscErrorCode SVDStoppingBasic(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)
89 : {
90 2080 : PetscFunctionBegin;
91 2080 : *reason = SVD_CONVERGED_ITERATING;
92 2080 : if (nconv >= nsv) {
93 123 : PetscCall(PetscInfo(svd,"Singular value solver finished successfully: %" PetscInt_FMT " singular triplets converged at iteration %" PetscInt_FMT "\n",nconv,its));
94 123 : *reason = SVD_CONVERGED_TOL;
95 1957 : } else if (its >= max_it) {
96 0 : if (svd->conv == SVD_CONV_MAXIT) *reason = SVD_CONVERGED_MAXIT;
97 : else {
98 0 : *reason = SVD_DIVERGED_ITS;
99 0 : PetscCall(PetscInfo(svd,"Singular value solver iteration reached maximum number of iterations (%" PetscInt_FMT ")\n",its));
100 : }
101 : }
102 2080 : PetscFunctionReturn(PETSC_SUCCESS);
103 : }
104 :
105 : /*@C
106 : SVDStoppingThreshold - Routine to determine whether the outer singular value
107 : solver iteration must be stopped, according to some threshold for the computed values.
108 :
109 : Collective
110 :
111 : Input Parameters:
112 : + svd - singular value solver context obtained from SVDCreate()
113 : . its - current number of iterations
114 : . max_it - maximum number of iterations
115 : . nconv - number of currently converged singular triplets (ignored here)
116 : . nsv - number of requested singular triplets (ignored here)
117 : - ctx - context containing additional data (SVDStoppingCtx)
118 :
119 : Output Parameter:
120 : . reason - result of the stopping test
121 :
122 : Notes:
123 : A positive value of reason indicates that the iteration has finished successfully
124 : (converged), and a negative value indicates an error condition (diverged). If
125 : the iteration needs to be continued, reason must be set to SVD_CONVERGED_ITERATING
126 : (zero).
127 :
128 : SVDStoppingThreshold() will stop when one of the computed singular values is not
129 : above/below the threshold given at SVDSetThreshold(). If a number of wanted singular
130 : values has been specified via SVDSetDimensions() then it is also taken into account,
131 : and the solver will stop when one of the two conditions (threshold or number of
132 : converged values) is met.
133 :
134 : Use SVDSetStoppingTest() to provide your own test instead of using this one.
135 :
136 : Level: advanced
137 :
138 : .seealso: SVDSetStoppingTest(), SVDStoppingBasic(), SVDSetThreshold(), SVDSetDimensions(), SVDConvergedReason, SVDGetConvergedReason()
139 : @*/
140 54 : PetscErrorCode SVDStoppingThreshold(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)
141 : {
142 54 : PetscReal thres,firstsv,lastsv;
143 54 : PetscBool rel;
144 54 : SVDWhich which;
145 :
146 54 : PetscFunctionBegin;
147 54 : *reason = SVD_CONVERGED_ITERATING;
148 54 : firstsv = ((SVDStoppingCtx)ctx)->firstsv;
149 54 : lastsv = ((SVDStoppingCtx)ctx)->lastsv;
150 54 : thres = ((SVDStoppingCtx)ctx)->thres;
151 54 : rel = ((SVDStoppingCtx)ctx)->threlative;
152 54 : which = ((SVDStoppingCtx)ctx)->which;
153 54 : if (nconv && ((which==SVD_LARGEST && ((rel && lastsv<thres*firstsv) || (!rel && lastsv<thres))) || (which==SVD_SMALLEST && lastsv>thres))) {
154 8 : if (which==SVD_SMALLEST) PetscCall(PetscInfo(svd,"Singular value solver finished successfully: singular value %g is above the threshold %g\n",(double)lastsv,(double)thres));
155 8 : else if (!rel) PetscCall(PetscInfo(svd,"Singular value solver finished successfully: singular value %g is below the threshold %g\n",(double)lastsv,(double)thres));
156 8 : else PetscCall(PetscInfo(svd,"Singular value solver finished successfully: the ratio %g/%g is below the threshold %g\n",(double)lastsv,(double)firstsv,(double)thres));
157 8 : *reason = SVD_CONVERGED_TOL;
158 46 : } else if (nsv && nconv >= nsv) {
159 0 : PetscCall(PetscInfo(svd,"Singular value solver finished successfully: %" PetscInt_FMT " singular triplets converged at iteration %" PetscInt_FMT "\n",nconv,its));
160 0 : *reason = SVD_CONVERGED_TOL;
161 46 : } else if (its >= max_it) {
162 0 : if (svd->conv == SVD_CONV_MAXIT) *reason = SVD_CONVERGED_MAXIT;
163 : else {
164 0 : *reason = SVD_DIVERGED_ITS;
165 0 : PetscCall(PetscInfo(svd,"Singular value solver iteration reached maximum number of iterations (%" PetscInt_FMT ")\n",its));
166 : }
167 : }
168 54 : PetscFunctionReturn(PETSC_SUCCESS);
169 : }
170 :
171 : /*@
172 : SVDSetWorkVecs - Sets a number of work vectors into an SVD object.
173 :
174 : Collective
175 :
176 : Input Parameters:
177 : + svd - singular value solver context
178 : . nleft - number of work vectors of dimension equal to left singular vector
179 : - nright - number of work vectors of dimension equal to right singular vector
180 :
181 : Developer Notes:
182 : This is SLEPC_EXTERN because it may be required by user plugin SVD
183 : implementations.
184 :
185 : Level: developer
186 :
187 : .seealso: SVDSetUp()
188 : @*/
189 1437 : PetscErrorCode SVDSetWorkVecs(SVD svd,PetscInt nleft,PetscInt nright)
190 : {
191 1437 : Vec t;
192 :
193 1437 : PetscFunctionBegin;
194 1437 : PetscValidHeaderSpecific(svd,SVD_CLASSID,1);
195 4311 : PetscValidLogicalCollectiveInt(svd,nleft,2);
196 4311 : PetscValidLogicalCollectiveInt(svd,nright,3);
197 1437 : PetscCheck(nleft>=0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nleft must be >= 0: nleft = %" PetscInt_FMT,nleft);
198 1437 : PetscCheck(nright>=0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nright must be >= 0: nright = %" PetscInt_FMT,nright);
199 1437 : PetscCheck(nleft>0 || nright>0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"nleft and nright cannot be both zero");
200 1437 : if (svd->nworkl < nleft) {
201 226 : PetscCall(VecDestroyVecs(svd->nworkl,&svd->workl));
202 226 : svd->nworkl = nleft;
203 226 : if (svd->isgeneralized) PetscCall(SVDCreateLeftTemplate(svd,&t));
204 148 : else PetscCall(MatCreateVecsEmpty(svd->OP,NULL,&t));
205 226 : PetscCall(VecDuplicateVecs(t,nleft,&svd->workl));
206 226 : PetscCall(VecDestroy(&t));
207 : }
208 1437 : if (svd->nworkr < nright) {
209 267 : PetscCall(VecDestroyVecs(svd->nworkr,&svd->workr));
210 267 : svd->nworkr = nright;
211 267 : PetscCall(MatCreateVecsEmpty(svd->OP,&t,NULL));
212 267 : PetscCall(VecDuplicateVecs(t,nright,&svd->workr));
213 267 : PetscCall(VecDestroy(&t));
214 : }
215 1437 : PetscFunctionReturn(PETSC_SUCCESS);
216 : }
|