Actual source code: svdopts.c

slepc-3.10.1 2018-10-23
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2018, 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
 95:    (use PETSC_DECIDE to compute an educated guess based on basis and matrix sizes)

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

100:    Level: intermediate

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

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

131:    Not Collective

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

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

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

143:    Level: intermediate

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

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

160:    Logically Collective on SVD

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

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

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

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

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

186:    Level: intermediate

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

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

219:    Not Collective

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

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

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

232:    Level: intermediate

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

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

250:     Logically Collective on SVD

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

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

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

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

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

268:     Level: intermediate

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

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

295:     Not Collective

297:     Input Parameter:
298: .   svd - singular value solver context obtained from SVDCreate()

300:     Output Parameter:
301: .   which - which singular triplets are to be sought

303:     Notes:
304:     See SVDSetWhichSingularTriplets() for possible values of which

306:     Level: intermediate

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

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

323:    Logically Collective on SVD

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

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

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

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

344:    Level: advanced

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

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

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

373:    Logically Collective on SVD

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

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

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

390:    Level: intermediate

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

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

417:    Not Collective

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

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

425:    Level: intermediate

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

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

442:    Logically Collective on SVD

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

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

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

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

467:    Level: advanced

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

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

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

495:    Logically Collective on SVD

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

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

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

510:    Level: advanced

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

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

536:    Not Collective

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

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

544:    Level: advanced

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

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

561:    Collective on SVD

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

571:    Level: developer

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

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

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

600:    Collective on SVD

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

609:    Level: developer

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

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

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

636:    Collective on SVD

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

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

644:    Level: beginner

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

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

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

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

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

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

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

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

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

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

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

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

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

757:    Logically Collective on SVD

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

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

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

772:    Level: developer

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

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

789:    Not Collective

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

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

797:    Level: developer

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


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

815:    Logically Collective on SVD

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

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

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

833:    Level: advanced

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

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

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

856:    Logically Collective on SVD

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

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

866:    Level: advanced

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

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

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

889:    Not Collective

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

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

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

901:    Level: advanced

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

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