Actual source code: epsbasic.c

slepc-3.16.0 2021-09-30
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 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: .  eps - 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: {
 51:   EPS            eps;

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

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

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

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

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

126:   PetscNewLog(eps,&eps->sc);
127:   *outeps = eps;
128:   return(0);
129: }

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

134:    Logically Collective on eps

136:    Input Parameters:
137: +  eps  - the eigensolver context
138: -  type - a known method

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

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

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

156:    Level: intermediate

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


169:   PetscObjectTypeCompare((PetscObject)eps,type,&match);
170:   if (match) return(0);

172:   PetscFunctionListFind(EPSList,type,&r);
173:   if (!r) SETERRQ1(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown EPS type given: %s",type);

175:   if (eps->ops->destroy) { (*eps->ops->destroy)(eps); }
176:   PetscMemzero(eps->ops,sizeof(struct _EPSOps));

178:   eps->state = EPS_STATE_INITIAL;
179:   PetscObjectChangeTypeName((PetscObject)eps,type);
180:   (*r)(eps);
181:   return(0);
182: }

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

187:    Not Collective

189:    Input Parameter:
190: .  eps - the eigensolver context

192:    Output Parameter:
193: .  name - name of EPS method

195:    Level: intermediate

197: .seealso: EPSSetType()
198: @*/
199: PetscErrorCode EPSGetType(EPS eps,EPSType *type)
200: {
204:   *type = ((PetscObject)eps)->type_name;
205:   return(0);
206: }

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

211:    Not Collective

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

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

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

225:    Then, your solver can be chosen with the procedural interface via
226: $     EPSSetType(eps,"my_solver")
227:    or at runtime via the option
228: $     -eps_type my_solver

230:    Level: advanced

232: .seealso: EPSRegisterAll()
233: @*/
234: PetscErrorCode EPSRegister(const char *name,PetscErrorCode (*function)(EPS))
235: {

239:   EPSInitializePackage();
240:   PetscFunctionListAdd(&EPSList,name,function);
241:   return(0);
242: }

244: /*@C
245:    EPSMonitorRegister - Adds EPS monitor routine.

247:    Not Collective

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

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

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

265:    Then, your monitor can be chosen with the procedural interface via
266: $      EPSMonitorSetFromOptions(eps,"-eps_monitor_my_monitor","my_monitor",NULL)
267:    or at runtime via the option
268: $      -eps_monitor_my_monitor

270:    Level: advanced

272: .seealso: EPSMonitorRegisterAll()
273: @*/
274: 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**))
275: {
276:   char           key[PETSC_MAX_PATH_LEN];

280:   EPSInitializePackage();
281:   SlepcMonitorMakeKey_Internal(name,vtype,format,key);
282:   PetscFunctionListAdd(&EPSMonitorList,key,monitor);
283:   if (create)  { PetscFunctionListAdd(&EPSMonitorCreateList,key,create); }
284:   if (destroy) { PetscFunctionListAdd(&EPSMonitorDestroyList,key,destroy); }
285:   return(0);
286: }

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

292:    Collective on eps

294:    Input Parameter:
295: .  eps - eigensolver context obtained from EPSCreate()

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

303:    Level: advanced

305: .seealso: EPSDestroy()
306: @*/
307: PetscErrorCode EPSReset(EPS eps)
308: {

313:   if (!eps) return(0);
314:   if (eps->ops->reset) { (eps->ops->reset)(eps); }
315:   if (eps->st) { STReset(eps->st); }
316:   VecDestroy(&eps->D);
317:   BVDestroy(&eps->V);
318:   BVDestroy(&eps->W);
319:   VecDestroyVecs(eps->nwork,&eps->work);
320:   eps->nwork = 0;
321:   eps->state = EPS_STATE_INITIAL;
322:   return(0);
323: }

325: /*@C
326:    EPSDestroy - Destroys the EPS context.

328:    Collective on eps

330:    Input Parameter:
331: .  eps - eigensolver context obtained from EPSCreate()

333:    Level: beginner

335: .seealso: EPSCreate(), EPSSetUp(), EPSSolve()
336: @*/
337: PetscErrorCode EPSDestroy(EPS *eps)
338: {

342:   if (!*eps) return(0);
344:   if (--((PetscObject)(*eps))->refct > 0) { *eps = 0; return(0); }
345:   EPSReset(*eps);
346:   if ((*eps)->ops->destroy) { (*(*eps)->ops->destroy)(*eps); }
347:   if ((*eps)->eigr) {
348:     PetscFree4((*eps)->eigr,(*eps)->eigi,(*eps)->errest,(*eps)->perm);
349:   }
350:   if ((*eps)->rr) {
351:     PetscFree2((*eps)->rr,(*eps)->ri);
352:   }
353:   STDestroy(&(*eps)->st);
354:   RGDestroy(&(*eps)->rg);
355:   DSDestroy(&(*eps)->ds);
356:   PetscFree((*eps)->sc);
357:   /* just in case the initial vectors have not been used */
358:   SlepcBasisDestroy_Private(&(*eps)->nds,&(*eps)->defl);
359:   SlepcBasisDestroy_Private(&(*eps)->nini,&(*eps)->IS);
360:   SlepcBasisDestroy_Private(&(*eps)->ninil,&(*eps)->ISL);
361:   if ((*eps)->convergeddestroy) {
362:     (*(*eps)->convergeddestroy)((*eps)->convergedctx);
363:   }
364:   EPSMonitorCancel(*eps);
365:   PetscHeaderDestroy(eps);
366:   return(0);
367: }

369: /*@
370:    EPSSetTarget - Sets the value of the target.

372:    Logically Collective on eps

374:    Input Parameters:
375: +  eps    - eigensolver context
376: -  target - the value of the target

378:    Options Database Key:
379: .  -eps_target <scalar> - the value of the target

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

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

389:    Level: intermediate

391: .seealso: EPSGetTarget(), EPSSetWhichEigenpairs()
392: @*/
393: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
394: {

400:   eps->target = target;
401:   if (!eps->st) { EPSGetST(eps,&eps->st); }
402:   STSetDefaultShift(eps->st,target);
403:   return(0);
404: }

406: /*@
407:    EPSGetTarget - Gets the value of the target.

409:    Not Collective

411:    Input Parameter:
412: .  eps - eigensolver context

414:    Output Parameter:
415: .  target - the value of the target

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

420:    Level: intermediate

422: .seealso: EPSSetTarget()
423: @*/
424: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
425: {
429:   *target = eps->target;
430:   return(0);
431: }

433: /*@
434:    EPSSetInterval - Defines the computational interval for spectrum slicing.

436:    Logically Collective on eps

438:    Input Parameters:
439: +  eps  - eigensolver context
440: .  inta - left end of the interval
441: -  intb - right end of the interval

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

446:    Notes:
447:    Spectrum slicing is a technique employed for computing all eigenvalues of
448:    symmetric eigenproblems in a given interval. This function provides the
449:    interval to be considered. It must be used in combination with EPS_ALL, see
450:    EPSSetWhichEigenpairs().

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

457:    Level: intermediate

459: .seealso: EPSGetInterval(), EPSSetWhichEigenpairs()
460: @*/
461: PetscErrorCode EPSSetInterval(EPS eps,PetscReal inta,PetscReal intb)
462: {
467:   if (inta>=intb) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
468:   if (eps->inta != inta || eps->intb != intb) {
469:     eps->inta = inta;
470:     eps->intb = intb;
471:     eps->state = EPS_STATE_INITIAL;
472:   }
473:   return(0);
474: }

476: /*@
477:    EPSGetInterval - Gets the computational interval for spectrum slicing.

479:    Not Collective

481:    Input Parameter:
482: .  eps - eigensolver context

484:    Output Parameters:
485: +  inta - left end of the interval
486: -  intb - right end of the interval

488:    Level: intermediate

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

493: .seealso: EPSSetInterval()
494: @*/
495: PetscErrorCode EPSGetInterval(EPS eps,PetscReal* inta,PetscReal* intb)
496: {
499:   if (inta) *inta = eps->inta;
500:   if (intb) *intb = eps->intb;
501:   return(0);
502: }

504: /*@
505:    EPSSetST - Associates a spectral transformation object to the eigensolver.

507:    Collective on eps

509:    Input Parameters:
510: +  eps - eigensolver context obtained from EPSCreate()
511: -  st   - the spectral transformation object

513:    Note:
514:    Use EPSGetST() to retrieve the spectral transformation context (for example,
515:    to free it at the end of the computations).

517:    Level: advanced

519: .seealso: EPSGetST()
520: @*/
521: PetscErrorCode EPSSetST(EPS eps,ST st)
522: {

529:   PetscObjectReference((PetscObject)st);
530:   STDestroy(&eps->st);
531:   eps->st = st;
532:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
533:   return(0);
534: }

536: /*@
537:    EPSGetST - Obtain the spectral transformation (ST) object associated
538:    to the eigensolver object.

540:    Not Collective

542:    Input Parameters:
543: .  eps - eigensolver context obtained from EPSCreate()

545:    Output Parameter:
546: .  st - spectral transformation context

548:    Level: intermediate

550: .seealso: EPSSetST()
551: @*/
552: PetscErrorCode EPSGetST(EPS eps,ST *st)
553: {

559:   if (!eps->st) {
560:     STCreate(PetscObjectComm((PetscObject)eps),&eps->st);
561:     PetscObjectIncrementTabLevel((PetscObject)eps->st,(PetscObject)eps,0);
562:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
563:     PetscObjectSetOptions((PetscObject)eps->st,((PetscObject)eps)->options);
564:   }
565:   *st = eps->st;
566:   return(0);
567: }

569: /*@
570:    EPSSetBV - Associates a basis vectors object to the eigensolver.

572:    Collective on eps

574:    Input Parameters:
575: +  eps - eigensolver context obtained from EPSCreate()
576: -  V   - the basis vectors object

578:    Level: advanced

580: .seealso: EPSGetBV()
581: @*/
582: PetscErrorCode EPSSetBV(EPS eps,BV V)
583: {

590:   PetscObjectReference((PetscObject)V);
591:   BVDestroy(&eps->V);
592:   eps->V = V;
593:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
594:   return(0);
595: }

597: /*@
598:    EPSGetBV - Obtain the basis vectors object associated to the eigensolver object.

600:    Not Collective

602:    Input Parameters:
603: .  eps - eigensolver context obtained from EPSCreate()

605:    Output Parameter:
606: .  V - basis vectors context

608:    Level: advanced

610: .seealso: EPSSetBV()
611: @*/
612: PetscErrorCode EPSGetBV(EPS eps,BV *V)
613: {

619:   if (!eps->V) {
620:     BVCreate(PetscObjectComm((PetscObject)eps),&eps->V);
621:     PetscObjectIncrementTabLevel((PetscObject)eps->V,(PetscObject)eps,0);
622:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
623:     PetscObjectSetOptions((PetscObject)eps->V,((PetscObject)eps)->options);
624:   }
625:   *V = eps->V;
626:   return(0);
627: }

629: /*@
630:    EPSSetRG - Associates a region object to the eigensolver.

632:    Collective on eps

634:    Input Parameters:
635: +  eps - eigensolver context obtained from EPSCreate()
636: -  rg  - the region object

638:    Note:
639:    Use EPSGetRG() to retrieve the region context (for example,
640:    to free it at the end of the computations).

642:    Level: advanced

644: .seealso: EPSGetRG()
645: @*/
646: PetscErrorCode EPSSetRG(EPS eps,RG rg)
647: {

652:   if (rg) {
655:   }
656:   PetscObjectReference((PetscObject)rg);
657:   RGDestroy(&eps->rg);
658:   eps->rg = rg;
659:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
660:   return(0);
661: }

663: /*@
664:    EPSGetRG - Obtain the region object associated to the eigensolver.

666:    Not Collective

668:    Input Parameters:
669: .  eps - eigensolver context obtained from EPSCreate()

671:    Output Parameter:
672: .  rg - region context

674:    Level: advanced

676: .seealso: EPSSetRG()
677: @*/
678: PetscErrorCode EPSGetRG(EPS eps,RG *rg)
679: {

685:   if (!eps->rg) {
686:     RGCreate(PetscObjectComm((PetscObject)eps),&eps->rg);
687:     PetscObjectIncrementTabLevel((PetscObject)eps->rg,(PetscObject)eps,0);
688:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
689:     PetscObjectSetOptions((PetscObject)eps->rg,((PetscObject)eps)->options);
690:   }
691:   *rg = eps->rg;
692:   return(0);
693: }

695: /*@
696:    EPSSetDS - Associates a direct solver object to the eigensolver.

698:    Collective on eps

700:    Input Parameters:
701: +  eps - eigensolver context obtained from EPSCreate()
702: -  ds  - the direct solver object

704:    Note:
705:    Use EPSGetDS() to retrieve the direct solver context (for example,
706:    to free it at the end of the computations).

708:    Level: advanced

710: .seealso: EPSGetDS()
711: @*/
712: PetscErrorCode EPSSetDS(EPS eps,DS ds)
713: {

720:   PetscObjectReference((PetscObject)ds);
721:   DSDestroy(&eps->ds);
722:   eps->ds = ds;
723:   PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
724:   return(0);
725: }

727: /*@
728:    EPSGetDS - Obtain the direct solver object associated to the eigensolver object.

730:    Not Collective

732:    Input Parameters:
733: .  eps - eigensolver context obtained from EPSCreate()

735:    Output Parameter:
736: .  ds - direct solver context

738:    Level: advanced

740: .seealso: EPSSetDS()
741: @*/
742: PetscErrorCode EPSGetDS(EPS eps,DS *ds)
743: {

749:   if (!eps->ds) {
750:     DSCreate(PetscObjectComm((PetscObject)eps),&eps->ds);
751:     PetscObjectIncrementTabLevel((PetscObject)eps->ds,(PetscObject)eps,0);
752:     PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
753:     PetscObjectSetOptions((PetscObject)eps->ds,((PetscObject)eps)->options);
754:   }
755:   *ds = eps->ds;
756:   return(0);
757: }

759: /*@
760:    EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
761:    eigenvalue problem.

763:    Not collective

765:    Input Parameter:
766: .  eps - the eigenproblem solver context

768:    Output Parameter:
769: .  is - the answer

771:    Level: intermediate

773: .seealso: EPSIsHermitian(), EPSIsPositive()
774: @*/
775: PetscErrorCode EPSIsGeneralized(EPS eps,PetscBool* is)
776: {
780:   *is = eps->isgeneralized;
781:   return(0);
782: }

784: /*@
785:    EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
786:    eigenvalue problem.

788:    Not collective

790:    Input Parameter:
791: .  eps - the eigenproblem solver context

793:    Output Parameter:
794: .  is - the answer

796:    Level: intermediate

798: .seealso: EPSIsGeneralized(), EPSIsPositive()
799: @*/
800: PetscErrorCode EPSIsHermitian(EPS eps,PetscBool* is)
801: {
805:   *is = eps->ishermitian;
806:   return(0);
807: }

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

813:    Not collective

815:    Input Parameter:
816: .  eps - the eigenproblem solver context

818:    Output Parameter:
819: .  is - the answer

821:    Level: intermediate

823: .seealso: EPSIsGeneralized(), EPSIsHermitian()
824: @*/
825: PetscErrorCode EPSIsPositive(EPS eps,PetscBool* is)
826: {
830:   *is = eps->ispositive;
831:   return(0);
832: }