Actual source code: svdopts.c

slepc-3.22.1 2024-10-28
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

 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: {
 46:   PetscFunctionBegin;
 49:   if (svd->impltrans!=impl) {
 50:     svd->impltrans = impl;
 51:     svd->state     = SVD_STATE_INITIAL;
 52:   }
 53:   PetscFunctionReturn(PETSC_SUCCESS);
 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: {
 74:   PetscFunctionBegin;
 76:   PetscAssertPointer(impl,2);
 77:   *impl = svd->impltrans;
 78:   PetscFunctionReturn(PETSC_SUCCESS);
 79: }

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

 85:    Logically Collective

 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_CURRENT to retain the current value of any of the parameters.
 98:    Use PETSC_DETERMINE for either argument to assign a default value computed
 99:    internally (may be different in each solver).
100:    For maxits use PETSC_UMLIMITED to indicate there is no upper bound on this value.

102:    Level: intermediate

104: .seealso: SVDGetTolerances()
105: @*/
106: PetscErrorCode SVDSetTolerances(SVD svd,PetscReal tol,PetscInt maxits)
107: {
108:   PetscFunctionBegin;
112:   if (tol == (PetscReal)PETSC_DETERMINE) {
113:     svd->tol   = PETSC_DETERMINE;
114:     svd->state = SVD_STATE_INITIAL;
115:   } else if (tol != (PetscReal)PETSC_CURRENT) {
116:     PetscCheck(tol>0.0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
117:     svd->tol = tol;
118:   }
119:   if (maxits == PETSC_DETERMINE) {
120:     svd->max_it = PETSC_DETERMINE;
121:     svd->state  = SVD_STATE_INITIAL;
122:   } else if (maxits == PETSC_UNLIMITED) {
123:     svd->max_it = PETSC_INT_MAX;
124:   } else if (maxits != PETSC_CURRENT) {
125:     PetscCheck(maxits>0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
126:     svd->max_it = maxits;
127:   }
128:   PetscFunctionReturn(PETSC_SUCCESS);
129: }

131: /*@
132:    SVDGetTolerances - Gets the tolerance and maximum
133:    iteration count used by the default SVD convergence tests.

135:    Not Collective

137:    Input Parameter:
138: .  svd - the singular value solver context

140:    Output Parameters:
141: +  tol - the convergence tolerance
142: -  maxits - maximum number of iterations

144:    Notes:
145:    The user can specify NULL for any parameter that is not needed.

147:    Level: intermediate

149: .seealso: SVDSetTolerances()
150: @*/
151: PetscErrorCode SVDGetTolerances(SVD svd,PetscReal *tol,PetscInt *maxits)
152: {
153:   PetscFunctionBegin;
155:   if (tol)    *tol    = svd->tol;
156:   if (maxits) *maxits = svd->max_it;
157:   PetscFunctionReturn(PETSC_SUCCESS);
158: }

160: /*@
161:    SVDSetDimensions - Sets the number of singular values to compute
162:    and the dimension of the subspace.

164:    Logically Collective

166:    Input Parameters:
167: +  svd - the singular value solver context
168: .  nsv - number of singular values to compute
169: .  ncv - the maximum dimension of the subspace to be used by the solver
170: -  mpd - the maximum dimension allowed for the projected problem

172:    Options Database Keys:
173: +  -svd_nsv <nsv> - Sets the number of singular values
174: .  -svd_ncv <ncv> - Sets the dimension of the subspace
175: -  -svd_mpd <mpd> - Sets the maximum projected dimension

177:    Notes:
178:    Use PETSC_DETERMINE for ncv and mpd to assign a reasonably good value, which is
179:    dependent on the solution method and the number of singular values required. For
180:    any of the arguments, use PETSC_CURRENT to preserve the current value.

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

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

191:    Level: intermediate

193: .seealso: SVDGetDimensions()
194: @*/
195: PetscErrorCode SVDSetDimensions(SVD svd,PetscInt nsv,PetscInt ncv,PetscInt mpd)
196: {
197:   PetscFunctionBegin;
202:   PetscCheck(nsv>0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
203:   if (nsv != PETSC_CURRENT) {
204:     PetscCheck(nsv>0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nsv. Must be > 0");
205:     svd->nsv = nsv;
206:   }
207:   if (ncv == PETSC_DETERMINE) {
208:     svd->ncv = PETSC_DETERMINE;
209:   } else if (ncv != PETSC_CURRENT) {
210:     PetscCheck(ncv>0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
211:     svd->ncv = ncv;
212:   }
213:   if (mpd == PETSC_DETERMINE) {
214:     svd->mpd = PETSC_DETERMINE;
215:   } else if (mpd != PETSC_CURRENT) {
216:     PetscCheck(mpd>0,PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
217:     svd->mpd = mpd;
218:   }
219:   svd->state = SVD_STATE_INITIAL;
220:   PetscFunctionReturn(PETSC_SUCCESS);
221: }

223: /*@
224:    SVDGetDimensions - Gets the number of singular values to compute
225:    and the dimension of the subspace.

227:    Not Collective

229:    Input Parameter:
230: .  svd - the singular value context

232:    Output Parameters:
233: +  nsv - number of singular values to compute
234: .  ncv - the maximum dimension of the subspace to be used by the solver
235: -  mpd - the maximum dimension allowed for the projected problem

237:    Notes:
238:    The user can specify NULL for any parameter that is not needed.

240:    Level: intermediate

242: .seealso: SVDSetDimensions()
243: @*/
244: PetscErrorCode SVDGetDimensions(SVD svd,PetscInt *nsv,PetscInt *ncv,PetscInt *mpd)
245: {
246:   PetscFunctionBegin;
248:   if (nsv) *nsv = svd->nsv;
249:   if (ncv) *ncv = svd->ncv;
250:   if (mpd) *mpd = svd->mpd;
251:   PetscFunctionReturn(PETSC_SUCCESS);
252: }

254: /*@
255:     SVDSetWhichSingularTriplets - Specifies which singular triplets are
256:     to be sought.

258:     Logically Collective

260:     Input Parameter:
261: .   svd - singular value solver context obtained from SVDCreate()

263:     Output Parameter:
264: .   which - which singular triplets are to be sought

266:     Options Database Keys:
267: +   -svd_largest  - Sets largest singular values
268: -   -svd_smallest - Sets smallest singular values

270:     Notes:
271:     The parameter 'which' can have one of these values

273: +     SVD_LARGEST  - largest singular values
274: -     SVD_SMALLEST - smallest singular values

276:     Level: intermediate

278: .seealso: SVDGetWhichSingularTriplets(), SVDWhich
279: @*/
280: PetscErrorCode SVDSetWhichSingularTriplets(SVD svd,SVDWhich which)
281: {
282:   PetscFunctionBegin;
285:   switch (which) {
286:     case SVD_LARGEST:
287:     case SVD_SMALLEST:
288:       if (svd->which != which) {
289:         svd->state = SVD_STATE_INITIAL;
290:         svd->which = which;
291:       }
292:       break;
293:   default:
294:     SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' parameter");
295:   }
296:   PetscFunctionReturn(PETSC_SUCCESS);
297: }

299: /*@
300:     SVDGetWhichSingularTriplets - Returns which singular triplets are
301:     to be sought.

303:     Not Collective

305:     Input Parameter:
306: .   svd - singular value solver context obtained from SVDCreate()

308:     Output Parameter:
309: .   which - which singular triplets are to be sought

311:     Notes:
312:     See SVDSetWhichSingularTriplets() for possible values of which

314:     Level: intermediate

316: .seealso: SVDSetWhichSingularTriplets(), SVDWhich
317: @*/
318: PetscErrorCode SVDGetWhichSingularTriplets(SVD svd,SVDWhich *which)
319: {
320:   PetscFunctionBegin;
322:   PetscAssertPointer(which,2);
323:   *which = svd->which;
324:   PetscFunctionReturn(PETSC_SUCCESS);
325: }

327: /*@C
328:    SVDSetConvergenceTestFunction - Sets a function to compute the error estimate
329:    used in the convergence test.

331:    Logically Collective

333:    Input Parameters:
334: +  svd     - singular value solver context obtained from SVDCreate()
335: .  conv    - the convergence test function, see SVDConvergenceTestFn for the calling sequence
336: .  ctx     - context for private data for the convergence routine (may be NULL)
337: -  destroy - a routine for destroying the context (may be NULL)

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,SVDConvergenceTestFn *conv,void* ctx,PetscErrorCode (*destroy)(void*))
348: {
349:   PetscFunctionBegin;
351:   if (svd->convergeddestroy) PetscCall((*svd->convergeddestroy)(svd->convergedctx));
352:   svd->convergeduser    = conv;
353:   svd->convergeddestroy = destroy;
354:   svd->convergedctx     = ctx;
355:   if (conv == SVDConvergedAbsolute) svd->conv = SVD_CONV_ABS;
356:   else if (conv == SVDConvergedRelative) svd->conv = SVD_CONV_REL;
357:   else if (conv == SVDConvergedNorm) svd->conv = SVD_CONV_NORM;
358:   else if (conv == SVDConvergedMaxIt) svd->conv = SVD_CONV_MAXIT;
359:   else {
360:     svd->conv      = SVD_CONV_USER;
361:     svd->converged = svd->convergeduser;
362:   }
363:   PetscFunctionReturn(PETSC_SUCCESS);
364: }

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

370:    Logically Collective

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

376:    Options Database Keys:
377: +  -svd_conv_abs   - Sets the absolute convergence test
378: .  -svd_conv_rel   - Sets the convergence test relative to the singular value
379: .  -svd_conv_norm  - Sets the convergence test relative to the matrix norm
380: .  -svd_conv_maxit - Forces the maximum number of iterations as set by -svd_max_it
381: -  -svd_conv_user  - Selects the user-defined convergence test

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

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

393:    Level: intermediate

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

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

422:    Not Collective

424:    Input Parameters:
425: .  svd   - singular value solver context obtained from SVDCreate()

427:    Output Parameters:
428: .  conv  - the type of convergence test

430:    Level: intermediate

432: .seealso: SVDSetConvergenceTest(), SVDConv
433: @*/
434: PetscErrorCode SVDGetConvergenceTest(SVD svd,SVDConv *conv)
435: {
436:   PetscFunctionBegin;
438:   PetscAssertPointer(conv,2);
439:   *conv = svd->conv;
440:   PetscFunctionReturn(PETSC_SUCCESS);
441: }

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

447:    Logically Collective

449:    Input Parameters:
450: +  svd     - singular value solver context obtained from SVDCreate()
451: .  stop    - the stopping test function, see SVDStoppingTestFn for the calling sequence
452: .  ctx     - context for private data for the stopping routine (may be NULL)
453: -  destroy - a routine for destroying the context (may be NULL)

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,SVDStoppingTestFn *stop,void* ctx,PetscErrorCode (*destroy)(void*))
466: {
467:   PetscFunctionBegin;
469:   if (svd->stoppingdestroy) PetscCall((*svd->stoppingdestroy)(svd->stoppingctx));
470:   svd->stoppinguser    = stop;
471:   svd->stoppingdestroy = destroy;
472:   svd->stoppingctx     = ctx;
473:   if (stop == SVDStoppingBasic) svd->stop = SVD_STOP_BASIC;
474:   else {
475:     svd->stop     = SVD_STOP_USER;
476:     svd->stopping = svd->stoppinguser;
477:   }
478:   PetscFunctionReturn(PETSC_SUCCESS);
479: }

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

485:    Logically Collective

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

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

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

500:    Level: advanced

502: .seealso: SVDGetStoppingTest(), SVDSetStoppingTestFunction(), SVDSetConvergenceTest(), SVDStop
503: @*/
504: PetscErrorCode SVDSetStoppingTest(SVD svd,SVDStop stop)
505: {
506:   PetscFunctionBegin;
509:   switch (stop) {
510:     case SVD_STOP_BASIC: svd->stopping = SVDStoppingBasic; break;
511:     case SVD_STOP_USER:
512:       PetscCheck(svd->stoppinguser,PetscObjectComm((PetscObject)svd),PETSC_ERR_ORDER,"Must call SVDSetStoppingTestFunction() first");
513:       svd->stopping = svd->stoppinguser;
514:       break;
515:     default:
516:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
517:   }
518:   svd->stop = stop;
519:   PetscFunctionReturn(PETSC_SUCCESS);
520: }

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

526:    Not Collective

528:    Input Parameters:
529: .  svd   - singular value solver context obtained from SVDCreate()

531:    Output Parameters:
532: .  stop  - the type of stopping test

534:    Level: advanced

536: .seealso: SVDSetStoppingTest(), SVDStop
537: @*/
538: PetscErrorCode SVDGetStoppingTest(SVD svd,SVDStop *stop)
539: {
540:   PetscFunctionBegin;
542:   PetscAssertPointer(stop,2);
543:   *stop = svd->stop;
544:   PetscFunctionReturn(PETSC_SUCCESS);
545: }

547: /*@C
548:    SVDMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
549:    indicated by the user.

551:    Collective

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

560:    Level: developer

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

576:   PetscFunctionBegin;
577:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)svd),((PetscObject)svd)->options,((PetscObject)svd)->prefix,opt,&viewer,&format,&flg));
578:   if (!flg) PetscFunctionReturn(PETSC_SUCCESS);

580:   PetscCall(PetscViewerGetType(viewer,&vtype));
581:   PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
582:   PetscCall(PetscFunctionListFind(SVDMonitorList,key,&mfunc));
583:   PetscCheck(mfunc,PetscObjectComm((PetscObject)svd),PETSC_ERR_SUP,"Specified viewer and format not supported");
584:   PetscCall(PetscFunctionListFind(SVDMonitorCreateList,key,&cfunc));
585:   PetscCall(PetscFunctionListFind(SVDMonitorDestroyList,key,&dfunc));
586:   if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
587:   if (!dfunc) dfunc = PetscViewerAndFormatDestroy;

589:   PetscCall((*cfunc)(viewer,format,ctx,&vf));
590:   PetscCall(PetscViewerDestroy(&viewer));
591:   PetscCall(SVDMonitorSet(svd,mfunc,vf,(PetscErrorCode(*)(void **))dfunc));
592:   if (trackall) PetscCall(SVDSetTrackAll(svd,PETSC_TRUE));
593:   PetscFunctionReturn(PETSC_SUCCESS);
594: }

596: /*@
597:    SVDSetFromOptions - Sets SVD options from the options database.
598:    This routine must be called before SVDSetUp() if the user is to be
599:    allowed to set the solver type.

601:    Collective

603:    Input Parameters:
604: .  svd - the singular value solver context

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

609:    Level: beginner

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

620:   PetscFunctionBegin;
622:   PetscCall(SVDRegisterAll());
623:   PetscObjectOptionsBegin((PetscObject)svd);
624:     PetscCall(PetscOptionsFList("-svd_type","SVD solver method","SVDSetType",SVDList,(char*)(((PetscObject)svd)->type_name?((PetscObject)svd)->type_name:SVDCROSS),type,sizeof(type),&flg));
625:     if (flg) PetscCall(SVDSetType(svd,type));
626:     else if (!((PetscObject)svd)->type_name) PetscCall(SVDSetType(svd,SVDCROSS));

628:     PetscCall(PetscOptionsBoolGroupBegin("-svd_standard","Singular value decomposition (SVD)","SVDSetProblemType",&flg));
629:     if (flg) PetscCall(SVDSetProblemType(svd,SVD_STANDARD));
630:     PetscCall(PetscOptionsBoolGroup("-svd_generalized","Generalized singular value decomposition (GSVD)","SVDSetProblemType",&flg));
631:     if (flg) PetscCall(SVDSetProblemType(svd,SVD_GENERALIZED));
632:     PetscCall(PetscOptionsBoolGroupEnd("-svd_hyperbolic","Hyperbolic singular value decomposition (HSVD)","SVDSetProblemType",&flg));
633:     if (flg) PetscCall(SVDSetProblemType(svd,SVD_HYPERBOLIC));

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

638:     i = svd->max_it;
639:     PetscCall(PetscOptionsInt("-svd_max_it","Maximum number of iterations","SVDSetTolerances",svd->max_it,&i,&flg1));
640:     r = svd->tol;
641:     PetscCall(PetscOptionsReal("-svd_tol","Tolerance","SVDSetTolerances",SlepcDefaultTol(svd->tol),&r,&flg2));
642:     if (flg1 || flg2) PetscCall(SVDSetTolerances(svd,r,i));

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

655:     PetscCall(PetscOptionsBoolGroupBegin("-svd_stop_basic","Stop iteration if all singular values converged or max_it reached","SVDSetStoppingTest",&flg));
656:     if (flg) PetscCall(SVDSetStoppingTest(svd,SVD_STOP_BASIC));
657:     PetscCall(PetscOptionsBoolGroupEnd("-svd_stop_user","User-defined stopping test","SVDSetStoppingTest",&flg));
658:     if (flg) PetscCall(SVDSetStoppingTest(svd,SVD_STOP_USER));

660:     i = svd->nsv;
661:     PetscCall(PetscOptionsInt("-svd_nsv","Number of singular values to compute","SVDSetDimensions",svd->nsv,&i,&flg1));
662:     j = svd->ncv;
663:     PetscCall(PetscOptionsInt("-svd_ncv","Number of basis vectors","SVDSetDimensions",svd->ncv,&j,&flg2));
664:     k = svd->mpd;
665:     PetscCall(PetscOptionsInt("-svd_mpd","Maximum dimension of projected problem","SVDSetDimensions",svd->mpd,&k,&flg3));
666:     if (flg1 || flg2 || flg3) PetscCall(SVDSetDimensions(svd,i,j,k));

668:     PetscCall(PetscOptionsBoolGroupBegin("-svd_largest","Compute largest singular values","SVDSetWhichSingularTriplets",&flg));
669:     if (flg) PetscCall(SVDSetWhichSingularTriplets(svd,SVD_LARGEST));
670:     PetscCall(PetscOptionsBoolGroupEnd("-svd_smallest","Compute smallest singular values","SVDSetWhichSingularTriplets",&flg));
671:     if (flg) PetscCall(SVDSetWhichSingularTriplets(svd,SVD_SMALLEST));

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

684:     /* -----------------------------------------------------------------------*/
685:     PetscCall(PetscOptionsName("-svd_view","Print detailed information on solver used","SVDView",&set));
686:     PetscCall(PetscOptionsName("-svd_view_vectors","View computed singular vectors","SVDVectorsView",&set));
687:     PetscCall(PetscOptionsName("-svd_view_values","View computed singular values","SVDValuesView",&set));
688:     PetscCall(PetscOptionsName("-svd_converged_reason","Print reason for convergence, and number of iterations","SVDConvergedReasonView",&set));
689:     PetscCall(PetscOptionsName("-svd_error_absolute","Print absolute errors of each singular triplet","SVDErrorView",&set));
690:     PetscCall(PetscOptionsName("-svd_error_relative","Print relative errors of each singular triplet","SVDErrorView",&set));
691:     PetscCall(PetscOptionsName("-svd_error_norm","Print errors relative to the matrix norms of each singular triplet","SVDErrorView",&set));

693:     PetscTryTypeMethod(svd,setfromoptions,PetscOptionsObject);
694:     PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)svd,PetscOptionsObject));
695:   PetscOptionsEnd();

697:   if (!svd->V) PetscCall(SVDGetBV(svd,&svd->V,NULL));
698:   PetscCall(BVSetFromOptions(svd->V));
699:   if (!svd->U) PetscCall(SVDGetBV(svd,NULL,&svd->U));
700:   PetscCall(BVSetFromOptions(svd->U));
701:   if (!svd->ds) PetscCall(SVDGetDS(svd,&svd->ds));
702:   PetscCall(SVDSetDSType(svd));
703:   PetscCall(DSSetFromOptions(svd->ds));
704:   PetscFunctionReturn(PETSC_SUCCESS);
705: }

707: /*@
708:    SVDSetProblemType - Specifies the type of the singular value problem.

710:    Logically Collective

712:    Input Parameters:
713: +  svd  - the singular value solver context
714: -  type - a known type of singular value problem

716:    Options Database Keys:
717: +  -svd_standard    - standard singular value decomposition (SVD)
718: .  -svd_generalized - generalized singular value problem (GSVD)
719: -  -svd_hyperbolic  - hyperbolic singular value problem (HSVD)

721:    Notes:
722:    The GSVD requires that two matrices have been passed via SVDSetOperators().
723:    The HSVD requires that a signature matrix has been passed via SVDSetSignature().

725:    Level: intermediate

727: .seealso: SVDSetOperators(), SVDSetSignature(), SVDSetType(), SVDGetProblemType(), SVDProblemType
728: @*/
729: PetscErrorCode SVDSetProblemType(SVD svd,SVDProblemType type)
730: {
731:   PetscFunctionBegin;
734:   if (type == svd->problem_type) PetscFunctionReturn(PETSC_SUCCESS);
735:   switch (type) {
736:     case SVD_STANDARD:
737:       svd->isgeneralized = PETSC_FALSE;
738:       svd->ishyperbolic  = PETSC_FALSE;
739:       break;
740:     case SVD_GENERALIZED:
741:       svd->isgeneralized = PETSC_TRUE;
742:       svd->ishyperbolic  = PETSC_FALSE;
743:       break;
744:     case SVD_HYPERBOLIC:
745:       svd->isgeneralized = PETSC_FALSE;
746:       svd->ishyperbolic  = PETSC_TRUE;
747:       break;
748:     default:
749:       SETERRQ(PetscObjectComm((PetscObject)svd),PETSC_ERR_ARG_WRONG,"Unknown singular value problem type");
750:   }
751:   svd->problem_type = type;
752:   svd->state = SVD_STATE_INITIAL;
753:   PetscFunctionReturn(PETSC_SUCCESS);
754: }

756: /*@
757:    SVDGetProblemType - Gets the problem type from the SVD object.

759:    Not Collective

761:    Input Parameter:
762: .  svd - the singular value solver context

764:    Output Parameter:
765: .  type - the problem type

767:    Level: intermediate

769: .seealso: SVDSetProblemType(), SVDProblemType
770: @*/
771: PetscErrorCode SVDGetProblemType(SVD svd,SVDProblemType *type)
772: {
773:   PetscFunctionBegin;
775:   PetscAssertPointer(type,2);
776:   *type = svd->problem_type;
777:   PetscFunctionReturn(PETSC_SUCCESS);
778: }

780: /*@
781:    SVDIsGeneralized - Ask if the SVD object corresponds to a generalized
782:    singular value problem.

784:    Not Collective

786:    Input Parameter:
787: .  svd - the singular value solver context

789:    Output Parameter:
790: .  is - the answer

792:    Level: intermediate

794: .seealso: SVDIsHyperbolic()
795: @*/
796: PetscErrorCode SVDIsGeneralized(SVD svd,PetscBool* is)
797: {
798:   PetscFunctionBegin;
800:   PetscAssertPointer(is,2);
801:   *is = svd->isgeneralized;
802:   PetscFunctionReturn(PETSC_SUCCESS);
803: }

805: /*@
806:    SVDIsHyperbolic - Ask if the SVD object corresponds to a hyperbolic
807:    singular value problem.

809:    Not Collective

811:    Input Parameter:
812: .  svd - the singular value solver context

814:    Output Parameter:
815: .  is - the answer

817:    Level: intermediate

819: .seealso: SVDIsGeneralized()
820: @*/
821: PetscErrorCode SVDIsHyperbolic(SVD svd,PetscBool* is)
822: {
823:   PetscFunctionBegin;
825:   PetscAssertPointer(is,2);
826:   *is = svd->ishyperbolic;
827:   PetscFunctionReturn(PETSC_SUCCESS);
828: }

830: /*@
831:    SVDSetTrackAll - Specifies if the solver must compute the residual norm of all
832:    approximate singular value or not.

834:    Logically Collective

836:    Input Parameters:
837: +  svd      - the singular value solver context
838: -  trackall - whether to compute all residuals or not

840:    Notes:
841:    If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
842:    the residual norm for each singular value approximation. Computing the residual is
843:    usually an expensive operation and solvers commonly compute only the residual
844:    associated to the first unconverged singular value.

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

848:    Level: developer

850: .seealso: SVDGetTrackAll()
851: @*/
852: PetscErrorCode SVDSetTrackAll(SVD svd,PetscBool trackall)
853: {
854:   PetscFunctionBegin;
857:   svd->trackall = trackall;
858:   PetscFunctionReturn(PETSC_SUCCESS);
859: }

861: /*@
862:    SVDGetTrackAll - Returns the flag indicating whether all residual norms must
863:    be computed or not.

865:    Not Collective

867:    Input Parameter:
868: .  svd - the singular value solver context

870:    Output Parameter:
871: .  trackall - the returned flag

873:    Level: developer

875: .seealso: SVDSetTrackAll()
876: @*/
877: PetscErrorCode SVDGetTrackAll(SVD svd,PetscBool *trackall)
878: {
879:   PetscFunctionBegin;
881:   PetscAssertPointer(trackall,2);
882:   *trackall = svd->trackall;
883:   PetscFunctionReturn(PETSC_SUCCESS);
884: }

886: /*@
887:    SVDSetOptionsPrefix - Sets the prefix used for searching for all
888:    SVD options in the database.

890:    Logically Collective

892:    Input Parameters:
893: +  svd - the singular value solver context
894: -  prefix - the prefix string to prepend to all SVD option requests

896:    Notes:
897:    A hyphen (-) must NOT be given at the beginning of the prefix name.
898:    The first character of all runtime options is AUTOMATICALLY the
899:    hyphen.

901:    For example, to distinguish between the runtime options for two
902:    different SVD contexts, one could call
903: .vb
904:       SVDSetOptionsPrefix(svd1,"svd1_")
905:       SVDSetOptionsPrefix(svd2,"svd2_")
906: .ve

908:    Level: advanced

910: .seealso: SVDAppendOptionsPrefix(), SVDGetOptionsPrefix()
911: @*/
912: PetscErrorCode SVDSetOptionsPrefix(SVD svd,const char *prefix)
913: {
914:   PetscFunctionBegin;
916:   if (!svd->V) PetscCall(SVDGetBV(svd,&svd->V,&svd->U));
917:   PetscCall(BVSetOptionsPrefix(svd->V,prefix));
918:   PetscCall(BVSetOptionsPrefix(svd->U,prefix));
919:   if (!svd->ds) PetscCall(SVDGetDS(svd,&svd->ds));
920:   PetscCall(DSSetOptionsPrefix(svd->ds,prefix));
921:   PetscCall(PetscObjectSetOptionsPrefix((PetscObject)svd,prefix));
922:   PetscFunctionReturn(PETSC_SUCCESS);
923: }

925: /*@
926:    SVDAppendOptionsPrefix - Appends to the prefix used for searching for all
927:    SVD options in the database.

929:    Logically Collective

931:    Input Parameters:
932: +  svd - the singular value solver context
933: -  prefix - the prefix string to prepend to all SVD option requests

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

939:    Level: advanced

941: .seealso: SVDSetOptionsPrefix(), SVDGetOptionsPrefix()
942: @*/
943: PetscErrorCode SVDAppendOptionsPrefix(SVD svd,const char *prefix)
944: {
945:   PetscFunctionBegin;
947:   if (!svd->V) PetscCall(SVDGetBV(svd,&svd->V,&svd->U));
948:   PetscCall(BVAppendOptionsPrefix(svd->V,prefix));
949:   PetscCall(BVAppendOptionsPrefix(svd->U,prefix));
950:   if (!svd->ds) PetscCall(SVDGetDS(svd,&svd->ds));
951:   PetscCall(DSAppendOptionsPrefix(svd->ds,prefix));
952:   PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)svd,prefix));
953:   PetscFunctionReturn(PETSC_SUCCESS);
954: }

956: /*@
957:    SVDGetOptionsPrefix - Gets the prefix used for searching for all
958:    SVD options in the database.

960:    Not Collective

962:    Input Parameters:
963: .  svd - the singular value solver context

965:    Output Parameters:
966: .  prefix - pointer to the prefix string used is returned

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

972:    Level: advanced

974: .seealso: SVDSetOptionsPrefix(), SVDAppendOptionsPrefix()
975: @*/
976: PetscErrorCode SVDGetOptionsPrefix(SVD svd,const char *prefix[])
977: {
978:   PetscFunctionBegin;
980:   PetscAssertPointer(prefix,2);
981:   PetscCall(PetscObjectGetOptionsPrefix((PetscObject)svd,prefix));
982:   PetscFunctionReturn(PETSC_SUCCESS);
983: }