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 \<type\> - sets the method; use `-help` for a list of available methods
144: Notes:
145: See `EPSType` for available methods. The default is `EPSKRYLOVSCHUR`.
147: Normally, it is best to use the `EPSSetFromOptions()` command and
148: then set the `EPS` 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 `EPSSetType()` 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: [](ch:eps), `STSetType()`, `EPSType`
158: @*/
159: PetscErrorCode EPSSetType(EPS eps,EPSType type)
160: {
161: PetscErrorCode (*r)(EPS);
162: PetscBool match;
164: PetscFunctionBegin;
166: PetscAssertPointer(type,2);
168: PetscCall(PetscObjectTypeCompare((PetscObject)eps,type,&match));
169: if (match) PetscFunctionReturn(PETSC_SUCCESS);
171: PetscCall(PetscFunctionListFind(EPSList,type,&r));
172: PetscCheck(r,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown EPS type given: %s",type);
174: PetscTryTypeMethod(eps,destroy);
175: PetscCall(PetscMemzero(eps->ops,sizeof(struct _EPSOps)));
177: eps->state = EPS_STATE_INITIAL;
178: PetscCall(PetscObjectChangeTypeName((PetscObject)eps,type));
179: PetscCall((*r)(eps));
180: PetscFunctionReturn(PETSC_SUCCESS);
181: }
183: /*@
184: EPSGetType - Gets the `EPS` type as a string from the `EPS` object.
186: Not Collective
188: Input Parameter:
189: . eps - the linear eigensolver context
191: Output Parameter:
192: . type - name of `EPS` method
194: Level: intermediate
196: .seealso: [](ch:eps), `EPSSetType()`
197: @*/
198: PetscErrorCode EPSGetType(EPS eps,EPSType *type)
199: {
200: PetscFunctionBegin;
202: PetscAssertPointer(type,2);
203: *type = ((PetscObject)eps)->type_name;
204: PetscFunctionReturn(PETSC_SUCCESS);
205: }
207: /*@C
208: EPSRegister - Adds a method to the 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: Note:
217: `EPSRegister()` may be called multiple times to add several user-defined solvers.
219: Example Usage:
220: .vb
221: EPSRegister("my_solver",MySolverCreate);
222: .ve
224: Then, your solver can be chosen with the procedural interface via
225: .vb
226: EPSSetType(eps,"my_solver")
227: .ve
228: or at runtime via the option `-eps_type my_solver`.
230: Level: advanced
232: .seealso: [](ch:eps), `EPSRegisterAll()`
233: @*/
234: PetscErrorCode EPSRegister(const char *name,PetscErrorCode (*function)(EPS))
235: {
236: PetscFunctionBegin;
237: PetscCall(EPSInitializePackage());
238: PetscCall(PetscFunctionListAdd(&EPSList,name,function));
239: PetscFunctionReturn(PETSC_SUCCESS);
240: }
242: /*@C
243: EPSMonitorRegister - Registers an `EPS` monitor routine that may be accessed with
244: `EPSMonitorSetFromOptions()`.
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, see `EPSMonitorRegisterFn`
253: . create - creation routine, or `NULL`
254: - destroy - destruction routine, or `NULL`
256: Notes:
257: `EPSMonitorRegister()` may be called multiple times to add several user-defined monitors.
259: The calling sequence for the given function matches the calling sequence of `EPSMonitorFn`
260: functions passed to `EPSMonitorSet()` with the additional requirement that its final argument
261: be a `PetscViewerAndFormat`.
263: Example Usage:
264: .vb
265: EPSMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
266: .ve
268: Then, your monitor can be chosen with the procedural interface via
269: .vb
270: EPSMonitorSetFromOptions(eps,"-eps_monitor_my_monitor","my_monitor",NULL)
271: .ve
272: or at runtime via the option `-eps_monitor_my_monitor`.
274: Level: advanced
276: .seealso: [](ch:eps), `EPSMonitorSet()`, `EPSMonitorRegisterAll()`, `EPSMonitorSetFromOptions()`
277: @*/
278: PetscErrorCode EPSMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,EPSMonitorRegisterFn *monitor,EPSMonitorRegisterCreateFn *create,EPSMonitorRegisterDestroyFn *destroy)
279: {
280: char key[PETSC_MAX_PATH_LEN];
282: PetscFunctionBegin;
283: PetscCall(EPSInitializePackage());
284: PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
285: PetscCall(PetscFunctionListAdd(&EPSMonitorList,key,monitor));
286: if (create) PetscCall(PetscFunctionListAdd(&EPSMonitorCreateList,key,create));
287: if (destroy) PetscCall(PetscFunctionListAdd(&EPSMonitorDestroyList,key,destroy));
288: PetscFunctionReturn(PETSC_SUCCESS);
289: }
291: /*@
292: EPSReset - Resets the `EPS` context to the initial state (prior to setup)
293: and destroys any allocated `Vec`s and `Mat`s.
295: Collective
297: Input Parameter:
298: . eps - the linear eigensolver context
300: Note:
301: This can be used when a problem of different matrix size wants to be solved.
302: All options that have previously been set are preserved, so in a next use
303: the solver configuration is the same, but new sizes for matrices and vectors
304: are allowed.
306: Level: advanced
308: .seealso: [](ch:eps), `EPSDestroy()`
309: @*/
310: PetscErrorCode EPSReset(EPS eps)
311: {
312: PetscFunctionBegin;
314: if (!eps) PetscFunctionReturn(PETSC_SUCCESS);
315: PetscTryTypeMethod(eps,reset);
316: if (eps->st) PetscCall(STReset(eps->st));
317: PetscCall(VecDestroy(&eps->D));
318: PetscCall(BVDestroy(&eps->V));
319: PetscCall(BVDestroy(&eps->W));
320: PetscCall(VecDestroyVecs(eps->nwork,&eps->work));
321: eps->nwork = 0;
322: eps->state = EPS_STATE_INITIAL;
323: PetscFunctionReturn(PETSC_SUCCESS);
324: }
326: /*@
327: EPSDestroy - Destroys the `EPS` context.
329: Collective
331: Input Parameter:
332: . eps - the linear eigensolver context
334: Level: beginner
336: .seealso: [](ch:eps), `EPSCreate()`, `EPSSetUp()`, `EPSSolve()`
337: @*/
338: PetscErrorCode EPSDestroy(EPS *eps)
339: {
340: PetscFunctionBegin;
341: if (!*eps) PetscFunctionReturn(PETSC_SUCCESS);
343: if (--((PetscObject)*eps)->refct > 0) { *eps = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
344: PetscCall(EPSReset(*eps));
345: PetscTryTypeMethod(*eps,destroy);
346: if ((*eps)->eigr) PetscCall(PetscFree4((*eps)->eigr,(*eps)->eigi,(*eps)->errest,(*eps)->perm));
347: if ((*eps)->rr) PetscCall(PetscFree2((*eps)->rr,(*eps)->ri));
348: PetscCall(STDestroy(&(*eps)->st));
349: PetscCall(RGDestroy(&(*eps)->rg));
350: PetscCall(DSDestroy(&(*eps)->ds));
351: PetscCall(PetscFree((*eps)->sc));
352: /* just in case the initial vectors have not been used */
353: PetscCall(SlepcBasisDestroy_Private(&(*eps)->nds,&(*eps)->defl));
354: PetscCall(SlepcBasisDestroy_Private(&(*eps)->nini,&(*eps)->IS));
355: PetscCall(SlepcBasisDestroy_Private(&(*eps)->ninil,&(*eps)->ISL));
356: if ((*eps)->convergeddestroy) PetscCall((*(*eps)->convergeddestroy)(&(*eps)->convergedctx));
357: if ((*eps)->stoppingdestroy) PetscCall((*(*eps)->stoppingdestroy)(&(*eps)->stoppingctx));
358: PetscCall(EPSMonitorCancel(*eps));
359: PetscCall(PetscHeaderDestroy(eps));
360: PetscFunctionReturn(PETSC_SUCCESS);
361: }
363: /*@
364: EPSSetTarget - Sets the value of the target.
366: Logically Collective
368: Input Parameters:
369: + eps - the linear eigensolver context
370: - target - the value of the target
372: Options Database Key:
373: . -eps_target \<target\> - the value of the target
375: Notes:
376: The target is a scalar value used to determine the portion of the spectrum
377: of interest. It is used in combination with `EPSSetWhichEigenpairs()`.
379: When PETSc is built with real scalars, it is not possible to specify a
380: complex target.
382: In the case of complex scalars, a complex value can be provided in the
383: command line with `[+/-][realnumber][+/-]realnumberi` with no spaces, e.g.
384: `-eps_target 1.0+2.0i`.
386: Level: intermediate
388: .seealso: [](ch:eps), `EPSGetTarget()`, `EPSSetWhichEigenpairs()`
389: @*/
390: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
391: {
392: PetscFunctionBegin;
395: eps->target = target;
396: if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
397: PetscCall(STSetDefaultShift(eps->st,target));
398: PetscFunctionReturn(PETSC_SUCCESS);
399: }
401: /*@
402: EPSGetTarget - Gets the value of the target.
404: Not Collective
406: Input Parameter:
407: . eps - the linear eigensolver context
409: Output Parameter:
410: . target - the value of the target
412: Note:
413: If the target was not set by the user, then zero is returned.
415: Level: intermediate
417: .seealso: [](ch:eps), `EPSSetTarget()`
418: @*/
419: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
420: {
421: PetscFunctionBegin;
423: PetscAssertPointer(target,2);
424: *target = eps->target;
425: PetscFunctionReturn(PETSC_SUCCESS);
426: }
428: /*@
429: EPSSetInterval - Defines the computational interval for spectrum slicing.
431: Logically Collective
433: Input Parameters:
434: + eps - the linear eigensolver context
435: . inta - left end of the interval
436: - intb - right end of the interval
438: Options Database Key:
439: . -eps_interval <a,b> - set $[a,b]$ as the interval of interest
441: Notes:
442: Spectrum slicing is a technique employed for computing all eigenvalues of
443: symmetric eigenproblems in a given interval, see section [](#sec:slice).
444: This function provides the interval to be considered. It must be used in
445: combination with `EPS_ALL`, see `EPSSetWhichEigenpairs()`.
447: A computational interval is also needed when using polynomial filters,
448: see `STFILTER` and section [](#sec:filter).
450: In the command-line option, two values must be provided. For an open interval,
451: one can give an infinite, e.g., `-eps_interval 1.0,inf` or `-eps_interval -inf,1.0`.
452: An open interval in the programmatic interface can be specified with
453: `PETSC_MAX_REAL` and -`PETSC_MAX_REAL`.
455: Level: intermediate
457: .seealso: [](ch:eps), [](#sec:slice), [](#sec:filter), `EPSGetInterval()`, `EPSSetWhichEigenpairs()`, `STFILTER`
458: @*/
459: PetscErrorCode EPSSetInterval(EPS eps,PetscReal inta,PetscReal intb)
460: {
461: PetscFunctionBegin;
465: PetscCheck(inta<intb,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
466: if (eps->inta != inta || eps->intb != intb) {
467: eps->inta = inta;
468: eps->intb = intb;
469: eps->state = EPS_STATE_INITIAL;
470: }
471: PetscFunctionReturn(PETSC_SUCCESS);
472: }
474: /*@
475: EPSGetInterval - Gets the computational interval for spectrum slicing.
477: Not Collective
479: Input Parameter:
480: . eps - the linear eigensolver context
482: Output Parameters:
483: + inta - left end of the interval
484: - intb - right end of the interval
486: Level: intermediate
488: Note:
489: If the interval was not set by the user, then zeros are returned.
491: .seealso: [](ch:eps), `EPSSetInterval()`
492: @*/
493: PetscErrorCode EPSGetInterval(EPS eps,PetscReal* inta,PetscReal* intb)
494: {
495: PetscFunctionBegin;
497: if (inta) *inta = eps->inta;
498: if (intb) *intb = eps->intb;
499: PetscFunctionReturn(PETSC_SUCCESS);
500: }
502: /*@
503: EPSSetST - Associates a spectral transformation object to the eigensolver.
505: Collective
507: Input Parameters:
508: + eps - the linear eigensolver context
509: - st - the spectral transformation object
511: Note:
512: Use `EPSGetST()` to retrieve the spectral transformation context at a later time
513: (for example, to free it at the end of the computations).
515: Level: advanced
517: .seealso: [](ch:eps), `EPSGetST()`
518: @*/
519: PetscErrorCode EPSSetST(EPS eps,ST st)
520: {
521: PetscFunctionBegin;
524: PetscCheckSameComm(eps,1,st,2);
525: PetscCall(PetscObjectReference((PetscObject)st));
526: PetscCall(STDestroy(&eps->st));
527: eps->st = st;
528: PetscFunctionReturn(PETSC_SUCCESS);
529: }
531: /*@
532: EPSGetST - Obtain the spectral transformation (`ST`) object associated
533: to the eigensolver object.
535: Not Collective
537: Input Parameter:
538: . eps - the linear eigensolver context
540: Output Parameter:
541: . st - spectral transformation context
543: Level: intermediate
545: .seealso: [](ch:eps), `EPSSetST()`
546: @*/
547: PetscErrorCode EPSGetST(EPS eps,ST *st)
548: {
549: PetscFunctionBegin;
551: PetscAssertPointer(st,2);
552: if (!eps->st) {
553: PetscCall(STCreate(PetscObjectComm((PetscObject)eps),&eps->st));
554: PetscCall(PetscObjectIncrementTabLevel((PetscObject)eps->st,(PetscObject)eps,0));
555: PetscCall(PetscObjectSetOptions((PetscObject)eps->st,((PetscObject)eps)->options));
556: }
557: *st = eps->st;
558: PetscFunctionReturn(PETSC_SUCCESS);
559: }
561: /*@
562: EPSSetBV - Associates a basis vectors object to the eigensolver.
564: Collective
566: Input Parameters:
567: + eps - the linear eigensolver context
568: - V - the basis vectors object
570: Level: advanced
572: .seealso: [](ch:eps), `EPSGetBV()`
573: @*/
574: PetscErrorCode EPSSetBV(EPS eps,BV V)
575: {
576: PetscFunctionBegin;
579: PetscCheckSameComm(eps,1,V,2);
580: PetscCall(PetscObjectReference((PetscObject)V));
581: PetscCall(BVDestroy(&eps->V));
582: eps->V = V;
583: PetscFunctionReturn(PETSC_SUCCESS);
584: }
586: /*@
587: EPSGetBV - Obtain the basis vectors object associated to the eigensolver object.
589: Not Collective
591: Input Parameter:
592: . eps - the linear eigensolver context
594: Output Parameter:
595: . V - basis vectors context
597: Level: advanced
599: .seealso: [](ch:eps), `EPSSetBV()`
600: @*/
601: PetscErrorCode EPSGetBV(EPS eps,BV *V)
602: {
603: PetscFunctionBegin;
605: PetscAssertPointer(V,2);
606: if (!eps->V) {
607: PetscCall(BVCreate(PetscObjectComm((PetscObject)eps),&eps->V));
608: PetscCall(PetscObjectIncrementTabLevel((PetscObject)eps->V,(PetscObject)eps,0));
609: PetscCall(PetscObjectSetOptions((PetscObject)eps->V,((PetscObject)eps)->options));
610: }
611: *V = eps->V;
612: PetscFunctionReturn(PETSC_SUCCESS);
613: }
615: /*@
616: EPSSetRG - Associates a region object to the eigensolver.
618: Collective
620: Input Parameters:
621: + eps - the linear eigensolver context
622: - rg - the region object
624: Note:
625: Use `EPSGetRG()` to retrieve the region context at a later time (for example,
626: to free it at the end of the computations).
628: Level: advanced
630: .seealso: [](ch:eps), `EPSGetRG()`
631: @*/
632: PetscErrorCode EPSSetRG(EPS eps,RG rg)
633: {
634: PetscFunctionBegin;
636: if (rg) {
638: PetscCheckSameComm(eps,1,rg,2);
639: }
640: PetscCall(PetscObjectReference((PetscObject)rg));
641: PetscCall(RGDestroy(&eps->rg));
642: eps->rg = rg;
643: PetscFunctionReturn(PETSC_SUCCESS);
644: }
646: /*@
647: EPSGetRG - Obtain the region object associated to the eigensolver.
649: Not Collective
651: Input Parameter:
652: . eps - the linear eigensolver context
654: Output Parameter:
655: . rg - region context
657: Level: advanced
659: .seealso: [](ch:eps), `EPSSetRG()`
660: @*/
661: PetscErrorCode EPSGetRG(EPS eps,RG *rg)
662: {
663: PetscFunctionBegin;
665: PetscAssertPointer(rg,2);
666: if (!eps->rg) {
667: PetscCall(RGCreate(PetscObjectComm((PetscObject)eps),&eps->rg));
668: PetscCall(PetscObjectIncrementTabLevel((PetscObject)eps->rg,(PetscObject)eps,0));
669: PetscCall(PetscObjectSetOptions((PetscObject)eps->rg,((PetscObject)eps)->options));
670: }
671: *rg = eps->rg;
672: PetscFunctionReturn(PETSC_SUCCESS);
673: }
675: /*@
676: EPSSetDS - Associates a direct solver object to the eigensolver.
678: Collective
680: Input Parameters:
681: + eps - the linear eigensolver context
682: - ds - the direct solver object
684: Note:
685: Use `EPSGetDS()` to retrieve the direct solver context at a later time (for example,
686: to free it at the end of the computations).
688: Level: advanced
690: .seealso: [](ch:eps), `EPSGetDS()`
691: @*/
692: PetscErrorCode EPSSetDS(EPS eps,DS ds)
693: {
694: PetscFunctionBegin;
697: PetscCheckSameComm(eps,1,ds,2);
698: PetscCall(PetscObjectReference((PetscObject)ds));
699: PetscCall(DSDestroy(&eps->ds));
700: eps->ds = ds;
701: PetscFunctionReturn(PETSC_SUCCESS);
702: }
704: /*@
705: EPSGetDS - Obtain the direct solver object associated to the eigensolver object.
707: Not Collective
709: Input Parameter:
710: . eps - the linear eigensolver context
712: Output Parameter:
713: . ds - direct solver context
715: Level: advanced
717: .seealso: [](ch:eps), `EPSSetDS()`
718: @*/
719: PetscErrorCode EPSGetDS(EPS eps,DS *ds)
720: {
721: PetscFunctionBegin;
723: PetscAssertPointer(ds,2);
724: if (!eps->ds) {
725: PetscCall(DSCreate(PetscObjectComm((PetscObject)eps),&eps->ds));
726: PetscCall(PetscObjectIncrementTabLevel((PetscObject)eps->ds,(PetscObject)eps,0));
727: PetscCall(PetscObjectSetOptions((PetscObject)eps->ds,((PetscObject)eps)->options));
728: }
729: *ds = eps->ds;
730: PetscFunctionReturn(PETSC_SUCCESS);
731: }
733: /*@
734: EPSIsGeneralized - Ask if the `EPS` object corresponds to a generalized
735: eigenvalue problem.
737: Not Collective
739: Input Parameter:
740: . eps - the linear eigensolver context
742: Output Parameter:
743: . is - `PETSC_TRUE` if the problem is generalized
745: Level: intermediate
747: .seealso: [](ch:eps), `EPSIsHermitian()`, `EPSIsPositive()`, `EPSIsStructured()`
748: @*/
749: PetscErrorCode EPSIsGeneralized(EPS eps,PetscBool* is)
750: {
751: PetscFunctionBegin;
753: PetscAssertPointer(is,2);
754: *is = eps->isgeneralized;
755: PetscFunctionReturn(PETSC_SUCCESS);
756: }
758: /*@
759: EPSIsHermitian - Ask if the `EPS` object corresponds to a Hermitian
760: eigenvalue problem.
762: Not Collective
764: Input Parameter:
765: . eps - the linear eigensolver context
767: Output Parameter:
768: . is - `PETSC_TRUE` if the problem is Hermitian
770: Level: intermediate
772: .seealso: [](ch:eps), `EPSIsGeneralized()`, `EPSIsPositive()`, `EPSIsStructured()`
773: @*/
774: PetscErrorCode EPSIsHermitian(EPS eps,PetscBool* is)
775: {
776: PetscFunctionBegin;
778: PetscAssertPointer(is,2);
779: *is = eps->ishermitian;
780: PetscFunctionReturn(PETSC_SUCCESS);
781: }
783: /*@
784: EPSIsPositive - Ask if the `EPS` object corresponds to an eigenvalue
785: problem type that requires a positive (semi-) definite matrix $B$.
787: Not Collective
789: Input Parameter:
790: . eps - the linear eigensolver context
792: Output Parameter:
793: . is - `PETSC_TRUE` if the problem is positive (semi-) definite
795: Level: intermediate
797: .seealso: [](ch:eps), `EPSIsGeneralized()`, `EPSIsHermitian()`, `EPSIsStructured()`
798: @*/
799: PetscErrorCode EPSIsPositive(EPS eps,PetscBool* is)
800: {
801: PetscFunctionBegin;
803: PetscAssertPointer(is,2);
804: *is = eps->ispositive;
805: PetscFunctionReturn(PETSC_SUCCESS);
806: }
808: /*@
809: EPSIsStructured - Ask if the `EPS` object corresponds to a structured
810: eigenvalue problem.
812: Not Collective
814: Input Parameter:
815: . eps - the linear eigensolver context
817: Output Parameter:
818: . is - `PETSC_TRUE` if the problem is structured
820: Note:
821: The result will be true if the problem type has been set to some
822: structured type such as `EPS_BSE`. This is independent of whether the
823: input matrix has been built with a certain structure with a helper function.
825: Level: intermediate
827: .seealso: [](ch:eps), `EPSIsGeneralized()`, `EPSIsHermitian()`, `EPSIsPositive()`, `EPSSetProblemType()`
828: @*/
829: PetscErrorCode EPSIsStructured(EPS eps,PetscBool* is)
830: {
831: PetscFunctionBegin;
833: PetscAssertPointer(is,2);
834: *is = eps->isstructured;
835: PetscFunctionReturn(PETSC_SUCCESS);
836: }