Actual source code: epsbasic.c
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 = NULL;
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 `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: [](ch:eps), `EPSSetUp()`, `EPSSolve()`, `EPSDestroy()`, `EPS`
47: @*/
48: PetscErrorCode EPSCreate(MPI_Comm comm,EPS *outeps)
49: {
50: EPS eps;
52: PetscFunctionBegin;
53: PetscAssertPointer(outeps,2);
54: PetscCall(EPSInitializePackage());
55: PetscCall(SlepcHeaderCreate(eps,EPS_CLASSID,"EPS","Eigenvalue Problem Solver","EPS",comm,EPSDestroy,EPSView));
57: eps->max_it = PETSC_DETERMINE;
58: eps->nev = 0;
59: eps->ncv = PETSC_DETERMINE;
60: eps->mpd = PETSC_DETERMINE;
61: eps->nini = 0;
62: eps->nds = 0;
63: eps->target = 0.0;
64: eps->tol = PETSC_DETERMINE;
65: eps->thres = PETSC_MIN_REAL;
66: eps->threlative = PETSC_FALSE;
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->isstructured = PETSC_FALSE;
125: eps->reason = EPS_CONVERGED_ITERATING;
127: PetscCall(PetscNew(&eps->sc));
128: *outeps = eps;
129: PetscFunctionReturn(PETSC_SUCCESS);
130: }
132: /*@
133: EPSSetType - Selects the particular solver to be used in the `EPS` object.
135: Logically Collective
137: Input Parameters:
138: + eps - the linear eigensolver context
139: - type - a known method
141: Options Database Key:
142: . -eps_type <method> - Sets the method; use `-help` for a list
143: of available methods
145: Notes:
146: See `EPSType` for available methods. The default 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: [](ch:eps), `STSetType()`, `EPSType`
159: @*/
160: PetscErrorCode EPSSetType(EPS eps,EPSType type)
161: {
162: PetscErrorCode (*r)(EPS);
163: PetscBool match;
165: PetscFunctionBegin;
167: PetscAssertPointer(type,2);
169: PetscCall(PetscObjectTypeCompare((PetscObject)eps,type,&match));
170: if (match) PetscFunctionReturn(PETSC_SUCCESS);
172: PetscCall(PetscFunctionListFind(EPSList,type,&r));
173: PetscCheck(r,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown EPS type given: %s",type);
175: PetscTryTypeMethod(eps,destroy);
176: PetscCall(PetscMemzero(eps->ops,sizeof(struct _EPSOps)));
178: eps->state = EPS_STATE_INITIAL;
179: PetscCall(PetscObjectChangeTypeName((PetscObject)eps,type));
180: PetscCall((*r)(eps));
181: PetscFunctionReturn(PETSC_SUCCESS);
182: }
184: /*@
185: EPSGetType - Gets the `EPS` type as a string from the `EPS` object.
187: Not Collective
189: Input Parameter:
190: . eps - the linear eigensolver context
192: Output Parameter:
193: . type - name of `EPS` method
195: Level: intermediate
197: .seealso: [](ch:eps), `EPSSetType()`
198: @*/
199: PetscErrorCode EPSGetType(EPS eps,EPSType *type)
200: {
201: PetscFunctionBegin;
203: PetscAssertPointer(type,2);
204: *type = ((PetscObject)eps)->type_name;
205: PetscFunctionReturn(PETSC_SUCCESS);
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: Note:
218: `EPSRegister()` may be called multiple times to add several user-defined solvers.
220: Example Usage:
221: .vb
222: EPSRegister("my_solver",MySolverCreate);
223: .ve
225: Then, your solver can be chosen with the procedural interface via
226: .vb
227: EPSSetType(eps,"my_solver")
228: .ve
229: or at runtime via the option `-eps_type my_solver`.
231: Level: advanced
233: .seealso: [](ch:eps), `EPSRegisterAll()`
234: @*/
235: PetscErrorCode EPSRegister(const char *name,PetscErrorCode (*function)(EPS))
236: {
237: PetscFunctionBegin;
238: PetscCall(EPSInitializePackage());
239: PetscCall(PetscFunctionListAdd(&EPSList,name,function));
240: PetscFunctionReturn(PETSC_SUCCESS);
241: }
243: /*@C
244: EPSMonitorRegister - Registers an `EPS` monitor routine that may be accessed with
245: `EPSMonitorSetFromOptions()`.
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, see `EPSMonitorRegisterFn`
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: The calling sequence for the given function matches the calling sequence of `EPSMonitorFn`
261: functions passed to `EPSMonitorSet()` with the additional requirement that its final argument
262: be a `PetscViewerAndFormat`.
264: Example Usage:
265: .vb
266: EPSMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
267: .ve
269: Then, your monitor can be chosen with the procedural interface via
270: .vb
271: EPSMonitorSetFromOptions(eps,"-eps_monitor_my_monitor","my_monitor",NULL)
272: .ve
273: or at runtime via the option `-eps_monitor_my_monitor`.
275: Level: advanced
277: .seealso: [](ch:eps), `EPSMonitorSet()`, `EPSMonitorRegisterAll()`, `EPSMonitorSetFromOptions()`
278: @*/
279: PetscErrorCode EPSMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,EPSMonitorRegisterFn *monitor,EPSMonitorRegisterCreateFn *create,EPSMonitorRegisterDestroyFn *destroy)
280: {
281: char key[PETSC_MAX_PATH_LEN];
283: PetscFunctionBegin;
284: PetscCall(EPSInitializePackage());
285: PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
286: PetscCall(PetscFunctionListAdd(&EPSMonitorList,key,monitor));
287: if (create) PetscCall(PetscFunctionListAdd(&EPSMonitorCreateList,key,create));
288: if (destroy) PetscCall(PetscFunctionListAdd(&EPSMonitorDestroyList,key,destroy));
289: PetscFunctionReturn(PETSC_SUCCESS);
290: }
292: /*@
293: EPSReset - Resets the EPS context to the initial state (prior to setup)
294: and destroys any allocated Vecs and Mats.
296: Collective
298: Input Parameter:
299: . eps - the linear eigensolver context
301: Note:
302: This can be used when a problem of different matrix size wants to be solved.
303: All options that have previously been set are preserved, so in a next use
304: the solver configuration is the same, but new sizes for matrices and vectors
305: are allowed.
307: Level: advanced
309: .seealso: [](ch:eps), `EPSDestroy()`
310: @*/
311: PetscErrorCode EPSReset(EPS eps)
312: {
313: PetscFunctionBegin;
315: if (!eps) PetscFunctionReturn(PETSC_SUCCESS);
316: PetscTryTypeMethod(eps,reset);
317: if (eps->st) PetscCall(STReset(eps->st));
318: PetscCall(VecDestroy(&eps->D));
319: PetscCall(BVDestroy(&eps->V));
320: PetscCall(BVDestroy(&eps->W));
321: PetscCall(VecDestroyVecs(eps->nwork,&eps->work));
322: eps->nwork = 0;
323: eps->state = EPS_STATE_INITIAL;
324: PetscFunctionReturn(PETSC_SUCCESS);
325: }
327: /*@
328: EPSDestroy - Destroys the `EPS` context.
330: Collective
332: Input Parameter:
333: . eps - the linear eigensolver context
335: Level: beginner
337: .seealso: [](ch:eps), `EPSCreate()`, `EPSSetUp()`, `EPSSolve()`
338: @*/
339: PetscErrorCode EPSDestroy(EPS *eps)
340: {
341: PetscFunctionBegin;
342: if (!*eps) PetscFunctionReturn(PETSC_SUCCESS);
344: if (--((PetscObject)*eps)->refct > 0) { *eps = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
345: PetscCall(EPSReset(*eps));
346: PetscTryTypeMethod(*eps,destroy);
347: if ((*eps)->eigr) PetscCall(PetscFree4((*eps)->eigr,(*eps)->eigi,(*eps)->errest,(*eps)->perm));
348: if ((*eps)->rr) PetscCall(PetscFree2((*eps)->rr,(*eps)->ri));
349: PetscCall(STDestroy(&(*eps)->st));
350: PetscCall(RGDestroy(&(*eps)->rg));
351: PetscCall(DSDestroy(&(*eps)->ds));
352: PetscCall(PetscFree((*eps)->sc));
353: /* just in case the initial vectors have not been used */
354: PetscCall(SlepcBasisDestroy_Private(&(*eps)->nds,&(*eps)->defl));
355: PetscCall(SlepcBasisDestroy_Private(&(*eps)->nini,&(*eps)->IS));
356: PetscCall(SlepcBasisDestroy_Private(&(*eps)->ninil,&(*eps)->ISL));
357: if ((*eps)->convergeddestroy) PetscCall((*(*eps)->convergeddestroy)(&(*eps)->convergedctx));
358: if ((*eps)->stoppingdestroy) PetscCall((*(*eps)->stoppingdestroy)(&(*eps)->stoppingctx));
359: PetscCall(EPSMonitorCancel(*eps));
360: PetscCall(PetscHeaderDestroy(eps));
361: PetscFunctionReturn(PETSC_SUCCESS);
362: }
364: /*@
365: EPSSetTarget - Sets the value of the target.
367: Logically Collective
369: Input Parameters:
370: + eps - the linear eigensolver context
371: - target - the value of the target
373: Options Database Key:
374: . -eps_target <scalar> - the value of the target
376: Notes:
377: The target is a scalar value used to determine the portion of the spectrum
378: of interest. It is used in combination with EPSSetWhichEigenpairs().
380: In the case of complex scalars, a complex value can be provided in the
381: command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
382: -eps_target 1.0+2.0i
384: Level: intermediate
386: .seealso: [](ch:eps), `EPSGetTarget()`, `EPSSetWhichEigenpairs()`
387: @*/
388: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
389: {
390: PetscFunctionBegin;
393: eps->target = target;
394: if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
395: PetscCall(STSetDefaultShift(eps->st,target));
396: PetscFunctionReturn(PETSC_SUCCESS);
397: }
399: /*@
400: EPSGetTarget - Gets the value of the target.
402: Not Collective
404: Input Parameter:
405: . eps - the linear eigensolver context
407: Output Parameter:
408: . target - the value of the target
410: Note:
411: If the target was not set by the user, then zero is returned.
413: Level: intermediate
415: .seealso: [](ch:eps), `EPSSetTarget()`
416: @*/
417: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
418: {
419: PetscFunctionBegin;
421: PetscAssertPointer(target,2);
422: *target = eps->target;
423: PetscFunctionReturn(PETSC_SUCCESS);
424: }
426: /*@
427: EPSSetInterval - Defines the computational interval for spectrum slicing.
429: Logically Collective
431: Input Parameters:
432: + eps - the linear eigensolver context
433: . inta - left end of the interval
434: - intb - right end of the interval
436: Options Database Key:
437: . -eps_interval <a,b> - set [a,b] as the interval of interest
439: Notes:
440: Spectrum slicing is a technique employed for computing all eigenvalues of
441: symmetric eigenproblems in a given interval. This function provides the
442: interval to be considered. It must be used in combination with EPS_ALL, see
443: EPSSetWhichEigenpairs().
445: In the command-line option, two values must be provided. For an open interval,
446: one can give an infinite, e.g., -eps_interval 1.0,inf or -eps_interval -inf,1.0.
447: An open interval in the programmatic interface can be specified with
448: PETSC_MAX_REAL and -PETSC_MAX_REAL.
450: Level: intermediate
452: .seealso: [](ch:eps), `EPSGetInterval()`, `EPSSetWhichEigenpairs()`
453: @*/
454: PetscErrorCode EPSSetInterval(EPS eps,PetscReal inta,PetscReal intb)
455: {
456: PetscFunctionBegin;
460: PetscCheck(inta<intb,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
461: if (eps->inta != inta || eps->intb != intb) {
462: eps->inta = inta;
463: eps->intb = intb;
464: eps->state = EPS_STATE_INITIAL;
465: }
466: PetscFunctionReturn(PETSC_SUCCESS);
467: }
469: /*@
470: EPSGetInterval - Gets the computational interval for spectrum slicing.
472: Not Collective
474: Input Parameter:
475: . eps - the linear eigensolver context
477: Output Parameters:
478: + inta - left end of the interval
479: - intb - right end of the interval
481: Level: intermediate
483: Note:
484: If the interval was not set by the user, then zeros are returned.
486: .seealso: [](ch:eps), `EPSSetInterval()`
487: @*/
488: PetscErrorCode EPSGetInterval(EPS eps,PetscReal* inta,PetscReal* intb)
489: {
490: PetscFunctionBegin;
492: if (inta) *inta = eps->inta;
493: if (intb) *intb = eps->intb;
494: PetscFunctionReturn(PETSC_SUCCESS);
495: }
497: /*@
498: EPSSetST - Associates a spectral transformation object to the eigensolver.
500: Collective
502: Input Parameters:
503: + eps - the linear eigensolver context
504: - st - the spectral transformation object
506: Note:
507: Use EPSGetST() to retrieve the spectral transformation context (for example,
508: to free it at the end of the computations).
510: Level: advanced
512: .seealso: [](ch:eps), `EPSGetST()`
513: @*/
514: PetscErrorCode EPSSetST(EPS eps,ST st)
515: {
516: PetscFunctionBegin;
519: PetscCheckSameComm(eps,1,st,2);
520: PetscCall(PetscObjectReference((PetscObject)st));
521: PetscCall(STDestroy(&eps->st));
522: eps->st = st;
523: PetscFunctionReturn(PETSC_SUCCESS);
524: }
526: /*@
527: EPSGetST - Obtain the spectral transformation (ST) object associated
528: to the eigensolver object.
530: Not Collective
532: Input Parameter:
533: . eps - the linear eigensolver context
535: Output Parameter:
536: . st - spectral transformation context
538: Level: intermediate
540: .seealso: [](ch:eps), `EPSSetST()`
541: @*/
542: PetscErrorCode EPSGetST(EPS eps,ST *st)
543: {
544: PetscFunctionBegin;
546: PetscAssertPointer(st,2);
547: if (!eps->st) {
548: PetscCall(STCreate(PetscObjectComm((PetscObject)eps),&eps->st));
549: PetscCall(PetscObjectIncrementTabLevel((PetscObject)eps->st,(PetscObject)eps,0));
550: PetscCall(PetscObjectSetOptions((PetscObject)eps->st,((PetscObject)eps)->options));
551: }
552: *st = eps->st;
553: PetscFunctionReturn(PETSC_SUCCESS);
554: }
556: /*@
557: EPSSetBV - Associates a basis vectors object to the eigensolver.
559: Collective
561: Input Parameters:
562: + eps - the linear eigensolver context
563: - V - the basis vectors object
565: Level: advanced
567: .seealso: [](ch:eps), `EPSGetBV()`
568: @*/
569: PetscErrorCode EPSSetBV(EPS eps,BV V)
570: {
571: PetscFunctionBegin;
574: PetscCheckSameComm(eps,1,V,2);
575: PetscCall(PetscObjectReference((PetscObject)V));
576: PetscCall(BVDestroy(&eps->V));
577: eps->V = V;
578: PetscFunctionReturn(PETSC_SUCCESS);
579: }
581: /*@
582: EPSGetBV - Obtain the basis vectors object associated to the eigensolver object.
584: Not Collective
586: Input Parameter:
587: . eps - the linear eigensolver context
589: Output Parameter:
590: . V - basis vectors context
592: Level: advanced
594: .seealso: [](ch:eps), `EPSSetBV()`
595: @*/
596: PetscErrorCode EPSGetBV(EPS eps,BV *V)
597: {
598: PetscFunctionBegin;
600: PetscAssertPointer(V,2);
601: if (!eps->V) {
602: PetscCall(BVCreate(PetscObjectComm((PetscObject)eps),&eps->V));
603: PetscCall(PetscObjectIncrementTabLevel((PetscObject)eps->V,(PetscObject)eps,0));
604: PetscCall(PetscObjectSetOptions((PetscObject)eps->V,((PetscObject)eps)->options));
605: }
606: *V = eps->V;
607: PetscFunctionReturn(PETSC_SUCCESS);
608: }
610: /*@
611: EPSSetRG - Associates a region object to the eigensolver.
613: Collective
615: Input Parameters:
616: + eps - the linear eigensolver context
617: - rg - the region object
619: Note:
620: Use EPSGetRG() to retrieve the region context (for example,
621: to free it at the end of the computations).
623: Level: advanced
625: .seealso: [](ch:eps), `EPSGetRG()`
626: @*/
627: PetscErrorCode EPSSetRG(EPS eps,RG rg)
628: {
629: PetscFunctionBegin;
631: if (rg) {
633: PetscCheckSameComm(eps,1,rg,2);
634: }
635: PetscCall(PetscObjectReference((PetscObject)rg));
636: PetscCall(RGDestroy(&eps->rg));
637: eps->rg = rg;
638: PetscFunctionReturn(PETSC_SUCCESS);
639: }
641: /*@
642: EPSGetRG - Obtain the region object associated to the eigensolver.
644: Not Collective
646: Input Parameter:
647: . eps - the linear eigensolver context
649: Output Parameter:
650: . rg - region context
652: Level: advanced
654: .seealso: [](ch:eps), `EPSSetRG()`
655: @*/
656: PetscErrorCode EPSGetRG(EPS eps,RG *rg)
657: {
658: PetscFunctionBegin;
660: PetscAssertPointer(rg,2);
661: if (!eps->rg) {
662: PetscCall(RGCreate(PetscObjectComm((PetscObject)eps),&eps->rg));
663: PetscCall(PetscObjectIncrementTabLevel((PetscObject)eps->rg,(PetscObject)eps,0));
664: PetscCall(PetscObjectSetOptions((PetscObject)eps->rg,((PetscObject)eps)->options));
665: }
666: *rg = eps->rg;
667: PetscFunctionReturn(PETSC_SUCCESS);
668: }
670: /*@
671: EPSSetDS - Associates a direct solver object to the eigensolver.
673: Collective
675: Input Parameters:
676: + eps - the linear eigensolver context
677: - ds - the direct solver object
679: Note:
680: Use EPSGetDS() to retrieve the direct solver context (for example,
681: to free it at the end of the computations).
683: Level: advanced
685: .seealso: [](ch:eps), `EPSGetDS()`
686: @*/
687: PetscErrorCode EPSSetDS(EPS eps,DS ds)
688: {
689: PetscFunctionBegin;
692: PetscCheckSameComm(eps,1,ds,2);
693: PetscCall(PetscObjectReference((PetscObject)ds));
694: PetscCall(DSDestroy(&eps->ds));
695: eps->ds = ds;
696: PetscFunctionReturn(PETSC_SUCCESS);
697: }
699: /*@
700: EPSGetDS - Obtain the direct solver object associated to the eigensolver object.
702: Not Collective
704: Input Parameter:
705: . eps - the linear eigensolver context
707: Output Parameter:
708: . ds - direct solver context
710: Level: advanced
712: .seealso: [](ch:eps), `EPSSetDS()`
713: @*/
714: PetscErrorCode EPSGetDS(EPS eps,DS *ds)
715: {
716: PetscFunctionBegin;
718: PetscAssertPointer(ds,2);
719: if (!eps->ds) {
720: PetscCall(DSCreate(PetscObjectComm((PetscObject)eps),&eps->ds));
721: PetscCall(PetscObjectIncrementTabLevel((PetscObject)eps->ds,(PetscObject)eps,0));
722: PetscCall(PetscObjectSetOptions((PetscObject)eps->ds,((PetscObject)eps)->options));
723: }
724: *ds = eps->ds;
725: PetscFunctionReturn(PETSC_SUCCESS);
726: }
728: /*@
729: EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
730: eigenvalue problem.
732: Not Collective
734: Input Parameter:
735: . eps - the linear eigensolver context
737: Output Parameter:
738: . is - the answer
740: Level: intermediate
742: .seealso: [](ch:eps), `EPSIsHermitian()`, `EPSIsPositive()`, `EPSIsStructured()`
743: @*/
744: PetscErrorCode EPSIsGeneralized(EPS eps,PetscBool* is)
745: {
746: PetscFunctionBegin;
748: PetscAssertPointer(is,2);
749: *is = eps->isgeneralized;
750: PetscFunctionReturn(PETSC_SUCCESS);
751: }
753: /*@
754: EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
755: eigenvalue problem.
757: Not Collective
759: Input Parameter:
760: . eps - the linear eigensolver context
762: Output Parameter:
763: . is - the answer
765: Level: intermediate
767: .seealso: [](ch:eps), `EPSIsGeneralized()`, `EPSIsPositive()`, `EPSIsStructured()`
768: @*/
769: PetscErrorCode EPSIsHermitian(EPS eps,PetscBool* is)
770: {
771: PetscFunctionBegin;
773: PetscAssertPointer(is,2);
774: *is = eps->ishermitian;
775: PetscFunctionReturn(PETSC_SUCCESS);
776: }
778: /*@
779: EPSIsPositive - Ask if the EPS object corresponds to an eigenvalue
780: problem type that requires a positive (semi-) definite matrix B.
782: Not Collective
784: Input Parameter:
785: . eps - the linear eigensolver context
787: Output Parameter:
788: . is - the answer
790: Level: intermediate
792: .seealso: [](ch:eps), `EPSIsGeneralized()`, `EPSIsHermitian()`, `EPSIsStructured()`
793: @*/
794: PetscErrorCode EPSIsPositive(EPS eps,PetscBool* is)
795: {
796: PetscFunctionBegin;
798: PetscAssertPointer(is,2);
799: *is = eps->ispositive;
800: PetscFunctionReturn(PETSC_SUCCESS);
801: }
803: /*@
804: EPSIsStructured - Ask if the EPS object corresponds to a structured
805: eigenvalue problem.
807: Not Collective
809: Input Parameter:
810: . eps - the linear eigensolver context
812: Output Parameter:
813: . is - the answer
815: Note:
816: The result will be true if the problem type has been set to some
817: structured type such as EPS_BSE. This is independent of whether the
818: input matrix has been built with a certain structure with a helper function.
820: Level: intermediate
822: .seealso: [](ch:eps), `EPSIsGeneralized()`, `EPSIsHermitian()`, `EPSIsPositive()`, `EPSSetProblemType()`
823: @*/
824: PetscErrorCode EPSIsStructured(EPS eps,PetscBool* is)
825: {
826: PetscFunctionBegin;
828: PetscAssertPointer(is,2);
829: *is = eps->isstructured;
830: PetscFunctionReturn(PETSC_SUCCESS);
831: }