Actual source code: epsbasic.c

slepc-3.17.2 2022-08-09
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:    Basic EPS routines
 12: */

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

 16: /* Logging support */
 17: PetscClassId      EPS_CLASSID = 0;
 18: PetscLogEvent     EPS_SetUp = 0,EPS_Solve = 0,EPS_CISS_SVD = 0;

 20: /* List of registered EPS routines */
 21: PetscFunctionList EPSList = 0;
 22: PetscBool         EPSRegisterAllCalled = PETSC_FALSE;

 24: /* List of registered EPS monitors */
 25: PetscFunctionList EPSMonitorList              = NULL;
 26: PetscFunctionList EPSMonitorCreateList        = NULL;
 27: PetscFunctionList EPSMonitorDestroyList       = NULL;
 28: PetscBool         EPSMonitorRegisterAllCalled = PETSC_FALSE;

 30: /*@
 31:    EPSCreate - Creates the default EPS context.

 33:    Collective

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

 38:    Output Parameter:
 39: .  outeps - location to put the EPS context

 41:    Note:
 42:    The default EPS type is EPSKRYLOVSCHUR

 44:    Level: beginner

 46: .seealso: EPSSetUp(), EPSSolve(), EPSDestroy(), EPS
 47: @*/
 48: PetscErrorCode EPSCreate(MPI_Comm comm,EPS *outeps)
 49: {
 50:   EPS            eps;

 53:   *outeps = 0;
 54:   EPSInitializePackage();
 55:   SlepcHeaderCreate(eps,EPS_CLASSID,"EPS","Eigenvalue Problem Solver","EPS",comm,EPSDestroy,EPSView);

 57:   eps->max_it          = PETSC_DEFAULT;
 58:   eps->nev             = 1;
 59:   eps->ncv             = PETSC_DEFAULT;
 60:   eps->mpd             = PETSC_DEFAULT;
 61:   eps->nini            = 0;
 62:   eps->nds             = 0;
 63:   eps->target          = 0.0;
 64:   eps->tol             = PETSC_DEFAULT;
 65:   eps->conv            = EPS_CONV_REL;
 66:   eps->stop            = EPS_STOP_BASIC;
 67:   eps->which           = (EPSWhich)0;
 68:   eps->inta            = 0.0;
 69:   eps->intb            = 0.0;
 70:   eps->problem_type    = (EPSProblemType)0;
 71:   eps->extraction      = EPS_RITZ;
 72:   eps->balance         = EPS_BALANCE_NONE;
 73:   eps->balance_its     = 5;
 74:   eps->balance_cutoff  = 1e-8;
 75:   eps->trueres         = PETSC_FALSE;
 76:   eps->trackall        = PETSC_FALSE;
 77:   eps->purify          = PETSC_TRUE;
 78:   eps->twosided        = PETSC_FALSE;

 80:   eps->converged       = EPSConvergedRelative;
 81:   eps->convergeduser   = NULL;
 82:   eps->convergeddestroy= NULL;
 83:   eps->stopping        = EPSStoppingBasic;
 84:   eps->stoppinguser    = NULL;
 85:   eps->stoppingdestroy = NULL;
 86:   eps->arbitrary       = NULL;
 87:   eps->convergedctx    = NULL;
 88:   eps->stoppingctx     = NULL;
 89:   eps->arbitraryctx    = NULL;
 90:   eps->numbermonitors  = 0;

 92:   eps->st              = NULL;
 93:   eps->ds              = NULL;
 94:   eps->V               = NULL;
 95:   eps->W               = NULL;
 96:   eps->rg              = NULL;
 97:   eps->D               = NULL;
 98:   eps->IS              = NULL;
 99:   eps->ISL             = NULL;
100:   eps->defl            = NULL;
101:   eps->eigr            = NULL;
102:   eps->eigi            = NULL;
103:   eps->errest          = NULL;
104:   eps->rr              = NULL;
105:   eps->ri              = NULL;
106:   eps->perm            = NULL;
107:   eps->nwork           = 0;
108:   eps->work            = NULL;
109:   eps->data            = NULL;

111:   eps->state           = EPS_STATE_INITIAL;
112:   eps->categ           = EPS_CATEGORY_KRYLOV;
113:   eps->nconv           = 0;
114:   eps->its             = 0;
115:   eps->nloc            = 0;
116:   eps->nrma            = 0.0;
117:   eps->nrmb            = 0.0;
118:   eps->useds           = PETSC_FALSE;
119:   eps->isgeneralized   = PETSC_FALSE;
120:   eps->ispositive      = PETSC_FALSE;
121:   eps->ishermitian     = PETSC_FALSE;
122:   eps->reason          = EPS_CONVERGED_ITERATING;

124:   PetscNewLog(eps,&eps->sc);
125:   *outeps = eps;
126:   PetscFunctionReturn(0);
127: }

129: /*@C
130:    EPSSetType - Selects the particular solver to be used in the EPS object.

132:    Logically Collective on eps

134:    Input Parameters:
135: +  eps  - the eigensolver context
136: -  type - a known method

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

142:    Notes:
143:    See "slepc/include/slepceps.h" for available methods. The default
144:    is EPSKRYLOVSCHUR.

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

154:    Level: intermediate

156: .seealso: STSetType(), EPSType
157: @*/
158: PetscErrorCode EPSSetType(EPS eps,EPSType type)
159: {
160:   PetscErrorCode (*r)(EPS);
161:   PetscBool      match;


166:   PetscObjectTypeCompare((PetscObject)eps,type,&match);
167:   if (match) PetscFunctionReturn(0);

169:   PetscFunctionListFind(EPSList,type,&r);

172:   if (eps->ops->destroy) (*eps->ops->destroy)(eps);
173:   PetscMemzero(eps->ops,sizeof(struct _EPSOps));

175:   eps->state = EPS_STATE_INITIAL;
176:   PetscObjectChangeTypeName((PetscObject)eps,type);
177:   (*r)(eps);
178:   PetscFunctionReturn(0);
179: }

181: /*@C
182:    EPSGetType - Gets the EPS type as a string from the EPS object.

184:    Not Collective

186:    Input Parameter:
187: .  eps - the eigensolver context

189:    Output Parameter:
190: .  type - name of EPS method

192:    Level: intermediate

194: .seealso: EPSSetType()
195: @*/
196: PetscErrorCode EPSGetType(EPS eps,EPSType *type)
197: {
200:   *type = ((PetscObject)eps)->type_name;
201:   PetscFunctionReturn(0);
202: }

204: /*@C
205:    EPSRegister - Adds a method to the eigenproblem solver package.

207:    Not Collective

209:    Input Parameters:
210: +  name - name of a new user-defined solver
211: -  function - routine to create the solver context

213:    Notes:
214:    EPSRegister() may be called multiple times to add several user-defined solvers.

216:    Sample usage:
217: .vb
218:     EPSRegister("my_solver",MySolverCreate);
219: .ve

221:    Then, your solver can be chosen with the procedural interface via
222: $     EPSSetType(eps,"my_solver")
223:    or at runtime via the option
224: $     -eps_type my_solver

226:    Level: advanced

228: .seealso: EPSRegisterAll()
229: @*/
230: PetscErrorCode EPSRegister(const char *name,PetscErrorCode (*function)(EPS))
231: {
232:   EPSInitializePackage();
233:   PetscFunctionListAdd(&EPSList,name,function);
234:   PetscFunctionReturn(0);
235: }

237: /*@C
238:    EPSMonitorRegister - Adds EPS monitor routine.

240:    Not Collective

242:    Input Parameters:
243: +  name    - name of a new monitor routine
244: .  vtype   - a PetscViewerType for the output
245: .  format  - a PetscViewerFormat for the output
246: .  monitor - monitor routine
247: .  create  - creation routine, or NULL
248: -  destroy - destruction routine, or NULL

250:    Notes:
251:    EPSMonitorRegister() may be called multiple times to add several user-defined monitors.

253:    Sample usage:
254: .vb
255:    EPSMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
256: .ve

258:    Then, your monitor can be chosen with the procedural interface via
259: $      EPSMonitorSetFromOptions(eps,"-eps_monitor_my_monitor","my_monitor",NULL)
260:    or at runtime via the option
261: $      -eps_monitor_my_monitor

263:    Level: advanced

265: .seealso: EPSMonitorRegisterAll()
266: @*/
267: PetscErrorCode EPSMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,PetscErrorCode (*monitor)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,PetscViewerAndFormat*),PetscErrorCode (*create)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**),PetscErrorCode (*destroy)(PetscViewerAndFormat**))
268: {
269:   char           key[PETSC_MAX_PATH_LEN];

271:   EPSInitializePackage();
272:   SlepcMonitorMakeKey_Internal(name,vtype,format,key);
273:   PetscFunctionListAdd(&EPSMonitorList,key,monitor);
274:   if (create)  PetscFunctionListAdd(&EPSMonitorCreateList,key,create);
275:   if (destroy) PetscFunctionListAdd(&EPSMonitorDestroyList,key,destroy);
276:   PetscFunctionReturn(0);
277: }

279: /*@
280:    EPSReset - Resets the EPS context to the initial state (prior to setup)
281:    and destroys any allocated Vecs and Mats.

283:    Collective on eps

285:    Input Parameter:
286: .  eps - eigensolver context obtained from EPSCreate()

288:    Note:
289:    This can be used when a problem of different matrix size wants to be solved.
290:    All options that have previously been set are preserved, so in a next use
291:    the solver configuration is the same, but new sizes for matrices and vectors
292:    are allowed.

294:    Level: advanced

296: .seealso: EPSDestroy()
297: @*/
298: PetscErrorCode EPSReset(EPS eps)
299: {
301:   if (!eps) PetscFunctionReturn(0);
302:   if (eps->ops->reset) (eps->ops->reset)(eps);
303:   if (eps->st) STReset(eps->st);
304:   VecDestroy(&eps->D);
305:   BVDestroy(&eps->V);
306:   BVDestroy(&eps->W);
307:   VecDestroyVecs(eps->nwork,&eps->work);
308:   eps->nwork = 0;
309:   eps->state = EPS_STATE_INITIAL;
310:   PetscFunctionReturn(0);
311: }

313: /*@C
314:    EPSDestroy - Destroys the EPS context.

316:    Collective on eps

318:    Input Parameter:
319: .  eps - eigensolver context obtained from EPSCreate()

321:    Level: beginner

323: .seealso: EPSCreate(), EPSSetUp(), EPSSolve()
324: @*/
325: PetscErrorCode EPSDestroy(EPS *eps)
326: {
327:   if (!*eps) PetscFunctionReturn(0);
329:   if (--((PetscObject)(*eps))->refct > 0) { *eps = 0; PetscFunctionReturn(0); }
330:   EPSReset(*eps);
331:   if ((*eps)->ops->destroy) (*(*eps)->ops->destroy)(*eps);
332:   if ((*eps)->eigr) PetscFree4((*eps)->eigr,(*eps)->eigi,(*eps)->errest,(*eps)->perm);
333:   if ((*eps)->rr) PetscFree2((*eps)->rr,(*eps)->ri);
334:   STDestroy(&(*eps)->st);
335:   RGDestroy(&(*eps)->rg);
336:   DSDestroy(&(*eps)->ds);
337:   PetscFree((*eps)->sc);
338:   /* just in case the initial vectors have not been used */
339:   SlepcBasisDestroy_Private(&(*eps)->nds,&(*eps)->defl);
340:   SlepcBasisDestroy_Private(&(*eps)->nini,&(*eps)->IS);
341:   SlepcBasisDestroy_Private(&(*eps)->ninil,&(*eps)->ISL);
342:   if ((*eps)->convergeddestroy) (*(*eps)->convergeddestroy)((*eps)->convergedctx);
343:   EPSMonitorCancel(*eps);
344:   PetscHeaderDestroy(eps);
345:   PetscFunctionReturn(0);
346: }

348: /*@
349:    EPSSetTarget - Sets the value of the target.

351:    Logically Collective on eps

353:    Input Parameters:
354: +  eps    - eigensolver context
355: -  target - the value of the target

357:    Options Database Key:
358: .  -eps_target <scalar> - the value of the target

360:    Notes:
361:    The target is a scalar value used to determine the portion of the spectrum
362:    of interest. It is used in combination with EPSSetWhichEigenpairs().

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

368:    Level: intermediate

370: .seealso: EPSGetTarget(), EPSSetWhichEigenpairs()
371: @*/
372: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
373: {
376:   eps->target = target;
377:   if (!eps->st) EPSGetST(eps,&eps->st);
378:   STSetDefaultShift(eps->st,target);
379:   PetscFunctionReturn(0);
380: }

382: /*@
383:    EPSGetTarget - Gets the value of the target.

385:    Not Collective

387:    Input Parameter:
388: .  eps - eigensolver context

390:    Output Parameter:
391: .  target - the value of the target

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

396:    Level: intermediate

398: .seealso: EPSSetTarget()
399: @*/
400: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
401: {
404:   *target = eps->target;
405:   PetscFunctionReturn(0);
406: }

408: /*@
409:    EPSSetInterval - Defines the computational interval for spectrum slicing.

411:    Logically Collective on eps

413:    Input Parameters:
414: +  eps  - eigensolver context
415: .  inta - left end of the interval
416: -  intb - right end of the interval

418:    Options Database Key:
419: .  -eps_interval <a,b> - set [a,b] as the interval of interest

421:    Notes:
422:    Spectrum slicing is a technique employed for computing all eigenvalues of
423:    symmetric eigenproblems in a given interval. This function provides the
424:    interval to be considered. It must be used in combination with EPS_ALL, see
425:    EPSSetWhichEigenpairs().

427:    In the command-line option, two values must be provided. For an open interval,
428:    one can give an infinite, e.g., -eps_interval 1.0,inf or -eps_interval -inf,1.0.
429:    An open interval in the programmatic interface can be specified with
430:    PETSC_MAX_REAL and -PETSC_MAX_REAL.

432:    Level: intermediate

434: .seealso: EPSGetInterval(), EPSSetWhichEigenpairs()
435: @*/
436: PetscErrorCode EPSSetInterval(EPS eps,PetscReal inta,PetscReal intb)
437: {
442:   if (eps->inta != inta || eps->intb != intb) {
443:     eps->inta = inta;
444:     eps->intb = intb;
445:     eps->state = EPS_STATE_INITIAL;
446:   }
447:   PetscFunctionReturn(0);
448: }

450: /*@
451:    EPSGetInterval - Gets the computational interval for spectrum slicing.

453:    Not Collective

455:    Input Parameter:
456: .  eps - eigensolver context

458:    Output Parameters:
459: +  inta - left end of the interval
460: -  intb - right end of the interval

462:    Level: intermediate

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

467: .seealso: EPSSetInterval()
468: @*/
469: PetscErrorCode EPSGetInterval(EPS eps,PetscReal* inta,PetscReal* intb)
470: {
472:   if (inta) *inta = eps->inta;
473:   if (intb) *intb = eps->intb;
474:   PetscFunctionReturn(0);
475: }

477: /*@
478:    EPSSetST - Associates a spectral transformation object to the eigensolver.

480:    Collective on eps

482:    Input Parameters:
483: +  eps - eigensolver context obtained from EPSCreate()
484: -  st   - the spectral transformation object

486:    Note:
487:    Use EPSGetST() to retrieve the spectral transformation context (for example,
488:    to free it at the end of the computations).

490:    Level: advanced

492: .seealso: EPSGetST()
493: @*/
494: PetscErrorCode EPSSetST(EPS eps,ST st)
495: {
499:   PetscObjectReference((PetscObject)st);
500:   STDestroy(&eps->st);
501:   eps->st = st;
502:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
503:   PetscFunctionReturn(0);
504: }

506: /*@
507:    EPSGetST - Obtain the spectral transformation (ST) object associated
508:    to the eigensolver object.

510:    Not Collective

512:    Input Parameters:
513: .  eps - eigensolver context obtained from EPSCreate()

515:    Output Parameter:
516: .  st - spectral transformation context

518:    Level: intermediate

520: .seealso: EPSSetST()
521: @*/
522: PetscErrorCode EPSGetST(EPS eps,ST *st)
523: {
526:   if (!eps->st) {
527:     STCreate(PetscObjectComm((PetscObject)eps),&eps->st);
528:     PetscObjectIncrementTabLevel((PetscObject)eps->st,(PetscObject)eps,0);
529:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
530:     PetscObjectSetOptions((PetscObject)eps->st,((PetscObject)eps)->options);
531:   }
532:   *st = eps->st;
533:   PetscFunctionReturn(0);
534: }

536: /*@
537:    EPSSetBV - Associates a basis vectors object to the eigensolver.

539:    Collective on eps

541:    Input Parameters:
542: +  eps - eigensolver context obtained from EPSCreate()
543: -  V   - the basis vectors object

545:    Level: advanced

547: .seealso: EPSGetBV()
548: @*/
549: PetscErrorCode EPSSetBV(EPS eps,BV V)
550: {
554:   PetscObjectReference((PetscObject)V);
555:   BVDestroy(&eps->V);
556:   eps->V = V;
557:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
558:   PetscFunctionReturn(0);
559: }

561: /*@
562:    EPSGetBV - Obtain the basis vectors object associated to the eigensolver object.

564:    Not Collective

566:    Input Parameters:
567: .  eps - eigensolver context obtained from EPSCreate()

569:    Output Parameter:
570: .  V - basis vectors context

572:    Level: advanced

574: .seealso: EPSSetBV()
575: @*/
576: PetscErrorCode EPSGetBV(EPS eps,BV *V)
577: {
580:   if (!eps->V) {
581:     BVCreate(PetscObjectComm((PetscObject)eps),&eps->V);
582:     PetscObjectIncrementTabLevel((PetscObject)eps->V,(PetscObject)eps,0);
583:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
584:     PetscObjectSetOptions((PetscObject)eps->V,((PetscObject)eps)->options);
585:   }
586:   *V = eps->V;
587:   PetscFunctionReturn(0);
588: }

590: /*@
591:    EPSSetRG - Associates a region object to the eigensolver.

593:    Collective on eps

595:    Input Parameters:
596: +  eps - eigensolver context obtained from EPSCreate()
597: -  rg  - the region object

599:    Note:
600:    Use EPSGetRG() to retrieve the region context (for example,
601:    to free it at the end of the computations).

603:    Level: advanced

605: .seealso: EPSGetRG()
606: @*/
607: PetscErrorCode EPSSetRG(EPS eps,RG rg)
608: {
610:   if (rg) {
613:   }
614:   PetscObjectReference((PetscObject)rg);
615:   RGDestroy(&eps->rg);
616:   eps->rg = rg;
617:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
618:   PetscFunctionReturn(0);
619: }

621: /*@
622:    EPSGetRG - Obtain the region object associated to the eigensolver.

624:    Not Collective

626:    Input Parameters:
627: .  eps - eigensolver context obtained from EPSCreate()

629:    Output Parameter:
630: .  rg - region context

632:    Level: advanced

634: .seealso: EPSSetRG()
635: @*/
636: PetscErrorCode EPSGetRG(EPS eps,RG *rg)
637: {
640:   if (!eps->rg) {
641:     RGCreate(PetscObjectComm((PetscObject)eps),&eps->rg);
642:     PetscObjectIncrementTabLevel((PetscObject)eps->rg,(PetscObject)eps,0);
643:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
644:     PetscObjectSetOptions((PetscObject)eps->rg,((PetscObject)eps)->options);
645:   }
646:   *rg = eps->rg;
647:   PetscFunctionReturn(0);
648: }

650: /*@
651:    EPSSetDS - Associates a direct solver object to the eigensolver.

653:    Collective on eps

655:    Input Parameters:
656: +  eps - eigensolver context obtained from EPSCreate()
657: -  ds  - the direct solver object

659:    Note:
660:    Use EPSGetDS() to retrieve the direct solver context (for example,
661:    to free it at the end of the computations).

663:    Level: advanced

665: .seealso: EPSGetDS()
666: @*/
667: PetscErrorCode EPSSetDS(EPS eps,DS ds)
668: {
672:   PetscObjectReference((PetscObject)ds);
673:   DSDestroy(&eps->ds);
674:   eps->ds = ds;
675:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
676:   PetscFunctionReturn(0);
677: }

679: /*@
680:    EPSGetDS - Obtain the direct solver object associated to the eigensolver object.

682:    Not Collective

684:    Input Parameters:
685: .  eps - eigensolver context obtained from EPSCreate()

687:    Output Parameter:
688: .  ds - direct solver context

690:    Level: advanced

692: .seealso: EPSSetDS()
693: @*/
694: PetscErrorCode EPSGetDS(EPS eps,DS *ds)
695: {
698:   if (!eps->ds) {
699:     DSCreate(PetscObjectComm((PetscObject)eps),&eps->ds);
700:     PetscObjectIncrementTabLevel((PetscObject)eps->ds,(PetscObject)eps,0);
701:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
702:     PetscObjectSetOptions((PetscObject)eps->ds,((PetscObject)eps)->options);
703:   }
704:   *ds = eps->ds;
705:   PetscFunctionReturn(0);
706: }

708: /*@
709:    EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
710:    eigenvalue problem.

712:    Not collective

714:    Input Parameter:
715: .  eps - the eigenproblem solver context

717:    Output Parameter:
718: .  is - the answer

720:    Level: intermediate

722: .seealso: EPSIsHermitian(), EPSIsPositive()
723: @*/
724: PetscErrorCode EPSIsGeneralized(EPS eps,PetscBool* is)
725: {
728:   *is = eps->isgeneralized;
729:   PetscFunctionReturn(0);
730: }

732: /*@
733:    EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
734:    eigenvalue problem.

736:    Not collective

738:    Input Parameter:
739: .  eps - the eigenproblem solver context

741:    Output Parameter:
742: .  is - the answer

744:    Level: intermediate

746: .seealso: EPSIsGeneralized(), EPSIsPositive()
747: @*/
748: PetscErrorCode EPSIsHermitian(EPS eps,PetscBool* is)
749: {
752:   *is = eps->ishermitian;
753:   PetscFunctionReturn(0);
754: }

756: /*@
757:    EPSIsPositive - Ask if the EPS object corresponds to an eigenvalue
758:    problem type that requires a positive (semi-) definite matrix B.

760:    Not collective

762:    Input Parameter:
763: .  eps - the eigenproblem solver context

765:    Output Parameter:
766: .  is - the answer

768:    Level: intermediate

770: .seealso: EPSIsGeneralized(), EPSIsHermitian()
771: @*/
772: PetscErrorCode EPSIsPositive(EPS eps,PetscBool* is)
773: {
776:   *is = eps->ispositive;
777:   PetscFunctionReturn(0);
778: }