1: /*
2: SVD routines for setting solver options.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2015, 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: SVDSetImplicitTranspose - Indicates 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: - impl - how to handle the transpose (implicitly or not)
38: Options Database Key:
39: . -svd_implicittranspose - Activate the implicit transpose mode.
41: Notes:
42: By default, the transpose of the matrix is explicitly built (if the matrix
43: has defined the MatTranspose operation).
45: If this flag is set to true, the solver does not build the transpose, but
46: handles it implicitly via MatMultTranspose() (or MatMultHermitianTranspose()
47: in the complex case) operations. This is likely to be more inefficient
48: than the default behaviour, both in sequential and in parallel, but
49: requires less storage.
51: Level: advanced
53: .seealso: SVDGetImplicitTranspose(), SVDSolve(), SVDSetOperator()
54: @*/
55: PetscErrorCode SVDSetImplicitTranspose(SVD svd,PetscBool impl) 56: {
60: if (svd->impltrans!=impl) {
61: svd->impltrans = impl;
62: svd->state = SVD_STATE_INITIAL;
63: }
64: return(0);
65: }
69: /*@
70: SVDGetImplicitTranspose - Gets the mode used to handle the transpose
71: of the matrix associated with the singular value problem.
73: Not Collective
75: Input Parameter:
76: . svd - the singular value solver context
78: Output paramter:
79: . impl - how to handle the transpose (implicitly or not)
81: Level: advanced
83: .seealso: SVDSetImplicitTranspose(), SVDSolve(), SVDSetOperator()
84: @*/
85: PetscErrorCode SVDGetImplicitTranspose(SVD svd,PetscBool *impl) 86: {
90: *impl = svd->impltrans;
91: return(0);
92: }
96: /*@
97: SVDSetTolerances - Sets the tolerance and maximum
98: iteration count used by the default SVD convergence testers.
100: Logically Collective on SVD102: Input Parameters:
103: + svd - the singular value solver context
104: . tol - the convergence tolerance
105: - maxits - maximum number of iterations to use
107: Options Database Keys:
108: + -svd_tol <tol> - Sets the convergence tolerance
109: - -svd_max_it <maxits> - Sets the maximum number of iterations allowed
110: (use PETSC_DECIDE to compute an educated guess based on basis and matrix sizes)
112: Note:
113: Use PETSC_DEFAULT for either argument to assign a reasonably good value.
115: Level: intermediate
117: .seealso: SVDGetTolerances()
118: @*/
119: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,PetscInt maxits)120: {
125: if (tol == PETSC_DEFAULT) {
126: svd->tol = PETSC_DEFAULT;
127: svd->state = SVD_STATE_INITIAL;
128: } else {
129: if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
130: svd->tol = tol;
131: }
132: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
133: svd->max_it = 0;
134: svd->state = SVD_STATE_INITIAL;
135: } else {
136: if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
137: svd->max_it = maxits;
138: }
139: return(0);
140: }
144: /*@
145: SVDGetTolerances - Gets the tolerance and maximum
146: iteration count used by the default SVD convergence tests.
148: Not Collective
150: Input Parameter:
151: . svd - the singular value solver context
153: Output Parameters:
154: + tol - the convergence tolerance
155: - maxits - maximum number of iterations
157: Notes:
158: The user can specify NULL for any parameter that is not needed.
160: Level: intermediate
162: .seealso: SVDSetTolerances()
163: @*/
164: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,PetscInt *maxits)165: {
168: if (tol) *tol = svd->tol;
169: if (maxits) *maxits = svd->max_it;
170: return(0);
171: }
175: /*@
176: SVDSetDimensions - Sets the number of singular values to compute
177: and the dimension of the subspace.
179: Logically Collective on SVD181: Input Parameters:
182: + svd - the singular value solver context
183: . nsv - number of singular values to compute
184: . ncv - the maximum dimension of the subspace to be used by the solver
185: - mpd - the maximum dimension allowed for the projected problem
187: Options Database Keys:
188: + -svd_nsv <nsv> - Sets the number of singular values
189: . -svd_ncv <ncv> - Sets the dimension of the subspace
190: - -svd_mpd <mpd> - Sets the maximum projected dimension
192: Notes:
193: Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
194: dependent on the solution method and the number of singular values required.
196: The parameters ncv and mpd are intimately related, so that the user is advised
197: to set one of them at most. Normal usage is that
198: (a) in cases where nsv is small, the user sets ncv (a reasonable default is 2*nsv); and
199: (b) in cases where nsv is large, the user sets mpd.
201: The value of ncv should always be between nsv and (nsv+mpd), typically
202: ncv=nsv+mpd. If nev is not too large, mpd=nsv is a reasonable choice, otherwise
203: a smaller value should be used.
205: Level: intermediate
207: .seealso: SVDGetDimensions()
208: @*/
209: PetscErrorCode SVDSetDimensions(SVD svd,PetscInt nsv,PetscInt ncv,PetscInt mpd)210: {
216: if (nsv<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
217: svd->nsv = nsv;
218: if (ncv == PETSC_DEFAULT || ncv == PETSC_DECIDE) {
219: svd->ncv = 0;
220: } else {
221: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
222: svd->ncv = ncv;
223: }
224: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
225: svd->mpd = 0;
226: } else {
227: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
228: svd->mpd = mpd;
229: }
230: svd->state = SVD_STATE_INITIAL;
231: return(0);
232: }
236: /*@
237: SVDGetDimensions - Gets the number of singular values to compute
238: and the dimension of the subspace.
240: Not Collective
242: Input Parameter:
243: . svd - the singular value context
245: Output Parameters:
246: + nsv - number of singular values to compute
247: . ncv - the maximum dimension of the subspace to be used by the solver
248: - mpd - the maximum dimension allowed for the projected problem
250: Notes:
251: The user can specify NULL for any parameter that is not needed.
253: Level: intermediate
255: .seealso: SVDSetDimensions()
256: @*/
257: PetscErrorCode SVDGetDimensions(SVD svd,PetscInt *nsv,PetscInt *ncv,PetscInt *mpd)258: {
261: if (nsv) *nsv = svd->nsv;
262: if (ncv) *ncv = svd->ncv;
263: if (mpd) *mpd = svd->mpd;
264: return(0);
265: }
269: /*@
270: SVDSetWhichSingularTriplets - Specifies which singular triplets are
271: to be sought.
273: Logically Collective on SVD275: Input Parameter:
276: . svd - singular value solver context obtained from SVDCreate()
278: Output Parameter:
279: . which - which singular triplets are to be sought
281: Possible values:
282: The parameter 'which' can have one of these values:
284: + SVD_LARGEST - largest singular values
285: - SVD_SMALLEST - smallest singular values
287: Options Database Keys:
288: + -svd_largest - Sets largest singular values
289: - -svd_smallest - Sets smallest singular values
291: Level: intermediate
293: .seealso: SVDGetWhichSingularTriplets(), SVDWhich294: @*/
295: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)296: {
300: switch (which) {
301: case SVD_LARGEST:
302: case SVD_SMALLEST:
303: if (svd->which != which) {
304: svd->state = SVD_STATE_INITIAL;
305: svd->which = which;
306: }
307: break;
308: default:309: SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
310: }
311: return(0);
312: }
316: /*@
317: SVDGetWhichSingularTriplets - Returns which singular triplets are
318: to be sought.
320: Not Collective
322: Input Parameter:
323: . svd - singular value solver context obtained from SVDCreate()
325: Output Parameter:
326: . which - which singular triplets are to be sought
328: Notes:
329: See SVDSetWhichSingularTriplets() for possible values of which
331: Level: intermediate
333: .seealso: SVDSetWhichSingularTriplets(), SVDWhich334: @*/
335: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)336: {
340: *which = svd->which;
341: return(0);
342: }
346: /*@
347: SVDSetFromOptions - Sets SVD options from the options database.
348: This routine must be called before SVDSetUp() if the user is to be
349: allowed to set the solver type.
351: Collective on SVD353: Input Parameters:
354: . svd - the singular value solver context
356: Notes:
357: To see all options, run your program with the -help option.
359: Level: beginner
361: .seealso:
362: @*/
363: PetscErrorCode SVDSetFromOptions(SVD svd)364: {
365: PetscErrorCode ierr;
366: char type[256],monfilename[PETSC_MAX_PATH_LEN];
367: PetscBool flg,val,flg1,flg2,flg3;
368: PetscInt i,j,k;
369: PetscReal r;
370: PetscViewer monviewer;
371: SlepcConvMonitor ctx;
375: SVDRegisterAll();
376: PetscObjectOptionsBegin((PetscObject)svd);
377: PetscOptionsFList("-svd_type","Singular Value Solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,256,&flg);
378: if (flg) {
379: SVDSetType(svd,type);
380: } else if (!((PetscObject)svd)->type_name) {
381: SVDSetType(svd,SVDCROSS);
382: }
384: PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",&flg);
385: PetscOptionsName("-svd_view_vectors","View computed singular vectors","SVDVectorsView",0);
386: PetscOptionsName("-svd_view_values","View computed singular values","SVDValuesView",0);
387: PetscOptionsName("-svd_converged_reason","Print reason for convergence, and number of iterations","SVDReasonView",0);
388: PetscOptionsName("-svd_error_absolute","Print absolute errors of each singular triplet","SVDErrorView",0);
389: PetscOptionsName("-svd_error_relative","Print relative errors of each singular triplet","SVDErrorView",0);
391: PetscOptionsBool("-svd_implicittranspose","Handle matrix transpose implicitly","SVDSetImplicitTranspose",svd->impltrans,&val,&flg);
392: if (flg) {
393: SVDSetImplicitTranspose(svd,val);
394: }
396: i = svd->max_it? svd->max_it: PETSC_DEFAULT;
397: PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,&flg1);
398: r = svd->tol;
399: PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",svd->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:svd->tol,&r,&flg2);
400: if (flg1 || flg2) {
401: SVDSetTolerances(svd,r,i);
402: }
404: i = svd->nsv;
405: PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,&flg1);
406: j = svd->ncv? svd->ncv: PETSC_DEFAULT;
407: PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,&flg2);
408: k = svd->mpd? svd->mpd: PETSC_DEFAULT;
409: PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,&flg3);
410: if (flg1 || flg2 || flg3) {
411: SVDSetDimensions(svd,i,j,k);
412: }
414: PetscOptionsBoolGroupBegin("-svd_largest","compute largest singular values","SVDSetWhichSingularTriplets",&flg);
415: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_LARGEST); }
416: PetscOptionsBoolGroupEnd("-svd_smallest","compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
417: if (flg) { SVDSetWhichSingularTriplets(svd,SVD_SMALLEST); }
419: /* -----------------------------------------------------------------------*/
420: /*
421: Cancels all monitors hardwired into code before call to SVDSetFromOptions()
422: */
423: flg = PETSC_FALSE;
424: PetscOptionsBool("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",flg,&flg,NULL);
425: if (flg) {
426: SVDMonitorCancel(svd);
427: }
429: PetscOptionsString("-svd_monitor_all","Monitor approximate singular values and error estimates","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
430: if (flg) {
431: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)svd),monfilename,&monviewer);
432: SVDMonitorSet(svd,SVDMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
433: SVDSetTrackAll(svd,PETSC_TRUE);
434: }
435: PetscOptionsString("-svd_monitor_conv","Monitor approximate singular values and error estimates as they converge","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
436: if (flg) {
437: PetscNew(&ctx);
438: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)svd),monfilename,&ctx->viewer);
439: SVDMonitorSet(svd,SVDMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
440: }
441: PetscOptionsString("-svd_monitor","Monitor first unconverged approximate singular value and error estimate","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
442: if (flg) {
443: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)svd),monfilename,&monviewer);
444: SVDMonitorSet(svd,SVDMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
445: }
446: flg = PETSC_FALSE;
447: PetscOptionsBool("-svd_monitor_lg","Monitor first unconverged approximate singular value and error estimate graphically","SVDMonitorSet",flg,&flg,NULL);
448: if (flg) {
449: SVDMonitorSet(svd,SVDMonitorLG,NULL,NULL);
450: }
451: flg = PETSC_FALSE;
452: PetscOptionsBool("-svd_monitor_lg_all","Monitor error estimates graphically","SVDMonitorSet",flg,&flg,NULL);
453: if (flg) {
454: SVDMonitorSet(svd,SVDMonitorLGAll,NULL,NULL);
455: SVDSetTrackAll(svd,PETSC_TRUE);
456: }
458: if (svd->ops->setfromoptions) {
459: (*svd->ops->setfromoptions)(PetscOptionsObject,svd);
460: }
461: PetscObjectProcessOptionsHandlers((PetscObject)svd);
462: PetscOptionsEnd();
464: if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
465: BVSetFromOptions(svd->V);
466: BVSetFromOptions(svd->U);
467: if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
468: DSSetFromOptions(svd->ds);
469: PetscRandomSetFromOptions(svd->rand);
470: return(0);
471: }
475: /*@
476: SVDSetTrackAll - Specifies if the solver must compute the residual norm of all
477: approximate singular value or not.
479: Logically Collective on SVD481: Input Parameters:
482: + svd - the singular value solver context
483: - trackall - whether to compute all residuals or not
485: Notes:
486: If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
487: the residual norm for each singular value approximation. Computing the residual is
488: usually an expensive operation and solvers commonly compute only the residual
489: associated to the first unconverged singular value.
491: The options '-svd_monitor_all' and '-svd_monitor_lg_all' automatically
492: activate this option.
494: Level: developer
496: .seealso: SVDGetTrackAll()
497: @*/
498: PetscErrorCode SVDSetTrackAll(SVD svd,PetscBool trackall)499: {
503: svd->trackall = trackall;
504: return(0);
505: }
509: /*@
510: SVDGetTrackAll - Returns the flag indicating whether all residual norms must
511: be computed or not.
513: Not Collective
515: Input Parameter:
516: . svd - the singular value solver context
518: Output Parameter:
519: . trackall - the returned flag
521: Level: developer
523: .seealso: SVDSetTrackAll()
524: @*/
525: PetscErrorCode SVDGetTrackAll(SVD svd,PetscBool *trackall)526: {
530: *trackall = svd->trackall;
531: return(0);
532: }
537: /*@C
538: SVDSetOptionsPrefix - Sets the prefix used for searching for all
539: SVD options in the database.
541: Logically Collective on SVD543: Input Parameters:
544: + svd - the singular value solver context
545: - prefix - the prefix string to prepend to all SVD option requests
547: Notes:
548: A hyphen (-) must NOT be given at the beginning of the prefix name.
549: The first character of all runtime options is AUTOMATICALLY the
550: hyphen.
552: For example, to distinguish between the runtime options for two
553: different SVD contexts, one could call
554: .vb
555: SVDSetOptionsPrefix(svd1,"svd1_")
556: SVDSetOptionsPrefix(svd2,"svd2_")
557: .ve
559: Level: advanced
561: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
562: @*/
563: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)564: {
566: PetscBool flg1,flg2;
567: EPS eps;
571: if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
572: BVSetOptionsPrefix(svd->V,prefix);
573: BVSetOptionsPrefix(svd->U,prefix);
574: if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
575: DSSetOptionsPrefix(svd->ds,prefix);
576: PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
577: PetscObjectTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
578: PetscObjectTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
579: if (flg1) {
580: SVDCrossGetEPS(svd,&eps);
581: } else if (flg2) {
582: SVDCyclicGetEPS(svd,&eps);
583: }
584: if (flg1 || flg2) {
585: EPSSetOptionsPrefix(eps,prefix);
586: EPSAppendOptionsPrefix(eps,"svd_");
587: }
588: return(0);
589: }
593: /*@C
594: SVDAppendOptionsPrefix - Appends to the prefix used for searching for all
595: SVD options in the database.
597: Logically Collective on SVD599: Input Parameters:
600: + svd - the singular value solver context
601: - prefix - the prefix string to prepend to all SVD option requests
603: Notes:
604: A hyphen (-) must NOT be given at the beginning of the prefix name.
605: The first character of all runtime options is AUTOMATICALLY the hyphen.
607: Level: advanced
609: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
610: @*/
611: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)612: {
614: PetscBool flg1,flg2;
615: EPS eps;
619: if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
620: BVSetOptionsPrefix(svd->V,prefix);
621: BVSetOptionsPrefix(svd->U,prefix);
622: if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
623: DSSetOptionsPrefix(svd->ds,prefix);
624: PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
625: PetscObjectTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
626: PetscObjectTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
627: if (flg1) {
628: SVDCrossGetEPS(svd,&eps);
629: } else if (flg2) {
630: SVDCyclicGetEPS(svd,&eps);
631: }
632: if (flg1 || flg2) {
633: EPSSetOptionsPrefix(eps,((PetscObject)svd)->prefix);
634: EPSAppendOptionsPrefix(eps,"svd_");
635: }
636: return(0);
637: }
641: /*@C
642: SVDGetOptionsPrefix - Gets the prefix used for searching for all
643: SVD options in the database.
645: Not Collective
647: Input Parameters:
648: . svd - the singular value solver context
650: Output Parameters:
651: . prefix - pointer to the prefix string used is returned
653: Notes: On the fortran side, the user should pass in a string 'prefix' of
654: sufficient length to hold the prefix.
656: Level: advanced
658: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
659: @*/
660: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])661: {
667: PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
668: return(0);
669: }