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-2009, 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 and Mat
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: switch (mode) {
65: case PETSC_DEFAULT:
66: mode = (SVDTransposeMode)PETSC_DECIDE;
67: case SVD_TRANSPOSE_EXPLICIT:
68: case SVD_TRANSPOSE_IMPLICIT:
69: case PETSC_DECIDE:
70: svd->transmode = mode;
71: svd->setupcalled = 0;
72: break;
73: default:
74: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid transpose mode");
75: }
76: return(0);
77: }
81: /*@C
82: SVDGetTransposeMode - Gets the mode use to compute the transpose
83: of the matrix associated with the singular value problem.
85: Not collective
87: Input Parameter:
88: . svd - the singular value solver context
90: Output paramter:
91: . mode - how to compute the transpose, one of SVD_TRANSPOSE_EXPLICIT
92: or SVD_TRANSPOSE_IMPLICIT
93:
94: Level: advanced
95:
96: .seealso: SVDSetTransposeMode(), SVDSolve(), SVDSetOperator(),
97: SVDGetOperator(), SVDTransposeMode
98: @*/
99: PetscErrorCode SVDGetTransposeMode(SVD svd,SVDTransposeMode *mode)
100: {
104: *mode = svd->transmode;
105: return(0);
106: }
110: /*@
111: SVDSetTolerances - Sets the tolerance and maximum
112: iteration count used by the default SVD convergence testers.
114: Collective on SVD
116: Input Parameters:
117: + svd - the singular value solver context
118: . tol - the convergence tolerance
119: - maxits - maximum number of iterations to use
121: Options Database Keys:
122: + -svd_tol <tol> - Sets the convergence tolerance
123: - -svd_max_it <maxits> - Sets the maximum number of iterations allowed
124: (use PETSC_DECIDE to compute an educated guess based on basis and matrix sizes)
126: Notes:
127: Use PETSC_IGNORE to retain the previous value of any parameter.
129: Level: intermediate
131: .seealso: SVDGetTolerances()
132: @*/
133: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,PetscInt maxits)
134: {
137: if (tol != PETSC_IGNORE) {
138: if (tol == PETSC_DEFAULT) {
139: tol = 1e-7;
140: } else {
141: if (tol < 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
142: svd->tol = tol;
143: }
144: }
145: if (maxits != PETSC_IGNORE) {
146: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
147: svd->max_it = 0;
148: svd->setupcalled = 0;
149: } else {
150: if (maxits < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
151: svd->max_it = maxits;
152: }
153: }
154: return(0);
155: }
159: /*@
160: SVDGetTolerances - Gets the tolerance and maximum
161: iteration count used by the default SVD convergence tests.
163: Not Collective
165: Input Parameter:
166: . svd - the singular value solver context
167:
168: Output Parameters:
169: + tol - the convergence tolerance
170: - maxits - maximum number of iterations
172: Notes:
173: The user can specify PETSC_NULL for any parameter that is not needed.
175: Level: intermediate
177: .seealso: SVDSetTolerances()
178: @*/
179: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,PetscInt *maxits)
180: {
183: if (tol) *tol = svd->tol;
184: if (maxits) *maxits = svd->max_it;
185: return(0);
186: }
190: /*@
191: SVDSetDimensions - Sets the number of singular values to compute
192: and the dimension of the subspace.
194: Collective on SVD
196: Input Parameters:
197: + svd - the singular value solver context
198: . nsv - number of singular values to compute
199: . ncv - the maximum dimension of the subspace to be used by the solver
200: - mpd - the maximum dimension allowed for the projected problem
202: Options Database Keys:
203: + -svd_nsv <nsv> - Sets the number of singular values
204: . -svd_ncv <ncv> - Sets the dimension of the subspace
205: - -svd_mpd <mpd> - Sets the maximum projected dimension
207: Notes:
208: Use PETSC_IGNORE to retain the previous value of any parameter.
210: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
211: dependent on the solution method and the number of singular values required.
213: The parameters ncv and mpd are intimately related, so that the user is advised
214: to set one of them at most. Normal usage is the following:
215: + - In cases where nsv is small, the user sets ncv (a reasonable default is 2*nsv).
216: - - In cases where nsv is large, the user sets mpd.
218: The value of ncv should always be between nsv and (nsv+mpd), typically
219: ncv=nsv+mpd. If nev is not too large, mpd=nsv is a reasonable choice, otherwise
220: a smaller value should be used.
222: Level: intermediate
224: .seealso: SVDGetDimensions()
225: @*/
226: PetscErrorCode SVDSetDimensions(SVD svd,PetscInt nsv,PetscInt ncv,PetscInt mpd)
227: {
231: if (nsv != PETSC_IGNORE) {
232: if (nsv<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
233: svd->nsv = nsv;
234: }
235: if (ncv != PETSC_IGNORE) {
236: if (ncv == PETSC_DEFAULT || ncv == PETSC_DECIDE) {
237: svd->ncv = 0;
238: } else {
239: if (ncv<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
240: svd->ncv = ncv;
241: }
242: svd->setupcalled = 0;
243: }
244: if( mpd != PETSC_IGNORE ) {
245: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
246: svd->mpd = 0;
247: } else {
248: if (mpd<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
249: svd->mpd = mpd;
250: }
251: }
252: return(0);
253: }
257: /*@
258: SVDGetDimensions - Gets the number of singular values to compute
259: and the dimension of the subspace.
261: Not Collective
263: Input Parameter:
264: . svd - the singular value context
265:
266: Output Parameters:
267: + nsv - number of singular values to compute
268: . ncv - the maximum dimension of the subspace to be used by the solver
269: - mpd - the maximum dimension allowed for the projected problem
271: Notes:
272: The user can specify PETSC_NULL for any parameter that is not needed.
274: Level: intermediate
276: .seealso: SVDSetDimensions()
277: @*/
278: PetscErrorCode SVDGetDimensions(SVD svd,PetscInt *nsv,PetscInt *ncv,PetscInt *mpd)
279: {
282: if (nsv) *nsv = svd->nsv;
283: if (ncv) *ncv = svd->ncv;
284: if (mpd) *mpd = svd->mpd;
285: return(0);
286: }
290: /*@
291: SVDSetWhichSingularTriplets - Specifies which singular triplets are
292: to be sought.
294: Collective on SVD
296: Input Parameter:
297: . svd - singular value solver context obtained from SVDCreate()
299: Output Parameter:
300: . which - which singular triplets are to be sought
302: Possible values:
303: The parameter 'which' can have one of these values:
304:
305: + SVD_LARGEST - largest singular values
306: - SVD_SMALLEST - smallest singular values
308: Options Database Keys:
309: + -svd_largest - Sets largest singular values
310: - -svd_smallest - Sets smallest singular values
311:
312: Level: intermediate
314: .seealso: SVDGetWhichSingularTriplets(), SVDWhich
315: @*/
316: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)
317: {
320: switch (which) {
321: case SVD_LARGEST:
322: case SVD_SMALLEST:
323: if (svd->which != which) {
324: svd->setupcalled = 0;
325: svd->which = which;
326: }
327: break;
328: default:
329: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
330: }
331: return(0);
332: }
336: /*@C
337: SVDGetWhichSingularTriplets - Returns which singular triplets are
338: to be sought.
340: Not Collective
342: Input Parameter:
343: . svd - singular value solver context obtained from SVDCreate()
345: Output Parameter:
346: . which - which singular triplets are to be sought
348: Notes:
349: See SVDSetWhichSingularTriplets() for possible values of which
351: Level: intermediate
353: .seealso: SVDSetWhichSingularTriplets(), SVDWhich
354: @*/
355: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)
356: {
360: *which = svd->which;
361: return(0);
362: }
366: /*@
367: SVDSetFromOptions - Sets SVD options from the options database.
368: This routine must be called before SVDSetUp() if the user is to be
369: allowed to set the solver type.
371: Collective on SVD
373: Input Parameters:
374: . svd - the singular value solver context
376: Notes:
377: To see all options, run your program with the -help option.
379: Level: beginner
381: .seealso:
382: @*/
383: PetscErrorCode SVDSetFromOptions(SVD svd)
384: {
386: char type[256],monfilename[PETSC_MAX_PATH_LEN];
387: PetscTruth flg;
388: const char *mode_list[2] = { "explicit", "implicit" };
389: PetscInt i,j,k;
390: PetscReal r;
391: PetscViewerASCIIMonitor monviewer;
395: svd->setupcalled = 0;
396: PetscOptionsBegin(((PetscObject)svd)->comm,((PetscObject)svd)->prefix,"Singular Value Solver (SVD) Options","SVD");
398: PetscOptionsList("-svd_type","Singular Value Solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,256,&flg);
399: if (flg) {
400: SVDSetType(svd,type);
401: } else if (!((PetscObject)svd)->type_name) {
402: SVDSetType(svd,SVDCROSS);
403: }
405: PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",0);
407: PetscOptionsEList("-svd_transpose_mode","Transpose SVD mode","SVDSetTransposeMode",mode_list,2,svd->transmode == PETSC_DECIDE ? "decide" : mode_list[svd->transmode],&i,&flg);
408: if (flg) {
409: SVDSetTransposeMode(svd,(SVDTransposeMode)i);
410: }
412: r = i = PETSC_IGNORE;
413: PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,PETSC_NULL);
414: PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",svd->tol,&r,PETSC_NULL);
415: SVDSetTolerances(svd,r,i);
417: i = j = k = PETSC_IGNORE;
418: PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,PETSC_NULL);
419: PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,PETSC_NULL);
420: PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,PETSC_NULL);
421: SVDSetDimensions(svd,i,j,k);
423: PetscOptionsTruthGroupBegin("-svd_largest","compute largest singular values","SVDSetWhichSingularTriplets",&flg);
424: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_LARGEST); }
425: PetscOptionsTruthGroupEnd("-svd_smallest","compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
426: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_SMALLEST); }
428: PetscOptionsName("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",&flg);
429: if (flg) {
430: SVDMonitorCancel(svd);
431: }
433: PetscOptionsString("-svd_monitor","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,SVDMonitorDefault,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
437: }
438: PetscOptionsName("-svd_monitor_draw","Monitor error estimates graphically","SVDMonitorSet",&flg);
439: if (flg) {
440: SVDMonitorSet(svd,SVDMonitorLG,PETSC_NULL,PETSC_NULL);
441: }
443: PetscOptionsEnd();
444: IPSetFromOptions(svd->ip);
445: if (svd->ops->setfromoptions) {
446: (*svd->ops->setfromoptions)(svd);
447: }
448: return(0);
449: }
453: /*@C
454: SVDSetOptionsPrefix - Sets the prefix used for searching for all
455: SVD options in the database.
457: Collective on SVD
459: Input Parameters:
460: + svd - the singular value solver context
461: - prefix - the prefix string to prepend to all SVD option requests
463: Notes:
464: A hyphen (-) must NOT be given at the beginning of the prefix name.
465: The first character of all runtime options is AUTOMATICALLY the
466: hyphen.
468: For example, to distinguish between the runtime options for two
469: different SVD contexts, one could call
470: .vb
471: SVDSetOptionsPrefix(svd1,"svd1_")
472: SVDSetOptionsPrefix(svd2,"svd2_")
473: .ve
475: Level: advanced
477: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
478: @*/
479: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)
480: {
482: PetscTruth flg1,flg2;
483: EPS eps;
487: PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
488: PetscTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
489: PetscTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
490: if (flg1) {
491: SVDCrossGetEPS(svd,&eps);
492: } else if (flg2) {
493: SVDCyclicGetEPS(svd,&eps);
494: }
495: if (flg1 || flg2) {
496: EPSSetOptionsPrefix(eps,prefix);
497: EPSAppendOptionsPrefix(eps,"svd_");
498: }
499: IPSetOptionsPrefix(svd->ip,prefix);
500: IPAppendOptionsPrefix(svd->ip,"svd_");
501: return(0);
502: }
506: /*@C
507: SVDAppendOptionsPrefix - Appends to the prefix used for searching for all
508: SVD options in the database.
510: Collective on SVD
512: Input Parameters:
513: + svd - the singular value solver context
514: - prefix - the prefix string to prepend to all SVD option requests
516: Notes:
517: A hyphen (-) must NOT be given at the beginning of the prefix name.
518: The first character of all runtime options is AUTOMATICALLY the hyphen.
520: Level: advanced
522: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
523: @*/
524: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
525: {
527: PetscTruth flg1,flg2;
528: EPS eps;
529:
532: PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
533: PetscTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
534: PetscTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
535: if (flg1) {
536: SVDCrossGetEPS(svd,&eps);
537: } else if (flg2) {
538: SVDCyclicGetEPS(svd,&eps);
539: }
540: if (flg1 || flg2) {
541: EPSSetOptionsPrefix(eps,((PetscObject)svd)->prefix);
542: EPSAppendOptionsPrefix(eps,"svd_");
543: }
544: IPSetOptionsPrefix(svd->ip,((PetscObject)svd)->prefix);
545: IPAppendOptionsPrefix(svd->ip,"svd_");
546: return(0);
547: }
551: /*@C
552: SVDGetOptionsPrefix - Gets the prefix used for searching for all
553: SVD options in the database.
555: Not Collective
557: Input Parameters:
558: . svd - the singular value solver context
560: Output Parameters:
561: . prefix - pointer to the prefix string used is returned
563: Notes: On the fortran side, the user should pass in a string 'prefix' of
564: sufficient length to hold the prefix.
566: Level: advanced
568: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
569: @*/
570: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
571: {
576: PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
577: return(0);
578: }