Actual source code: svdopts.c

slepc-3.18.2 2023-01-26
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, 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>
 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(), SVDSetOperators()
 43: @*/
 44: PetscErrorCode SVDSetImplicitTranspose(SVD svd,PetscBool impl)
 45: {
 48:   if (svd->impltrans!=impl) {
 49:     svd->impltrans = impl;
 50:     svd->state     = SVD_STATE_INITIAL;
 51:   }
 52:   return 0;
 53: }

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

 59:    Not Collective

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

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

 67:    Level: advanced

 69: .seealso: SVDSetImplicitTranspose(), SVDSolve(), SVDSetOperators()
 70: @*/
 71: PetscErrorCode SVDGetImplicitTranspose(SVD svd,PetscBool *impl)
 72: {
 75:   *impl = svd->impltrans;
 76:   return 0;
 77: }

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

 83:    Logically Collective on svd

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

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

 94:    Note:
 95:    Use PETSC_DEFAULT for either argument to assign a reasonably good value.

 97:    Level: intermediate

 99: .seealso: SVDGetTolerances()
100: @*/
101: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,PetscInt maxits)
102: {
106:   if (tol == PETSC_DEFAULT) {
107:     svd->tol   = PETSC_DEFAULT;
108:     svd->state = SVD_STATE_INITIAL;
109:   } else {
111:     svd->tol = tol;
112:   }
113:   if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
114:     svd->max_it = PETSC_DEFAULT;
115:     svd->state  = SVD_STATE_INITIAL;
116:   } else {
118:     svd->max_it = maxits;
119:   }
120:   return 0;
121: }

123: /*@C
124:    SVDGetTolerances - Gets the tolerance and maximum
125:    iteration count used by the default SVD convergence tests.

127:    Not Collective

129:    Input Parameter:
130: .  svd - the singular value solver context

132:    Output Parameters:
133: +  tol - the convergence tolerance
134: -  maxits - maximum number of iterations

136:    Notes:
137:    The user can specify NULL for any parameter that is not needed.

139:    Level: intermediate

141: .seealso: SVDSetTolerances()
142: @*/
143: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,PetscInt *maxits)
144: {
146:   if (tol)    *tol    = svd->tol;
147:   if (maxits) *maxits = svd->max_it;
148:   return 0;
149: }

151: /*@
152:    SVDSetDimensions - Sets the number of singular values to compute
153:    and the dimension of the subspace.

155:    Logically Collective on svd

157:    Input Parameters:
158: +  svd - the singular value solver context
159: .  nsv - number of singular values to compute
160: .  ncv - the maximum dimension of the subspace to be used by the solver
161: -  mpd - the maximum dimension allowed for the projected problem

163:    Options Database Keys:
164: +  -svd_nsv <nsv> - Sets the number of singular values
165: .  -svd_ncv <ncv> - Sets the dimension of the subspace
166: -  -svd_mpd <mpd> - Sets the maximum projected dimension

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

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

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

181:    Level: intermediate

183: .seealso: SVDGetDimensions()
184: @*/
185: PetscErrorCode SVDSetDimensions(SVD svd,PetscInt nsv,PetscInt ncv,PetscInt mpd)
186: {
192:   svd->nsv = nsv;
193:   if (ncv == PETSC_DEFAULT || ncv == PETSC_DECIDE) {
194:     svd->ncv = PETSC_DEFAULT;
195:   } else {
197:     svd->ncv = ncv;
198:   }
199:   if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
200:     svd->mpd = PETSC_DEFAULT;
201:   } else {
203:     svd->mpd = mpd;
204:   }
205:   svd->state = SVD_STATE_INITIAL;
206:   return 0;
207: }

209: /*@C
210:    SVDGetDimensions - Gets the number of singular values to compute
211:    and the dimension of the subspace.

213:    Not Collective

215:    Input Parameter:
216: .  svd - the singular value context

218:    Output Parameters:
219: +  nsv - number of singular values to compute
220: .  ncv - the maximum dimension of the subspace to be used by the solver
221: -  mpd - the maximum dimension allowed for the projected problem

223:    Notes:
224:    The user can specify NULL for any parameter that is not needed.

226:    Level: intermediate

228: .seealso: SVDSetDimensions()
229: @*/
230: PetscErrorCode SVDGetDimensions(SVD svd,PetscInt *nsv,PetscInt *ncv,PetscInt *mpd)
231: {
233:   if (nsv) *nsv = svd->nsv;
234:   if (ncv) *ncv = svd->ncv;
235:   if (mpd) *mpd = svd->mpd;
236:   return 0;
237: }

239: /*@
240:     SVDSetWhichSingularTriplets - Specifies which singular triplets are
241:     to be sought.

243:     Logically Collective on svd

245:     Input Parameter:
246: .   svd - singular value solver context obtained from SVDCreate()

248:     Output Parameter:
249: .   which - which singular triplets are to be sought

251:     Possible values:
252:     The parameter 'which' can have one of these values

254: +     SVD_LARGEST  - largest singular values
255: -     SVD_SMALLEST - smallest singular values

257:     Options Database Keys:
258: +   -svd_largest  - Sets largest singular values
259: -   -svd_smallest - Sets smallest singular values

261:     Level: intermediate

263: .seealso: SVDGetWhichSingularTriplets(), SVDWhich
264: @*/
265: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)
266: {
269:   switch (which) {
270:     case SVD_LARGEST:
271:     case SVD_SMALLEST:
272:       if (svd->which != which) {
273:         svd->state = SVD_STATE_INITIAL;
274:         svd->which = which;
275:       }
276:       break;
277:   default:
278:     SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
279:   }
280:   return 0;
281: }

283: /*@
284:     SVDGetWhichSingularTriplets - Returns which singular triplets are
285:     to be sought.

287:     Not Collective

289:     Input Parameter:
290: .   svd - singular value solver context obtained from SVDCreate()

292:     Output Parameter:
293: .   which - which singular triplets are to be sought

295:     Notes:
296:     See SVDSetWhichSingularTriplets() for possible values of which

298:     Level: intermediate

300: .seealso: SVDSetWhichSingularTriplets(), SVDWhich
301: @*/
302: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)
303: {
306:   *which = svd->which;
307:   return 0;
308: }

310: /*@C
311:    SVDSetConvergenceTestFunction - Sets a function to compute the error estimate
312:    used in the convergence test.

314:    Logically Collective on svd

316:    Input Parameters:
317: +  svd     - singular value solver context obtained from SVDCreate()
318: .  func    - a pointer to the convergence test function
319: .  ctx     - context for private data for the convergence routine (may be null)
320: -  destroy - a routine for destroying the context (may be null)

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

325: +   svd    - singular value solver context obtained from SVDCreate()
326: .   sigma  - computed singular value
327: .   res    - residual norm associated to the singular triplet
328: .   errest - (output) computed error estimate
329: -   ctx    - optional context, as set by SVDSetConvergenceTestFunction()

331:    Note:
332:    If the error estimate returned by the convergence test function is less than
333:    the tolerance, then the singular value is accepted as converged.

335:    Level: advanced

337: .seealso: SVDSetConvergenceTest(), SVDSetTolerances()
338: @*/
339: PetscErrorCode SVDSetConvergenceTestFunction(SVD svd,PetscErrorCode (*func)(SVD,PetscReal,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
340: {
342:   if (svd->convergeddestroy) (*svd->convergeddestroy)(svd->convergedctx);
343:   svd->convergeduser    = func;
344:   svd->convergeddestroy = destroy;
345:   svd->convergedctx     = ctx;
346:   if (func == SVDConvergedAbsolute) svd->conv = SVD_CONV_ABS;
347:   else if (func == SVDConvergedRelative) svd->conv = SVD_CONV_REL;
348:   else if (func == SVDConvergedNorm) svd->conv = SVD_CONV_NORM;
349:   else if (func == SVDConvergedMaxIt) svd->conv = SVD_CONV_MAXIT;
350:   else {
351:     svd->conv      = SVD_CONV_USER;
352:     svd->converged = svd->convergeduser;
353:   }
354:   return 0;
355: }

357: /*@
358:    SVDSetConvergenceTest - Specifies how to compute the error estimate
359:    used in the convergence test.

361:    Logically Collective on svd

363:    Input Parameters:
364: +  svd  - singular value solver context obtained from SVDCreate()
365: -  conv - the type of convergence test

367:    Options Database Keys:
368: +  -svd_conv_abs   - Sets the absolute convergence test
369: .  -svd_conv_rel   - Sets the convergence test relative to the singular value
370: .  -svd_conv_norm  - Sets the convergence test relative to the matrix norm
371: .  -svd_conv_maxit - Forces the maximum number of iterations as set by -svd_max_it
372: -  -svd_conv_user  - Selects the user-defined convergence test

374:    Notes:
375:    The parameter 'conv' can have one of these values
376: +     SVD_CONV_ABS   - absolute error ||r||
377: .     SVD_CONV_REL   - error relative to the singular value sigma, ||r||/sigma
378: .     SVD_CONV_NORM  - error relative to the matrix norms, ||r||/||Z||, with Z=A or Z=[A;B]
379: .     SVD_CONV_MAXIT - no convergence until maximum number of iterations has been reached
380: -     SVD_CONV_USER  - function set by SVDSetConvergenceTestFunction()

382:    The default in standard SVD is SVD_CONV_REL, while in GSVD the default is SVD_CONV_NORM.

384:    Level: intermediate

386: .seealso: SVDGetConvergenceTest(), SVDSetConvergenceTestFunction(), SVDSetStoppingTest(), SVDConv
387: @*/
388: PetscErrorCode SVDSetConvergenceTest(SVD svd,SVDConv conv)
389: {
392:   switch (conv) {
393:     case SVD_CONV_ABS:   svd->converged = SVDConvergedAbsolute; break;
394:     case SVD_CONV_REL:   svd->converged = SVDConvergedRelative; break;
395:     case SVD_CONV_NORM:  svd->converged = SVDConvergedNorm; break;
396:     case SVD_CONV_MAXIT: svd->converged = SVDConvergedMaxIt; break;
397:     case SVD_CONV_USER:
399:       svd->converged = svd->convergeduser;
400:       break;
401:     default:
402:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
403:   }
404:   svd->conv = conv;
405:   return 0;
406: }

408: /*@
409:    SVDGetConvergenceTest - Gets the method used to compute the error estimate
410:    used in the convergence test.

412:    Not Collective

414:    Input Parameters:
415: .  svd   - singular value solver context obtained from SVDCreate()

417:    Output Parameters:
418: .  conv  - the type of convergence test

420:    Level: intermediate

422: .seealso: SVDSetConvergenceTest(), SVDConv
423: @*/
424: PetscErrorCode SVDGetConvergenceTest(SVD svd,SVDConv *conv)
425: {
428:   *conv = svd->conv;
429:   return 0;
430: }

432: /*@C
433:    SVDSetStoppingTestFunction - Sets a function to decide when to stop the outer
434:    iteration of the singular value solver.

436:    Logically Collective on svd

438:    Input Parameters:
439: +  svd     - singular value solver context obtained from SVDCreate()
440: .  func    - pointer to the stopping test function
441: .  ctx     - context for private data for the stopping routine (may be null)
442: -  destroy - a routine for destroying the context (may be null)

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

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

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

461:    Level: advanced

463: .seealso: SVDSetStoppingTest(), SVDStoppingBasic()
464: @*/
465: PetscErrorCode SVDSetStoppingTestFunction(SVD svd,PetscErrorCode (*func)(SVD,PetscInt,PetscInt,PetscInt,PetscInt,SVDConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
466: {
468:   if (svd->stoppingdestroy) (*svd->stoppingdestroy)(svd->stoppingctx);
469:   svd->stoppinguser    = func;
470:   svd->stoppingdestroy = destroy;
471:   svd->stoppingctx     = ctx;
472:   if (func == SVDStoppingBasic) svd->stop = SVD_STOP_BASIC;
473:   else {
474:     svd->stop     = SVD_STOP_USER;
475:     svd->stopping = svd->stoppinguser;
476:   }
477:   return 0;
478: }

480: /*@
481:    SVDSetStoppingTest - Specifies how to decide the termination of the outer
482:    loop of the singular value solver.

484:    Logically Collective on svd

486:    Input Parameters:
487: +  svd  - singular value solver context obtained from SVDCreate()
488: -  stop - the type of stopping test

490:    Options Database Keys:
491: +  -svd_stop_basic - Sets the default stopping test
492: -  -svd_stop_user  - Selects the user-defined stopping test

494:    Note:
495:    The parameter 'stop' can have one of these values
496: +     SVD_STOP_BASIC - default stopping test
497: -     SVD_STOP_USER  - function set by SVDSetStoppingTestFunction()

499:    Level: advanced

501: .seealso: SVDGetStoppingTest(), SVDSetStoppingTestFunction(), SVDSetConvergenceTest(), SVDStop
502: @*/
503: PetscErrorCode SVDSetStoppingTest(SVD svd,SVDStop stop)
504: {
507:   switch (stop) {
508:     case SVD_STOP_BASIC: svd->stopping = SVDStoppingBasic; break;
509:     case SVD_STOP_USER:
511:       svd->stopping = svd->stoppinguser;
512:       break;
513:     default:
514:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
515:   }
516:   svd->stop = stop;
517:   return 0;
518: }

520: /*@
521:    SVDGetStoppingTest - Gets the method used to decide the termination of the outer
522:    loop of the singular value solver.

524:    Not Collective

526:    Input Parameters:
527: .  svd   - singular value solver context obtained from SVDCreate()

529:    Output Parameters:
530: .  stop  - the type of stopping test

532:    Level: advanced

534: .seealso: SVDSetStoppingTest(), SVDStop
535: @*/
536: PetscErrorCode SVDGetStoppingTest(SVD svd,SVDStop *stop)
537: {
540:   *stop = svd->stop;
541:   return 0;
542: }

544: /*@C
545:    SVDMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
546:    indicated by the user.

548:    Collective on svd

550:    Input Parameters:
551: +  svd      - the singular value solver context
552: .  opt      - the command line option for this monitor
553: .  name     - the monitor type one is seeking
554: .  ctx      - an optional user context for the monitor, or NULL
555: -  trackall - whether this monitor tracks all singular values or not

557:    Level: developer

559: .seealso: SVDMonitorSet(), SVDSetTrackAll()
560: @*/
561: PetscErrorCode SVDMonitorSetFromOptions(SVD svd,const char opt[],const char name[],void *ctx,PetscBool trackall)
562: {
563:   PetscErrorCode       (*mfunc)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*);
564:   PetscErrorCode       (*cfunc)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**);
565:   PetscErrorCode       (*dfunc)(PetscViewerAndFormat**);
566:   PetscViewerAndFormat *vf;
567:   PetscViewer          viewer;
568:   PetscViewerFormat    format;
569:   PetscViewerType      vtype;
570:   char                 key[PETSC_MAX_PATH_LEN];
571:   PetscBool            flg;

573:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,opt,&viewer,&format,&flg);
574:   if (!flg) return 0;

576:   PetscViewerGetType(viewer,&vtype);
577:   SlepcMonitorMakeKey_Internal(name,vtype,format,key);
578:   PetscFunctionListFind(SVDMonitorList,key,&mfunc);
580:   PetscFunctionListFind(SVDMonitorCreateList,key,&cfunc);
581:   PetscFunctionListFind(SVDMonitorDestroyList,key,&dfunc);
582:   if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
583:   if (!dfunc) dfunc = PetscViewerAndFormatDestroy;

585:   (*cfunc)(viewer,format,ctx,&vf);
586:   PetscObjectDereference((PetscObject)viewer);
587:   SVDMonitorSet(svd,mfunc,vf,(PetscErrorCode(*)(void **))dfunc);
588:   if (trackall) SVDSetTrackAll(svd,PETSC_TRUE);
589:   return 0;
590: }

592: /*@
593:    SVDSetFromOptions - Sets SVD options from the options database.
594:    This routine must be called before SVDSetUp() if the user is to be
595:    allowed to set the solver type.

597:    Collective on svd

599:    Input Parameters:
600: .  svd - the singular value solver context

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

605:    Level: beginner

607: .seealso: SVDSetOptionsPrefix()
608: @*/
609: PetscErrorCode SVDSetFromOptions(SVD svd)
610: {
611:   char           type[256];
612:   PetscBool      set,flg,val,flg1,flg2,flg3;
613:   PetscInt       i,j,k;
614:   PetscReal      r;

617:   SVDRegisterAll();
618:   PetscObjectOptionsBegin((PetscObject)svd);
619:     PetscOptionsFList("-svd_type","SVD solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,sizeof(type),&flg);
620:     if (flg) SVDSetType(svd,type);
621:     else if (!((PetscObject)svd)->type_name) SVDSetType(svd,SVDCROSS);

623:     PetscOptionsBoolGroupBegin("-svd_standard","Singular value decomposition (SVD)","SVDSetProblemType",&flg);
624:     if (flg) SVDSetProblemType(svd,SVD_STANDARD);
625:     PetscOptionsBoolGroup("-svd_generalized","Generalized singular value decomposition (GSVD)","SVDSetProblemType",&flg);
626:     if (flg) SVDSetProblemType(svd,SVD_GENERALIZED);
627:     PetscOptionsBoolGroupEnd("-svd_hyperbolic","Hyperbolic singular value decomposition (HSVD)","SVDSetProblemType",&flg);
628:     if (flg) SVDSetProblemType(svd,SVD_HYPERBOLIC);

630:     PetscOptionsBool("-svd_implicittranspose","Handle matrix transpose implicitly","SVDSetImplicitTranspose",svd->impltrans,&val,&flg);
631:     if (flg) SVDSetImplicitTranspose(svd,val);

633:     i = svd->max_it;
634:     PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,&flg1);
635:     r = svd->tol;
636:     PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",SlepcDefaultTol(svd->tol),&r,&flg2);
637:     if (flg1 || flg2) SVDSetTolerances(svd,r,i);

639:     PetscOptionsBoolGroupBegin("-svd_conv_abs","Absolute error convergence test","SVDSetConvergenceTest",&flg);
640:     if (flg) SVDSetConvergenceTest(svd,SVD_CONV_ABS);
641:     PetscOptionsBoolGroup("-svd_conv_rel","Relative error convergence test","SVDSetConvergenceTest",&flg);
642:     if (flg) SVDSetConvergenceTest(svd,SVD_CONV_REL);
643:     PetscOptionsBoolGroup("-svd_conv_norm","Convergence test relative to the matrix norms","SVDSetConvergenceTest",&flg);
644:     if (flg) SVDSetConvergenceTest(svd,SVD_CONV_NORM);
645:     PetscOptionsBoolGroup("-svd_conv_maxit","Maximum iterations convergence test","SVDSetConvergenceTest",&flg);
646:     if (flg) SVDSetConvergenceTest(svd,SVD_CONV_MAXIT);
647:     PetscOptionsBoolGroupEnd("-svd_conv_user","User-defined convergence test","SVDSetConvergenceTest",&flg);
648:     if (flg) SVDSetConvergenceTest(svd,SVD_CONV_USER);

650:     PetscOptionsBoolGroupBegin("-svd_stop_basic","Stop iteration if all singular values converged or max_it reached","SVDSetStoppingTest",&flg);
651:     if (flg) SVDSetStoppingTest(svd,SVD_STOP_BASIC);
652:     PetscOptionsBoolGroupEnd("-svd_stop_user","User-defined stopping test","SVDSetStoppingTest",&flg);
653:     if (flg) SVDSetStoppingTest(svd,SVD_STOP_USER);

655:     i = svd->nsv;
656:     PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,&flg1);
657:     j = svd->ncv;
658:     PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,&flg2);
659:     k = svd->mpd;
660:     PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,&flg3);
661:     if (flg1 || flg2 || flg3) SVDSetDimensions(svd,i,j,k);

663:     PetscOptionsBoolGroupBegin("-svd_largest","Compute largest singular values","SVDSetWhichSingularTriplets",&flg);
664:     if (flg) SVDSetWhichSingularTriplets(svd,SVD_LARGEST);
665:     PetscOptionsBoolGroupEnd("-svd_smallest","Compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
666:     if (flg) SVDSetWhichSingularTriplets(svd,SVD_SMALLEST);

668:     /* -----------------------------------------------------------------------*/
669:     /*
670:       Cancels all monitors hardwired into code before call to SVDSetFromOptions()
671:     */
672:     PetscOptionsBool("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",PETSC_FALSE,&flg,&set);
673:     if (set && flg) SVDMonitorCancel(svd);
674:     SVDMonitorSetFromOptions(svd,"-svd_monitor","first_approximation",NULL,PETSC_FALSE);
675:     SVDMonitorSetFromOptions(svd,"-svd_monitor_all","all_approximations",NULL,PETSC_TRUE);
676:     SVDMonitorSetFromOptions(svd,"-svd_monitor_conv","convergence_history",NULL,PETSC_FALSE);
677:     SVDMonitorSetFromOptions(svd,"-svd_monitor_conditioning","conditioning",NULL,PETSC_FALSE);

679:     /* -----------------------------------------------------------------------*/
680:     PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",NULL);
681:     PetscOptionsName("-svd_view_vectors","View computed singular vectors","SVDVectorsView",NULL);
682:     PetscOptionsName("-svd_view_values","View computed singular values","SVDValuesView",NULL);
683:     PetscOptionsName("-svd_converged_reason","Print reason for convergence, and number of iterations","SVDConvergedReasonView",NULL);
684:     PetscOptionsName("-svd_error_absolute","Print absolute errors of each singular triplet","SVDErrorView",NULL);
685:     PetscOptionsName("-svd_error_relative","Print relative errors of each singular triplet","SVDErrorView",NULL);
686:     PetscOptionsName("-svd_error_norm","Print errors relative to the matrix norms of each singular triplet","SVDErrorView",NULL);

688:     PetscTryTypeMethod(svd,setfromoptions,PetscOptionsObject);
689:     PetscObjectProcessOptionsHandlers((PetscObject)svd,PetscOptionsObject);
690:   PetscOptionsEnd();

692:   if (!svd->V) SVDGetBV(svd,&svd->V,NULL);
693:   BVSetFromOptions(svd->V);
694:   if (!svd->U) SVDGetBV(svd,NULL,&svd->U);
695:   BVSetFromOptions(svd->U);
696:   if (!svd->ds) SVDGetDS(svd,&svd->ds);
697:   DSSetFromOptions(svd->ds);
698:   return 0;
699: }

701: /*@
702:    SVDSetProblemType - Specifies the type of the singular value problem.

704:    Logically Collective on svd

706:    Input Parameters:
707: +  svd  - the singular value solver context
708: -  type - a known type of singular value problem

710:    Options Database Keys:
711: +  -svd_standard    - standard singular value decomposition (SVD)
712: .  -svd_generalized - generalized singular value problem (GSVD)
713: -  -svd_hyperbolic  - hyperbolic singular value problem (HSVD)

715:    Notes:
716:    The GSVD requires that two matrices have been passed via SVDSetOperators().
717:    The HSVD requires that a signature matrix has been passed via SVDSetSignature().

719:    Level: intermediate

721: .seealso: SVDSetOperators(), SVDSetSignature(), SVDSetType(), SVDGetProblemType(), SVDProblemType
722: @*/
723: PetscErrorCode SVDSetProblemType(SVD svd,SVDProblemType type)
724: {
727:   if (type == svd->problem_type) return 0;
728:   switch (type) {
729:     case SVD_STANDARD:
730:       svd->isgeneralized = PETSC_FALSE;
731:       svd->ishyperbolic  = PETSC_FALSE;
732:       break;
733:     case SVD_GENERALIZED:
734:       svd->isgeneralized = PETSC_TRUE;
735:       svd->ishyperbolic  = PETSC_FALSE;
736:       break;
737:     case SVD_HYPERBOLIC:
738:       svd->isgeneralized = PETSC_FALSE;
739:       svd->ishyperbolic  = PETSC_TRUE;
740:       break;
741:     default:
742:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_WRONG,"Unknown singular value problem type");
743:   }
744:   svd->problem_type = type;
745:   svd->state = SVD_STATE_INITIAL;
746:   return 0;
747: }

749: /*@
750:    SVDGetProblemType - Gets the problem type from the SVD object.

752:    Not Collective

754:    Input Parameter:
755: .  svd - the singular value solver context

757:    Output Parameter:
758: .  type - the problem type

760:    Level: intermediate

762: .seealso: SVDSetProblemType(), SVDProblemType
763: @*/
764: PetscErrorCode SVDGetProblemType(SVD svd,SVDProblemType *type)
765: {
768:   *type = svd->problem_type;
769:   return 0;
770: }

772: /*@
773:    SVDIsGeneralized - Ask if the SVD object corresponds to a generalized
774:    singular value problem.

776:    Not collective

778:    Input Parameter:
779: .  svd - the singular value solver context

781:    Output Parameter:
782: .  is - the answer

784:    Level: intermediate

786: .seealso: SVDIsHyperbolic()
787: @*/
788: PetscErrorCode SVDIsGeneralized(SVD svd,PetscBool* is)
789: {
792:   *is = svd->isgeneralized;
793:   return 0;
794: }

796: /*@
797:    SVDIsHyperbolic - Ask if the SVD object corresponds to a hyperbolic
798:    singular value problem.

800:    Not collective

802:    Input Parameter:
803: .  svd - the singular value solver context

805:    Output Parameter:
806: .  is - the answer

808:    Level: intermediate

810: .seealso: SVDIsGeneralized()
811: @*/
812: PetscErrorCode SVDIsHyperbolic(SVD svd,PetscBool* is)
813: {
816:   *is = svd->ishyperbolic;
817:   return 0;
818: }

820: /*@
821:    SVDSetTrackAll - Specifies if the solver must compute the residual norm of all
822:    approximate singular value or not.

824:    Logically Collective on svd

826:    Input Parameters:
827: +  svd      - the singular value solver context
828: -  trackall - whether to compute all residuals or not

830:    Notes:
831:    If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
832:    the residual norm for each singular value approximation. Computing the residual is
833:    usually an expensive operation and solvers commonly compute only the residual
834:    associated to the first unconverged singular value.

836:    The option '-svd_monitor_all' automatically activates this option.

838:    Level: developer

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

850: /*@
851:    SVDGetTrackAll - Returns the flag indicating whether all residual norms must
852:    be computed or not.

854:    Not Collective

856:    Input Parameter:
857: .  svd - the singular value solver context

859:    Output Parameter:
860: .  trackall - the returned flag

862:    Level: developer

864: .seealso: SVDSetTrackAll()
865: @*/
866: PetscErrorCode SVDGetTrackAll(SVD svd,PetscBool *trackall)
867: {
870:   *trackall = svd->trackall;
871:   return 0;
872: }

874: /*@C
875:    SVDSetOptionsPrefix - Sets the prefix used for searching for all
876:    SVD options in the database.

878:    Logically Collective on svd

880:    Input Parameters:
881: +  svd - the singular value solver context
882: -  prefix - the prefix string to prepend to all SVD option requests

884:    Notes:
885:    A hyphen (-) must NOT be given at the beginning of the prefix name.
886:    The first character of all runtime options is AUTOMATICALLY the
887:    hyphen.

889:    For example, to distinguish between the runtime options for two
890:    different SVD contexts, one could call
891: .vb
892:       SVDSetOptionsPrefix(svd1,"svd1_")
893:       SVDSetOptionsPrefix(svd2,"svd2_")
894: .ve

896:    Level: advanced

898: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
899: @*/
900: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)
901: {
903:   if (!svd->V) SVDGetBV(svd,&svd->V,&svd->U);
904:   BVSetOptionsPrefix(svd->V,prefix);
905:   BVSetOptionsPrefix(svd->U,prefix);
906:   if (!svd->ds) SVDGetDS(svd,&svd->ds);
907:   DSSetOptionsPrefix(svd->ds,prefix);
908:   PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
909:   return 0;
910: }

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

916:    Logically Collective on svd

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

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

926:    Level: advanced

928: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
929: @*/
930: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
931: {
933:   if (!svd->V) SVDGetBV(svd,&svd->V,&svd->U);
934:   BVAppendOptionsPrefix(svd->V,prefix);
935:   BVAppendOptionsPrefix(svd->U,prefix);
936:   if (!svd->ds) SVDGetDS(svd,&svd->ds);
937:   DSAppendOptionsPrefix(svd->ds,prefix);
938:   PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
939:   return 0;
940: }

942: /*@C
943:    SVDGetOptionsPrefix - Gets the prefix used for searching for all
944:    SVD options in the database.

946:    Not Collective

948:    Input Parameters:
949: .  svd - the singular value solver context

951:    Output Parameters:
952: .  prefix - pointer to the prefix string used is returned

954:    Note:
955:    On the Fortran side, the user should pass in a string 'prefix' of
956:    sufficient length to hold the prefix.

958:    Level: advanced

960: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
961: @*/
962: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
963: {
966:   PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
967:   return 0;
968: }