Actual source code: svdopts.c
1: /*
2: SVD routines for setting solver options.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2013, 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: SVDSetTransposeMode - Sets how to handle the transpose of the matrix
30: associated with the singular value problem.
32: Logically Collective on SVD
34: Input Parameters:
35: + svd - the singular value solver context
36: - mode - how to compute the transpose, one of SVD_TRANSPOSE_EXPLICIT
37: or SVD_TRANSPOSE_IMPLICIT (see notes below)
39: Options Database Key:
40: . -svd_transpose_mode <mode> - Indicates the mode flag, where <mode>
41: is one of 'explicit' or 'implicit'.
43: Notes:
44: In the SVD_TRANSPOSE_EXPLICIT mode, the transpose of the matrix is
45: explicitly built.
47: The option SVD_TRANSPOSE_IMPLICIT does not build the transpose, but
48: handles it implicitly via MatMultTranspose() (or MatMultHermitianTranspose()
49: in the complex case) operations. This is
50: likely to be more inefficient than SVD_TRANSPOSE_EXPLICIT, both in
51: sequential and in parallel, but requires less storage.
53: The default is SVD_TRANSPOSE_EXPLICIT if the matrix has defined the
54: MatTranspose operation, and SVD_TRANSPOSE_IMPLICIT otherwise.
56: Level: advanced
58: .seealso: SVDGetTransposeMode(), SVDSolve(), SVDSetOperator(),
59: SVDGetOperator(), SVDTransposeMode
60: @*/
61: PetscErrorCode SVDSetTransposeMode(SVD svd,SVDTransposeMode mode)
62: {
66: if (mode == PETSC_DEFAULT || mode == PETSC_DECIDE) mode = (SVDTransposeMode)PETSC_DECIDE;
67: else switch (mode) {
68: case SVD_TRANSPOSE_EXPLICIT:
69: case SVD_TRANSPOSE_IMPLICIT:
70: if (svd->transmode!=mode) {
71: svd->transmode = mode;
72: svd->setupcalled = 0;
73: }
74: break;
75: default:
76: SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid transpose mode");
77: }
78: return(0);
79: }
83: /*@C
84: SVDGetTransposeMode - Gets the mode used to compute the transpose
85: of the matrix associated with the singular value problem.
87: Not Collective
89: Input Parameter:
90: . svd - the singular value solver context
92: Output paramter:
93: . mode - how to compute the transpose, one of SVD_TRANSPOSE_EXPLICIT
94: or SVD_TRANSPOSE_IMPLICIT
96: Level: advanced
98: .seealso: SVDSetTransposeMode(), SVDSolve(), SVDSetOperator(),
99: SVDGetOperator(), SVDTransposeMode
100: @*/
101: PetscErrorCode SVDGetTransposeMode(SVD svd,SVDTransposeMode *mode)
102: {
106: *mode = svd->transmode;
107: return(0);
108: }
112: /*@
113: SVDSetTolerances - Sets the tolerance and maximum
114: iteration count used by the default SVD convergence testers.
116: Logically Collective on SVD
118: Input Parameters:
119: + svd - the singular value solver context
120: . tol - the convergence tolerance
121: - maxits - maximum number of iterations to use
123: Options Database Keys:
124: + -svd_tol <tol> - Sets the convergence tolerance
125: - -svd_max_it <maxits> - Sets the maximum number of iterations allowed
126: (use PETSC_DECIDE to compute an educated guess based on basis and matrix sizes)
128: Notes:
129: Pass 0 to retain the previous value of any parameter.
131: Level: intermediate
133: .seealso: SVDGetTolerances()
134: @*/
135: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,PetscInt maxits)
136: {
141: if (tol) {
142: if (tol == PETSC_DEFAULT) {
143: tol = PETSC_DEFAULT;
144: } else {
145: if (tol < 0.0) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
146: svd->tol = tol;
147: }
148: }
149: if (maxits) {
150: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
151: svd->max_it = 0;
152: svd->setupcalled = 0;
153: } else {
154: if (maxits < 0) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
155: svd->max_it = maxits;
156: }
157: }
158: return(0);
159: }
163: /*@
164: SVDGetTolerances - Gets the tolerance and maximum
165: iteration count used by the default SVD convergence tests.
167: Not Collective
169: Input Parameter:
170: . svd - the singular value solver context
172: Output Parameters:
173: + tol - the convergence tolerance
174: - maxits - maximum number of iterations
176: Notes:
177: The user can specify NULL for any parameter that is not needed.
179: Level: intermediate
181: .seealso: SVDSetTolerances()
182: @*/
183: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,PetscInt *maxits)
184: {
187: if (tol) *tol = svd->tol;
188: if (maxits) *maxits = svd->max_it;
189: return(0);
190: }
194: /*@
195: SVDSetDimensions - Sets the number of singular values to compute
196: and the dimension of the subspace.
198: Logically Collective on SVD
200: Input Parameters:
201: + svd - the singular value solver context
202: . nsv - number of singular values to compute
203: . ncv - the maximum dimension of the subspace to be used by the solver
204: - mpd - the maximum dimension allowed for the projected problem
206: Options Database Keys:
207: + -svd_nsv <nsv> - Sets the number of singular values
208: . -svd_ncv <ncv> - Sets the dimension of the subspace
209: - -svd_mpd <mpd> - Sets the maximum projected dimension
211: Notes:
212: Pass 0 to retain the previous value of any parameter.
214: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
215: dependent on the solution method and the number of singular values required.
217: The parameters ncv and mpd are intimately related, so that the user is advised
218: to set one of them at most. Normal usage is that
219: (a) in cases where nsv is small, the user sets ncv (a reasonable default is 2*nsv); and
220: (b) in cases where nsv is large, the user sets mpd.
222: The value of ncv should always be between nsv and (nsv+mpd), typically
223: ncv=nsv+mpd. If nev is not too large, mpd=nsv is a reasonable choice, otherwise
224: a smaller value should be used.
226: Level: intermediate
228: .seealso: SVDGetDimensions()
229: @*/
230: PetscErrorCode SVDSetDimensions(SVD svd,PetscInt nsv,PetscInt ncv,PetscInt mpd)
231: {
237: if (nsv) {
238: if (nsv<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
239: svd->nsv = nsv;
240: }
241: if (ncv) {
242: if (ncv == PETSC_DEFAULT || ncv == PETSC_DECIDE) {
243: svd->ncv = 0;
244: } else {
245: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
246: svd->ncv = ncv;
247: }
248: svd->setupcalled = 0;
249: }
250: if (mpd) {
251: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
252: svd->mpd = 0;
253: } else {
254: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
255: svd->mpd = mpd;
256: }
257: }
258: return(0);
259: }
263: /*@
264: SVDGetDimensions - Gets the number of singular values to compute
265: and the dimension of the subspace.
267: Not Collective
269: Input Parameter:
270: . svd - the singular value context
272: Output Parameters:
273: + nsv - number of singular values to compute
274: . ncv - the maximum dimension of the subspace to be used by the solver
275: - mpd - the maximum dimension allowed for the projected problem
277: Notes:
278: The user can specify NULL for any parameter that is not needed.
280: Level: intermediate
282: .seealso: SVDSetDimensions()
283: @*/
284: PetscErrorCode SVDGetDimensions(SVD svd,PetscInt *nsv,PetscInt *ncv,PetscInt *mpd)
285: {
288: if (nsv) *nsv = svd->nsv;
289: if (ncv) *ncv = svd->ncv;
290: if (mpd) *mpd = svd->mpd;
291: return(0);
292: }
296: /*@
297: SVDSetWhichSingularTriplets - Specifies which singular triplets are
298: to be sought.
300: Logically Collective on SVD
302: Input Parameter:
303: . svd - singular value solver context obtained from SVDCreate()
305: Output Parameter:
306: . which - which singular triplets are to be sought
308: Possible values:
309: The parameter 'which' can have one of these values:
311: + SVD_LARGEST - largest singular values
312: - SVD_SMALLEST - smallest singular values
314: Options Database Keys:
315: + -svd_largest - Sets largest singular values
316: - -svd_smallest - Sets smallest singular values
318: Level: intermediate
320: .seealso: SVDGetWhichSingularTriplets(), SVDWhich
321: @*/
322: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)
323: {
327: switch (which) {
328: case SVD_LARGEST:
329: case SVD_SMALLEST:
330: if (svd->which != which) {
331: svd->setupcalled = 0;
332: svd->which = which;
333: }
334: break;
335: default:
336: SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
337: }
338: return(0);
339: }
343: /*@C
344: SVDGetWhichSingularTriplets - Returns which singular triplets are
345: to be sought.
347: Not Collective
349: Input Parameter:
350: . svd - singular value solver context obtained from SVDCreate()
352: Output Parameter:
353: . which - which singular triplets are to be sought
355: Notes:
356: See SVDSetWhichSingularTriplets() for possible values of which
358: Level: intermediate
360: .seealso: SVDSetWhichSingularTriplets(), SVDWhich
361: @*/
362: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)
363: {
367: *which = svd->which;
368: return(0);
369: }
373: /*@
374: SVDSetFromOptions - Sets SVD options from the options database.
375: This routine must be called before SVDSetUp() if the user is to be
376: allowed to set the solver type.
378: Collective on SVD
380: Input Parameters:
381: . svd - the singular value solver context
383: Notes:
384: To see all options, run your program with the -help option.
386: Level: beginner
388: .seealso:
389: @*/
390: PetscErrorCode SVDSetFromOptions(SVD svd)
391: {
392: PetscErrorCode ierr;
393: char type[256],monfilename[PETSC_MAX_PATH_LEN];
394: PetscBool flg;
395: const char *mode_list[2] = {"explicit","implicit"};
396: PetscInt i,j,k;
397: PetscReal r;
398: PetscViewer monviewer;
399: SlepcConvMonitor ctx;
403: svd->setupcalled = 0;
404: if (!SVDRegisterAllCalled) { SVDRegisterAll(); }
405: PetscObjectOptionsBegin((PetscObject)svd);
407: PetscOptionsList("-svd_type","Singular Value Solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,256,&flg);
408: if (flg) {
409: SVDSetType(svd,type);
410: } else if (!((PetscObject)svd)->type_name) {
411: SVDSetType(svd,SVDCROSS);
412: }
414: PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",&flg);
416: PetscOptionsEList("-svd_transpose_mode","Transpose SVD mode","SVDSetTransposeMode",mode_list,2,svd->transmode == PETSC_DECIDE ? "decide" : mode_list[svd->transmode],&i,&flg);
417: if (flg) {
418: SVDSetTransposeMode(svd,(SVDTransposeMode)i);
419: }
421: r = i = 0;
422: PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,NULL);
423: PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",svd->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:svd->tol,&r,NULL);
424: SVDSetTolerances(svd,r,i);
426: i = j = k = 0;
427: PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,NULL);
428: PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,NULL);
429: PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,NULL);
430: SVDSetDimensions(svd,i,j,k);
432: PetscOptionsBoolGroupBegin("-svd_largest","compute largest singular values","SVDSetWhichSingularTriplets",&flg);
433: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_LARGEST); }
434: PetscOptionsBoolGroupEnd("-svd_smallest","compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
435: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_SMALLEST); }
437: flg = PETSC_FALSE;
438: PetscOptionsBool("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",flg,&flg,NULL);
439: if (flg) {
440: SVDMonitorCancel(svd);
441: }
443: PetscOptionsString("-svd_monitor_all","Monitor approximate singular values and error estimates","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
444: if (flg) {
445: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)svd),monfilename,&monviewer);
446: SVDMonitorSet(svd,SVDMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
447: SVDSetTrackAll(svd,PETSC_TRUE);
448: }
449: PetscOptionsString("-svd_monitor_conv","Monitor approximate singular values and error estimates as they converge","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
450: if (flg) {
451: PetscNew(struct _n_SlepcConvMonitor,&ctx);
452: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)svd),monfilename,&ctx->viewer);
453: SVDMonitorSet(svd,SVDMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
454: }
455: PetscOptionsString("-svd_monitor","Monitor first unconverged approximate singular value and error estimate","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
456: if (flg) {
457: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)svd),monfilename,&monviewer);
458: SVDMonitorSet(svd,SVDMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
459: }
460: flg = PETSC_FALSE;
461: PetscOptionsBool("-svd_monitor_lg","Monitor first unconverged approximate singular value and error estimate graphically","SVDMonitorSet",flg,&flg,NULL);
462: if (flg) {
463: SVDMonitorSet(svd,SVDMonitorLG,NULL,NULL);
464: }
465: flg = PETSC_FALSE;
466: PetscOptionsBool("-svd_monitor_lg_all","Monitor error estimates graphically","SVDMonitorSet",flg,&flg,NULL);
467: if (flg) {
468: SVDMonitorSet(svd,SVDMonitorLGAll,NULL,NULL);
469: SVDSetTrackAll(svd,PETSC_TRUE);
470: }
471: if (svd->ops->setfromoptions) {
472: (*svd->ops->setfromoptions)(svd);
473: }
475: PetscObjectProcessOptionsHandlers((PetscObject)svd);
476: PetscOptionsEnd();
478: if (!svd->ip) { SVDGetIP(svd,&svd->ip); }
479: IPSetFromOptions(svd->ip);
480: if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
481: DSSetFromOptions(svd->ds);
482: PetscRandomSetFromOptions(svd->rand);
483: return(0);
484: }
488: /*@
489: SVDSetTrackAll - Specifies if the solver must compute the residual norm of all
490: approximate singular value or not.
492: Logically Collective on SVD
494: Input Parameters:
495: + svd - the singular value solver context
496: - trackall - whether to compute all residuals or not
498: Notes:
499: If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
500: the residual norm for each singular value approximation. Computing the residual is
501: usually an expensive operation and solvers commonly compute only the residual
502: associated to the first unconverged singular value.
504: The options '-svd_monitor_all' and '-svd_monitor_lg_all' automatically
505: activate this option.
507: Level: intermediate
509: .seealso: SVDGetTrackAll()
510: @*/
511: PetscErrorCode SVDSetTrackAll(SVD svd,PetscBool trackall)
512: {
516: svd->trackall = trackall;
517: return(0);
518: }
522: /*@
523: SVDGetTrackAll - Returns the flag indicating whether all residual norms must
524: be computed or not.
526: Not Collective
528: Input Parameter:
529: . svd - the singular value solver context
531: Output Parameter:
532: . trackall - the returned flag
534: Level: intermediate
536: .seealso: SVDSetTrackAll()
537: @*/
538: PetscErrorCode SVDGetTrackAll(SVD svd,PetscBool *trackall)
539: {
543: *trackall = svd->trackall;
544: return(0);
545: }
550: /*@C
551: SVDSetOptionsPrefix - Sets the prefix used for searching for all
552: SVD options in the database.
554: Logically Collective on SVD
556: Input Parameters:
557: + svd - the singular value solver context
558: - prefix - the prefix string to prepend to all SVD option requests
560: Notes:
561: A hyphen (-) must NOT be given at the beginning of the prefix name.
562: The first character of all runtime options is AUTOMATICALLY the
563: hyphen.
565: For example, to distinguish between the runtime options for two
566: different SVD contexts, one could call
567: .vb
568: SVDSetOptionsPrefix(svd1,"svd1_")
569: SVDSetOptionsPrefix(svd2,"svd2_")
570: .ve
572: Level: advanced
574: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
575: @*/
576: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)
577: {
579: PetscBool flg1,flg2;
580: EPS eps;
584: if (!svd->ip) { SVDGetIP(svd,&svd->ip); }
585: IPSetOptionsPrefix(svd->ip,prefix);
586: if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
587: DSSetOptionsPrefix(svd->ds,prefix);
588: PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
589: PetscObjectTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
590: PetscObjectTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
591: if (flg1) {
592: SVDCrossGetEPS(svd,&eps);
593: } else if (flg2) {
594: SVDCyclicGetEPS(svd,&eps);
595: }
596: if (flg1 || flg2) {
597: EPSSetOptionsPrefix(eps,prefix);
598: EPSAppendOptionsPrefix(eps,"svd_");
599: }
600: return(0);
601: }
605: /*@C
606: SVDAppendOptionsPrefix - Appends to the prefix used for searching for all
607: SVD options in the database.
609: Logically Collective on SVD
611: Input Parameters:
612: + svd - the singular value solver context
613: - prefix - the prefix string to prepend to all SVD option requests
615: Notes:
616: A hyphen (-) must NOT be given at the beginning of the prefix name.
617: The first character of all runtime options is AUTOMATICALLY the hyphen.
619: Level: advanced
621: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
622: @*/
623: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
624: {
626: PetscBool flg1,flg2;
627: EPS eps;
631: if (!svd->ip) { SVDGetIP(svd,&svd->ip); }
632: IPSetOptionsPrefix(svd->ip,prefix);
633: if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
634: DSSetOptionsPrefix(svd->ds,prefix);
635: PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
636: PetscObjectTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
637: PetscObjectTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
638: if (flg1) {
639: SVDCrossGetEPS(svd,&eps);
640: } else if (flg2) {
641: SVDCyclicGetEPS(svd,&eps);
642: }
643: if (flg1 || flg2) {
644: EPSSetOptionsPrefix(eps,((PetscObject)svd)->prefix);
645: EPSAppendOptionsPrefix(eps,"svd_");
646: }
647: return(0);
648: }
652: /*@C
653: SVDGetOptionsPrefix - Gets the prefix used for searching for all
654: SVD options in the database.
656: Not Collective
658: Input Parameters:
659: . svd - the singular value solver context
661: Output Parameters:
662: . prefix - pointer to the prefix string used is returned
664: Notes: On the fortran side, the user should pass in a string 'prefix' of
665: sufficient length to hold the prefix.
667: Level: advanced
669: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
670: @*/
671: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
672: {
678: PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
679: return(0);
680: }