Actual source code: pepbasic.c
slepc-3.5.2 2014-10-10
1: /*
2: The basic PEP routines, Create, View, etc. are here.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2014, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include <slepc-private/pepimpl.h> /*I "slepcpep.h" I*/
26: PetscFunctionList PEPList = 0;
27: PetscBool PEPRegisterAllCalled = PETSC_FALSE;
28: PetscClassId PEP_CLASSID = 0;
29: PetscLogEvent PEP_SetUp = 0,PEP_Solve = 0,PEP_Refine = 0;
33: /*@C
34: PEPView - Prints the PEP data structure.
36: Collective on PEP
38: Input Parameters:
39: + pep - the polynomial eigenproblem solver context
40: - viewer - optional visualization context
42: Options Database Key:
43: . -pep_view - Calls PEPView() at end of PEPSolve()
45: Note:
46: The available visualization contexts include
47: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
48: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
49: output where only the first processor opens
50: the file. All other processors send their
51: data to the first processor to print.
53: The user can open an alternative visualization context with
54: PetscViewerASCIIOpen() - output to a specified file.
56: Level: beginner
58: .seealso: PetscViewerASCIIOpen()
59: @*/
60: PetscErrorCode PEPView(PEP pep,PetscViewer viewer)
61: {
63: const char *type;
64: char str[50];
65: PetscBool isascii,islinear,istrivial;
66: PetscInt i;
70: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
74: #if defined(PETSC_USE_COMPLEX)
75: #define HERM "hermitian"
76: #else
77: #define HERM "symmetric"
78: #endif
79: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
80: if (isascii) {
81: PetscObjectPrintClassNamePrefixType((PetscObject)pep,viewer);
82: if (pep->ops->view) {
83: PetscViewerASCIIPushTab(viewer);
84: (*pep->ops->view)(pep,viewer);
85: PetscViewerASCIIPopTab(viewer);
86: }
87: if (pep->problem_type) {
88: switch (pep->problem_type) {
89: case PEP_GENERAL: type = "general polynomial eigenvalue problem"; break;
90: case PEP_HERMITIAN: type = HERM " polynomial eigenvalue problem"; break;
91: case PEP_GYROSCOPIC: type = "gyroscopic polynomial eigenvalue problem"; break;
92: default: SETERRQ(PetscObjectComm((PetscObject)pep),1,"Wrong value of pep->problem_type");
93: }
94: } else type = "not yet set";
95: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
96: PetscViewerASCIIPrintf(viewer," polynomial represented in %s basis\n",PEPBasisTypes[pep->basis]);
97: switch (pep->scale) {
98: case PEP_SCALE_NONE:
99: break;
100: case PEP_SCALE_SCALAR:
101: PetscViewerASCIIPrintf(viewer," scalar balancing enabled, with scaling factor=%g\n",(double)pep->sfactor);
102: break;
103: case PEP_SCALE_DIAGONAL:
104: PetscViewerASCIIPrintf(viewer," diagonal balancing enabled, with its=%D and lambda=%g\n",pep->sits,(double)pep->slambda);
105: break;
106: case PEP_SCALE_BOTH:
107: PetscViewerASCIIPrintf(viewer," scalar & diagonal balancing enabled, with scaling factor=%g, its=%D and lambda=%g\n",(double)pep->sfactor,pep->sits,(double)pep->slambda);
108: break;
109: }
110: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",PEPExtractTypes[pep->extract]);
111: PetscViewerASCIIPrintf(viewer," iterative refinement: %s%s\n",PEPRefineTypes[pep->refine],pep->schur?", with a Schur complement approach":"");
112: if (pep->refine) {
113: PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%D\n",(double)pep->rtol,pep->rits);
114: if (pep->npart>1) {
115: PetscViewerASCIIPrintf(viewer," splitting communicator in %D partitions for refinement\n",pep->npart);
116: }
117: }
118: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
119: SlepcSNPrintfScalar(str,50,pep->target,PETSC_FALSE);
120: if (!pep->which) {
121: PetscViewerASCIIPrintf(viewer,"not yet set\n");
122: } else switch (pep->which) {
123: case PEP_TARGET_MAGNITUDE:
124: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
125: break;
126: case PEP_TARGET_REAL:
127: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
128: break;
129: case PEP_TARGET_IMAGINARY:
130: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
131: break;
132: case PEP_LARGEST_MAGNITUDE:
133: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
134: break;
135: case PEP_SMALLEST_MAGNITUDE:
136: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
137: break;
138: case PEP_LARGEST_REAL:
139: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
140: break;
141: case PEP_SMALLEST_REAL:
142: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
143: break;
144: case PEP_LARGEST_IMAGINARY:
145: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
146: break;
147: case PEP_SMALLEST_IMAGINARY:
148: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
149: break;
150: default: SETERRQ(PetscObjectComm((PetscObject)pep),1,"Wrong value of pep->which");
151: }
152: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",pep->nev);
153: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",pep->ncv);
154: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",pep->mpd);
155: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",pep->max_it);
156: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)pep->tol);
157: PetscViewerASCIIPrintf(viewer," convergence test: ");
158: switch (pep->conv) {
159: case PEP_CONV_ABS:
160: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
161: case PEP_CONV_EIG:
162: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
163: case PEP_CONV_NORM:
164: PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
165: if (pep->nrma) {
166: PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)pep->nrma[0]);
167: for (i=1;i<pep->nmat;i++) {
168: PetscViewerASCIIPrintf(viewer,", %g",(double)pep->nrma[i]);
169: }
170: PetscViewerASCIIPrintf(viewer,"\n");
171: }
172: break;
173: case PEP_CONV_USER:
174: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
175: }
176: if (pep->nini) {
177: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(pep->nini));
178: }
179: } else {
180: if (pep->ops->view) {
181: (*pep->ops->view)(pep,viewer);
182: }
183: }
184: PetscObjectTypeCompare((PetscObject)pep,PEPLINEAR,&islinear);
185: if (!islinear) {
186: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
187: if (!pep->V) { PEPGetBV(pep,&pep->V); }
188: BVView(pep->V,viewer);
189: if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
190: RGIsTrivial(pep->rg,&istrivial);
191: if (!istrivial) { RGView(pep->rg,viewer); }
192: if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
193: DSView(pep->ds,viewer);
194: PetscViewerPopFormat(viewer);
195: if (!pep->st) { PEPGetST(pep,&pep->st); }
196: STView(pep->st,viewer);
197: }
198: return(0);
199: }
203: /*@
204: PEPPrintSolution - Prints the computed eigenvalues.
206: Collective on PEP
208: Input Parameters:
209: + pep - the eigensolver context
210: - viewer - optional visualization context
212: Options Database Key:
213: . -pep_terse - print only minimal information
215: Note:
216: By default, this function prints a table with eigenvalues and associated
217: relative errors. With -pep_terse only the eigenvalues are printed.
219: Level: intermediate
221: .seealso: PetscViewerASCIIOpen()
222: @*/
223: PetscErrorCode PEPPrintSolution(PEP pep,PetscViewer viewer)
224: {
225: PetscBool terse,errok,isascii;
226: PetscReal error,re,im;
227: PetscScalar kr,ki;
228: PetscInt i,j;
233: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
236: PEPCheckSolved(pep,1);
237: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
238: if (!isascii) return(0);
240: PetscOptionsHasName(NULL,"-pep_terse",&terse);
241: if (terse) {
242: if (pep->nconv<pep->nev) {
243: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",pep->nev);
244: } else {
245: errok = PETSC_TRUE;
246: for (i=0;i<pep->nev;i++) {
247: PEPComputeRelativeError(pep,i,&error);
248: errok = (errok && error<5.0*pep->tol)? PETSC_TRUE: PETSC_FALSE;
249: }
250: if (errok) {
251: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
252: for (i=0;i<=(pep->nev-1)/8;i++) {
253: PetscViewerASCIIPrintf(viewer,"\n ");
254: for (j=0;j<PetscMin(8,pep->nev-8*i);j++) {
255: PEPGetEigenpair(pep,8*i+j,&kr,&ki,NULL,NULL);
256: #if defined(PETSC_USE_COMPLEX)
257: re = PetscRealPart(kr);
258: im = PetscImaginaryPart(kr);
259: #else
260: re = kr;
261: im = ki;
262: #endif
263: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
264: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
265: if (im!=0.0) {
266: PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
267: } else {
268: PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
269: }
270: if (8*i+j+1<pep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
271: }
272: }
273: PetscViewerASCIIPrintf(viewer,"\n\n");
274: } else {
275: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",pep->nev);
276: }
277: }
278: } else {
279: PetscViewerASCIIPrintf(viewer," Number of converged approximate eigenpairs: %D\n\n",pep->nconv);
280: if (pep->nconv>0) {
281: PetscViewerASCIIPrintf(viewer,
282: " k ||P(k)x||/||kx||\n"
283: " ----------------- -------------------------\n");
284: for (i=0;i<pep->nconv;i++) {
285: PEPGetEigenpair(pep,i,&kr,&ki,NULL,NULL);
286: PEPComputeRelativeError(pep,i,&error);
287: #if defined(PETSC_USE_COMPLEX)
288: re = PetscRealPart(kr);
289: im = PetscImaginaryPart(kr);
290: #else
291: re = kr;
292: im = ki;
293: #endif
294: if (im!=0.0) {
295: PetscViewerASCIIPrintf(viewer," % 9f%+9f i %12g\n",(double)re,(double)im,(double)error);
296: } else {
297: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
298: }
299: }
300: PetscViewerASCIIPrintf(viewer,"\n");
301: }
302: }
303: return(0);
304: }
308: /*@C
309: PEPCreate - Creates the default PEP context.
311: Collective on MPI_Comm
313: Input Parameter:
314: . comm - MPI communicator
316: Output Parameter:
317: . pep - location to put the PEP context
319: Note:
320: The default PEP type is PEPTOAR
322: Level: beginner
324: .seealso: PEPSetUp(), PEPSolve(), PEPDestroy(), PEP
325: @*/
326: PetscErrorCode PEPCreate(MPI_Comm comm,PEP *outpep)
327: {
329: PEP pep;
333: *outpep = 0;
334: PEPInitializePackage();
335: SlepcHeaderCreate(pep,_p_PEP,struct _PEPOps,PEP_CLASSID,"PEP","Polynomial Eigenvalue Problem","PEP",comm,PEPDestroy,PEPView);
337: pep->max_it = 0;
338: pep->nev = 1;
339: pep->ncv = 0;
340: pep->mpd = 0;
341: pep->nini = 0;
342: pep->target = 0.0;
343: pep->tol = PETSC_DEFAULT;
344: pep->conv = PEP_CONV_NORM;
345: pep->which = (PEPWhich)0;
346: pep->basis = PEP_BASIS_MONOMIAL;
347: pep->problem_type = (PEPProblemType)0;
348: pep->scale = PEP_SCALE_NONE;
349: pep->sfactor = 1.0;
350: pep->sits = 5;
351: pep->slambda = 1.0;
352: pep->refine = PEP_REFINE_NONE;
353: pep->npart = 1;
354: pep->rtol = PETSC_DEFAULT;
355: pep->rits = PETSC_DEFAULT;
356: pep->schur = PETSC_FALSE;
357: pep->extract = PEP_EXTRACT_NORM;
358: pep->trackall = PETSC_FALSE;
360: pep->converged = PEPConvergedNormRelative;
361: pep->convergeddestroy= NULL;
362: pep->convergedctx = NULL;
363: pep->numbermonitors = 0;
365: pep->st = NULL;
366: pep->ds = NULL;
367: pep->V = NULL;
368: pep->rg = NULL;
369: pep->rand = NULL;
370: pep->A = NULL;
371: pep->nmat = 0;
372: pep->Dl = NULL;
373: pep->Dr = NULL;
374: pep->IS = NULL;
375: pep->eigr = NULL;
376: pep->eigi = NULL;
377: pep->errest = NULL;
378: pep->perm = NULL;
379: pep->pbc = NULL;
380: pep->solvematcoeffs = NULL;
381: pep->nwork = 0;
382: pep->work = NULL;
383: pep->data = NULL;
385: pep->state = PEP_STATE_INITIAL;
386: pep->nconv = 0;
387: pep->its = 0;
388: pep->n = 0;
389: pep->nloc = 0;
390: pep->nrma = NULL;
391: pep->sfactor_set = PETSC_FALSE;
392: pep->reason = PEP_CONVERGED_ITERATING;
394: PetscNewLog(pep,&pep->sc);
395: PetscRandomCreate(comm,&pep->rand);
396: PetscRandomSetSeed(pep->rand,0x12345678);
397: PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->rand);
398: *outpep = pep;
399: return(0);
400: }
404: /*@C
405: PEPSetType - Selects the particular solver to be used in the PEP object.
407: Logically Collective on PEP
409: Input Parameters:
410: + pep - the polynomial eigensolver context
411: - type - a known method
413: Options Database Key:
414: . -pep_type <method> - Sets the method; use -help for a list
415: of available methods
417: Notes:
418: See "slepc/include/slepcpep.h" for available methods. The default
419: is PEPTOAR.
421: Normally, it is best to use the PEPSetFromOptions() command and
422: then set the PEP type from the options database rather than by using
423: this routine. Using the options database provides the user with
424: maximum flexibility in evaluating the different available methods.
425: The PEPSetType() routine is provided for those situations where it
426: is necessary to set the iterative solver independently of the command
427: line or options database.
429: Level: intermediate
431: .seealso: PEPType
432: @*/
433: PetscErrorCode PEPSetType(PEP pep,PEPType type)
434: {
435: PetscErrorCode ierr,(*r)(PEP);
436: PetscBool match;
442: PetscObjectTypeCompare((PetscObject)pep,type,&match);
443: if (match) return(0);
445: PetscFunctionListFind(PEPList,type,&r);
446: if (!r) SETERRQ1(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown PEP type given: %s",type);
448: if (pep->ops->destroy) { (*pep->ops->destroy)(pep); }
449: PetscMemzero(pep->ops,sizeof(struct _PEPOps));
451: pep->state = PEP_STATE_INITIAL;
452: PetscObjectChangeTypeName((PetscObject)pep,type);
453: (*r)(pep);
454: return(0);
455: }
459: /*@C
460: PEPGetType - Gets the PEP type as a string from the PEP object.
462: Not Collective
464: Input Parameter:
465: . pep - the eigensolver context
467: Output Parameter:
468: . name - name of PEP method
470: Level: intermediate
472: .seealso: PEPSetType()
473: @*/
474: PetscErrorCode PEPGetType(PEP pep,PEPType *type)
475: {
479: *type = ((PetscObject)pep)->type_name;
480: return(0);
481: }
485: /*@C
486: PEPRegister - Adds a method to the polynomial eigenproblem solver package.
488: Not Collective
490: Input Parameters:
491: + name - name of a new user-defined solver
492: - function - routine to create the solver context
494: Notes:
495: PEPRegister() may be called multiple times to add several user-defined solvers.
497: Sample usage:
498: .vb
499: PEPRegister("my_solver",MySolverCreate);
500: .ve
502: Then, your solver can be chosen with the procedural interface via
503: $ PEPSetType(pep,"my_solver")
504: or at runtime via the option
505: $ -pep_type my_solver
507: Level: advanced
509: .seealso: PEPRegisterAll()
510: @*/
511: PetscErrorCode PEPRegister(const char *name,PetscErrorCode (*function)(PEP))
512: {
516: PetscFunctionListAdd(&PEPList,name,function);
517: return(0);
518: }
522: /*@
523: PEPReset - Resets the PEP context to the initial state and removes any
524: allocated objects.
526: Collective on PEP
528: Input Parameter:
529: . pep - eigensolver context obtained from PEPCreate()
531: Level: advanced
533: .seealso: PEPDestroy()
534: @*/
535: PetscErrorCode PEPReset(PEP pep)
536: {
538: PetscInt ncols;
542: if (pep->ops->reset) { (pep->ops->reset)(pep); }
543: if (pep->st) { STReset(pep->st); }
544: if (pep->ds) { DSReset(pep->ds); }
545: if (pep->nmat) {
546: MatDestroyMatrices(pep->nmat,&pep->A);
547: PetscFree3(pep->pbc,pep->solvematcoeffs,pep->nrma);
548: pep->nmat = 0;
549: }
550: VecDestroy(&pep->Dl);
551: VecDestroy(&pep->Dr);
552: BVGetSizes(pep->V,NULL,NULL,&ncols);
553: if (ncols) {
554: PetscFree4(pep->eigr,pep->eigi,pep->errest,pep->perm);
555: }
556: BVDestroy(&pep->V);
557: VecDestroyVecs(pep->nwork,&pep->work);
558: pep->nwork = 0;
559: pep->state = PEP_STATE_INITIAL;
560: return(0);
561: }
565: /*@C
566: PEPDestroy - Destroys the PEP context.
568: Collective on PEP
570: Input Parameter:
571: . pep - eigensolver context obtained from PEPCreate()
573: Level: beginner
575: .seealso: PEPCreate(), PEPSetUp(), PEPSolve()
576: @*/
577: PetscErrorCode PEPDestroy(PEP *pep)
578: {
582: if (!*pep) return(0);
584: if (--((PetscObject)(*pep))->refct > 0) { *pep = 0; return(0); }
585: PEPReset(*pep);
586: if ((*pep)->ops->destroy) { (*(*pep)->ops->destroy)(*pep); }
587: STDestroy(&(*pep)->st);
588: RGDestroy(&(*pep)->rg);
589: DSDestroy(&(*pep)->ds);
590: PetscRandomDestroy(&(*pep)->rand);
591: PetscFree((*pep)->sc);
592: /* just in case the initial vectors have not been used */
593: SlepcBasisDestroy_Private(&(*pep)->nini,&(*pep)->IS);
594: if ((*pep)->convergeddestroy) {
595: (*(*pep)->convergeddestroy)((*pep)->convergedctx);
596: }
597: PEPMonitorCancel(*pep);
598: PetscHeaderDestroy(pep);
599: return(0);
600: }
604: /*@
605: PEPSetBV - Associates a basis vectors object to the polynomial eigensolver.
607: Collective on PEP
609: Input Parameters:
610: + pep - eigensolver context obtained from PEPCreate()
611: - bv - the basis vectors object
613: Note:
614: Use PEPGetBV() to retrieve the basis vectors context (for example,
615: to free it at the end of the computations).
617: Level: advanced
619: .seealso: PEPGetBV()
620: @*/
621: PetscErrorCode PEPSetBV(PEP pep,BV bv)
622: {
629: PetscObjectReference((PetscObject)bv);
630: BVDestroy(&pep->V);
631: pep->V = bv;
632: PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->V);
633: return(0);
634: }
638: /*@C
639: PEPGetBV - Obtain the basis vectors object associated to the polynomial
640: eigensolver object.
642: Not Collective
644: Input Parameters:
645: . pep - eigensolver context obtained from PEPCreate()
647: Output Parameter:
648: . bv - basis vectors context
650: Level: advanced
652: .seealso: PEPSetBV()
653: @*/
654: PetscErrorCode PEPGetBV(PEP pep,BV *bv)
655: {
661: if (!pep->V) {
662: BVCreate(PetscObjectComm((PetscObject)pep),&pep->V);
663: PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->V);
664: }
665: *bv = pep->V;
666: return(0);
667: }
671: /*@
672: PEPSetRG - Associates a region object to the polynomial eigensolver.
674: Collective on PEP
676: Input Parameters:
677: + pep - eigensolver context obtained from PEPCreate()
678: - rg - the region object
680: Note:
681: Use PEPGetRG() to retrieve the region context (for example,
682: to free it at the end of the computations).
684: Level: advanced
686: .seealso: PEPGetRG()
687: @*/
688: PetscErrorCode PEPSetRG(PEP pep,RG rg)
689: {
696: PetscObjectReference((PetscObject)rg);
697: RGDestroy(&pep->rg);
698: pep->rg = rg;
699: PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->rg);
700: return(0);
701: }
705: /*@C
706: PEPGetRG - Obtain the region object associated to the
707: polynomial eigensolver object.
709: Not Collective
711: Input Parameters:
712: . pep - eigensolver context obtained from PEPCreate()
714: Output Parameter:
715: . rg - region context
717: Level: advanced
719: .seealso: PEPSetRG()
720: @*/
721: PetscErrorCode PEPGetRG(PEP pep,RG *rg)
722: {
728: if (!pep->rg) {
729: RGCreate(PetscObjectComm((PetscObject)pep),&pep->rg);
730: PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->rg);
731: }
732: *rg = pep->rg;
733: return(0);
734: }
738: /*@
739: PEPSetDS - Associates a direct solver object to the polynomial eigensolver.
741: Collective on PEP
743: Input Parameters:
744: + pep - eigensolver context obtained from PEPCreate()
745: - ds - the direct solver object
747: Note:
748: Use PEPGetDS() to retrieve the direct solver context (for example,
749: to free it at the end of the computations).
751: Level: advanced
753: .seealso: PEPGetDS()
754: @*/
755: PetscErrorCode PEPSetDS(PEP pep,DS ds)
756: {
763: PetscObjectReference((PetscObject)ds);
764: DSDestroy(&pep->ds);
765: pep->ds = ds;
766: PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->ds);
767: return(0);
768: }
772: /*@C
773: PEPGetDS - Obtain the direct solver object associated to the
774: polynomial eigensolver object.
776: Not Collective
778: Input Parameters:
779: . pep - eigensolver context obtained from PEPCreate()
781: Output Parameter:
782: . ds - direct solver context
784: Level: advanced
786: .seealso: PEPSetDS()
787: @*/
788: PetscErrorCode PEPGetDS(PEP pep,DS *ds)
789: {
795: if (!pep->ds) {
796: DSCreate(PetscObjectComm((PetscObject)pep),&pep->ds);
797: PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->ds);
798: }
799: *ds = pep->ds;
800: return(0);
801: }
805: /*@
806: PEPSetST - Associates a spectral transformation object to the eigensolver.
808: Collective on PEP
810: Input Parameters:
811: + pep - eigensolver context obtained from PEPCreate()
812: - st - the spectral transformation object
814: Note:
815: Use PEPGetST() to retrieve the spectral transformation context (for example,
816: to free it at the end of the computations).
818: Level: developer
820: .seealso: PEPGetST()
821: @*/
822: PetscErrorCode PEPSetST(PEP pep,ST st)
823: {
830: PetscObjectReference((PetscObject)st);
831: STDestroy(&pep->st);
832: pep->st = st;
833: PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->st);
834: return(0);
835: }
839: /*@C
840: PEPGetST - Obtain the spectral transformation (ST) object associated
841: to the eigensolver object.
843: Not Collective
845: Input Parameters:
846: . pep - eigensolver context obtained from PEPCreate()
848: Output Parameter:
849: . st - spectral transformation context
851: Level: beginner
853: .seealso: PEPSetST()
854: @*/
855: PetscErrorCode PEPGetST(PEP pep,ST *st)
856: {
862: if (!pep->st) {
863: STCreate(PetscObjectComm((PetscObject)pep),&pep->st);
864: PetscLogObjectParent((PetscObject)pep,(PetscObject)pep->st);
865: }
866: *st = pep->st;
867: return(0);
868: }
872: /*@
873: PEPSetTarget - Sets the value of the target.
875: Logically Collective on PEP
877: Input Parameters:
878: + pep - eigensolver context
879: - target - the value of the target
881: Options Database Key:
882: . -pep_target <scalar> - the value of the target
884: Notes:
885: The target is a scalar value used to determine the portion of the spectrum
886: of interest. It is used in combination with PEPSetWhichEigenpairs().
888: In the case of complex scalars, a complex value can be provided in the
889: command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
890: -pep_target 1.0+2.0i
892: Level: beginner
894: .seealso: PEPGetTarget(), PEPSetWhichEigenpairs()
895: @*/
896: PetscErrorCode PEPSetTarget(PEP pep,PetscScalar target)
897: {
903: pep->target = target;
904: if (!pep->st) { PEPGetST(pep,&pep->st); }
905: STSetDefaultShift(pep->st,target);
906: return(0);
907: }
911: /*@
912: PEPGetTarget - Gets the value of the target.
914: Not Collective
916: Input Parameter:
917: . pep - eigensolver context
919: Output Parameter:
920: . target - the value of the target
922: Note:
923: If the target was not set by the user, then zero is returned.
925: Level: beginner
927: .seealso: PEPSetTarget()
928: @*/
929: PetscErrorCode PEPGetTarget(PEP pep,PetscScalar* target)
930: {
934: *target = pep->target;
935: return(0);
936: }