Actual source code: svdopts.c

slepc-3.16.3 2022-04-11
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2021, 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: {
 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(), SVDSetOperators()
 71: @*/
 72: PetscErrorCode SVDGetImplicitTranspose(SVD svd,PetscBool *impl)
 73: {
 77:   *impl = svd->impltrans;
 78:   return(0);
 79: }

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

 85:    Logically Collective on svd

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

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

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

 99:    Level: intermediate

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

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

130:    Not Collective

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

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

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

142:    Level: intermediate

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

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

159:    Logically Collective on svd

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

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

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

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

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

185:    Level: intermediate

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

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

218:    Not Collective

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

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

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

231:    Level: intermediate

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

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

249:     Logically Collective on svd

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

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

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

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

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

267:     Level: intermediate

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

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

294:     Not Collective

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

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

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

305:     Level: intermediate

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

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

322:    Logically Collective on svd

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

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

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

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

343:    Level: advanced

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

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

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

374:    Logically Collective on svd

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

380:    Options Database Keys:
381: +  -svd_conv_abs   - Sets the absolute convergence test
382: .  -svd_conv_rel   - Sets the convergence test relative to the singular value
383: .  -svd_conv_norm  - Sets the convergence test relative to the matrix norms
384: .  -svd_conv_maxit - Forces the maximum number of iterations as set by -svd_max_it
385: -  -svd_conv_user  - Selects the user-defined convergence test

387:    Note:
388:    The parameter 'conv' can have one of these values
389: +     SVD_CONV_ABS   - absolute error ||r||
390: .     SVD_CONV_REL   - error relative to the singular value l, ||r||/sigma
391: .     SVD_CONV_NORM  - error relative to the matrix norms, ||r||/||A||
392: .     SVD_CONV_MAXIT - no convergence until maximum number of iterations has been reached
393: -     SVD_CONV_USER  - function set by SVDSetConvergenceTestFunction()

395:    Level: intermediate

397: .seealso: SVDGetConvergenceTest(), SVDSetConvergenceTestFunction(), SVDSetStoppingTest(), SVDConv
398: @*/
399: PetscErrorCode SVDSetConvergenceTest(SVD svd,SVDConv conv)
400: {
404:   switch (conv) {
405:     case SVD_CONV_ABS:   svd->converged = SVDConvergedAbsolute; break;
406:     case SVD_CONV_REL:   svd->converged = SVDConvergedRelative; break;
407:     case SVD_CONV_NORM:  svd->converged = SVDConvergedNorm; break;
408:     case SVD_CONV_MAXIT: svd->converged = SVDConvergedMaxIt; break;
409:     case SVD_CONV_USER:
410:       if (!svd->convergeduser) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ORDER,"Must call SVDSetConvergenceTestFunction() first");
411:       svd->converged = svd->convergeduser;
412:       break;
413:     default:
414:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
415:   }
416:   svd->conv = conv;
417:   return(0);
418: }

420: /*@
421:    SVDGetConvergenceTest - Gets the method used to compute the error estimate
422:    used in the convergence test.

424:    Not Collective

426:    Input Parameters:
427: .  svd   - singular value solver context obtained from SVDCreate()

429:    Output Parameters:
430: .  conv  - the type of convergence test

432:    Level: intermediate

434: .seealso: SVDSetConvergenceTest(), SVDConv
435: @*/
436: PetscErrorCode SVDGetConvergenceTest(SVD svd,SVDConv *conv)
437: {
441:   *conv = svd->conv;
442:   return(0);
443: }

445: /*@C
446:    SVDSetStoppingTestFunction - Sets a function to decide when to stop the outer
447:    iteration of the singular value solver.

449:    Logically Collective on svd

451:    Input Parameters:
452: +  svd     - singular value solver context obtained from SVDCreate()
453: .  func    - pointer to the stopping test function
454: .  ctx     - context for private data for the stopping routine (may be null)
455: -  destroy - a routine for destroying the context (may be null)

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

460: +   svd    - singular value solver context obtained from SVDCreate()
461: .   its    - current number of iterations
462: .   max_it - maximum number of iterations
463: .   nconv  - number of currently converged singular triplets
464: .   nsv    - number of requested singular triplets
465: .   reason - (output) result of the stopping test
466: -   ctx    - optional context, as set by SVDSetStoppingTestFunction()

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

474:    Level: advanced

476: .seealso: SVDSetStoppingTest(), SVDStoppingBasic()
477: @*/
478: PetscErrorCode SVDSetStoppingTestFunction(SVD svd,PetscErrorCode (*func)(SVD,PetscInt,PetscInt,PetscInt,PetscInt,SVDConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
479: {

484:   if (svd->stoppingdestroy) {
485:     (*svd->stoppingdestroy)(svd->stoppingctx);
486:   }
487:   svd->stoppinguser    = func;
488:   svd->stoppingdestroy = destroy;
489:   svd->stoppingctx     = ctx;
490:   if (func == SVDStoppingBasic) svd->stop = SVD_STOP_BASIC;
491:   else {
492:     svd->stop     = SVD_STOP_USER;
493:     svd->stopping = svd->stoppinguser;
494:   }
495:   return(0);
496: }

498: /*@
499:    SVDSetStoppingTest - Specifies how to decide the termination of the outer
500:    loop of the singular value solver.

502:    Logically Collective on svd

504:    Input Parameters:
505: +  svd  - singular value solver context obtained from SVDCreate()
506: -  stop - the type of stopping test

508:    Options Database Keys:
509: +  -svd_stop_basic - Sets the default stopping test
510: -  -svd_stop_user  - Selects the user-defined stopping test

512:    Note:
513:    The parameter 'stop' can have one of these values
514: +     SVD_STOP_BASIC - default stopping test
515: -     SVD_STOP_USER  - function set by SVDSetStoppingTestFunction()

517:    Level: advanced

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

539: /*@
540:    SVDGetStoppingTest - Gets the method used to decide the termination of the outer
541:    loop of the singular value solver.

543:    Not Collective

545:    Input Parameters:
546: .  svd   - singular value solver context obtained from SVDCreate()

548:    Output Parameters:
549: .  stop  - the type of stopping test

551:    Level: advanced

553: .seealso: SVDSetStoppingTest(), SVDStop
554: @*/
555: PetscErrorCode SVDGetStoppingTest(SVD svd,SVDStop *stop)
556: {
560:   *stop = svd->stop;
561:   return(0);
562: }

564: /*@C
565:    SVDMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
566:    indicated by the user.

568:    Collective on svd

570:    Input Parameters:
571: +  svd      - the singular value solver context
572: .  opt      - the command line option for this monitor
573: .  name     - the monitor type one is seeking
574: .  ctx      - an optional user context for the monitor, or NULL
575: -  trackall - whether this monitor tracks all singular values or not

577:    Level: developer

579: .seealso: SVDMonitorSet(), SVDSetTrackAll()
580: @*/
581: PetscErrorCode SVDMonitorSetFromOptions(SVD svd,const char opt[],const char name[],void *ctx,PetscBool trackall)
582: {
583:   PetscErrorCode       (*mfunc)(SVD,PetscInt,PetscInt,PetscReal*,PetscReal*,PetscInt,void*);
584:   PetscErrorCode       (*cfunc)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**);
585:   PetscErrorCode       (*dfunc)(PetscViewerAndFormat**);
586:   PetscViewerAndFormat *vf;
587:   PetscViewer          viewer;
588:   PetscViewerFormat    format;
589:   PetscViewerType      vtype;
590:   char                 key[PETSC_MAX_PATH_LEN];
591:   PetscBool            flg;
592:   PetscErrorCode       ierr;

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

598:   PetscViewerGetType(viewer,&vtype);
599:   SlepcMonitorMakeKey_Internal(name,vtype,format,key);
600:   PetscFunctionListFind(SVDMonitorList,key,&mfunc);
601:   if (!mfunc) SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"Specified viewer and format not supported");
602:   PetscFunctionListFind(SVDMonitorCreateList,key,&cfunc);
603:   PetscFunctionListFind(SVDMonitorDestroyList,key,&dfunc);
604:   if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
605:   if (!dfunc) dfunc = PetscViewerAndFormatDestroy;

607:   (*cfunc)(viewer,format,ctx,&vf);
608:   PetscObjectDereference((PetscObject)viewer);
609:   SVDMonitorSet(svd,mfunc,vf,(PetscErrorCode(*)(void **))dfunc);
610:   if (trackall) {
611:     SVDSetTrackAll(svd,PETSC_TRUE);
612:   }
613:   return(0);
614: }

616: /*@
617:    SVDSetFromOptions - Sets SVD options from the options database.
618:    This routine must be called before SVDSetUp() if the user is to be
619:    allowed to set the solver type.

621:    Collective on svd

623:    Input Parameters:
624: .  svd - the singular value solver context

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

629:    Level: beginner

631: .seealso:
632: @*/
633: PetscErrorCode SVDSetFromOptions(SVD svd)
634: {
636:   char           type[256];
637:   PetscBool      set,flg,val,flg1,flg2,flg3;
638:   PetscInt       i,j,k;
639:   PetscReal      r;

643:   SVDRegisterAll();
644:   PetscObjectOptionsBegin((PetscObject)svd);
645:     PetscOptionsFList("-svd_type","SVD solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,sizeof(type),&flg);
646:     if (flg) {
647:       SVDSetType(svd,type);
648:     } else if (!((PetscObject)svd)->type_name) {
649:       SVDSetType(svd,SVDCROSS);
650:     }

652:     PetscOptionsBoolGroupBegin("-svd_standard","Singular value decomposition (SVD)","SVDSetProblemType",&flg);
653:     if (flg) { SVDSetProblemType(svd,SVD_STANDARD); }
654:     PetscOptionsBoolGroupEnd("-svd_generalized","Generalized singular value decomposition (GSVD)","SVDSetProblemType",&flg);
655:     if (flg) { SVDSetProblemType(svd,SVD_GENERALIZED); }

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

660:     i = svd->max_it;
661:     PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,&flg1);
662:     r = svd->tol;
663:     PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",SlepcDefaultTol(svd->tol),&r,&flg2);
664:     if (flg1 || flg2) { SVDSetTolerances(svd,r,i); }

666:     PetscOptionsBoolGroupBegin("-svd_conv_abs","Absolute error convergence test","SVDSetConvergenceTest",&flg);
667:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_ABS); }
668:     PetscOptionsBoolGroup("-svd_conv_rel","Relative error convergence test","SVDSetConvergenceTest",&flg);
669:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_REL); }
670:     PetscOptionsBoolGroup("-svd_conv_norm","Convergence test relative to the matrix norms","SVDSetConvergenceTest",&flg);
671:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_NORM); }
672:     PetscOptionsBoolGroup("-svd_conv_maxit","Maximum iterations convergence test","SVDSetConvergenceTest",&flg);
673:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_MAXIT); }
674:     PetscOptionsBoolGroupEnd("-svd_conv_user","User-defined convergence test","SVDSetConvergenceTest",&flg);
675:     if (flg) { SVDSetConvergenceTest(svd,SVD_CONV_USER); }

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

682:     i = svd->nsv;
683:     PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,&flg1);
684:     j = svd->ncv;
685:     PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,&flg2);
686:     k = svd->mpd;
687:     PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,&flg3);
688:     if (flg1 || flg2 || flg3) { SVDSetDimensions(svd,i,j,k); }

690:     PetscOptionsBoolGroupBegin("-svd_largest","Compute largest singular values","SVDSetWhichSingularTriplets",&flg);
691:     if (flg) { SVDSetWhichSingularTriplets(svd,SVD_LARGEST); }
692:     PetscOptionsBoolGroupEnd("-svd_smallest","Compute smallest singular values","SVDSetWhichSingularTriplets",&flg);
693:     if (flg) { SVDSetWhichSingularTriplets(svd,SVD_SMALLEST); }

695:     /* -----------------------------------------------------------------------*/
696:     /*
697:       Cancels all monitors hardwired into code before call to SVDSetFromOptions()
698:     */
699:     PetscOptionsBool("-svd_monitor_cancel","Remove any hardwired monitor routines","SVDMonitorCancel",PETSC_FALSE,&flg,&set);
700:     if (set && flg) { SVDMonitorCancel(svd); }
701:     SVDMonitorSetFromOptions(svd,"-svd_monitor","first_approximation",NULL,PETSC_FALSE);
702:     SVDMonitorSetFromOptions(svd,"-svd_monitor_all","all_approximations",NULL,PETSC_TRUE);
703:     SVDMonitorSetFromOptions(svd,"-svd_monitor_conv","convergence_history",NULL,PETSC_FALSE);

705:     /* -----------------------------------------------------------------------*/
706:     PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",NULL);
707:     PetscOptionsName("-svd_view_vectors","View computed singular vectors","SVDVectorsView",NULL);
708:     PetscOptionsName("-svd_view_values","View computed singular values","SVDValuesView",NULL);
709:     PetscOptionsName("-svd_converged_reason","Print reason for convergence, and number of iterations","SVDConvergedReasonView",NULL);
710:     PetscOptionsName("-svd_error_absolute","Print absolute errors of each singular triplet","SVDErrorView",NULL);
711:     PetscOptionsName("-svd_error_relative","Print relative errors of each singular triplet","SVDErrorView",NULL);

713:     if (svd->ops->setfromoptions) {
714:       (*svd->ops->setfromoptions)(PetscOptionsObject,svd);
715:     }
716:     PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)svd);
717:   PetscOptionsEnd();

719:   if (!svd->V) { SVDGetBV(svd,&svd->V,NULL); }
720:   BVSetFromOptions(svd->V);
721:   if (!svd->U) { SVDGetBV(svd,NULL,&svd->U); }
722:   BVSetFromOptions(svd->U);
723:   if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
724:   DSSetFromOptions(svd->ds);
725:   return(0);
726: }

728: /*@
729:    SVDSetProblemType - Specifies the type of the singular value problem.

731:    Logically Collective on svd

733:    Input Parameters:
734: +  svd  - the singular value solver context
735: -  type - a known type of singular value problem

737:    Options Database Keys:
738: +  -svd_standard    - standard singular value decomposition (SVD)
739: -  -svd_generalized - generalized singular value problem (GSVD)

741:    Notes:
742:    The GSVD requires that two matrices have been passed via SVDSetOperators().

744:    Level: intermediate

746: .seealso: SVDSetOperators(), SVDSetType(), SVDGetProblemType(), SVDProblemType
747: @*/
748: PetscErrorCode SVDSetProblemType(SVD svd,SVDProblemType type)
749: {
753:   if (type == svd->problem_type) return(0);
754:   switch (type) {
755:     case SVD_STANDARD:
756:       svd->isgeneralized = PETSC_FALSE;
757:       break;
758:     case SVD_GENERALIZED:
759:       svd->isgeneralized = PETSC_TRUE;
760:       break;
761:     default:
762:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_WRONG,"Unknown singular value problem type");
763:   }
764:   svd->problem_type = type;
765:   svd->state = SVD_STATE_INITIAL;
766:   return(0);
767: }

769: /*@
770:    SVDGetProblemType - Gets the problem type from the SVD object.

772:    Not Collective

774:    Input Parameter:
775: .  svd - the singular value solver context

777:    Output Parameter:
778: .  type - the problem type

780:    Level: intermediate

782: .seealso: SVDSetProblemType(), SVDProblemType
783: @*/
784: PetscErrorCode SVDGetProblemType(SVD svd,SVDProblemType *type)
785: {
789:   *type = svd->problem_type;
790:   return(0);
791: }

793: /*@
794:    SVDIsGeneralized - Ask if the SVD object corresponds to a generalized
795:    singular value problem.

797:    Not collective

799:    Input Parameter:
800: .  svd - the singular value solver context

802:    Output Parameter:
803: .  is - the answer

805:    Level: intermediate

807: .seealso: SVDIsHermitian(), SVDIsPositive()
808: @*/
809: PetscErrorCode SVDIsGeneralized(SVD svd,PetscBool* is)
810: {
814:   *is = svd->isgeneralized;
815:   return(0);
816: }

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

822:    Logically Collective on svd

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

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

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

836:    Level: developer

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

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

853:    Not Collective

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

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

861:    Level: developer

863: .seealso: SVDSetTrackAll()
864: @*/
865: PetscErrorCode SVDGetTrackAll(SVD svd,PetscBool *trackall)
866: {
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: {

906:   if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
907:   BVSetOptionsPrefix(svd->V,prefix);
908:   BVSetOptionsPrefix(svd->U,prefix);
909:   if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
910:   DSSetOptionsPrefix(svd->ds,prefix);
911:   PetscObjectSetOptionsPrefix((PetscObject)svd,prefix);
912:   return(0);
913: }

915: /*@C
916:    SVDAppendOptionsPrefix - Appends to the prefix used for searching for all
917:    SVD options in the database.

919:    Logically Collective on svd

921:    Input Parameters:
922: +  svd - the singular value solver context
923: -  prefix - the prefix string to prepend to all SVD option requests

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

929:    Level: advanced

931: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
932: @*/
933: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
934: {

939:   if (!svd->V) { SVDGetBV(svd,&svd->V,&svd->U); }
940:   BVAppendOptionsPrefix(svd->V,prefix);
941:   BVAppendOptionsPrefix(svd->U,prefix);
942:   if (!svd->ds) { SVDGetDS(svd,&svd->ds); }
943:   DSAppendOptionsPrefix(svd->ds,prefix);
944:   PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix);
945:   return(0);
946: }

948: /*@C
949:    SVDGetOptionsPrefix - Gets the prefix used for searching for all
950:    SVD options in the database.

952:    Not Collective

954:    Input Parameters:
955: .  svd - the singular value solver context

957:    Output Parameters:
958: .  prefix - pointer to the prefix string used is returned

960:    Note:
961:    On the Fortran side, the user should pass in a string 'prefix' of
962:    sufficient length to hold the prefix.

964:    Level: advanced

966: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
967: @*/
968: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
969: {

975:   PetscObjectGetOptionsPrefix((PetscObject)svd,prefix);
976:   return(0);
977: }