Actual source code: svdopts.c

slepc-3.11.2 2019-07-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    SVD routines for setting solver options
 12: */

 14: #include <slepc/private/svdimpl.h>      /*I "slepcsvd.h" I*/
 15: #include <petscdraw.h>

 17: /*@
 18:    SVDSetImplicitTranspose - Indicates how to handle the transpose of the matrix
 19:    associated with the singular value problem.

 21:    Logically Collective on SVD

 23:    Input Parameters:
 24: +  svd  - the singular value solver context
 25: -  impl - how to handle the transpose (implicitly or not)

 27:    Options Database Key:
 28: .  -svd_implicittranspose - Activate the implicit transpose mode.

 30:    Notes:
 31:    By default, the transpose of the matrix is explicitly built (if the matrix
 32:    has defined the MatTranspose operation).

 34:    If this flag is set to true, the solver does not build the transpose, but
 35:    handles it implicitly via MatMultTranspose() (or MatMultHermitianTranspose()
 36:    in the complex case) operations. This is likely to be more inefficient
 37:    than the default behaviour, both in sequential and in parallel, but
 38:    requires less storage.

 40:    Level: advanced

 42: .seealso: SVDGetImplicitTranspose(), SVDSolve(), SVDSetOperator()
 43: @*/
 44: PetscErrorCode SVDSetImplicitTranspose(SVD svd,PetscBool impl)
 45: {
 49:   if (svd->impltrans!=impl) {
 50:     svd->impltrans = impl;
 51:     svd->state     = SVD_STATE_INITIAL;
 52:   }
 53:   return(0);
 54: }

 56: /*@
 57:    SVDGetImplicitTranspose - Gets the mode used to handle the transpose
 58:    of the matrix associated with the singular value problem.

 60:    Not Collective

 62:    Input Parameter:
 63: .  svd  - the singular value solver context

 65:    Output Parameter:
 66: .  impl - how to handle the transpose (implicitly or not)

 68:    Level: advanced

 70: .seealso: SVDSetImplicitTranspose(), SVDSolve(), SVDSetOperator()
 71: @*/
 72: PetscErrorCode SVDGetImplicitTranspose(SVD svd,PetscBool *impl)
 73: {
 77:   *impl = svd->impltrans;
 78:   return(0);
 79: }

 81: /*@
 82:    SVDSetTolerances - Sets the tolerance and maximum
 83:    iteration count used by the default SVD convergence testers.

 85:    Logically Collective on SVD

 87:    Input Parameters:
 88: +  svd - the singular value solver context
 89: .  tol - the convergence tolerance
 90: -  maxits - maximum number of iterations to use

 92:    Options Database Keys:
 93: +  -svd_tol <tol> - Sets the convergence tolerance
 94: -  -svd_max_it <maxits> - Sets the maximum number of iterations allowed

 96:    Note:
 97:    Use PETSC_DEFAULT for either argument to assign a reasonably good value.

 99:    Level: intermediate

101: .seealso: SVDGetTolerances()
102: @*/
103: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,PetscInt maxits)
104: {
109:   if (tol == PETSC_DEFAULT) {
110:     svd->tol   = PETSC_DEFAULT;
111:     svd->state = SVD_STATE_INITIAL;
112:   } else {
113:     if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
114:     svd->tol = tol;
115:   }
116:   if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
117:     svd->max_it = 0;
118:     svd->state  = SVD_STATE_INITIAL;
119:   } else {
120:     if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
121:     svd->max_it = maxits;
122:   }
123:   return(0);
124: }

126: /*@C
127:    SVDGetTolerances - Gets the tolerance and maximum
128:    iteration count used by the default SVD convergence tests.

130:    Not Collective

132:    Input Parameter:
133: .  svd - the singular value solver context

135:    Output Parameters:
136: +  tol - the convergence tolerance
137: -  maxits - maximum number of iterations

139:    Notes:
140:    The user can specify NULL for any parameter that is not needed.

142:    Level: intermediate

144: .seealso: SVDSetTolerances()
145: @*/
146: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,PetscInt *maxits)
147: {
150:   if (tol)    *tol    = svd->tol;
151:   if (maxits) *maxits = svd->max_it;
152:   return(0);
153: }

155: /*@
156:    SVDSetDimensions - Sets the number of singular values to compute
157:    and the dimension of the subspace.

159:    Logically Collective on SVD

161:    Input Parameters:
162: +  svd - the singular value solver context
163: .  nsv - number of singular values to compute
164: .  ncv - the maximum dimension of the subspace to be used by the solver
165: -  mpd - the maximum dimension allowed for the projected problem

167:    Options Database Keys:
168: +  -svd_nsv <nsv> - Sets the number of singular values
169: .  -svd_ncv <ncv> - Sets the dimension of the subspace
170: -  -svd_mpd <mpd> - Sets the maximum projected dimension

172:    Notes:
173:    Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
174:    dependent on the solution method and the number of singular values required.

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

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

185:    Level: intermediate

187: .seealso: SVDGetDimensions()
188: @*/
189: PetscErrorCode SVDSetDimensions(SVD svd,PetscInt nsv,PetscInt ncv,PetscInt mpd)
190: {
196:   if (nsv<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
197:   svd->nsv = nsv;
198:   if (ncv == PETSC_DEFAULT || ncv == PETSC_DECIDE) {
199:     svd->ncv = 0;
200:   } else {
201:     if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
202:     svd->ncv = ncv;
203:   }
204:   if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
205:     svd->mpd = 0;
206:   } else {
207:     if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
208:     svd->mpd = mpd;
209:   }
210:   svd->state = SVD_STATE_INITIAL;
211:   return(0);
212: }

214: /*@C
215:    SVDGetDimensions - Gets the number of singular values to compute
216:    and the dimension of the subspace.

218:    Not Collective

220:    Input Parameter:
221: .  svd - the singular value context

223:    Output Parameters:
224: +  nsv - number of singular values to compute
225: .  ncv - the maximum dimension of the subspace to be used by the solver
226: -  mpd - the maximum dimension allowed for the projected problem

228:    Notes:
229:    The user can specify NULL for any parameter that is not needed.

231:    Level: intermediate

233: .seealso: SVDSetDimensions()
234: @*/
235: PetscErrorCode SVDGetDimensions(SVD svd,PetscInt *nsv,PetscInt *ncv,PetscInt *mpd)
236: {
239:   if (nsv) *nsv = svd->nsv;
240:   if (ncv) *ncv = svd->ncv;
241:   if (mpd) *mpd = svd->mpd;
242:   return(0);
243: }

245: /*@
246:     SVDSetWhichSingularTriplets - Specifies which singular triplets are
247:     to be sought.

249:     Logically Collective on SVD

251:     Input Parameter:
252: .   svd - singular value solver context obtained from SVDCreate()

254:     Output Parameter:
255: .   which - which singular triplets are to be sought

257:     Possible values:
258:     The parameter 'which' can have one of these values

260: +     SVD_LARGEST  - largest singular values
261: -     SVD_SMALLEST - smallest singular values

263:     Options Database Keys:
264: +   -svd_largest  - Sets largest singular values
265: -   -svd_smallest - Sets smallest singular values

267:     Level: intermediate

269: .seealso: SVDGetWhichSingularTriplets(), SVDWhich
270: @*/
271: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)
272: {
276:   switch (which) {
277:     case SVD_LARGEST:
278:     case SVD_SMALLEST:
279:       if (svd->which != which) {
280:         svd->state = SVD_STATE_INITIAL;
281:         svd->which = which;
282:       }
283:       break;
284:   default:
285:     SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
286:   }
287:   return(0);
288: }

290: /*@
291:     SVDGetWhichSingularTriplets - Returns which singular triplets are
292:     to be sought.

294:     Not Collective

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:     Notes:
303:     See SVDSetWhichSingularTriplets() for possible values of which

305:     Level: intermediate

307: .seealso: SVDSetWhichSingularTriplets(), SVDWhich
308: @*/
309: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)
310: {
314:   *which = svd->which;
315:   return(0);
316: }

318: /*@C
319:    SVDSetConvergenceTestFunction - Sets a function to compute the error estimate
320:    used in the convergence test.

322:    Logically Collective on SVD

324:    Input Parameters:
325: +  svd     - singular value solver context obtained from SVDCreate()
326: .  func    - a pointer to the convergence test function
327: .  ctx     - context for private data for the convergence routine (may be null)
328: -  destroy - a routine for destroying the context (may be null)

330:    Calling Sequence of func:
331: $   func(SVD svd,PetscReal sigma,PetscReal res,PetscReal *errest,void *ctx)

333: +   svd    - singular value solver context obtained from SVDCreate()
334: .   sigma  - computed singular value
335: .   res    - residual norm associated to the singular triplet
336: .   errest - (output) computed error estimate
337: -   ctx    - optional context, as set by SVDSetConvergenceTestFunction()

339:    Note:
340:    If the error estimate returned by the convergence test function is less than
341:    the tolerance, then the singular value is accepted as converged.

343:    Level: advanced

345: .seealso: SVDSetConvergenceTest(), SVDSetTolerances()
346: @*/
347: PetscErrorCode SVDSetConvergenceTestFunction(SVD svd,PetscErrorCode (*func)(SVD,PetscReal,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
348: {

353:   if (svd->convergeddestroy) {
354:     (*svd->convergeddestroy)(svd->convergedctx);
355:   }
356:   svd->convergeduser    = func;
357:   svd->convergeddestroy = destroy;
358:   svd->convergedctx     = ctx;
359:   if (func == SVDConvergedRelative) svd->conv = SVD_CONV_REL;
360:   else if (func == SVDConvergedAbsolute) svd->conv = SVD_CONV_ABS;
361:   else {
362:     svd->conv      = SVD_CONV_USER;
363:     svd->converged = svd->convergeduser;
364:   }
365:   return(0);
366: }

368: /*@
369:    SVDSetConvergenceTest - Specifies how to compute the error estimate
370:    used in the convergence test.

372:    Logically Collective on SVD

374:    Input Parameters:
375: +  svd  - singular value solver context obtained from SVDCreate()
376: -  conv - the type of convergence test

378:    Options Database Keys:
379: +  -svd_conv_abs  - Sets the absolute convergence test
380: .  -svd_conv_rel  - Sets the convergence test relative to the singular value
381: -  -svd_conv_user - Selects the user-defined convergence test

383:    Note:
384:    The parameter 'conv' can have one of these values
385: +     SVD_CONV_ABS  - absolute error ||r||
386: .     SVD_CONV_REL  - error relative to the singular value l, ||r||/sigma
387: -     SVD_CONV_USER - function set by SVDSetConvergenceTestFunction()

389:    Level: intermediate

391: .seealso: SVDGetConvergenceTest(), SVDSetConvergenceTestFunction(), SVDSetStoppingTest(), SVDConv
392: @*/
393: PetscErrorCode SVDSetConvergenceTest(SVD svd,SVDConv conv)
394: {
398:   switch (conv) {
399:     case SVD_CONV_ABS:  svd->converged = SVDConvergedAbsolute; break;
400:     case SVD_CONV_REL:  svd->converged = SVDConvergedRelative; break;
401:     case SVD_CONV_USER:
402:       if (!svd->convergeduser) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ORDER,"Must call SVDSetConvergenceTestFunction() first");
403:       svd->converged = svd->convergeduser;
404:       break;
405:     default:
406:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
407:   }
408:   svd->conv = conv;
409:   return(0);
410: }

412: /*@
413:    SVDGetConvergenceTest - Gets the method used to compute the error estimate
414:    used in the convergence test.

416:    Not Collective

418:    Input Parameters:
419: .  svd   - singular value solver context obtained from SVDCreate()

421:    Output Parameters:
422: .  conv  - the type of convergence test

424:    Level: intermediate

426: .seealso: SVDSetConvergenceTest(), SVDConv
427: @*/
428: PetscErrorCode SVDGetConvergenceTest(SVD svd,SVDConv *conv)
429: {
433:   *conv = svd->conv;
434:   return(0);
435: }

437: /*@C
438:    SVDSetStoppingTestFunction - Sets a function to decide when to stop the outer
439:    iteration of the singular value solver.

441:    Logically Collective on SVD

443:    Input Parameters:
444: +  svd     - singular value solver context obtained from SVDCreate()
445: .  func    - pointer to the stopping test function
446: .  ctx     - context for private data for the stopping routine (may be null)
447: -  destroy - a routine for destroying the context (may be null)

449:    Calling Sequence of func:
450: $   func(SVD svd,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nsv,SVDConvergedReason *reason,void *ctx)

452: +   svd    - singular value solver context obtained from SVDCreate()
453: .   its    - current number of iterations
454: .   max_it - maximum number of iterations
455: .   nconv  - number of currently converged singular triplets
456: .   nsv    - number of requested singular triplets
457: .   reason - (output) result of the stopping test
458: -   ctx    - optional context, as set by SVDSetStoppingTestFunction()

460:    Note:
461:    Normal usage is to first call the default routine SVDStoppingBasic() and then
462:    set reason to SVD_CONVERGED_USER if some user-defined conditions have been
463:    met. To let the singular value solver continue iterating, the result must be
464:    left as SVD_CONVERGED_ITERATING.

466:    Level: advanced

468: .seealso: SVDSetStoppingTest(), SVDStoppingBasic()
469: @*/
470: PetscErrorCode SVDSetStoppingTestFunction(SVD svd,PetscErrorCode (*func)(SVD,PetscInt,PetscInt,PetscInt,PetscInt,SVDConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
471: {

476:   if (svd->stoppingdestroy) {
477:     (*svd->stoppingdestroy)(svd->stoppingctx);
478:   }
479:   svd->stoppinguser    = func;
480:   svd->stoppingdestroy = destroy;
481:   svd->stoppingctx     = ctx;
482:   if (func == SVDStoppingBasic) svd->stop = SVD_STOP_BASIC;
483:   else {
484:     svd->stop     = SVD_STOP_USER;
485:     svd->stopping = svd->stoppinguser;
486:   }
487:   return(0);
488: }

490: /*@
491:    SVDSetStoppingTest - Specifies how to decide the termination of the outer
492:    loop of the singular value solver.

494:    Logically Collective on SVD

496:    Input Parameters:
497: +  svd  - singular value solver context obtained from SVDCreate()
498: -  stop - the type of stopping test

500:    Options Database Keys:
501: +  -svd_stop_basic - Sets the default stopping test
502: -  -svd_stop_user  - Selects the user-defined stopping test

504:    Note:
505:    The parameter 'stop' can have one of these values
506: +     SVD_STOP_BASIC - default stopping test
507: -     SVD_STOP_USER  - function set by SVDSetStoppingTestFunction()

509:    Level: advanced

511: .seealso: SVDGetStoppingTest(), SVDSetStoppingTestFunction(), SVDSetConvergenceTest(), SVDStop
512: @*/
513: PetscErrorCode SVDSetStoppingTest(SVD svd,SVDStop stop)
514: {
518:   switch (stop) {
519:     case SVD_STOP_BASIC: svd->stopping = SVDStoppingBasic; break;
520:     case SVD_STOP_USER:
521:       if (!svd->stoppinguser) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ORDER,"Must call SVDSetStoppingTestFunction() first");
522:       svd->stopping = svd->stoppinguser;
523:       break;
524:     default:
525:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
526:   }
527:   svd->stop = stop;
528:   return(0);
529: }

531: /*@
532:    SVDGetStoppingTest - Gets the method used to decide the termination of the outer
533:    loop of the singular value solver.

535:    Not Collective

537:    Input Parameters:
538: .  svd   - singular value solver context obtained from SVDCreate()

540:    Output Parameters:
541: .  stop  - the type of stopping test

543:    Level: advanced

545: .seealso: SVDSetStoppingTest(), SVDStop
546: @*/
547: PetscErrorCode SVDGetStoppingTest(SVD svd,SVDStop *stop)
548: {
552:   *stop = svd->stop;
553:   return(0);
554: }

556: /*@C
557:    SVDMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
558:    indicated by the user.

560:    Collective on SVD

562:    Input Parameters:
563: +  svd      - the singular value solver context
564: .  name     - the monitor option name
565: .  help     - message indicating what monitoring is done
566: .  manual   - manual page for the monitor
567: .  monitor  - the monitor function, whose context is a PetscViewerAndFormat
568: -  trackall - whether this monitor tracks all singular values or not

570:    Level: developer

572: .seealso: SVDMonitorSet(), SVDSetTrackAll(), SVDConvMonitorSetFromOptions()
573: @*/
574: PetscErrorCode SVDMonitorSetFromOptions(SVD svd,const char name[],const char help[],const char manual[],PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,PetscViewerAndFormat*),PetscBool trackall)
575: {
576:   PetscErrorCode       ierr;
577:   PetscBool            flg;
578:   PetscViewer          viewer;
579:   PetscViewerFormat    format;
580:   PetscViewerAndFormat *vf;

583:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,name,&viewer,&format,&flg);
584:   if (flg) {
585:     PetscViewerAndFormatCreate(viewer,format,&vf);
586:     PetscObjectDereference((PetscObject)viewer);
587:     SVDMonitorSet(svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))monitor,vf,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
588:     if (trackall) {
589:       SVDSetTrackAll(svd,PETSC_TRUE);
590:     }
591:   }
592:   return(0);
593: }

595: /*@C
596:    SVDConvMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
597:    indicated by the user (for monitors that only show iteration numbers of convergence).

599:    Collective on SVD

601:    Input Parameters:
602: +  svd      - the singular value solver context
603: .  name     - the monitor option name
604: .  help     - message indicating what monitoring is done
605: .  manual   - manual page for the monitor
606: -  monitor  - the monitor function, whose context is a SlepcConvMonitor

608:    Level: developer

610: .seealso: SVDMonitorSet(), SVDMonitorSetFromOptions()
611: @*/
612: PetscErrorCode SVDConvMonitorSetFromOptions(SVD svd,const char name[],const char help[],const char manual[],PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,SlepcConvMonitor))
613: {
614:   PetscErrorCode    ierr;
615:   PetscBool         flg;
616:   PetscViewer       viewer;
617:   PetscViewerFormat format;
618:   SlepcConvMonitor  ctx;

621:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,name,&viewer,&format,&flg);
622:   if (flg) {
623:     SlepcConvMonitorCreate(viewer,format,&ctx);
624:     PetscObjectDereference((PetscObject)viewer);
625:     SVDMonitorSet(svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))monitor,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
626:   }
627:   return(0);
628: }

630: /*@
631:    SVDSetFromOptions - Sets SVD options from the options database.
632:    This routine must be called before SVDSetUp() if the user is to be
633:    allowed to set the solver type.

635:    Collective on SVD

637:    Input Parameters:
638: .  svd - the singular value solver context

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

643:    Level: beginner

645: .seealso:
646: @*/
647: PetscErrorCode SVDSetFromOptions(SVD svd)
648: {
650:   char           type[256];
651:   PetscBool      set,flg,val,flg1,flg2,flg3;
652:   PetscInt       i,j,k;
653:   PetscReal      r;
654:   PetscDrawLG    lg;

658:   SVDRegisterAll();
659:   PetscObjectOptionsBegin((PetscObject)svd);
660:     PetscOptionsFList("-svd_type","SVD solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,256,&flg);
661:     if (flg) {
662:       SVDSetType(svd,type);
663:     } else if (!((PetscObject)svd)->type_name) {
664:       SVDSetType(svd,SVDCROSS);
665:     }

667:     PetscOptionsBool("-svd_implicittranspose","Handle matrix transpose implicitly","SVDSetImplicitTranspose",svd->impltrans,&val,&flg);
668:     if (flg) { SVDSetImplicitTranspose(svd,val); }

670:     i = svd->max_it? svd->max_it: PETSC_DEFAULT;
671:     PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,&flg1);
672:     r = svd->tol;
673:     PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",svd->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:svd->tol,&r,&flg2);
674:     if (flg1 || flg2) { SVDSetTolerances(svd,r,i); }

676:     PetscOptionsBoolGroupBegin("-svd_conv_rel","Relative error convergence test","SVDSetConvergenceTest",&flg);
677:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_REL); }
678:     PetscOptionsBoolGroup("-svd_conv_abs","Absolute error convergence test","SVDSetConvergenceTest",&flg);
679:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_ABS); }
680:     PetscOptionsBoolGroupEnd("-svd_conv_user","User-defined convergence test","SVDSetConvergenceTest",&flg);
681:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_USER); }

683:     PetscOptionsBoolGroupBegin("-svd_stop_basic","Stop iteration if all singular values converged or max_it reached","SVDSetStoppingTest",&flg);
684:     if (flg) { SVDSetStoppingTest(svd,SVD_STOP_BASIC); }
685:     PetscOptionsBoolGroupEnd("-svd_stop_user","User-defined stopping test","SVDSetStoppingTest",&flg);
686:     if (flg) { SVDSetStoppingTest(svd,SVD_STOP_USER); }

688:     i = svd->nsv;
689:     PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,&flg1);
690:     j = svd->ncv? svd->ncv: PETSC_DEFAULT;
691:     PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,&flg2);
692:     k = svd->mpd? svd->mpd: PETSC_DEFAULT;
693:     PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,&flg3);
694:     if (flg1 || flg2 || flg3) { SVDSetDimensions(svd,i,j,k); }

696:     PetscOptionsBoolGroupBegin("-svd_largest","Compute largest singular values","SVDSetWhichSingularTriplets",&flg);
697:     if (flg) { SVDSetWhichSingularTriplets(svd,SVD_LARGEST); }
698:     PetscOptionsBoolGroupEnd("-svd_smallest","Compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
699:     if (flg) { SVDSetWhichSingularTriplets(svd,SVD_SMALLEST); }

701:     /* -----------------------------------------------------------------------*/
702:     /*
703:       Cancels all monitors hardwired into code before call to SVDSetFromOptions()
704:     */
705:     PetscOptionsBool("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",PETSC_FALSE,&flg,&set);
706:     if (set && flg) {
707:       SVDMonitorCancel(svd);
708:     }
709:     /*
710:       Text monitors
711:     */
712:     SVDMonitorSetFromOptions(svd,"-svd_monitor","Monitor first unconverged approximate singular value and error estimate","SVDMonitorFirst",SVDMonitorFirst,PETSC_FALSE);
713:     SVDConvMonitorSetFromOptions(svd,"-svd_monitor_conv","Monitor approximate singular values and error estimates as they converge","SVDMonitorConverged",SVDMonitorConverged);
714:     SVDMonitorSetFromOptions(svd,"-svd_monitor_all","Monitor approximate singular values and error estimates","SVDMonitorAll",SVDMonitorAll,PETSC_TRUE);
715:     /*
716:       Line graph monitors
717:     */
718:     PetscOptionsBool("-svd_monitor_lg","Monitor first unconverged approximate singular value and error estimate graphically","SVDMonitorSet",PETSC_FALSE,&flg,&set);
719:     if (set && flg) {
720:       SVDMonitorLGCreate(PetscObjectComm((PetscObject)svd),NULL,"Error estimates",PETSC_DECIDE,PETSC_DECIDE,300,300,&lg);
721:       SVDMonitorSet(svd,SVDMonitorLG,lg,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
722:     }
723:     PetscOptionsBool("-svd_monitor_lg_all","Monitor error estimates graphically","SVDMonitorSet",PETSC_FALSE,&flg,&set);
724:     if (set && flg) {
725:       SVDMonitorLGCreate(PetscObjectComm((PetscObject)svd),NULL,"Error estimates",PETSC_DECIDE,PETSC_DECIDE,300,300,&lg);
726:       SVDMonitorSet(svd,SVDMonitorLGAll,lg,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
727:       SVDSetTrackAll(svd,PETSC_TRUE);
728:     }

730:     /* -----------------------------------------------------------------------*/
731:     PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",NULL);
732:     PetscOptionsName("-svd_view_vectors","View computed singular vectors","SVDVectorsView",NULL);
733:     PetscOptionsName("-svd_view_values","View computed singular values","SVDValuesView",NULL);
734:     PetscOptionsName("-svd_converged_reason","Print reason for convergence, and number of iterations","SVDReasonView",NULL);
735:     PetscOptionsName("-svd_error_absolute","Print absolute errors of each singular triplet","SVDErrorView",NULL);
736:     PetscOptionsName("-svd_error_relative","Print relative errors of each singular triplet","SVDErrorView",NULL);

738:     if (svd->ops->setfromoptions) {
739:       (*svd->ops->setfromoptions)(PetscOptionsObject,svd);
740:     }
741:     PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)svd);
742:   PetscOptionsEnd();

744:   if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
745:   BVSetFromOptions(svd->V);
746:   BVSetFromOptions(svd->U);
747:   if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
748:   DSSetFromOptions(svd->ds);
749:   return(0);
750: }

752: /*@
753:    SVDSetTrackAll - Specifies if the solver must compute the residual norm of all
754:    approximate singular value or not.

756:    Logically Collective on SVD

758:    Input Parameters:
759: +  svd      - the singular value solver context
760: -  trackall - whether to compute all residuals or not

762:    Notes:
763:    If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
764:    the residual norm for each singular value approximation. Computing the residual is
765:    usually an expensive operation and solvers commonly compute only the residual
766:    associated to the first unconverged singular value.

768:    The options '-svd_monitor_all' and '-svd_monitor_lg_all' automatically
769:    activate this option.

771:    Level: developer

773: .seealso: SVDGetTrackAll()
774: @*/
775: PetscErrorCode SVDSetTrackAll(SVD svd,PetscBool trackall)
776: {
780:   svd->trackall = trackall;
781:   return(0);
782: }

784: /*@
785:    SVDGetTrackAll - Returns the flag indicating whether all residual norms must
786:    be computed or not.

788:    Not Collective

790:    Input Parameter:
791: .  svd - the singular value solver context

793:    Output Parameter:
794: .  trackall - the returned flag

796:    Level: developer

798: .seealso: SVDSetTrackAll()
799: @*/
800: PetscErrorCode SVDGetTrackAll(SVD svd,PetscBool *trackall)
801: {
805:   *trackall = svd->trackall;
806:   return(0);
807: }


810: /*@C
811:    SVDSetOptionsPrefix - Sets the prefix used for searching for all
812:    SVD options in the database.

814:    Logically Collective on SVD

816:    Input Parameters:
817: +  svd - the singular value solver context
818: -  prefix - the prefix string to prepend to all SVD option requests

820:    Notes:
821:    A hyphen (-) must NOT be given at the beginning of the prefix name.
822:    The first character of all runtime options is AUTOMATICALLY the
823:    hyphen.

825:    For example, to distinguish between the runtime options for two
826:    different SVD contexts, one could call
827: .vb
828:       SVDSetOptionsPrefix(svd1,"svd1_")
829:       SVDSetOptionsPrefix(svd2,"svd2_")
830: .ve

832:    Level: advanced

834: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
835: @*/
836: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)
837: {

842:   if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
843:   BVSetOptionsPrefix(svd->V,prefix);
844:   BVSetOptionsPrefix(svd->U,prefix);
845:   if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
846:   DSSetOptionsPrefix(svd->ds,prefix);
847:   PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
848:   return(0);
849: }

851: /*@C
852:    SVDAppendOptionsPrefix - Appends to the prefix used for searching for all
853:    SVD options in the database.

855:    Logically Collective on SVD

857:    Input Parameters:
858: +  svd - the singular value solver context
859: -  prefix - the prefix string to prepend to all SVD option requests

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

865:    Level: advanced

867: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
868: @*/
869: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
870: {

875:   if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
876:   BVAppendOptionsPrefix(svd->V,prefix);
877:   BVAppendOptionsPrefix(svd->U,prefix);
878:   if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
879:   DSAppendOptionsPrefix(svd->ds,prefix);
880:   PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
881:   return(0);
882: }

884: /*@C
885:    SVDGetOptionsPrefix - Gets the prefix used for searching for all
886:    SVD options in the database.

888:    Not Collective

890:    Input Parameters:
891: .  svd - the singular value solver context

893:    Output Parameters:
894: .  prefix - pointer to the prefix string used is returned

896:    Note:
897:    On the Fortran side, the user should pass in a string 'prefix' of
898:    sufficient length to hold the prefix.

900:    Level: advanced

902: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
903: @*/
904: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
905: {

911:   PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
912:   return(0);
913: }