Actual source code: pepbasic.c

slepc-3.16.1 2021-11-17
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:    Basic PEP routines
 12: */

 14: #include <slepc/private/pepimpl.h>

 16: /* Logging support */
 17: PetscClassId      PEP_CLASSID = 0;
 18: PetscLogEvent     PEP_SetUp = 0,PEP_Solve = 0,PEP_Refine = 0,PEP_CISS_SVD = 0;

 20: /* List of registered PEP routines */
 21: PetscFunctionList PEPList = 0;
 22: PetscBool         PEPRegisterAllCalled = PETSC_FALSE;

 24: /* List of registered PEP monitors */
 25: PetscFunctionList PEPMonitorList              = NULL;
 26: PetscFunctionList PEPMonitorCreateList        = NULL;
 27: PetscFunctionList PEPMonitorDestroyList       = NULL;
 28: PetscBool         PEPMonitorRegisterAllCalled = PETSC_FALSE;

 30: /*@
 31:    PEPCreate - Creates the default PEP context.

 33:    Collective

 35:    Input Parameter:
 36: .  comm - MPI communicator

 38:    Output Parameter:
 39: .  pep - location to put the PEP context

 41:    Note:
 42:    The default PEP type is PEPTOAR

 44:    Level: beginner

 46: .seealso: PEPSetUp(), PEPSolve(), PEPDestroy(), PEP
 47: @*/
 48: PetscErrorCode PEPCreate(MPI_Comm comm,PEP *outpep)
 49: {
 51:   PEP            pep;

 55:   *outpep = 0;
 56:   PEPInitializePackage();
 57:   SlepcHeaderCreate(pep,PEP_CLASSID,"PEP","Polynomial Eigenvalue Problem","PEP",comm,PEPDestroy,PEPView);

 59:   pep->max_it          = PETSC_DEFAULT;
 60:   pep->nev             = 1;
 61:   pep->ncv             = PETSC_DEFAULT;
 62:   pep->mpd             = PETSC_DEFAULT;
 63:   pep->nini            = 0;
 64:   pep->target          = 0.0;
 65:   pep->tol             = PETSC_DEFAULT;
 66:   pep->conv            = PEP_CONV_REL;
 67:   pep->stop            = PEP_STOP_BASIC;
 68:   pep->which           = (PEPWhich)0;
 69:   pep->basis           = PEP_BASIS_MONOMIAL;
 70:   pep->problem_type    = (PEPProblemType)0;
 71:   pep->scale           = PEP_SCALE_NONE;
 72:   pep->sfactor         = 1.0;
 73:   pep->dsfactor        = 1.0;
 74:   pep->sits            = 5;
 75:   pep->slambda         = 1.0;
 76:   pep->refine          = PEP_REFINE_NONE;
 77:   pep->npart           = 1;
 78:   pep->rtol            = PETSC_DEFAULT;
 79:   pep->rits            = PETSC_DEFAULT;
 80:   pep->scheme          = (PEPRefineScheme)0;
 81:   pep->extract         = (PEPExtract)0;
 82:   pep->trackall        = PETSC_FALSE;

 84:   pep->converged       = PEPConvergedRelative;
 85:   pep->convergeduser   = NULL;
 86:   pep->convergeddestroy= NULL;
 87:   pep->stopping        = PEPStoppingBasic;
 88:   pep->stoppinguser    = NULL;
 89:   pep->stoppingdestroy = NULL;
 90:   pep->convergedctx    = NULL;
 91:   pep->stoppingctx     = NULL;
 92:   pep->numbermonitors  = 0;

 94:   pep->st              = NULL;
 95:   pep->ds              = NULL;
 96:   pep->V               = NULL;
 97:   pep->rg              = NULL;
 98:   pep->A               = NULL;
 99:   pep->nmat            = 0;
100:   pep->Dl              = NULL;
101:   pep->Dr              = NULL;
102:   pep->IS              = NULL;
103:   pep->eigr            = NULL;
104:   pep->eigi            = NULL;
105:   pep->errest          = NULL;
106:   pep->perm            = NULL;
107:   pep->pbc             = NULL;
108:   pep->solvematcoeffs  = NULL;
109:   pep->nwork           = 0;
110:   pep->work            = NULL;
111:   pep->refineksp       = NULL;
112:   pep->refinesubc      = NULL;
113:   pep->data            = NULL;

115:   pep->state           = PEP_STATE_INITIAL;
116:   pep->nconv           = 0;
117:   pep->its             = 0;
118:   pep->n               = 0;
119:   pep->nloc            = 0;
120:   pep->nrma            = NULL;
121:   pep->sfactor_set     = PETSC_FALSE;
122:   pep->lineariz        = PETSC_FALSE;
123:   pep->reason          = PEP_CONVERGED_ITERATING;

125:   PetscNewLog(pep,&pep->sc);
126:   *outpep = pep;
127:   return(0);
128: }

130: /*@C
131:    PEPSetType - Selects the particular solver to be used in the PEP object.

133:    Logically Collective on pep

135:    Input Parameters:
136: +  pep      - the polynomial eigensolver context
137: -  type     - a known method

139:    Options Database Key:
140: .  -pep_type <method> - Sets the method; use -help for a list
141:     of available methods

143:    Notes:
144:    See "slepc/include/slepcpep.h" for available methods. The default
145:    is PEPTOAR.

147:    Normally, it is best to use the PEPSetFromOptions() command and
148:    then set the PEP type from the options database rather than by using
149:    this routine.  Using the options database provides the user with
150:    maximum flexibility in evaluating the different available methods.
151:    The PEPSetType() routine is provided for those situations where it
152:    is necessary to set the iterative solver independently of the command
153:    line or options database.

155:    Level: intermediate

157: .seealso: PEPType
158: @*/
159: PetscErrorCode PEPSetType(PEP pep,PEPType type)
160: {
161:   PetscErrorCode ierr,(*r)(PEP);
162:   PetscBool      match;


168:   PetscObjectTypeCompare((PetscObject)pep,type,&match);
169:   if (match) return(0);

171:   PetscFunctionListFind(PEPList,type,&r);
172:   if (!r) SETERRQ1(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown PEP type given: %s",type);

174:   if (pep->ops->destroy) { (*pep->ops->destroy)(pep); }
175:   PetscMemzero(pep->ops,sizeof(struct _PEPOps));

177:   pep->state = PEP_STATE_INITIAL;
178:   PetscObjectChangeTypeName((PetscObject)pep,type);
179:   (*r)(pep);
180:   return(0);
181: }

183: /*@C
184:    PEPGetType - Gets the PEP type as a string from the PEP object.

186:    Not Collective

188:    Input Parameter:
189: .  pep - the eigensolver context

191:    Output Parameter:
192: .  name - name of PEP method

194:    Level: intermediate

196: .seealso: PEPSetType()
197: @*/
198: PetscErrorCode PEPGetType(PEP pep,PEPType *type)
199: {
203:   *type = ((PetscObject)pep)->type_name;
204:   return(0);
205: }

207: /*@C
208:    PEPRegister - Adds a method to the polynomial eigenproblem solver package.

210:    Not Collective

212:    Input Parameters:
213: +  name - name of a new user-defined solver
214: -  function - routine to create the solver context

216:    Notes:
217:    PEPRegister() may be called multiple times to add several user-defined solvers.

219:    Sample usage:
220: .vb
221:     PEPRegister("my_solver",MySolverCreate);
222: .ve

224:    Then, your solver can be chosen with the procedural interface via
225: $     PEPSetType(pep,"my_solver")
226:    or at runtime via the option
227: $     -pep_type my_solver

229:    Level: advanced

231: .seealso: PEPRegisterAll()
232: @*/
233: PetscErrorCode PEPRegister(const char *name,PetscErrorCode (*function)(PEP))
234: {

238:   PEPInitializePackage();
239:   PetscFunctionListAdd(&PEPList,name,function);
240:   return(0);
241: }

243: /*@C
244:    PEPMonitorRegister - Adds PEP monitor routine.

246:    Not Collective

248:    Input Parameters:
249: +  name    - name of a new monitor routine
250: .  vtype   - a PetscViewerType for the output
251: .  format  - a PetscViewerFormat for the output
252: .  monitor - monitor routine
253: .  create  - creation routine, or NULL
254: -  destroy - destruction routine, or NULL

256:    Notes:
257:    PEPMonitorRegister() may be called multiple times to add several user-defined monitors.

259:    Sample usage:
260: .vb
261:    PEPMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
262: .ve

264:    Then, your monitor can be chosen with the procedural interface via
265: $      PEPMonitorSetFromOptions(pep,"-pep_monitor_my_monitor","my_monitor",NULL)
266:    or at runtime via the option
267: $      -pep_monitor_my_monitor

269:    Level: advanced

271: .seealso: PEPMonitorRegisterAll()
272: @*/
273: PetscErrorCode PEPMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,PetscErrorCode (*monitor)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,PetscViewerAndFormat*),PetscErrorCode (*create)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**),PetscErrorCode (*destroy)(PetscViewerAndFormat**))
274: {
275:   char           key[PETSC_MAX_PATH_LEN];

279:   PEPInitializePackage();
280:   SlepcMonitorMakeKey_Internal(name,vtype,format,key);
281:   PetscFunctionListAdd(&PEPMonitorList,key,monitor);
282:   if (create)  { PetscFunctionListAdd(&PEPMonitorCreateList,key,create); }
283:   if (destroy) { PetscFunctionListAdd(&PEPMonitorDestroyList,key,destroy); }
284:   return(0);
285: }

287: /*@
288:    PEPReset - Resets the PEP context to the initial state (prior to setup)
289:    and destroys any allocated Vecs and Mats.

291:    Collective on pep

293:    Input Parameter:
294: .  pep - eigensolver context obtained from PEPCreate()

296:    Level: advanced

298: .seealso: PEPDestroy()
299: @*/
300: PetscErrorCode PEPReset(PEP pep)
301: {

306:   if (!pep) return(0);
307:   if (pep->ops->reset) { (pep->ops->reset)(pep); }
308:   if (pep->st) { STReset(pep->st); }
309:   if (pep->refineksp) { KSPReset(pep->refineksp); }
310:   if (pep->nmat) {
311:     MatDestroyMatrices(pep->nmat,&pep->A);
312:     PetscFree2(pep->pbc,pep->nrma);
313:     PetscFree(pep->solvematcoeffs);
314:     pep->nmat = 0;
315:   }
316:   VecDestroy(&pep->Dl);
317:   VecDestroy(&pep->Dr);
318:   BVDestroy(&pep->V);
319:   VecDestroyVecs(pep->nwork,&pep->work);
320:   pep->nwork = 0;
321:   pep->state = PEP_STATE_INITIAL;
322:   return(0);
323: }

325: /*@C
326:    PEPDestroy - Destroys the PEP context.

328:    Collective on pep

330:    Input Parameter:
331: .  pep - eigensolver context obtained from PEPCreate()

333:    Level: beginner

335: .seealso: PEPCreate(), PEPSetUp(), PEPSolve()
336: @*/
337: PetscErrorCode PEPDestroy(PEP *pep)
338: {

342:   if (!*pep) return(0);
344:   if (--((PetscObject)(*pep))->refct > 0) { *pep = 0; return(0); }
345:   PEPReset(*pep);
346:   if ((*pep)->ops->destroy) { (*(*pep)->ops->destroy)(*pep); }
347:   if ((*pep)->eigr) {
348:     PetscFree4((*pep)->eigr,(*pep)->eigi,(*pep)->errest,(*pep)->perm);
349:   }
350:   STDestroy(&(*pep)->st);
351:   RGDestroy(&(*pep)->rg);
352:   DSDestroy(&(*pep)->ds);
353:   KSPDestroy(&(*pep)->refineksp);
354:   PetscSubcommDestroy(&(*pep)->refinesubc);
355:   PetscFree((*pep)->sc);
356:   /* just in case the initial vectors have not been used */
357:   SlepcBasisDestroy_Private(&(*pep)->nini,&(*pep)->IS);
358:   if ((*pep)->convergeddestroy) {
359:     (*(*pep)->convergeddestroy)((*pep)->convergedctx);
360:   }
361:   PEPMonitorCancel(*pep);
362:   PetscHeaderDestroy(pep);
363:   return(0);
364: }

366: /*@
367:    PEPSetBV - Associates a basis vectors object to the polynomial eigensolver.

369:    Collective on pep

371:    Input Parameters:
372: +  pep - eigensolver context obtained from PEPCreate()
373: -  bv  - the basis vectors object

375:    Note:
376:    Use PEPGetBV() to retrieve the basis vectors context (for example,
377:    to free it at the end of the computations).

379:    Level: advanced

381: .seealso: PEPGetBV()
382: @*/
383: PetscErrorCode PEPSetBV(PEP pep,BV bv)
384: {

391:   PetscObjectReference((PetscObject)bv);
392:   BVDestroy(&pep->V);
393:   pep->V = bv;
394:   PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->V);
395:   return(0);
396: }

398: /*@
399:    PEPGetBV - Obtain the basis vectors object associated to the polynomial
400:    eigensolver object.

402:    Not Collective

404:    Input Parameters:
405: .  pep - eigensolver context obtained from PEPCreate()

407:    Output Parameter:
408: .  bv - basis vectors context

410:    Level: advanced

412: .seealso: PEPSetBV()
413: @*/
414: PetscErrorCode PEPGetBV(PEP pep,BV *bv)
415: {

421:   if (!pep->V) {
422:     BVCreate(PetscObjectComm((PetscObject)pep),&pep->V);
423:     PetscObjectIncrementTabLevel((PetscObject)pep->V,(PetscObject)pep,0);
424:     PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->V);
425:     PetscObjectSetOptions((PetscObject)pep->V,((PetscObject)pep)->options);
426:   }
427:   *bv = pep->V;
428:   return(0);
429: }

431: /*@
432:    PEPSetRG - Associates a region object to the polynomial eigensolver.

434:    Collective on pep

436:    Input Parameters:
437: +  pep - eigensolver context obtained from PEPCreate()
438: -  rg  - the region object

440:    Note:
441:    Use PEPGetRG() to retrieve the region context (for example,
442:    to free it at the end of the computations).

444:    Level: advanced

446: .seealso: PEPGetRG()
447: @*/
448: PetscErrorCode PEPSetRG(PEP pep,RG rg)
449: {

454:   if (rg) {
457:   }
458:   PetscObjectReference((PetscObject)rg);
459:   RGDestroy(&pep->rg);
460:   pep->rg = rg;
461:   PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->rg);
462:   return(0);
463: }

465: /*@
466:    PEPGetRG - Obtain the region object associated to the
467:    polynomial eigensolver object.

469:    Not Collective

471:    Input Parameters:
472: .  pep - eigensolver context obtained from PEPCreate()

474:    Output Parameter:
475: .  rg - region context

477:    Level: advanced

479: .seealso: PEPSetRG()
480: @*/
481: PetscErrorCode PEPGetRG(PEP pep,RG *rg)
482: {

488:   if (!pep->rg) {
489:     RGCreate(PetscObjectComm((PetscObject)pep),&pep->rg);
490:     PetscObjectIncrementTabLevel((PetscObject)pep->rg,(PetscObject)pep,0);
491:     PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->rg);
492:     PetscObjectSetOptions((PetscObject)pep->rg,((PetscObject)pep)->options);
493:   }
494:   *rg = pep->rg;
495:   return(0);
496: }

498: /*@
499:    PEPSetDS - Associates a direct solver object to the polynomial eigensolver.

501:    Collective on pep

503:    Input Parameters:
504: +  pep - eigensolver context obtained from PEPCreate()
505: -  ds  - the direct solver object

507:    Note:
508:    Use PEPGetDS() to retrieve the direct solver context (for example,
509:    to free it at the end of the computations).

511:    Level: advanced

513: .seealso: PEPGetDS()
514: @*/
515: PetscErrorCode PEPSetDS(PEP pep,DS ds)
516: {

523:   PetscObjectReference((PetscObject)ds);
524:   DSDestroy(&pep->ds);
525:   pep->ds = ds;
526:   PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->ds);
527:   return(0);
528: }

530: /*@
531:    PEPGetDS - Obtain the direct solver object associated to the
532:    polynomial eigensolver object.

534:    Not Collective

536:    Input Parameters:
537: .  pep - eigensolver context obtained from PEPCreate()

539:    Output Parameter:
540: .  ds - direct solver context

542:    Level: advanced

544: .seealso: PEPSetDS()
545: @*/
546: PetscErrorCode PEPGetDS(PEP pep,DS *ds)
547: {

553:   if (!pep->ds) {
554:     DSCreate(PetscObjectComm((PetscObject)pep),&pep->ds);
555:     PetscObjectIncrementTabLevel((PetscObject)pep->ds,(PetscObject)pep,0);
556:     PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->ds);
557:     PetscObjectSetOptions((PetscObject)pep->ds,((PetscObject)pep)->options);
558:   }
559:   *ds = pep->ds;
560:   return(0);
561: }

563: /*@
564:    PEPSetST - Associates a spectral transformation object to the eigensolver.

566:    Collective on pep

568:    Input Parameters:
569: +  pep - eigensolver context obtained from PEPCreate()
570: -  st   - the spectral transformation object

572:    Note:
573:    Use PEPGetST() to retrieve the spectral transformation context (for example,
574:    to free it at the end of the computations).

576:    Level: advanced

578: .seealso: PEPGetST()
579: @*/
580: PetscErrorCode PEPSetST(PEP pep,ST st)
581: {

588:   PetscObjectReference((PetscObject)st);
589:   STDestroy(&pep->st);
590:   pep->st = st;
591:   PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->st);
592:   return(0);
593: }

595: /*@
596:    PEPGetST - Obtain the spectral transformation (ST) object associated
597:    to the eigensolver object.

599:    Not Collective

601:    Input Parameters:
602: .  pep - eigensolver context obtained from PEPCreate()

604:    Output Parameter:
605: .  st - spectral transformation context

607:    Level: intermediate

609: .seealso: PEPSetST()
610: @*/
611: PetscErrorCode PEPGetST(PEP pep,ST *st)
612: {

618:   if (!pep->st) {
619:     STCreate(PetscObjectComm((PetscObject)pep),&pep->st);
620:     PetscObjectIncrementTabLevel((PetscObject)pep->st,(PetscObject)pep,0);
621:     PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->st);
622:     PetscObjectSetOptions((PetscObject)pep->st,((PetscObject)pep)->options);
623:   }
624:   *st = pep->st;
625:   return(0);
626: }

628: /*@
629:    PEPRefineGetKSP - Obtain the ksp object used by the eigensolver
630:    object in the refinement phase.

632:    Not Collective

634:    Input Parameters:
635: .  pep - eigensolver context obtained from PEPCreate()

637:    Output Parameter:
638: .  ksp - ksp context

640:    Level: advanced

642: .seealso: PEPSetRefine()
643: @*/
644: PetscErrorCode PEPRefineGetKSP(PEP pep,KSP *ksp)
645: {

651:   if (!pep->refineksp) {
652:     if (pep->npart>1) {
653:       /* Split in subcomunicators */
654:       PetscSubcommCreate(PetscObjectComm((PetscObject)pep),&pep->refinesubc);
655:       PetscSubcommSetNumber(pep->refinesubc,pep->npart);
656:       PetscSubcommSetType(pep->refinesubc,PETSC_SUBCOMM_CONTIGUOUS);
657:       PetscLogObjectMemory((PetscObject)pep,sizeof(PetscSubcomm));
658:     }
659:     KSPCreate((pep->npart==1)?PetscObjectComm((PetscObject)pep):PetscSubcommChild(pep->refinesubc),&pep->refineksp);
660:     PetscObjectIncrementTabLevel((PetscObject)pep->refineksp,(PetscObject)pep,0);
661:     PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->refineksp);
662:     PetscObjectSetOptions((PetscObject)pep->refineksp,((PetscObject)pep)->options);
663:     KSPSetOptionsPrefix(*ksp,((PetscObject)pep)->prefix);
664:     KSPAppendOptionsPrefix(*ksp,"pep_refine_");
665:     KSPSetTolerances(pep->refineksp,SlepcDefaultTol(pep->rtol),PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
666:   }
667:   *ksp = pep->refineksp;
668:   return(0);
669: }

671: /*@
672:    PEPSetTarget - Sets the value of the target.

674:    Logically Collective on pep

676:    Input Parameters:
677: +  pep    - eigensolver context
678: -  target - the value of the target

680:    Options Database Key:
681: .  -pep_target <scalar> - the value of the target

683:    Notes:
684:    The target is a scalar value used to determine the portion of the spectrum
685:    of interest. It is used in combination with PEPSetWhichEigenpairs().

687:    In the case of complex scalars, a complex value can be provided in the
688:    command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
689:    -pep_target 1.0+2.0i

691:    Level: intermediate

693: .seealso: PEPGetTarget(), PEPSetWhichEigenpairs()
694: @*/
695: PetscErrorCode PEPSetTarget(PEP pep,PetscScalar target)
696: {

702:   pep->target = target;
703:   if (!pep->st) { PEPGetST(pep,&pep->st); }
704:   STSetDefaultShift(pep->st,target);
705:   return(0);
706: }

708: /*@
709:    PEPGetTarget - Gets the value of the target.

711:    Not Collective

713:    Input Parameter:
714: .  pep - eigensolver context

716:    Output Parameter:
717: .  target - the value of the target

719:    Note:
720:    If the target was not set by the user, then zero is returned.

722:    Level: intermediate

724: .seealso: PEPSetTarget()
725: @*/
726: PetscErrorCode PEPGetTarget(PEP pep,PetscScalar* target)
727: {
731:   *target = pep->target;
732:   return(0);
733: }

735: /*@
736:    PEPSetInterval - Defines the computational interval for spectrum slicing.

738:    Logically Collective on pep

740:    Input Parameters:
741: +  pep  - eigensolver context
742: .  inta - left end of the interval
743: -  intb - right end of the interval

745:    Options Database Key:
746: .  -pep_interval <a,b> - set [a,b] as the interval of interest

748:    Notes:
749:    Spectrum slicing is a technique employed for computing all eigenvalues of
750:    symmetric eigenproblems in a given interval. This function provides the
751:    interval to be considered. It must be used in combination with PEP_ALL, see
752:    PEPSetWhichEigenpairs().

754:    In the command-line option, two values must be provided. For an open interval,
755:    one can give an infinite, e.g., -pep_interval 1.0,inf or -pep_interval -inf,1.0.
756:    An open interval in the programmatic interface can be specified with
757:    PETSC_MAX_REAL and -PETSC_MAX_REAL.

759:    Level: intermediate

761: .seealso: PEPGetInterval(), PEPSetWhichEigenpairs()
762: @*/
763: PetscErrorCode PEPSetInterval(PEP pep,PetscReal inta,PetscReal intb)
764: {
769:   if (inta>=intb) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
770:   if (pep->inta != inta || pep->intb != intb) {
771:     pep->inta = inta;
772:     pep->intb = intb;
773:     pep->state = PEP_STATE_INITIAL;
774:   }
775:   return(0);
776: }

778: /*@
779:    PEPGetInterval - Gets the computational interval for spectrum slicing.

781:    Not Collective

783:    Input Parameter:
784: .  pep - eigensolver context

786:    Output Parameters:
787: +  inta - left end of the interval
788: -  intb - right end of the interval

790:    Level: intermediate

792:    Note:
793:    If the interval was not set by the user, then zeros are returned.

795: .seealso: PEPSetInterval()
796: @*/
797: PetscErrorCode PEPGetInterval(PEP pep,PetscReal* inta,PetscReal* intb)
798: {
801:   if (inta) *inta = pep->inta;
802:   if (intb) *intb = pep->intb;
803:   return(0);
804: }