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-2011, Universitat 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>      /*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() 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: {
 65:   if (mode == PETSC_DEFAULT || mode == PETSC_DECIDE) mode = (SVDTransposeMode)PETSC_DECIDE;
 66:   else switch (mode) {
 67:     case SVD_TRANSPOSE_EXPLICIT:
 68:     case SVD_TRANSPOSE_IMPLICIT:
 69:       if (svd->transmode!=mode) {
 70:         svd->transmode = mode;
 71:         svd->setupcalled = 0;
 72:       }
 73:       break;
 74:     default:
 75:       SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid transpose mode");
 76:   }
 77:   return(0);
 78: }

 82: /*@C
 83:    SVDGetTransposeMode - Gets the mode used to compute the transpose 
 84:    of the matrix associated with the singular value problem.

 86:    Not Collective

 88:    Input Parameter:
 89: .  svd  - the singular value solver context

 91:    Output paramter:
 92: .  mode - how to compute the transpose, one of SVD_TRANSPOSE_EXPLICIT
 93:           or SVD_TRANSPOSE_IMPLICIT
 94:    
 95:    Level: advanced
 96:    
 97: .seealso: SVDSetTransposeMode(), SVDSolve(), SVDSetOperator(), 
 98:    SVDGetOperator(), SVDTransposeMode
 99: @*/
100: PetscErrorCode SVDGetTransposeMode(SVD svd,SVDTransposeMode *mode)
101: {
105:   *mode = svd->transmode;
106:   return(0);
107: }

111: /*@
112:    SVDSetTolerances - Sets the tolerance and maximum
113:    iteration count used by the default SVD convergence testers. 

115:    Logically Collective on SVD

117:    Input Parameters:
118: +  svd - the singular value solver context
119: .  tol - the convergence tolerance
120: -  maxits - maximum number of iterations to use

122:    Options Database Keys:
123: +  -svd_tol <tol> - Sets the convergence tolerance
124: -  -svd_max_it <maxits> - Sets the maximum number of iterations allowed 
125:    (use PETSC_DECIDE to compute an educated guess based on basis and matrix sizes)

127:    Notes:
128:    Use PETSC_IGNORE to retain the previous value of any parameter. 

130:    Level: intermediate

132: .seealso: SVDGetTolerances()
133: @*/
134: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,PetscInt maxits)
135: {
140:   if (tol != PETSC_IGNORE) {
141:     if (tol == PETSC_DEFAULT) {
142:       tol = PETSC_DEFAULT;
143:     } else {
144:       if (tol < 0.0) SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
145:       svd->tol = tol;
146:     }
147:   }
148:   if (maxits != PETSC_IGNORE) {
149:     if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
150:       svd->max_it = 0;
151:       svd->setupcalled = 0;
152:     } else {
153:       if (maxits < 0) SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
154:       svd->max_it = maxits;
155:     }
156:   }
157:   return(0);
158: }

162: /*@
163:    SVDGetTolerances - Gets the tolerance and maximum
164:    iteration count used by the default SVD convergence tests. 

166:    Not Collective

168:    Input Parameter:
169: .  svd - the singular value solver context
170:   
171:    Output Parameters:
172: +  tol - the convergence tolerance
173: -  maxits - maximum number of iterations

175:    Notes:
176:    The user can specify PETSC_NULL for any parameter that is not needed.

178:    Level: intermediate

180: .seealso: SVDSetTolerances()
181: @*/
182: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,PetscInt *maxits)
183: {
186:   if (tol)    *tol    = svd->tol;
187:   if (maxits) *maxits = svd->max_it;
188:   return(0);
189: }

193: /*@
194:    SVDSetDimensions - Sets the number of singular values to compute
195:    and the dimension of the subspace.

197:    Logically Collective on SVD

199:    Input Parameters:
200: +  svd - the singular value solver context
201: .  nsv - number of singular values to compute
202: .  ncv - the maximum dimension of the subspace to be used by the solver
203: -  mpd - the maximum dimension allowed for the projected problem

205:    Options Database Keys:
206: +  -svd_nsv <nsv> - Sets the number of singular values
207: .  -svd_ncv <ncv> - Sets the dimension of the subspace
208: -  -svd_mpd <mpd> - Sets the maximum projected dimension

210:    Notes:
211:    Use PETSC_IGNORE to retain the previous value of any parameter.

213:    Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is 
214:    dependent on the solution method and the number of singular values required.

216:    The parameters ncv and mpd are intimately related, so that the user is advised
217:    to set one of them at most. Normal usage is that
218:    (a) in cases where nsv is small, the user sets ncv (a reasonable default is 2*nsv); and
219:    (b) in cases where nsv is large, the user sets mpd.

221:    The value of ncv should always be between nsv and (nsv+mpd), typically
222:    ncv=nsv+mpd. If nev is not too large, mpd=nsv is a reasonable choice, otherwise
223:    a smaller value should be used.

225:    Level: intermediate

227: .seealso: SVDGetDimensions()
228: @*/
229: PetscErrorCode SVDSetDimensions(SVD svd,PetscInt nsv,PetscInt ncv,PetscInt mpd)
230: {
236:   if (nsv != PETSC_IGNORE) {
237:     if (nsv<1) SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
238:     svd->nsv = nsv;
239:   }
240:   if (ncv != PETSC_IGNORE) {
241:     if (ncv == PETSC_DEFAULT || ncv == PETSC_DECIDE) {
242:       svd->ncv = 0;
243:     } else {
244:       if (ncv<1) SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
245:       svd->ncv = ncv;
246:     }
247:     svd->setupcalled = 0;
248:   }
249:   if (mpd != PETSC_IGNORE) {
250:     if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
251:       svd->mpd = 0;
252:     } else {
253:       if (mpd<1) SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
254:       svd->mpd = mpd;
255:     }
256:   }
257:   return(0);
258: }

262: /*@
263:    SVDGetDimensions - Gets the number of singular values to compute
264:    and the dimension of the subspace.

266:    Not Collective

268:    Input Parameter:
269: .  svd - the singular value context
270:   
271:    Output Parameters:
272: +  nsv - number of singular values to compute
273: .  ncv - the maximum dimension of the subspace to be used by the solver
274: -  mpd - the maximum dimension allowed for the projected problem

276:    Notes:
277:    The user can specify PETSC_NULL for any parameter that is not needed.

279:    Level: intermediate

281: .seealso: SVDSetDimensions()
282: @*/
283: PetscErrorCode SVDGetDimensions(SVD svd,PetscInt *nsv,PetscInt *ncv,PetscInt *mpd)
284: {
287:   if (nsv) *nsv = svd->nsv;
288:   if (ncv) *ncv = svd->ncv;
289:   if (mpd) *mpd = svd->mpd;
290:   return(0);
291: }

295: /*@
296:     SVDSetWhichSingularTriplets - Specifies which singular triplets are 
297:     to be sought.

299:     Logically Collective on SVD

301:     Input Parameter:
302: .   svd - singular value solver context obtained from SVDCreate()

304:     Output Parameter:
305: .   which - which singular triplets are to be sought

307:     Possible values:
308:     The parameter 'which' can have one of these values:
309:     
310: +     SVD_LARGEST  - largest singular values
311: -     SVD_SMALLEST - smallest singular values

313:     Options Database Keys:
314: +   -svd_largest  - Sets largest singular values
315: -   -svd_smallest - Sets smallest singular values
316:     
317:     Level: intermediate

319: .seealso: SVDGetWhichSingularTriplets(), SVDWhich
320: @*/
321: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)
322: {
326:   switch (which) {
327:     case SVD_LARGEST:
328:     case SVD_SMALLEST:
329:       if (svd->which != which) {
330:         svd->setupcalled = 0;
331:         svd->which = which;
332:       }
333:       break;
334:   default:
335:     SETERRQ(((PetscObject)svd)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
336:   }
337:   return(0);
338: }

342: /*@C
343:     SVDGetWhichSingularTriplets - Returns which singular triplets are
344:     to be sought.

346:     Not Collective

348:     Input Parameter:
349: .   svd - singular value solver context obtained from SVDCreate()

351:     Output Parameter:
352: .   which - which singular triplets are to be sought

354:     Notes:
355:     See SVDSetWhichSingularTriplets() for possible values of which

357:     Level: intermediate

359: .seealso: SVDSetWhichSingularTriplets(), SVDWhich
360: @*/
361: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)
362: {
366:   *which = svd->which;
367:   return(0);
368: }

372: /*@
373:    SVDSetFromOptions - Sets SVD options from the options database.
374:    This routine must be called before SVDSetUp() if the user is to be 
375:    allowed to set the solver type. 

377:    Collective on SVD

379:    Input Parameters:
380: .  svd - the singular value solver context

382:    Notes:  
383:    To see all options, run your program with the -help option.

385:    Level: beginner

387: .seealso: 
388: @*/
389: PetscErrorCode SVDSetFromOptions(SVD svd)
390: {
391:   PetscErrorCode   ierr;
392:   char             type[256],monfilename[PETSC_MAX_PATH_LEN];
393:   PetscBool        flg;
394:   const char       *mode_list[2] = {"explicit","implicit"};
395:   PetscInt         i,j,k;
396:   PetscReal        r;
397:   PetscViewer      monviewer;
398:   SlepcConvMonitor ctx;

402:   svd->setupcalled = 0;
403:   if (!SVDRegisterAllCalled) { SVDRegisterAll(PETSC_NULL); }
404:   PetscOptionsBegin(((PetscObject)svd)->comm,((PetscObject)svd)->prefix,"Singular Value Solver (SVD) Options","SVD");

406:   PetscOptionsList("-svd_type","Singular Value Solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,256,&flg);
407:   if (flg) {
408:     SVDSetType(svd,type);
409:   } else if (!((PetscObject)svd)->type_name) {
410:     SVDSetType(svd,SVDCROSS);
411:   }

413:   PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",&flg);

415:   PetscOptionsEList("-svd_transpose_mode","Transpose SVD mode","SVDSetTransposeMode",mode_list,2,svd->transmode == PETSC_DECIDE ? "decide" : mode_list[svd->transmode],&i,&flg);
416:   if (flg) {
417:     SVDSetTransposeMode(svd,(SVDTransposeMode)i);
418:   }

420:   r = i = PETSC_IGNORE;
421:   PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,PETSC_NULL);
422:   PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",svd->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:svd->tol,&r,PETSC_NULL);
423:   SVDSetTolerances(svd,r,i);

425:   i = j = k = PETSC_IGNORE;
426:   PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,PETSC_NULL);
427:   PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,PETSC_NULL);
428:   PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,PETSC_NULL);
429:   SVDSetDimensions(svd,i,j,k);

431:   PetscOptionsBoolGroupBegin("-svd_largest","compute largest singular values","SVDSetWhichSingularTriplets",&flg);
432:   if (flg) { SVDSetWhichSingularTriplets(svd,SVD_LARGEST); }
433:   PetscOptionsBoolGroupEnd("-svd_smallest","compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
434:   if (flg) { SVDSetWhichSingularTriplets(svd,SVD_SMALLEST); }

436:   flg = PETSC_FALSE;
437:   PetscOptionsBool("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",flg,&flg,PETSC_NULL);
438:   if (flg) {
439:     SVDMonitorCancel(svd);
440:   }

442:   PetscOptionsString("-svd_monitor_all","Monitor approximate singular values and error estimates","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
443:   if (flg) {
444:     PetscViewerASCIIOpen(((PetscObject)svd)->comm,monfilename,&monviewer);
445:     SVDMonitorSet(svd,SVDMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
446:     SVDSetTrackAll(svd,PETSC_TRUE);
447:   }
448:   PetscOptionsString("-svd_monitor_conv","Monitor approximate singular values and error estimates as they converge","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
449:   if (flg) {
450:       PetscNew(struct _n_SlepcConvMonitor,&ctx);
451:       PetscViewerASCIIOpen(((PetscObject)svd)->comm,monfilename,&ctx->viewer);
452:       SVDMonitorSet(svd,SVDMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
453:   }
454:   PetscOptionsString("-svd_monitor","Monitor first unconverged approximate singular value and error estimate","SVDMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
455:   if (flg) {
456:     PetscViewerASCIIOpen(((PetscObject)svd)->comm,monfilename,&monviewer);
457:     SVDMonitorSet(svd,SVDMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
458:   }
459:   flg = PETSC_FALSE;
460:   PetscOptionsBool("-svd_monitor_draw","Monitor first unconverged approximate singular value and error estimate graphically","SVDMonitorSet",flg,&flg,PETSC_NULL);
461:   if (flg) {
462:     SVDMonitorSet(svd,SVDMonitorLG,PETSC_NULL,PETSC_NULL);
463:   }
464:   flg = PETSC_FALSE;
465:   PetscOptionsBool("-svd_monitor_draw_all","Monitor error estimates graphically","SVDMonitorSet",flg,&flg,PETSC_NULL);
466:   if (flg) {
467:     SVDMonitorSet(svd,SVDMonitorLGAll,PETSC_NULL,PETSC_NULL);
468:     SVDSetTrackAll(svd,PETSC_TRUE);
469:   }
470:   if (svd->ops->setfromoptions) {
471:     (*svd->ops->setfromoptions)(svd);
472:   }

474:   PetscObjectProcessOptionsHandlers((PetscObject)svd);
475:   PetscOptionsEnd();

477:   if (!svd->ip) { SVDGetIP(svd,&svd->ip); }
478:   IPSetFromOptions(svd->ip);
479:   PetscRandomSetFromOptions(svd->rand);
480:   return(0);
481: }

485: /*@
486:    SVDSetTrackAll - Specifies if the solver must compute the residual norm of all
487:    approximate singular value or not.

489:    Logically Collective on SVD

491:    Input Parameters:
492: +  svd      - the singular value solver context
493: -  trackall - whether to compute all residuals or not

495:    Notes:
496:    If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
497:    the residual norm for each singular value approximation. Computing the residual is
498:    usually an expensive operation and solvers commonly compute only the residual 
499:    associated to the first unconverged singular value.

501:    The options '-svd_monitor_all' and '-svd_monitor_draw_all' automatically
502:    activate this option.

504:    Level: intermediate

506: .seealso: SVDGetTrackAll()
507: @*/
508: PetscErrorCode SVDSetTrackAll(SVD svd,PetscBool trackall)
509: {
513:   svd->trackall = trackall;
514:   return(0);
515: }

519: /*@
520:    SVDGetTrackAll - Returns the flag indicating whether all residual norms must
521:    be computed or not.

523:    Not Collective

525:    Input Parameter:
526: .  svd - the singular value solver context

528:    Output Parameter:
529: .  trackall - the returned flag

531:    Level: intermediate

533: .seealso: SVDSetTrackAll()
534: @*/
535: PetscErrorCode SVDGetTrackAll(SVD svd,PetscBool *trackall)
536: {
540:   *trackall = svd->trackall;
541:   return(0);
542: }


547: /*@C
548:    SVDSetOptionsPrefix - Sets the prefix used for searching for all 
549:    SVD options in the database.

551:    Logically Collective on SVD

553:    Input Parameters:
554: +  svd - the singular value solver context
555: -  prefix - the prefix string to prepend to all SVD option requests

557:    Notes:
558:    A hyphen (-) must NOT be given at the beginning of the prefix name.
559:    The first character of all runtime options is AUTOMATICALLY the
560:    hyphen.

562:    For example, to distinguish between the runtime options for two
563:    different SVD contexts, one could call
564: .vb
565:       SVDSetOptionsPrefix(svd1,"svd1_")
566:       SVDSetOptionsPrefix(svd2,"svd2_")
567: .ve

569:    Level: advanced

571: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
572: @*/
573: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)
574: {
576:   PetscBool      flg1,flg2;
577:   EPS            eps;

581:   if (!svd->ip) { SVDGetIP(svd,&svd->ip); }
582:   IPSetOptionsPrefix(svd->ip,prefix);
583:   PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
584:   PetscTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
585:   PetscTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
586:   if (flg1) {
587:     SVDCrossGetEPS(svd,&eps);
588:   } else if (flg2) {
589:       SVDCyclicGetEPS(svd,&eps);
590:   }
591:   if (flg1 || flg2) {
592:     EPSSetOptionsPrefix(eps,prefix);
593:     EPSAppendOptionsPrefix(eps,"svd_");
594:   }
595:   return(0);
596: }

600: /*@C
601:    SVDAppendOptionsPrefix - Appends to the prefix used for searching for all 
602:    SVD options in the database.

604:    Logically Collective on SVD

606:    Input Parameters:
607: +  svd - the singular value solver context
608: -  prefix - the prefix string to prepend to all SVD option requests

610:    Notes:
611:    A hyphen (-) must NOT be given at the beginning of the prefix name.
612:    The first character of all runtime options is AUTOMATICALLY the hyphen.

614:    Level: advanced

616: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
617: @*/
618: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
619: {
621:   PetscBool      flg1,flg2;
622:   EPS            eps;
623: 
626:   if (!svd->ip) { SVDGetIP(svd,&svd->ip); }
627:   IPSetOptionsPrefix(svd->ip,prefix);
628:   PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
629:   PetscTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
630:   PetscTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
631:   if (flg1) {
632:     SVDCrossGetEPS(svd,&eps);
633:   } else if (flg2) {
634:     SVDCyclicGetEPS(svd,&eps);
635:   }
636:   if (flg1 || flg2) {
637:     EPSSetOptionsPrefix(eps,((PetscObject)svd)->prefix);
638:     EPSAppendOptionsPrefix(eps,"svd_");
639:   }
640:   return(0);
641: }

645: /*@C
646:    SVDGetOptionsPrefix - Gets the prefix used for searching for all 
647:    SVD options in the database.

649:    Not Collective

651:    Input Parameters:
652: .  svd - the singular value solver context

654:    Output Parameters:
655: .  prefix - pointer to the prefix string used is returned

657:    Notes: On the fortran side, the user should pass in a string 'prefix' of
658:    sufficient length to hold the prefix.

660:    Level: advanced

662: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
663: @*/
664: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
665: {

671:   PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
672:   return(0);
673: }