Actual source code: svdopts.c

slepc-3.7.0 2016-05-16
Report Typos and Errors
  1: /*
  2:      SVD routines for setting solver options.

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2016, 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*/
 25: #include <petscdraw.h>

 29: /*@
 30:    SVDSetImplicitTranspose - Indicates how to handle the transpose of the matrix
 31:    associated with the singular value problem.

 33:    Logically Collective on SVD

 35:    Input Parameters:
 36: +  svd  - the singular value solver context
 37: -  impl - how to handle the transpose (implicitly or not)

 39:    Options Database Key:
 40: .  -svd_implicittranspose - Activate the implicit transpose mode.

 42:    Notes:
 43:    By default, the transpose of the matrix is explicitly built (if the matrix
 44:    has defined the MatTranspose operation).

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

 52:    Level: advanced

 54: .seealso: SVDGetImplicitTranspose(), SVDSolve(), SVDSetOperator()
 55: @*/
 56: PetscErrorCode SVDSetImplicitTranspose(SVD svd,PetscBool impl)
 57: {
 61:   if (svd->impltrans!=impl) {
 62:     svd->impltrans = impl;
 63:     svd->state     = SVD_STATE_INITIAL;
 64:   }
 65:   return(0);
 66: }

 70: /*@
 71:    SVDGetImplicitTranspose - Gets the mode used to handle the transpose
 72:    of the matrix associated with the singular value problem.

 74:    Not Collective

 76:    Input Parameter:
 77: .  svd  - the singular value solver context

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

 82:    Level: advanced

 84: .seealso: SVDSetImplicitTranspose(), SVDSolve(), SVDSetOperator()
 85: @*/
 86: PetscErrorCode SVDGetImplicitTranspose(SVD svd,PetscBool *impl)
 87: {
 91:   *impl = svd->impltrans;
 92:   return(0);
 93: }

 97: /*@
 98:    SVDSetTolerances - Sets the tolerance and maximum
 99:    iteration count used by the default SVD convergence testers.

101:    Logically Collective on SVD

103:    Input Parameters:
104: +  svd - the singular value solver context
105: .  tol - the convergence tolerance
106: -  maxits - maximum number of iterations to use

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

113:    Note:
114:    Use PETSC_DEFAULT for either argument to assign a reasonably good value.

116:    Level: intermediate

118: .seealso: SVDGetTolerances()
119: @*/
120: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,PetscInt maxits)
121: {
126:   if (tol == PETSC_DEFAULT) {
127:     svd->tol   = PETSC_DEFAULT;
128:     svd->state = SVD_STATE_INITIAL;
129:   } else {
130:     if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
131:     svd->tol = tol;
132:   }
133:   if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
134:     svd->max_it = 0;
135:     svd->state  = SVD_STATE_INITIAL;
136:   } else {
137:     if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
138:     svd->max_it = maxits;
139:   }
140:   return(0);
141: }

145: /*@
146:    SVDGetTolerances - Gets the tolerance and maximum
147:    iteration count used by the default SVD convergence tests.

149:    Not Collective

151:    Input Parameter:
152: .  svd - the singular value solver context

154:    Output Parameters:
155: +  tol - the convergence tolerance
156: -  maxits - maximum number of iterations

158:    Notes:
159:    The user can specify NULL for any parameter that is not needed.

161:    Level: intermediate

163: .seealso: SVDSetTolerances()
164: @*/
165: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,PetscInt *maxits)
166: {
169:   if (tol)    *tol    = svd->tol;
170:   if (maxits) *maxits = svd->max_it;
171:   return(0);
172: }

176: /*@
177:    SVDSetDimensions - Sets the number of singular values to compute
178:    and the dimension of the subspace.

180:    Logically Collective on SVD

182:    Input Parameters:
183: +  svd - the singular value solver context
184: .  nsv - number of singular values to compute
185: .  ncv - the maximum dimension of the subspace to be used by the solver
186: -  mpd - the maximum dimension allowed for the projected problem

188:    Options Database Keys:
189: +  -svd_nsv <nsv> - Sets the number of singular values
190: .  -svd_ncv <ncv> - Sets the dimension of the subspace
191: -  -svd_mpd <mpd> - Sets the maximum projected dimension

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

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

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

206:    Level: intermediate

208: .seealso: SVDGetDimensions()
209: @*/
210: PetscErrorCode SVDSetDimensions(SVD svd,PetscInt nsv,PetscInt ncv,PetscInt mpd)
211: {
217:   if (nsv<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
218:   svd->nsv = nsv;
219:   if (ncv == PETSC_DEFAULT || ncv == PETSC_DECIDE) {
220:     svd->ncv = 0;
221:   } else {
222:     if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
223:     svd->ncv = ncv;
224:   }
225:   if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
226:     svd->mpd = 0;
227:   } else {
228:     if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
229:     svd->mpd = mpd;
230:   }
231:   svd->state = SVD_STATE_INITIAL;
232:   return(0);
233: }

237: /*@
238:    SVDGetDimensions - Gets the number of singular values to compute
239:    and the dimension of the subspace.

241:    Not Collective

243:    Input Parameter:
244: .  svd - the singular value context

246:    Output Parameters:
247: +  nsv - number of singular values to compute
248: .  ncv - the maximum dimension of the subspace to be used by the solver
249: -  mpd - the maximum dimension allowed for the projected problem

251:    Notes:
252:    The user can specify NULL for any parameter that is not needed.

254:    Level: intermediate

256: .seealso: SVDSetDimensions()
257: @*/
258: PetscErrorCode SVDGetDimensions(SVD svd,PetscInt *nsv,PetscInt *ncv,PetscInt *mpd)
259: {
262:   if (nsv) *nsv = svd->nsv;
263:   if (ncv) *ncv = svd->ncv;
264:   if (mpd) *mpd = svd->mpd;
265:   return(0);
266: }

270: /*@
271:     SVDSetWhichSingularTriplets - Specifies which singular triplets are
272:     to be sought.

274:     Logically Collective on SVD

276:     Input Parameter:
277: .   svd - singular value solver context obtained from SVDCreate()

279:     Output Parameter:
280: .   which - which singular triplets are to be sought

282:     Possible values:
283:     The parameter 'which' can have one of these values:

285: +     SVD_LARGEST  - largest singular values
286: -     SVD_SMALLEST - smallest singular values

288:     Options Database Keys:
289: +   -svd_largest  - Sets largest singular values
290: -   -svd_smallest - Sets smallest singular values

292:     Level: intermediate

294: .seealso: SVDGetWhichSingularTriplets(), SVDWhich
295: @*/
296: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)
297: {
301:   switch (which) {
302:     case SVD_LARGEST:
303:     case SVD_SMALLEST:
304:       if (svd->which != which) {
305:         svd->state = SVD_STATE_INITIAL;
306:         svd->which = which;
307:       }
308:       break;
309:   default:
310:     SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
311:   }
312:   return(0);
313: }

317: /*@
318:     SVDGetWhichSingularTriplets - Returns which singular triplets are
319:     to be sought.

321:     Not Collective

323:     Input Parameter:
324: .   svd - singular value solver context obtained from SVDCreate()

326:     Output Parameter:
327: .   which - which singular triplets are to be sought

329:     Notes:
330:     See SVDSetWhichSingularTriplets() for possible values of which

332:     Level: intermediate

334: .seealso: SVDSetWhichSingularTriplets(), SVDWhich
335: @*/
336: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)
337: {
341:   *which = svd->which;
342:   return(0);
343: }

347: /*@C
348:    SVDSetConvergenceTestFunction - Sets a function to compute the error estimate
349:    used in the convergence test.

351:    Logically Collective on SVD

353:    Input Parameters:
354: +  svd     - singular value solver context obtained from SVDCreate()
355: .  func    - a pointer to the convergence test function
356: .  ctx     - context for private data for the convergence routine (may be null)
357: -  destroy - a routine for destroying the context (may be null)

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

362: +   svd    - singular value solver context obtained from SVDCreate()
363: .   sigma  - computed singular value
364: .   res    - residual norm associated to the singular triplet
365: .   errest - (output) computed error estimate
366: -   ctx    - optional context, as set by SVDSetConvergenceTestFunction()

368:    Note:
369:    If the error estimate returned by the convergence test function is less than
370:    the tolerance, then the singular value is accepted as converged.

372:    Level: advanced

374: .seealso: SVDSetConvergenceTest(), SVDSetTolerances()
375: @*/
376: PetscErrorCode SVDSetConvergenceTestFunction(SVD svd,PetscErrorCode (*func)(SVD,PetscReal,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
377: {

382:   if (svd->convergeddestroy) {
383:     (*svd->convergeddestroy)(svd->convergedctx);
384:   }
385:   svd->converged        = func;
386:   svd->convergeddestroy = destroy;
387:   svd->convergedctx     = ctx;
388:   if (func == SVDConvergedRelative) svd->conv = SVD_CONV_REL;
389:   else if (func == SVDConvergedAbsolute) svd->conv = SVD_CONV_ABS;
390:   else svd->conv = SVD_CONV_USER;
391:   return(0);
392: }

396: /*@
397:    SVDSetConvergenceTest - Specifies how to compute the error estimate
398:    used in the convergence test.

400:    Logically Collective on SVD

402:    Input Parameters:
403: +  svd  - singular value solver context obtained from SVDCreate()
404: -  conv - the type of convergence test

406:    Options Database Keys:
407: +  -svd_conv_abs  - Sets the absolute convergence test
408: .  -svd_conv_rel  - Sets the convergence test relative to the singular value
409: -  -svd_conv_user - Selects the user-defined convergence test

411:    Note:
412:    The parameter 'conv' can have one of these values
413: +     SVD_CONV_ABS  - absolute error ||r||
414: .     SVD_CONV_REL  - error relative to the singular value l, ||r||/sigma
415: -     SVD_CONV_USER - function set by SVDSetConvergenceTestFunction()

417:    Level: intermediate

419: .seealso: SVDGetConvergenceTest(), SVDSetConvergenceTestFunction(), SVDSetStoppingTest(), SVDConv
420: @*/
421: PetscErrorCode SVDSetConvergenceTest(SVD svd,SVDConv conv)
422: {
426:   switch (conv) {
427:     case SVD_CONV_ABS:  svd->converged = SVDConvergedAbsolute; break;
428:     case SVD_CONV_REL:  svd->converged = SVDConvergedRelative; break;
429:     case SVD_CONV_USER: break;
430:     default:
431:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
432:   }
433:   svd->conv = conv;
434:   return(0);
435: }

439: /*@
440:    SVDGetConvergenceTest - Gets the method used to compute the error estimate
441:    used in the convergence test.

443:    Not Collective

445:    Input Parameters:
446: .  svd   - singular value solver context obtained from SVDCreate()

448:    Output Parameters:
449: .  conv  - the type of convergence test

451:    Level: intermediate

453: .seealso: SVDSetConvergenceTest(), SVDConv
454: @*/
455: PetscErrorCode SVDGetConvergenceTest(SVD svd,SVDConv *conv)
456: {
460:   *conv = svd->conv;
461:   return(0);
462: }

466: /*@C
467:    SVDSetStoppingTestFunction - Sets a function to decide when to stop the outer
468:    iteration of the singular value solver.

470:    Logically Collective on SVD

472:    Input Parameters:
473: +  svd     - singular value solver context obtained from SVDCreate()
474: .  func    - pointer to the stopping test function
475: .  ctx     - context for private data for the stopping routine (may be null)
476: -  destroy - a routine for destroying the context (may be null)

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

481: +   svd    - singular value solver context obtained from SVDCreate()
482: .   its    - current number of iterations
483: .   max_it - maximum number of iterations
484: .   nconv  - number of currently converged singular triplets
485: .   nsv    - number of requested singular triplets
486: .   reason - (output) result of the stopping test
487: -   ctx    - optional context, as set by SVDSetStoppingTestFunction()

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

495:    Level: advanced

497: .seealso: SVDSetStoppingTest(), SVDStoppingBasic()
498: @*/
499: PetscErrorCode SVDSetStoppingTestFunction(SVD svd,PetscErrorCode (*func)(SVD,PetscInt,PetscInt,PetscInt,PetscInt,SVDConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
500: {

505:   if (svd->stoppingdestroy) {
506:     (*svd->stoppingdestroy)(svd->stoppingctx);
507:   }
508:   svd->stopping        = func;
509:   svd->stoppingdestroy = destroy;
510:   svd->stoppingctx     = ctx;
511:   if (func == SVDStoppingBasic) svd->stop = SVD_STOP_BASIC;
512:   else svd->stop = SVD_STOP_USER;
513:   return(0);
514: }

518: /*@
519:    SVDSetStoppingTest - Specifies how to decide the termination of the outer
520:    loop of the singular value solver.

522:    Logically Collective on SVD

524:    Input Parameters:
525: +  svd  - singular value solver context obtained from SVDCreate()
526: -  stop - the type of stopping test

528:    Options Database Keys:
529: +  -svd_stop_basic - Sets the default stopping test
530: -  -svd_stop_user  - Selects the user-defined stopping test

532:    Note:
533:    The parameter 'stop' can have one of these values
534: +     SVD_STOP_BASIC - default stopping test
535: -     SVD_STOP_USER  - function set by SVDSetStoppingTestFunction()

537:    Level: advanced

539: .seealso: SVDGetStoppingTest(), SVDSetStoppingTestFunction(), SVDSetConvergenceTest(), SVDStop
540: @*/
541: PetscErrorCode SVDSetStoppingTest(SVD svd,SVDStop stop)
542: {
546:   switch (stop) {
547:     case SVD_STOP_BASIC: svd->stopping = SVDStoppingBasic; break;
548:     case SVD_STOP_USER:  break;
549:     default:
550:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
551:   }
552:   svd->stop = stop;
553:   return(0);
554: }

558: /*@
559:    SVDGetStoppingTest - Gets the method used to decide the termination of the outer
560:    loop of the singular value solver.

562:    Not Collective

564:    Input Parameters:
565: .  svd   - singular value solver context obtained from SVDCreate()

567:    Output Parameters:
568: .  stop  - the type of stopping test

570:    Level: advanced

572: .seealso: SVDSetStoppingTest(), SVDStop
573: @*/
574: PetscErrorCode SVDGetStoppingTest(SVD svd,SVDStop *stop)
575: {
579:   *stop = svd->stop;
580:   return(0);
581: }

585: /*@C
586:    SVDMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
587:    indicated by the user.

589:    Collective on SVD

591:    Input Parameters:
592: +  svd      - the singular value solver context
593: .  name     - the monitor option name
594: .  help     - message indicating what monitoring is done
595: .  manual   - manual page for the monitor
596: .  monitor  - the monitor function, whose context is a PetscViewerAndFormat
597: -  trackall - whether this monitor tracks all singular values or not

599:    Level: developer

601: .seealso: SVDMonitorSet(), SVDSetTrackAll(), SVDConvMonitorSetFromOptions()
602: @*/
603: PetscErrorCode SVDMonitorSetFromOptions(SVD svd,const char name[],const char help[],const char manual[],PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,PetscViewerAndFormat*),PetscBool trackall)
604: {
605:   PetscErrorCode       ierr;
606:   PetscBool            flg;
607:   PetscViewer          viewer;
608:   PetscViewerFormat    format;
609:   PetscViewerAndFormat *vf;

612:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->prefix,name,&viewer,&format,&flg);
613:   if (flg) {
614:     PetscViewerAndFormatCreate(viewer,format,&vf);
615:     PetscObjectDereference((PetscObject)viewer);
616:     SVDMonitorSet(svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))monitor,vf,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
617:     if (trackall) {
618:       SVDSetTrackAll(svd,PETSC_TRUE);
619:     }
620:   }
621:   return(0);
622: }

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

630:    Collective on SVD

632:    Input Parameters:
633: +  svd      - the singular value solver context
634: .  name     - the monitor option name
635: .  help     - message indicating what monitoring is done
636: .  manual   - manual page for the monitor
637: -  monitor  - the monitor function, whose context is a SlepcConvMonitor

639:    Level: developer

641: .seealso: SVDMonitorSet(), SVDMonitorSetFromOptions()
642: @*/
643: PetscErrorCode SVDConvMonitorSetFromOptions(SVD svd,const char name[],const char help[],const char manual[],PetscErrorCode (*monitor)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,SlepcConvMonitor))
644: {
645:   PetscErrorCode    ierr;
646:   PetscBool         flg;
647:   PetscViewer       viewer;
648:   PetscViewerFormat format;
649:   SlepcConvMonitor  ctx;

652:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->prefix,name,&viewer,&format,&flg);
653:   if (flg) {
654:     SlepcConvMonitorCreate(viewer,format,&ctx);
655:     PetscObjectDereference((PetscObject)viewer);
656:     SVDMonitorSet(svd,(PetscErrorCode (*)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*))monitor,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
657:   }
658:   return(0);
659: }

663: /*@
664:    SVDSetFromOptions - Sets SVD options from the options database.
665:    This routine must be called before SVDSetUp() if the user is to be
666:    allowed to set the solver type.

668:    Collective on SVD

670:    Input Parameters:
671: .  svd - the singular value solver context

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

676:    Level: beginner

678: .seealso:
679: @*/
680: PetscErrorCode SVDSetFromOptions(SVD svd)
681: {
683:   char           type[256];
684:   PetscBool      set,flg,val,flg1,flg2,flg3;
685:   PetscInt       i,j,k;
686:   PetscReal      r;
687:   PetscDrawLG    lg;

691:   SVDRegisterAll();
692:   PetscObjectOptionsBegin((PetscObject)svd);
693:     PetscOptionsFList("-svd_type","Singular Value Solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,256,&flg);
694:     if (flg) {
695:       SVDSetType(svd,type);
696:     } else if (!((PetscObject)svd)->type_name) {
697:       SVDSetType(svd,SVDCROSS);
698:     }

700:     PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",NULL);
701:     PetscOptionsName("-svd_view_vectors","View computed singular vectors","SVDVectorsView",NULL);
702:     PetscOptionsName("-svd_view_values","View computed singular values","SVDValuesView",NULL);
703:     PetscOptionsName("-svd_converged_reason","Print reason for convergence, and number of iterations","SVDReasonView",NULL);
704:     PetscOptionsName("-svd_error_absolute","Print absolute errors of each singular triplet","SVDErrorView",NULL);
705:     PetscOptionsName("-svd_error_relative","Print relative errors of each singular triplet","SVDErrorView",NULL);

707:     PetscOptionsBool("-svd_implicittranspose","Handle matrix transpose implicitly","SVDSetImplicitTranspose",svd->impltrans,&val,&flg);
708:     if (flg) {
709:       SVDSetImplicitTranspose(svd,val);
710:     }

712:     i = svd->max_it? svd->max_it: PETSC_DEFAULT;
713:     PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,&flg1);
714:     r = svd->tol;
715:     PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",svd->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:svd->tol,&r,&flg2);
716:     if (flg1 || flg2) {
717:       SVDSetTolerances(svd,r,i);
718:     }

720:     PetscOptionsBoolGroupBegin("-svd_conv_rel","Relative error convergence test","SVDSetConvergenceTest",&flg);
721:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_REL); }
722:     PetscOptionsBoolGroup("-svd_conv_abs","Absolute error convergence test","SVDSetConvergenceTest",&flg);
723:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_ABS); }
724:     PetscOptionsBoolGroupEnd("-svd_conv_user","User-defined convergence test","SVDSetConvergenceTest",&flg);
725:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_USER); }

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

732:     i = svd->nsv;
733:     PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,&flg1);
734:     j = svd->ncv? svd->ncv: PETSC_DEFAULT;
735:     PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,&flg2);
736:     k = svd->mpd? svd->mpd: PETSC_DEFAULT;
737:     PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,&flg3);
738:     if (flg1 || flg2 || flg3) {
739:       SVDSetDimensions(svd,i,j,k);
740:     }

742:     PetscOptionsBoolGroupBegin("-svd_largest","compute largest singular values","SVDSetWhichSingularTriplets",&flg);
743:     if (flg) { SVDSetWhichSingularTriplets(svd,SVD_LARGEST); }
744:     PetscOptionsBoolGroupEnd("-svd_smallest","compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
745:     if (flg) { SVDSetWhichSingularTriplets(svd,SVD_SMALLEST); }

747:     /* -----------------------------------------------------------------------*/
748:     /*
749:       Cancels all monitors hardwired into code before call to SVDSetFromOptions()
750:     */
751:     PetscOptionsBool("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",PETSC_FALSE,&flg,&set);
752:     if (set && flg) {
753:       SVDMonitorCancel(svd);
754:     }
755:     /*
756:       Text monitors
757:     */
758:     SVDMonitorSetFromOptions(svd,"-svd_monitor","Monitor first unconverged approximate singular value and error estimate","SVDMonitorFirst",SVDMonitorFirst,PETSC_FALSE);
759:     SVDConvMonitorSetFromOptions(svd,"-svd_monitor_conv","Monitor approximate singular values and error estimates as they converge","SVDMonitorConverged",SVDMonitorConverged);
760:     SVDMonitorSetFromOptions(svd,"-svd_monitor_all","Monitor approximate singular values and error estimates","SVDMonitorAll",SVDMonitorAll,PETSC_TRUE);
761:     /*
762:       Line graph monitors
763:     */
764:     PetscOptionsBool("-svd_monitor_lg","Monitor first unconverged approximate singular value and error estimate graphically","SVDMonitorSet",PETSC_FALSE,&flg,&set);
765:     if (set && flg) {
766:       SVDMonitorLGCreate(PetscObjectComm((PetscObject)svd),NULL,"Error estimates",PETSC_DECIDE,PETSC_DECIDE,300,300,&lg);
767:       SVDMonitorSet(svd,SVDMonitorLG,lg,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
768:     }
769:     PetscOptionsBool("-svd_monitor_lg_all","Monitor error estimates graphically","SVDMonitorSet",PETSC_FALSE,&flg,&set);
770:     if (set && flg) {
771:       SVDMonitorLGCreate(PetscObjectComm((PetscObject)svd),NULL,"Error estimates",PETSC_DECIDE,PETSC_DECIDE,300,300,&lg);
772:       SVDMonitorSet(svd,SVDMonitorLGAll,lg,(PetscErrorCode (*)(void**))PetscDrawLGDestroy);
773:       SVDSetTrackAll(svd,PETSC_TRUE);
774:     }

776:     if (svd->ops->setfromoptions) {
777:       (*svd->ops->setfromoptions)(PetscOptionsObject,svd);
778:     }
779:     PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)svd);
780:   PetscOptionsEnd();

782:   if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
783:   BVSetFromOptions(svd->V);
784:   BVSetFromOptions(svd->U);
785:   if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
786:   DSSetFromOptions(svd->ds);
787:   return(0);
788: }

792: /*@
793:    SVDSetTrackAll - Specifies if the solver must compute the residual norm of all
794:    approximate singular value or not.

796:    Logically Collective on SVD

798:    Input Parameters:
799: +  svd      - the singular value solver context
800: -  trackall - whether to compute all residuals or not

802:    Notes:
803:    If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
804:    the residual norm for each singular value approximation. Computing the residual is
805:    usually an expensive operation and solvers commonly compute only the residual
806:    associated to the first unconverged singular value.

808:    The options '-svd_monitor_all' and '-svd_monitor_lg_all' automatically
809:    activate this option.

811:    Level: developer

813: .seealso: SVDGetTrackAll()
814: @*/
815: PetscErrorCode SVDSetTrackAll(SVD svd,PetscBool trackall)
816: {
820:   svd->trackall = trackall;
821:   return(0);
822: }

826: /*@
827:    SVDGetTrackAll - Returns the flag indicating whether all residual norms must
828:    be computed or not.

830:    Not Collective

832:    Input Parameter:
833: .  svd - the singular value solver context

835:    Output Parameter:
836: .  trackall - the returned flag

838:    Level: developer

840: .seealso: SVDSetTrackAll()
841: @*/
842: PetscErrorCode SVDGetTrackAll(SVD svd,PetscBool *trackall)
843: {
847:   *trackall = svd->trackall;
848:   return(0);
849: }


854: /*@C
855:    SVDSetOptionsPrefix - Sets the prefix used for searching for all
856:    SVD options in the database.

858:    Logically Collective on SVD

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

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

869:    For example, to distinguish between the runtime options for two
870:    different SVD contexts, one could call
871: .vb
872:       SVDSetOptionsPrefix(svd1,"svd1_")
873:       SVDSetOptionsPrefix(svd2,"svd2_")
874: .ve

876:    Level: advanced

878: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
879: @*/
880: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)
881: {
883:   PetscBool      flg1,flg2;
884:   EPS            eps;

888:   if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
889:   BVSetOptionsPrefix(svd->V,prefix);
890:   BVSetOptionsPrefix(svd->U,prefix);
891:   if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
892:   DSSetOptionsPrefix(svd->ds,prefix);
893:   PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
894:   PetscObjectTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
895:   PetscObjectTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
896:   if (flg1) {
897:     SVDCrossGetEPS(svd,&eps);
898:   } else if (flg2) {
899:       SVDCyclicGetEPS(svd,&eps);
900:   }
901:   if (flg1 || flg2) {
902:     EPSSetOptionsPrefix(eps,prefix);
903:     EPSAppendOptionsPrefix(eps,"svd_");
904:   }
905:   return(0);
906: }

910: /*@C
911:    SVDAppendOptionsPrefix - Appends to the prefix used for searching for all
912:    SVD options in the database.

914:    Logically Collective on SVD

916:    Input Parameters:
917: +  svd - the singular value solver context
918: -  prefix - the prefix string to prepend to all SVD option requests

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

924:    Level: advanced

926: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
927: @*/
928: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
929: {
931:   PetscBool      flg1,flg2;
932:   EPS            eps;

936:   if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
937:   BVSetOptionsPrefix(svd->V,prefix);
938:   BVSetOptionsPrefix(svd->U,prefix);
939:   if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
940:   DSSetOptionsPrefix(svd->ds,prefix);
941:   PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
942:   PetscObjectTypeCompare((PetscObject)svd,SVDCROSS,&flg1);
943:   PetscObjectTypeCompare((PetscObject)svd,SVDCYCLIC,&flg2);
944:   if (flg1) {
945:     SVDCrossGetEPS(svd,&eps);
946:   } else if (flg2) {
947:     SVDCyclicGetEPS(svd,&eps);
948:   }
949:   if (flg1 || flg2) {
950:     EPSSetOptionsPrefix(eps,((PetscObject)svd)->prefix);
951:     EPSAppendOptionsPrefix(eps,"svd_");
952:   }
953:   return(0);
954: }

958: /*@C
959:    SVDGetOptionsPrefix - Gets the prefix used for searching for all
960:    SVD options in the database.

962:    Not Collective

964:    Input Parameters:
965: .  svd - the singular value solver context

967:    Output Parameters:
968: .  prefix - pointer to the prefix string used is returned

970:    Note:
971:    On the Fortran side, the user should pass in a string 'prefix' of
972:    sufficient length to hold the prefix.

974:    Level: advanced

976: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
977: @*/
978: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
979: {

985:   PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
986:   return(0);
987: }