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-2007, Universidad Politecnica de Valencia, Spain
8: This file is part of SLEPc. See the README file for conditions of use
9: and additional information.
10: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
11: */
13: #include src/svd/svdimpl.h
17: /*@
18: SVDSetTransposeMode - Sets how to handle the transpose of the matrix
19: associated with the singular value problem.
21: Collective on SVD and Mat
23: Input Parameters:
24: + svd - the singular value solver context
25: - mode - how to compute the transpose, one of SVD_TRANSPOSE_EXPLICIT
26: or SVD_TRANSPOSE_MATMULT (see notes below)
28: Options Database Key:
29: . -svd_transpose_mode <mode> - Indicates the mode flag, where <mode>
30: is one of 'explicit' or 'matmult'.
32: Notes:
33: In the SVD_TRANSPOSE_EXPLICIT mode, the transpose of the matrix is
34: explicitly built.
36: The option SVD_TRANSPOSE_MATMULT does not build the transpose, but
37: handles it implicitly via MatMultTranspose() operations. This is
38: likely to be more inefficient than SVD_TRANSPOSE_EXPLICIT, both in
39: sequential and in parallel, but requires less storage.
41: The default is SVD_TRANSPOSE_EXPLICIT if the matrix has defined the
42: MatTranspose operation, and SVD_TRANSPOSE_MATMULT otherwise.
43:
44: Level: advanced
45:
46: .seealso: SVDGetTransposeMode(), SVDSolve(), SVDSetOperator(),
47: SVDGetOperator(), SVDTransposeMode
48: @*/
49: PetscErrorCode SVDSetTransposeMode(SVD svd,SVDTransposeMode mode)
50: {
53: switch (mode) {
54: case PETSC_DEFAULT:
55: mode = (SVDTransposeMode)PETSC_DECIDE;
56: case SVD_TRANSPOSE_EXPLICIT:
57: case SVD_TRANSPOSE_IMPLICIT:
58: case PETSC_DECIDE:
59: svd->transmode = mode;
60: svd->setupcalled = 0;
61: break;
62: default:
63: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid transpose mode");
64: }
65: return(0);
66: }
70: /*@C
71: SVDGetTransposeMode - Gets the mode use to compute the transpose
72: of the matrix associated with the singular value problem.
74: Not collective
76: Input Parameter:
77: + svd - the singular value solver context
79: Output paramter:
80: + mode - how to compute the transpose, one of SVD_TRANSPOSE_EXPLICIT
81: or SVD_TRANSPOSE_MATMULT
82:
83: Level: advanced
84:
85: .seealso: SVDSetTransposeMode(), SVDSolve(), SVDSetOperator(),
86: SVDGetOperator(), SVDTransposeMode
87: @*/
88: PetscErrorCode SVDGetTransposeMode(SVD svd,SVDTransposeMode *mode)
89: {
93: *mode = svd->transmode;
94: return(0);
95: }
99: /*@
100: SVDSetTolerances - Sets the tolerance and maximum
101: iteration count used by the default SVD convergence testers.
103: Collective on SVD
105: Input Parameters:
106: + svd - the singluar value solver context
107: . tol - the convergence tolerance
108: - maxits - maximum number of iterations to use
110: Options Database Keys:
111: + -svd_tol <tol> - Sets the convergence tolerance
112: - -svd_max_it <maxits> - Sets the maximum number of iterations allowed
113: (use PETSC_DECIDE to compute an educated guess based on basis and matrix sizes)
115: Notes:
116: Use PETSC_IGNORE to retain the previous value of any parameter.
118: Level: intermediate
120: .seealso: SVDGetTolerances()
121: @*/
122: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,int maxits)
123: {
126: if (tol != PETSC_IGNORE) {
127: if (tol == PETSC_DEFAULT) {
128: tol = 1e-7;
129: } else {
130: if (tol < 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
131: svd->tol = tol;
132: }
133: }
134: if (maxits != PETSC_IGNORE) {
135: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
136: svd->max_it = PETSC_DECIDE;
137: svd->setupcalled = 0;
138: } else {
139: if (maxits < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
140: svd->max_it = maxits;
141: }
142: }
143: return(0);
144: }
148: /*@
149: SVDGetTolerances - Gets the tolerance and maximum
150: iteration count used by the default SVD convergence tests.
152: Not Collective
154: Input Parameter:
155: . svd - the singular value solver context
156:
157: Output Parameters:
158: + tol - the convergence tolerance
159: - maxits - maximum number of iterations
161: Notes:
162: The user can specify PETSC_NULL for any parameter that is not needed.
164: Level: intermediate
166: .seealso: SVDSetTolerances()
167: @*/
168: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,int *maxits)
169: {
172: if (tol) *tol = svd->tol;
173: if (maxits) *maxits = svd->max_it;
174: return(0);
175: }
179: /*@
180: SVDSetDimensions - Sets the number of singular values to compute
181: and the dimension of the subspace.
183: Collective on SVD
185: Input Parameters:
186: + svd - the singular value solver context
187: . nsv - number of singular values to compute
188: - ncv - the maximum dimension of the subspace to be used by the solver
190: Options Database Keys:
191: + -svd_nsv <nsv> - Sets the number of singular values
192: - -svd_ncv <ncv> - Sets the dimension of the subspace
194: Notes:
195: Use PETSC_IGNORE to retain the previous value of any parameter.
197: Use PETSC_DECIDE for ncv to assign a reasonably good value, which is
198: dependent on the solution method and the number of singular values required.
200: Level: intermediate
202: .seealso: SVDGetDimensions()
203: @*/
204: PetscErrorCode SVDSetDimensions(SVD svd,int nsv,int ncv)
205: {
209: if (nsv != PETSC_IGNORE) {
210: if (nsv<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
211: svd->nsv = nsv;
212: }
213: if (ncv != PETSC_IGNORE) {
214: if (ncv == PETSC_DEFAULT || ncv == PETSC_DECIDE) {
215: svd->ncv = PETSC_DECIDE;
216: } else {
217: if (ncv<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
218: svd->ncv = ncv;
219: }
220: svd->setupcalled = 0;
221: }
222: return(0);
223: }
227: /*@
228: SVDGetDimensions - Gets the number of singular values to compute
229: and the dimension of the subspace.
231: Not Collective
233: Input Parameter:
234: . svd - the singular value context
235:
236: Output Parameters:
237: + nsv - number of singular values to compute
238: - ncv - the maximum dimension of the subspace to be used by the solver
240: Notes:
241: The user can specify PETSC_NULL for any parameter that is not needed.
243: Level: intermediate
245: .seealso: SVDSetDimensions()
246: @*/
247: PetscErrorCode SVDGetDimensions(SVD svd,int *nsv,int *ncv)
248: {
251: if (nsv) *nsv = svd->nsv;
252: if (ncv) *ncv = svd->ncv;
253: return(0);
254: }
258: /*@
259: SVDSetWhichSingularTriplets - Specifies which singular triplets are
260: to be sought.
262: Collective on SVD
264: Input Parameter:
265: . svd - singular value solver context obtained from SVDCreate()
267: Output Parameter:
268: . which - which singular triplets are to be sought
270: Possible values:
271: The parameter 'which' can have one of these values:
272:
273: + SVD_LARGEST - largest singular values
274: - SVD_SMALLEST - smallest singular values
276: Options Database Keys:
277: + -svd_largest - Sets largest singular values
278: - -svd_smallest - Sets smallest singular values
279:
280: Level: intermediate
282: .seealso: SVDGetWhichSingularTriplets(), SVDWhich
283: @*/
284: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)
285: {
288: switch (which) {
289: case SVD_LARGEST:
290: case SVD_SMALLEST:
291: svd->which = which;
292: break;
293: default:
294: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
295: }
296: return(0);
297: }
301: /*@C
302: SVDGetWhichSingularTriplets - Returns which singular triplets are
303: to be sought.
305: Not Collective
307: Input Parameter:
308: . svd - singular value solver context obtained from SVDCreate()
310: Output Parameter:
311: . which - which singular triplets are to be sought
313: Notes:
314: See SVDSetWhichSingularTriplets() for possible values of which
316: Level: intermediate
318: .seealso: SVDSetWhichSingularTriplets(), SVDWhich
319: @*/
320: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)
321: {
325: *which = svd->which;
326: return(0);
327: }
331: /*@
332: SVDSetFromOptions - Sets SVD options from the options database.
333: This routine must be called before SVDSetUp() if the user is to be
334: allowed to set the solver type.
336: Collective on SVD
338: Input Parameters:
339: . svd - the singular value solver context
341: Notes:
342: To see all options, run your program with the -help option.
344: Level: beginner
346: .seealso:
347: @*/
348: PetscErrorCode SVDSetFromOptions(SVD svd)
349: {
351: char type[256],monfilename[PETSC_MAX_PATH_LEN];
352: PetscTruth flg;
353: const char *mode_list[2] = { "explicit", "implicit" };
354: PetscInt i,j;
355: PetscReal r;
356: PetscViewer monviewer;
360: svd->setupcalled = 0;
361: PetscOptionsBegin(svd->comm,svd->prefix,"Singular Value Solver (SVD) Options","SVD");
363: PetscOptionsList("-svd_type","Singular Value Solver method","SVDSetType",SVDList,(char*)(svd->type_name?svd->type_name:SVDCROSS),type,256,&flg);
364: if (flg) {
365: SVDSetType(svd,type);
366: } else if (!svd->type_name) {
367: SVDSetType(svd,SVDCROSS);
368: }
370: PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",0);
372: PetscOptionsEList("-svd_transpose_mode","Transpose SVD mode","SVDSetTransposeMode",mode_list,2,svd->transmode == PETSC_DECIDE ? "decide" : mode_list[svd->transmode],&i,&flg);
373: if (flg) {
374: SVDSetTransposeMode(svd,(SVDTransposeMode)i);
375: }
377: r = i = PETSC_IGNORE;
378: PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,PETSC_NULL);
379: PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",svd->tol,&r,PETSC_NULL);
380: SVDSetTolerances(svd,r,i);
382: i = j = PETSC_IGNORE;
383: PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->ncv,&i,PETSC_NULL);
384: PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,PETSC_NULL);
385: SVDSetDimensions(svd,i,j);
387: PetscOptionsTruthGroupBegin("-svd_largest","compute largest singular values","SVDSetWhichSingularTriplets",&flg);
388: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_LARGEST); }
389: PetscOptionsTruthGroupEnd("-svd_smallest","compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
390: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_SMALLEST); }
392: PetscOptionsName("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",&flg);
393: if (flg) {
394: SVDMonitorCancel(svd);
395: }
397: PetscOptionsString("-svd_monitor","Monitor approximate singular values and error estimates","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
398: if (flg) {
399: PetscViewerASCIIOpen(svd->comm,monfilename,&monviewer);
400: SVDMonitorSet(svd,SVDMonitorDefault,monviewer,(PetscErrorCode (*)(void*))PetscViewerDestroy);
401: }
402: PetscOptionsName("-svd_monitor_draw","Monitor error estimates graphically","SVDMonitorSet",&flg);
403: if (flg) {
404: SVDMonitorSet(svd,SVDMonitorLG,PETSC_NULL,PETSC_NULL);
405: }
407: PetscOptionsEnd();
408: IPSetFromOptions(svd->ip);
409: if (svd->ops->setfromoptions) {
410: (*svd->ops->setfromoptions)(svd);
411: }
412: return(0);
413: }
417: /*@C
418: SVDSetOptionsPrefix - Sets the prefix used for searching for all
419: SVD options in the database.
421: Collective on SVD
423: Input Parameters:
424: + svd - the singular value solver context
425: - prefix - the prefix string to prepend to all SVD option requests
427: Notes:
428: A hyphen (-) must NOT be given at the beginning of the prefix name.
429: The first character of all runtime options is AUTOMATICALLY the
430: hyphen.
432: For example, to distinguish between the runtime options for two
433: different SVD contexts, one could call
434: .vb
435: SVDSetOptionsPrefix(svd1,"svd1_")
436: SVDSetOptionsPrefix(svd2,"svd2_")
437: .ve
439: Level: advanced
441: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
442: @*/
443: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)
444: {
446: PetscTruth flg1,flg2;
447: EPS eps;
451: PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
452: PetscTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
453: PetscTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
454: if (flg1) {
455: SVDCrossGetEPS(svd,&eps);
456: } else if (flg2) {
457: SVDCyclicGetEPS(svd,&eps);
458: }
459: if (flg1 || flg2) {
460: EPSSetOptionsPrefix(eps,prefix);
461: EPSAppendOptionsPrefix(eps,"svd_");
462: }
463: IPSetOptionsPrefix(svd->ip,prefix);
464: IPAppendOptionsPrefix(svd->ip,"svd_");
465: return(0);
466: }
470: /*@C
471: SVDAppendOptionsPrefix - Appends to the prefix used for searching for all
472: SVD options in the database.
474: Collective on SVD
476: Input Parameters:
477: + svd - the singular value solver context
478: - prefix - the prefix string to prepend to all SVD option requests
480: Notes:
481: A hyphen (-) must NOT be given at the beginning of the prefix name.
482: The first character of all runtime options is AUTOMATICALLY the hyphen.
484: Level: advanced
486: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
487: @*/
488: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
489: {
491: PetscTruth flg1,flg2;
492: EPS eps;
493:
496: PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
497: PetscTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
498: PetscTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
499: if (flg1) {
500: SVDCrossGetEPS(svd,&eps);
501: } else if (flg2) {
502: SVDCyclicGetEPS(svd,&eps);
503: }
504: if (flg1 || flg2) {
505: EPSSetOptionsPrefix(eps,svd->prefix);
506: EPSAppendOptionsPrefix(eps,"svd_");
507: }
508: IPSetOptionsPrefix(svd->ip,svd->prefix);
509: IPAppendOptionsPrefix(svd->ip,"svd_");
510: return(0);
511: }
515: /*@C
516: SVDGetOptionsPrefix - Gets the prefix used for searching for all
517: SVD options in the database.
519: Not Collective
521: Input Parameters:
522: . svd - the singular value solver context
524: Output Parameters:
525: . prefix - pointer to the prefix string used is returned
527: Notes: On the fortran side, the user should pass in a string 'prefix' of
528: sufficient length to hold the prefix.
530: Level: advanced
532: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
533: @*/
534: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
535: {
540: PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
541: return(0);
542: }