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-2010, Universidad Politecnica de Valencia, Spain
8: This file is part of SLEPc.
9:
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 private/svdimpl.h
28: /*@
29: SVDSetTransposeMode - Sets how to handle the transpose of the matrix
30: associated with the singular value problem.
32: 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() operations. This is
49: likely to be more inefficient than SVD_TRANSPOSE_EXPLICIT, both in
50: sequential and in parallel, but requires less storage.
52: The default is SVD_TRANSPOSE_EXPLICIT if the matrix has defined the
53: MatTranspose operation, and SVD_TRANSPOSE_IMPLICIT otherwise.
54:
55: Level: advanced
56:
57: .seealso: SVDGetTransposeMode(), SVDSolve(), SVDSetOperator(),
58: SVDGetOperator(), SVDTransposeMode
59: @*/
60: PetscErrorCode SVDSetTransposeMode(SVD svd,SVDTransposeMode mode)
61: {
64: if (mode == PETSC_DEFAULT || mode == PETSC_DECIDE) mode = (SVDTransposeMode)PETSC_DECIDE;
65: else switch (mode) {
66: case SVD_TRANSPOSE_EXPLICIT:
67: case SVD_TRANSPOSE_IMPLICIT:
68: svd->transmode = mode;
69: svd->setupcalled = 0;
70: break;
71: default:
72: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid transpose mode");
73: }
74: return(0);
75: }
79: /*@C
80: SVDGetTransposeMode - Gets the mode used to compute the transpose
81: of the matrix associated with the singular value problem.
83: Not collective
85: Input Parameter:
86: . svd - the singular value solver context
88: Output paramter:
89: . mode - how to compute the transpose, one of SVD_TRANSPOSE_EXPLICIT
90: or SVD_TRANSPOSE_IMPLICIT
91:
92: Level: advanced
93:
94: .seealso: SVDSetTransposeMode(), SVDSolve(), SVDSetOperator(),
95: SVDGetOperator(), SVDTransposeMode
96: @*/
97: PetscErrorCode SVDGetTransposeMode(SVD svd,SVDTransposeMode *mode)
98: {
102: *mode = svd->transmode;
103: return(0);
104: }
108: /*@
109: SVDSetTolerances - Sets the tolerance and maximum
110: iteration count used by the default SVD convergence testers.
112: Collective on SVD
114: Input Parameters:
115: + svd - the singular value solver context
116: . tol - the convergence tolerance
117: - maxits - maximum number of iterations to use
119: Options Database Keys:
120: + -svd_tol <tol> - Sets the convergence tolerance
121: - -svd_max_it <maxits> - Sets the maximum number of iterations allowed
122: (use PETSC_DECIDE to compute an educated guess based on basis and matrix sizes)
124: Notes:
125: Use PETSC_IGNORE to retain the previous value of any parameter.
127: Level: intermediate
129: .seealso: SVDGetTolerances()
130: @*/
131: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,PetscInt maxits)
132: {
135: if (tol != PETSC_IGNORE) {
136: if (tol == PETSC_DEFAULT) {
137: tol = 1e-7;
138: } else {
139: if (tol < 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
140: svd->tol = tol;
141: }
142: }
143: if (maxits != PETSC_IGNORE) {
144: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
145: svd->max_it = 0;
146: svd->setupcalled = 0;
147: } else {
148: if (maxits < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
149: svd->max_it = maxits;
150: }
151: }
152: return(0);
153: }
157: /*@
158: SVDGetTolerances - Gets the tolerance and maximum
159: iteration count used by the default SVD convergence tests.
161: Not Collective
163: Input Parameter:
164: . svd - the singular value solver context
165:
166: Output Parameters:
167: + tol - the convergence tolerance
168: - maxits - maximum number of iterations
170: Notes:
171: The user can specify PETSC_NULL for any parameter that is not needed.
173: Level: intermediate
175: .seealso: SVDSetTolerances()
176: @*/
177: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,PetscInt *maxits)
178: {
181: if (tol) *tol = svd->tol;
182: if (maxits) *maxits = svd->max_it;
183: return(0);
184: }
188: /*@
189: SVDSetDimensions - Sets the number of singular values to compute
190: and the dimension of the subspace.
192: Collective on SVD
194: Input Parameters:
195: + svd - the singular value solver context
196: . nsv - number of singular values to compute
197: . ncv - the maximum dimension of the subspace to be used by the solver
198: - mpd - the maximum dimension allowed for the projected problem
200: Options Database Keys:
201: + -svd_nsv <nsv> - Sets the number of singular values
202: . -svd_ncv <ncv> - Sets the dimension of the subspace
203: - -svd_mpd <mpd> - Sets the maximum projected dimension
205: Notes:
206: Use PETSC_IGNORE to retain the previous value of any parameter.
208: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
209: dependent on the solution method and the number of singular values required.
211: The parameters ncv and mpd are intimately related, so that the user is advised
212: to set one of them at most. Normal usage is the following:
213: + - In cases where nsv is small, the user sets ncv (a reasonable default is 2*nsv).
214: - - In cases where nsv is large, the user sets mpd.
216: The value of ncv should always be between nsv and (nsv+mpd), typically
217: ncv=nsv+mpd. If nev is not too large, mpd=nsv is a reasonable choice, otherwise
218: a smaller value should be used.
220: Level: intermediate
222: .seealso: SVDGetDimensions()
223: @*/
224: PetscErrorCode SVDSetDimensions(SVD svd,PetscInt nsv,PetscInt ncv,PetscInt mpd)
225: {
229: if (nsv != PETSC_IGNORE) {
230: if (nsv<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
231: svd->nsv = nsv;
232: }
233: if (ncv != PETSC_IGNORE) {
234: if (ncv == PETSC_DEFAULT || ncv == PETSC_DECIDE) {
235: svd->ncv = 0;
236: } else {
237: if (ncv<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
238: svd->ncv = ncv;
239: }
240: svd->setupcalled = 0;
241: }
242: if( mpd != PETSC_IGNORE ) {
243: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
244: svd->mpd = 0;
245: } else {
246: if (mpd<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
247: svd->mpd = mpd;
248: }
249: }
250: return(0);
251: }
255: /*@
256: SVDGetDimensions - Gets the number of singular values to compute
257: and the dimension of the subspace.
259: Not Collective
261: Input Parameter:
262: . svd - the singular value context
263:
264: Output Parameters:
265: + nsv - number of singular values to compute
266: . ncv - the maximum dimension of the subspace to be used by the solver
267: - mpd - the maximum dimension allowed for the projected problem
269: Notes:
270: The user can specify PETSC_NULL for any parameter that is not needed.
272: Level: intermediate
274: .seealso: SVDSetDimensions()
275: @*/
276: PetscErrorCode SVDGetDimensions(SVD svd,PetscInt *nsv,PetscInt *ncv,PetscInt *mpd)
277: {
280: if (nsv) *nsv = svd->nsv;
281: if (ncv) *ncv = svd->ncv;
282: if (mpd) *mpd = svd->mpd;
283: return(0);
284: }
288: /*@
289: SVDSetWhichSingularTriplets - Specifies which singular triplets are
290: to be sought.
292: Collective on SVD
294: Input Parameter:
295: . svd - singular value solver context obtained from SVDCreate()
297: Output Parameter:
298: . which - which singular triplets are to be sought
300: Possible values:
301: The parameter 'which' can have one of these values:
302:
303: + SVD_LARGEST - largest singular values
304: - SVD_SMALLEST - smallest singular values
306: Options Database Keys:
307: + -svd_largest - Sets largest singular values
308: - -svd_smallest - Sets smallest singular values
309:
310: Level: intermediate
312: .seealso: SVDGetWhichSingularTriplets(), SVDWhich
313: @*/
314: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)
315: {
318: switch (which) {
319: case SVD_LARGEST:
320: case SVD_SMALLEST:
321: if (svd->which != which) {
322: svd->setupcalled = 0;
323: svd->which = which;
324: }
325: break;
326: default:
327: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
328: }
329: return(0);
330: }
334: /*@C
335: SVDGetWhichSingularTriplets - Returns which singular triplets are
336: to be sought.
338: Not Collective
340: Input Parameter:
341: . svd - singular value solver context obtained from SVDCreate()
343: Output Parameter:
344: . which - which singular triplets are to be sought
346: Notes:
347: See SVDSetWhichSingularTriplets() for possible values of which
349: Level: intermediate
351: .seealso: SVDSetWhichSingularTriplets(), SVDWhich
352: @*/
353: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)
354: {
358: *which = svd->which;
359: return(0);
360: }
364: /*@
365: SVDSetFromOptions - Sets SVD options from the options database.
366: This routine must be called before SVDSetUp() if the user is to be
367: allowed to set the solver type.
369: Collective on SVD
371: Input Parameters:
372: . svd - the singular value solver context
374: Notes:
375: To see all options, run your program with the -help option.
377: Level: beginner
379: .seealso:
380: @*/
381: PetscErrorCode SVDSetFromOptions(SVD svd)
382: {
384: char type[256],monfilename[PETSC_MAX_PATH_LEN];
385: PetscTruth flg;
386: const char *mode_list[2] = { "explicit", "implicit" };
387: PetscInt i,j,k;
388: PetscReal r;
389: PetscViewerASCIIMonitor monviewer;
390: SVDMONITOR_CONV *ctx;
394: svd->setupcalled = 0;
395: PetscOptionsBegin(((PetscObject)svd)->comm,((PetscObject)svd)->prefix,"Singular Value Solver (SVD) Options","SVD");
397: PetscOptionsList("-svd_type","Singular Value Solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,256,&flg);
398: if (flg) {
399: SVDSetType(svd,type);
400: } else if (!((PetscObject)svd)->type_name) {
401: SVDSetType(svd,SVDCROSS);
402: }
404: PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",0);
406: PetscOptionsEList("-svd_transpose_mode","Transpose SVD mode","SVDSetTransposeMode",mode_list,2,svd->transmode == PETSC_DECIDE ? "decide" : mode_list[svd->transmode],&i,&flg);
407: if (flg) {
408: SVDSetTransposeMode(svd,(SVDTransposeMode)i);
409: }
411: r = i = PETSC_IGNORE;
412: PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,PETSC_NULL);
413: PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",svd->tol,&r,PETSC_NULL);
414: SVDSetTolerances(svd,r,i);
416: i = j = k = PETSC_IGNORE;
417: PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,PETSC_NULL);
418: PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,PETSC_NULL);
419: PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,PETSC_NULL);
420: SVDSetDimensions(svd,i,j,k);
422: PetscOptionsTruthGroupBegin("-svd_largest","compute largest singular values","SVDSetWhichSingularTriplets",&flg);
423: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_LARGEST); }
424: PetscOptionsTruthGroupEnd("-svd_smallest","compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
425: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_SMALLEST); }
427: flg = PETSC_FALSE;
428: PetscOptionsTruth("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",flg,&flg,PETSC_NULL);
429: if (flg) {
430: SVDMonitorCancel(svd);
431: }
433: PetscOptionsString("-svd_monitor_all","Monitor approximate singular values and error estimates","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
434: if (flg) {
435: PetscViewerASCIIMonitorCreate(((PetscObject)svd)->comm,monfilename,((PetscObject)svd)->tablevel,&monviewer);
436: SVDMonitorSet(svd,SVDMonitorAll,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
437: SVDSetTrackAll(svd,PETSC_TRUE);
438: }
439: PetscOptionsString("-svd_monitor_conv","Monitor approximate singular values and error estimates as they converge","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
440: if (flg) {
441: PetscNew(SVDMONITOR_CONV,&ctx);
442: PetscViewerASCIIMonitorCreate(((PetscObject)svd)->comm,monfilename,((PetscObject)svd)->tablevel,&ctx->viewer);
443: SVDMonitorSet(svd,SVDMonitorConverged,ctx,(PetscErrorCode (*)(void*))SVDMonitorDestroy_Converged);
444: }
445: PetscOptionsString("-svd_monitor","Monitor first unconverged approximate singular value and error estimate","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
446: if (flg) {
447: PetscViewerASCIIMonitorCreate(((PetscObject)svd)->comm,monfilename,((PetscObject)svd)->tablevel,&monviewer);
448: SVDMonitorSet(svd,SVDMonitorFirst,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
449: }
450: flg = PETSC_FALSE;
451: PetscOptionsTruth("-svd_monitor_draw","Monitor first unconverged approximate singular value and error estimate graphically","SVDMonitorSet",flg,&flg,PETSC_NULL);
452: if (flg) {
453: SVDMonitorSet(svd,SVDMonitorLG,PETSC_NULL,PETSC_NULL);
454: }
455: flg = PETSC_FALSE;
456: PetscOptionsTruth("-svd_monitor_draw_all","Monitor error estimates graphically","SVDMonitorSet",flg,&flg,PETSC_NULL);
457: if (flg) {
458: SVDMonitorSet(svd,SVDMonitorLGAll,PETSC_NULL,PETSC_NULL);
459: SVDSetTrackAll(svd,PETSC_TRUE);
460: }
462: PetscOptionsEnd();
463: IPSetFromOptions(svd->ip);
464: if (svd->ops->setfromoptions) {
465: (*svd->ops->setfromoptions)(svd);
466: }
467: return(0);
468: }
472: /*@
473: SVDSetTrackAll - Specifies if the solver must compute the residual norm of all
474: approximate singular value or not.
476: Collective on SVD
478: Input Parameters:
479: + svd - the singular value solver context
480: - trackall - whether to compute all residuals or not
482: Notes:
483: If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
484: the residual norm for each singular value approximation. Computing the residual is
485: usually an expensive operation and solvers commonly compute only the residual
486: associated to the first unconverged singular value.
488: The options '-svd_monitor_all' and '-svd_monitor_draw_all' automatically
489: activate this option.
491: Level: intermediate
493: .seealso: SVDGetTrackAll()
494: @*/
495: PetscErrorCode SVDSetTrackAll(SVD svd,PetscTruth trackall)
496: {
499: svd->trackall = trackall;
500: return(0);
501: }
505: /*@
506: SVDGetTrackAll - Returns the flag indicating whether all residual norms must
507: be computed or not.
509: Not Collective
511: Input Parameter:
512: . svd - the singular value solver context
514: Output Parameter:
515: . trackall - the returned flag
517: Level: intermediate
519: .seealso: SVDSetTrackAll()
520: @*/
521: PetscErrorCode SVDGetTrackAll(SVD svd,PetscTruth *trackall)
522: {
526: *trackall = svd->trackall;
527: return(0);
528: }
533: /*@C
534: SVDSetOptionsPrefix - Sets the prefix used for searching for all
535: SVD options in the database.
537: Collective on SVD
539: Input Parameters:
540: + svd - the singular value solver context
541: - prefix - the prefix string to prepend to all SVD option requests
543: Notes:
544: A hyphen (-) must NOT be given at the beginning of the prefix name.
545: The first character of all runtime options is AUTOMATICALLY the
546: hyphen.
548: For example, to distinguish between the runtime options for two
549: different SVD contexts, one could call
550: .vb
551: SVDSetOptionsPrefix(svd1,"svd1_")
552: SVDSetOptionsPrefix(svd2,"svd2_")
553: .ve
555: Level: advanced
557: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
558: @*/
559: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)
560: {
562: PetscTruth flg1,flg2;
563: EPS eps;
567: PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
568: PetscTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
569: PetscTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
570: if (flg1) {
571: SVDCrossGetEPS(svd,&eps);
572: } else if (flg2) {
573: SVDCyclicGetEPS(svd,&eps);
574: }
575: if (flg1 || flg2) {
576: EPSSetOptionsPrefix(eps,prefix);
577: EPSAppendOptionsPrefix(eps,"svd_");
578: }
579: IPSetOptionsPrefix(svd->ip,prefix);
580: IPAppendOptionsPrefix(svd->ip,"svd_");
581: return(0);
582: }
586: /*@C
587: SVDAppendOptionsPrefix - Appends to the prefix used for searching for all
588: SVD options in the database.
590: Collective on SVD
592: Input Parameters:
593: + svd - the singular value solver context
594: - prefix - the prefix string to prepend to all SVD option requests
596: Notes:
597: A hyphen (-) must NOT be given at the beginning of the prefix name.
598: The first character of all runtime options is AUTOMATICALLY the hyphen.
600: Level: advanced
602: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
603: @*/
604: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
605: {
607: PetscTruth flg1,flg2;
608: EPS eps;
609:
612: PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
613: PetscTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
614: PetscTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
615: if (flg1) {
616: SVDCrossGetEPS(svd,&eps);
617: } else if (flg2) {
618: SVDCyclicGetEPS(svd,&eps);
619: }
620: if (flg1 || flg2) {
621: EPSSetOptionsPrefix(eps,((PetscObject)svd)->prefix);
622: EPSAppendOptionsPrefix(eps,"svd_");
623: }
624: IPSetOptionsPrefix(svd->ip,((PetscObject)svd)->prefix);
625: IPAppendOptionsPrefix(svd->ip,"svd_");
626: return(0);
627: }
631: /*@C
632: SVDGetOptionsPrefix - Gets the prefix used for searching for all
633: SVD options in the database.
635: Not Collective
637: Input Parameters:
638: . svd - the singular value solver context
640: Output Parameters:
641: . prefix - pointer to the prefix string used is returned
643: Notes: On the fortran side, the user should pass in a string 'prefix' of
644: sufficient length to hold the prefix.
646: Level: advanced
648: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
649: @*/
650: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
651: {
656: PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
657: return(0);
658: }