Actual source code: pepbasic.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 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 = NULL;
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 `PEP` context.
33: Collective
35: Input Parameter:
36: . comm - MPI communicator
38: Output Parameter:
39: . outpep - location to put the `PEP` context
41: Note:
42: The default `PEP` type is `PEPTOAR`.
44: Level: beginner
46: .seealso: [](ch:pep), `PEPSetUp()`, `PEPSolve()`, `PEPDestroy()`, `PEP`
47: @*/
48: PetscErrorCode PEPCreate(MPI_Comm comm,PEP *outpep)
49: {
50: PEP pep;
52: PetscFunctionBegin;
53: PetscAssertPointer(outpep,2);
54: PetscCall(PEPInitializePackage());
55: PetscCall(SlepcHeaderCreate(pep,PEP_CLASSID,"PEP","Polynomial Eigenvalue Problem","PEP",comm,PEPDestroy,PEPView));
57: pep->max_it = PETSC_DETERMINE;
58: pep->nev = 1;
59: pep->ncv = PETSC_DETERMINE;
60: pep->mpd = PETSC_DETERMINE;
61: pep->nini = 0;
62: pep->target = 0.0;
63: pep->tol = PETSC_DETERMINE;
64: pep->conv = PEP_CONV_REL;
65: pep->stop = PEP_STOP_BASIC;
66: pep->which = (PEPWhich)0;
67: pep->basis = PEP_BASIS_MONOMIAL;
68: pep->problem_type = (PEPProblemType)0;
69: pep->scale = PEP_SCALE_NONE;
70: pep->sfactor = 1.0;
71: pep->dsfactor = 1.0;
72: pep->sits = 5;
73: pep->slambda = 1.0;
74: pep->refine = PEP_REFINE_NONE;
75: pep->npart = 1;
76: pep->rtol = PETSC_DETERMINE;
77: pep->rits = PETSC_DETERMINE;
78: pep->scheme = (PEPRefineScheme)0;
79: pep->extract = (PEPExtract)0;
80: pep->trackall = PETSC_FALSE;
82: pep->converged = PEPConvergedRelative;
83: pep->convergeduser = NULL;
84: pep->convergeddestroy= NULL;
85: pep->stopping = PEPStoppingBasic;
86: pep->stoppinguser = NULL;
87: pep->stoppingdestroy = NULL;
88: pep->convergedctx = NULL;
89: pep->stoppingctx = NULL;
90: pep->numbermonitors = 0;
92: pep->st = NULL;
93: pep->ds = NULL;
94: pep->V = NULL;
95: pep->rg = NULL;
96: pep->A = NULL;
97: pep->nmat = 0;
98: pep->Dl = NULL;
99: pep->Dr = NULL;
100: pep->IS = NULL;
101: pep->eigr = NULL;
102: pep->eigi = NULL;
103: pep->errest = NULL;
104: pep->perm = NULL;
105: pep->pbc = NULL;
106: pep->solvematcoeffs = NULL;
107: pep->nwork = 0;
108: pep->work = NULL;
109: pep->refineksp = NULL;
110: pep->refinesubc = NULL;
111: pep->data = NULL;
113: pep->state = PEP_STATE_INITIAL;
114: pep->nconv = 0;
115: pep->its = 0;
116: pep->n = 0;
117: pep->nloc = 0;
118: pep->nrma = NULL;
119: pep->sfactor_set = PETSC_FALSE;
120: pep->lineariz = PETSC_FALSE;
121: pep->reason = PEP_CONVERGED_ITERATING;
123: PetscCall(PetscNew(&pep->sc));
124: *outpep = pep;
125: PetscFunctionReturn(PETSC_SUCCESS);
126: }
128: /*@
129: PEPSetType - Selects the particular solver to be used in the `PEP` object.
131: Logically Collective
133: Input Parameters:
134: + pep - the polynomial eigensolver context
135: - type - a known method
137: Options Database Key:
138: . -pep_type type - sets the method; use `-help` for a list of available methods
140: Notes:
141: See `PEPType` for available methods. The default is `PEPTOAR`.
143: Normally, it is best to use the `PEPSetFromOptions()` command and
144: then set the `PEP` type from the options database rather than by using
145: this routine. Using the options database provides the user with
146: maximum flexibility in evaluating the different available methods.
147: The `PEPSetType()` routine is provided for those situations where it
148: is necessary to set the iterative solver independently of the command
149: line or options database.
151: Level: intermediate
153: .seealso: [](ch:pep), `PEPType`
154: @*/
155: PetscErrorCode PEPSetType(PEP pep,PEPType type)
156: {
157: PetscErrorCode (*r)(PEP);
158: PetscBool match;
160: PetscFunctionBegin;
162: PetscAssertPointer(type,2);
164: PetscCall(PetscObjectTypeCompare((PetscObject)pep,type,&match));
165: if (match) PetscFunctionReturn(PETSC_SUCCESS);
167: PetscCall(PetscFunctionListFind(PEPList,type,&r));
168: PetscCheck(r,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown PEP type given: %s",type);
170: PetscTryTypeMethod(pep,destroy);
171: PetscCall(PetscMemzero(pep->ops,sizeof(struct _PEPOps)));
173: pep->state = PEP_STATE_INITIAL;
174: PetscCall(PetscObjectChangeTypeName((PetscObject)pep,type));
175: PetscCall((*r)(pep));
176: PetscFunctionReturn(PETSC_SUCCESS);
177: }
179: /*@
180: PEPGetType - Gets the `PEP` type as a string from the `PEP` object.
182: Not Collective
184: Input Parameter:
185: . pep - the polynomial eigensolver context
187: Output Parameter:
188: . type - name of `PEP` method
190: Note:
191: `type` should not be retained for later use as it will be an invalid pointer
192: if the `PEPType` of `pep` is changed.
194: Level: intermediate
196: .seealso: [](ch:pep), `PEPSetType()`, `PetscObjectTypeCompare()`, `PetscObjectTypeCompareAny()`
197: @*/
198: PetscErrorCode PEPGetType(PEP pep,PEPType *type)
199: {
200: PetscFunctionBegin;
202: PetscAssertPointer(type,2);
203: *type = ((PetscObject)pep)->type_name;
204: PetscFunctionReturn(PETSC_SUCCESS);
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: Note:
217: `PEPRegister()` may be called multiple times to add several user-defined solvers.
219: Example Usage:
220: .vb
221: PEPRegister("my_solver",MySolverCreate);
222: .ve
224: Then, your solver can be chosen with the procedural interface via
225: .vb
226: PEPSetType(pep,"my_solver")
227: .ve
228: or at runtime via the option `-pep_type my_solver`.
230: Level: advanced
232: .seealso: [](ch:pep), `PEPRegisterAll()`
233: @*/
234: PetscErrorCode PEPRegister(const char *name,PetscErrorCode (*function)(PEP))
235: {
236: PetscFunctionBegin;
237: PetscCall(PEPInitializePackage());
238: PetscCall(PetscFunctionListAdd(&PEPList,name,function));
239: PetscFunctionReturn(PETSC_SUCCESS);
240: }
242: /*@C
243: PEPMonitorRegister - Registers a `PEP` monitor routine that may be accessed with
244: `PEPMonitorSetFromOptions()`.
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 `PEPMonitorRegisterFn`
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: The calling sequence for the given function matches the calling sequence of `PEPMonitorFn`
260: functions passed to `PEPMonitorSet()` with the additional requirement that its final argument
261: be a `PetscViewerAndFormat`.
263: Example Usage:
264: .vb
265: PEPMonitorRegister("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: PEPMonitorSetFromOptions(pep,"-pep_monitor_my_monitor","my_monitor",NULL);
271: .ve
272: or at runtime via the option `-pep_monitor_my_monitor`.
274: Level: advanced
276: .seealso: [](ch:pep), `PEPMonitorSet()`, `PEPMonitorRegisterAll()`, `PEPMonitorSetFromOptions()`
277: @*/
278: PetscErrorCode PEPMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,PEPMonitorRegisterFn *monitor,PEPMonitorRegisterCreateFn *create,PEPMonitorRegisterDestroyFn *destroy)
279: {
280: char key[PETSC_MAX_PATH_LEN];
282: PetscFunctionBegin;
283: PetscCall(PEPInitializePackage());
284: PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
285: PetscCall(PetscFunctionListAdd(&PEPMonitorList,key,monitor));
286: if (create) PetscCall(PetscFunctionListAdd(&PEPMonitorCreateList,key,create));
287: if (destroy) PetscCall(PetscFunctionListAdd(&PEPMonitorDestroyList,key,destroy));
288: PetscFunctionReturn(PETSC_SUCCESS);
289: }
291: /*@
292: PEPReset - Resets the `PEP` context to the initial state (prior to setup)
293: and destroys any allocated `Vec`s and `Mat`s.
295: Collective
297: Input Parameter:
298: . pep - the polynomial eigensolver context
300: Level: advanced
302: .seealso: [](ch:pep), `PEPDestroy()`
303: @*/
304: PetscErrorCode PEPReset(PEP pep)
305: {
306: PetscFunctionBegin;
308: if (!pep) PetscFunctionReturn(PETSC_SUCCESS);
309: PetscTryTypeMethod(pep,reset);
310: if (pep->st) PetscCall(STReset(pep->st));
311: if (pep->refineksp) PetscCall(KSPReset(pep->refineksp));
312: if (pep->nmat) {
313: PetscCall(MatDestroyMatrices(pep->nmat,&pep->A));
314: PetscCall(PetscFree2(pep->pbc,pep->nrma));
315: PetscCall(PetscFree(pep->solvematcoeffs));
316: pep->nmat = 0;
317: }
318: PetscCall(VecDestroy(&pep->Dl));
319: PetscCall(VecDestroy(&pep->Dr));
320: PetscCall(BVDestroy(&pep->V));
321: PetscCall(VecDestroyVecs(pep->nwork,&pep->work));
322: pep->nwork = 0;
323: pep->state = PEP_STATE_INITIAL;
324: PetscFunctionReturn(PETSC_SUCCESS);
325: }
327: /*@
328: PEPDestroy - Destroys the `PEP` context.
330: Collective
332: Input Parameter:
333: . pep - the polynomial eigensolver context
335: Level: beginner
337: .seealso: [](ch:pep), `PEPCreate()`, `PEPSetUp()`, `PEPSolve()`
338: @*/
339: PetscErrorCode PEPDestroy(PEP *pep)
340: {
341: PetscFunctionBegin;
342: if (!*pep) PetscFunctionReturn(PETSC_SUCCESS);
344: if (--((PetscObject)*pep)->refct > 0) { *pep = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
345: PetscCall(PEPReset(*pep));
346: PetscTryTypeMethod(*pep,destroy);
347: if ((*pep)->eigr) PetscCall(PetscFree4((*pep)->eigr,(*pep)->eigi,(*pep)->errest,(*pep)->perm));
348: PetscCall(STDestroy(&(*pep)->st));
349: PetscCall(RGDestroy(&(*pep)->rg));
350: PetscCall(DSDestroy(&(*pep)->ds));
351: PetscCall(KSPDestroy(&(*pep)->refineksp));
352: PetscCall(PetscSubcommDestroy(&(*pep)->refinesubc));
353: PetscCall(PetscFree((*pep)->sc));
354: /* just in case the initial vectors have not been used */
355: PetscCall(SlepcBasisDestroy_Private(&(*pep)->nini,&(*pep)->IS));
356: if ((*pep)->convergeddestroy) PetscCall((*(*pep)->convergeddestroy)(&(*pep)->convergedctx));
357: if ((*pep)->stoppingdestroy) PetscCall((*(*pep)->stoppingdestroy)(&(*pep)->stoppingctx));
358: PetscCall(PEPMonitorCancel(*pep));
359: PetscCall(PetscHeaderDestroy(pep));
360: PetscFunctionReturn(PETSC_SUCCESS);
361: }
363: /*@
364: PEPSetBV - Associates a basis vectors object to the polynomial eigensolver.
366: Collective
368: Input Parameters:
369: + pep - the polynomial eigensolver context
370: - bv - the basis vectors object
372: Note:
373: Use `PEPGetBV()` to retrieve the basis vectors context (for example,
374: to free it at the end of the computations).
376: Level: advanced
378: .seealso: [](ch:pep), `PEPGetBV()`
379: @*/
380: PetscErrorCode PEPSetBV(PEP pep,BV bv)
381: {
382: PetscFunctionBegin;
385: PetscCheckSameComm(pep,1,bv,2);
386: PetscCall(PetscObjectReference((PetscObject)bv));
387: PetscCall(BVDestroy(&pep->V));
388: pep->V = bv;
389: PetscFunctionReturn(PETSC_SUCCESS);
390: }
392: /*@
393: PEPGetBV - Obtain the basis vectors object associated to the polynomial
394: eigensolver object.
396: Not Collective
398: Input Parameter:
399: . pep - the polynomial eigensolver context
401: Output Parameter:
402: . bv - basis vectors context
404: Level: advanced
406: .seealso: [](ch:pep), `PEPSetBV()`
407: @*/
408: PetscErrorCode PEPGetBV(PEP pep,BV *bv)
409: {
410: PetscFunctionBegin;
412: PetscAssertPointer(bv,2);
413: if (!pep->V) {
414: PetscCall(BVCreate(PetscObjectComm((PetscObject)pep),&pep->V));
415: PetscCall(PetscObjectIncrementTabLevel((PetscObject)pep->V,(PetscObject)pep,0));
416: PetscCall(PetscObjectSetOptions((PetscObject)pep->V,((PetscObject)pep)->options));
417: }
418: *bv = pep->V;
419: PetscFunctionReturn(PETSC_SUCCESS);
420: }
422: /*@
423: PEPSetRG - Associates a region object to the polynomial eigensolver.
425: Collective
427: Input Parameters:
428: + pep - the polynomial eigensolver context
429: - rg - the region object
431: Note:
432: Use `PEPGetRG()` to retrieve the region context (for example,
433: to free it at the end of the computations).
435: Level: advanced
437: .seealso: [](ch:pep), `PEPGetRG()`
438: @*/
439: PetscErrorCode PEPSetRG(PEP pep,RG rg)
440: {
441: PetscFunctionBegin;
443: if (rg) {
445: PetscCheckSameComm(pep,1,rg,2);
446: }
447: PetscCall(PetscObjectReference((PetscObject)rg));
448: PetscCall(RGDestroy(&pep->rg));
449: pep->rg = rg;
450: PetscFunctionReturn(PETSC_SUCCESS);
451: }
453: /*@
454: PEPGetRG - Obtain the region object associated to the
455: polynomial eigensolver object.
457: Not Collective
459: Input Parameter:
460: . pep - the polynomial eigensolver context
462: Output Parameter:
463: . rg - region context
465: Level: advanced
467: .seealso: [](ch:pep), `PEPSetRG()`
468: @*/
469: PetscErrorCode PEPGetRG(PEP pep,RG *rg)
470: {
471: PetscFunctionBegin;
473: PetscAssertPointer(rg,2);
474: if (!pep->rg) {
475: PetscCall(RGCreate(PetscObjectComm((PetscObject)pep),&pep->rg));
476: PetscCall(PetscObjectIncrementTabLevel((PetscObject)pep->rg,(PetscObject)pep,0));
477: PetscCall(PetscObjectSetOptions((PetscObject)pep->rg,((PetscObject)pep)->options));
478: }
479: *rg = pep->rg;
480: PetscFunctionReturn(PETSC_SUCCESS);
481: }
483: /*@
484: PEPSetDS - Associates a direct solver object to the polynomial eigensolver.
486: Collective
488: Input Parameters:
489: + pep - the polynomial eigensolver context
490: - ds - the direct solver object
492: Note:
493: Use `PEPGetDS()` to retrieve the direct solver context (for example,
494: to free it at the end of the computations).
496: Level: advanced
498: .seealso: [](ch:pep), `PEPGetDS()`
499: @*/
500: PetscErrorCode PEPSetDS(PEP pep,DS ds)
501: {
502: PetscFunctionBegin;
505: PetscCheckSameComm(pep,1,ds,2);
506: PetscCall(PetscObjectReference((PetscObject)ds));
507: PetscCall(DSDestroy(&pep->ds));
508: pep->ds = ds;
509: PetscFunctionReturn(PETSC_SUCCESS);
510: }
512: /*@
513: PEPGetDS - Obtain the direct solver object associated to the
514: polynomial eigensolver object.
516: Not Collective
518: Input Parameter:
519: . pep - the polynomial eigensolver context
521: Output Parameter:
522: . ds - direct solver context
524: Level: advanced
526: .seealso: [](ch:pep), `PEPSetDS()`
527: @*/
528: PetscErrorCode PEPGetDS(PEP pep,DS *ds)
529: {
530: PetscFunctionBegin;
532: PetscAssertPointer(ds,2);
533: if (!pep->ds) {
534: PetscCall(DSCreate(PetscObjectComm((PetscObject)pep),&pep->ds));
535: PetscCall(PetscObjectIncrementTabLevel((PetscObject)pep->ds,(PetscObject)pep,0));
536: PetscCall(PetscObjectSetOptions((PetscObject)pep->ds,((PetscObject)pep)->options));
537: }
538: *ds = pep->ds;
539: PetscFunctionReturn(PETSC_SUCCESS);
540: }
542: /*@
543: PEPSetST - Associates a spectral transformation object to the eigensolver.
545: Collective
547: Input Parameters:
548: + pep - the polynomial eigensolver context
549: - st - the spectral transformation object
551: Note:
552: Use `PEPGetST()` to retrieve the spectral transformation context (for example,
553: to free it at the end of the computations).
555: Level: advanced
557: .seealso: [](ch:pep), `PEPGetST()`
558: @*/
559: PetscErrorCode PEPSetST(PEP pep,ST st)
560: {
561: PetscFunctionBegin;
564: PetscCheckSameComm(pep,1,st,2);
565: PetscCall(PetscObjectReference((PetscObject)st));
566: PetscCall(STDestroy(&pep->st));
567: pep->st = st;
568: PetscFunctionReturn(PETSC_SUCCESS);
569: }
571: /*@
572: PEPGetST - Obtain the spectral transformation (`ST`) object associated
573: to the eigensolver object.
575: Not Collective
577: Input Parameter:
578: . pep - the polynomial eigensolver context
580: Output Parameter:
581: . st - spectral transformation context
583: Level: intermediate
585: .seealso: [](ch:pep), `PEPSetST()`
586: @*/
587: PetscErrorCode PEPGetST(PEP pep,ST *st)
588: {
589: PetscFunctionBegin;
591: PetscAssertPointer(st,2);
592: if (!pep->st) {
593: PetscCall(STCreate(PetscObjectComm((PetscObject)pep),&pep->st));
594: PetscCall(PetscObjectIncrementTabLevel((PetscObject)pep->st,(PetscObject)pep,0));
595: PetscCall(PetscObjectSetOptions((PetscObject)pep->st,((PetscObject)pep)->options));
596: }
597: *st = pep->st;
598: PetscFunctionReturn(PETSC_SUCCESS);
599: }
601: /*@
602: PEPRefineGetKSP - Obtain the `KSP` object used by the eigensolver
603: object in the refinement phase.
605: Collective
607: Input Parameter:
608: . pep - the polynomial eigensolver context
610: Output Parameter:
611: . ksp - the linear solver context
613: Level: advanced
615: .seealso: [](ch:pep), `PEPSetRefine()`
616: @*/
617: PetscErrorCode PEPRefineGetKSP(PEP pep,KSP *ksp)
618: {
619: MPI_Comm comm;
621: PetscFunctionBegin;
623: PetscAssertPointer(ksp,2);
624: if (!pep->refineksp) {
625: if (pep->npart>1) {
626: /* Split in subcomunicators */
627: PetscCall(PetscSubcommCreate(PetscObjectComm((PetscObject)pep),&pep->refinesubc));
628: PetscCall(PetscSubcommSetNumber(pep->refinesubc,pep->npart));
629: PetscCall(PetscSubcommSetType(pep->refinesubc,PETSC_SUBCOMM_CONTIGUOUS));
630: PetscCall(PetscSubcommGetChild(pep->refinesubc,&comm));
631: } else PetscCall(PetscObjectGetComm((PetscObject)pep,&comm));
632: PetscCall(KSPCreate(comm,&pep->refineksp));
633: PetscCall(PetscObjectIncrementTabLevel((PetscObject)pep->refineksp,(PetscObject)pep,0));
634: PetscCall(PetscObjectSetOptions((PetscObject)pep->refineksp,((PetscObject)pep)->options));
635: PetscCall(KSPSetOptionsPrefix(*ksp,((PetscObject)pep)->prefix));
636: PetscCall(KSPAppendOptionsPrefix(*ksp,"pep_refine_"));
637: PetscCall(KSPSetTolerances(pep->refineksp,SlepcDefaultTol(pep->rtol),PETSC_CURRENT,PETSC_CURRENT,PETSC_CURRENT));
638: }
639: *ksp = pep->refineksp;
640: PetscFunctionReturn(PETSC_SUCCESS);
641: }
643: /*@
644: PEPSetTarget - Sets the value of the target.
646: Logically Collective
648: Input Parameters:
649: + pep - the polynomial eigensolver context
650: - target - the value of the target
652: Options Database Key:
653: . -pep_target target - the value of the target
655: Notes:
656: The target is a scalar value used to determine the portion of the spectrum
657: of interest. It is used in combination with `PEPSetWhichEigenpairs()`.
659: When PETSc is built with real scalars, it is not possible to specify a
660: complex target.
662: In the case of complex scalars, a complex value can be provided in the
663: command line with `[+/-][realnumber][+/-]realnumberi` with no spaces, e.g.
664: `-pep_target 1.0+2.0i`.
666: Level: intermediate
668: .seealso: [](ch:pep), `PEPGetTarget()`, `PEPSetWhichEigenpairs()`
669: @*/
670: PetscErrorCode PEPSetTarget(PEP pep,PetscScalar target)
671: {
672: PetscFunctionBegin;
675: pep->target = target;
676: if (!pep->st) PetscCall(PEPGetST(pep,&pep->st));
677: PetscCall(STSetDefaultShift(pep->st,target));
678: PetscFunctionReturn(PETSC_SUCCESS);
679: }
681: /*@
682: PEPGetTarget - Gets the value of the target.
684: Not Collective
686: Input Parameter:
687: . pep - the polynomial eigensolver context
689: Output Parameter:
690: . target - the value of the target
692: Note:
693: If the target was not set by the user, then zero is returned.
695: Level: intermediate
697: .seealso: [](ch:pep), `PEPSetTarget()`
698: @*/
699: PetscErrorCode PEPGetTarget(PEP pep,PetscScalar* target)
700: {
701: PetscFunctionBegin;
703: PetscAssertPointer(target,2);
704: *target = pep->target;
705: PetscFunctionReturn(PETSC_SUCCESS);
706: }
708: /*@
709: PEPSetInterval - Defines the computational interval for spectrum slicing.
711: Logically Collective
713: Input Parameters:
714: + pep - the polynomial eigensolver context
715: . inta - left end of the interval
716: - intb - right end of the interval
718: Options Database Key:
719: . -pep_interval a,b - set $[a,b]$ as the interval of interest
721: Notes:
722: Spectrum slicing is a technique employed for computing all eigenvalues of
723: symmetric quadratic eigenproblems in a given interval, see section [](#sec:qslice).
724: This function provides the interval to be considered. It must be used in
725: combination with `PEP_ALL`, see `PEPSetWhichEigenpairs()`. Note that in
726: polynomial eigenproblems spectrum slicing is implemented in `PEPSTOAR` only.
728: In the command-line option, two values must be provided. For an open interval,
729: one can give an infinite, e.g., `-pep_interval 1.0,inf` or `-pep_interval -inf,1.0`.
730: An open interval in the programmatic interface can be specified with
731: `PETSC_MAX_REAL` and -`PETSC_MAX_REAL`.
733: Level: intermediate
735: .seealso: [](ch:pep), [](#sec:qslice), `PEPGetInterval()`, `PEPSetWhichEigenpairs()`, `PEPSTOAR`
736: @*/
737: PetscErrorCode PEPSetInterval(PEP pep,PetscReal inta,PetscReal intb)
738: {
739: PetscFunctionBegin;
743: PetscCheck(inta<intb,PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
744: if (pep->inta != inta || pep->intb != intb) {
745: pep->inta = inta;
746: pep->intb = intb;
747: pep->state = PEP_STATE_INITIAL;
748: }
749: PetscFunctionReturn(PETSC_SUCCESS);
750: }
752: /*@
753: PEPGetInterval - Gets the computational interval for spectrum slicing.
755: Not Collective
757: Input Parameter:
758: . pep - the polynomial eigensolver context
760: Output Parameters:
761: + inta - left end of the interval
762: - intb - right end of the interval
764: Level: intermediate
766: Note:
767: If the interval was not set by the user, then zeros are returned.
769: .seealso: [](ch:pep), `PEPSetInterval()`
770: @*/
771: PetscErrorCode PEPGetInterval(PEP pep,PetscReal* inta,PetscReal* intb)
772: {
773: PetscFunctionBegin;
775: if (inta) *inta = pep->inta;
776: if (intb) *intb = pep->intb;
777: PetscFunctionReturn(PETSC_SUCCESS);
778: }