Actual source code: qepbasic.c
1: /*
2: The basic QEP routines, Create, View, etc. are here.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2013, 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/qepimpl.h> /*I "slepcqep.h" I*/
26: PetscFunctionList QEPList = 0;
27: PetscBool QEPRegisterAllCalled = PETSC_FALSE;
28: PetscClassId QEP_CLASSID = 0;
29: PetscLogEvent QEP_SetUp = 0,QEP_Solve = 0,QEP_Dense = 0;
30: static PetscBool QEPPackageInitialized = PETSC_FALSE;
34: /*@C
35: QEPFinalizePackage - This function destroys everything in the Slepc interface
36: to the QEP package. It is called from SlepcFinalize().
38: Level: developer
40: .seealso: SlepcFinalize()
41: @*/
42: PetscErrorCode QEPFinalizePackage(void)
43: {
47: PetscFunctionListDestroy(&QEPList);
48: QEPPackageInitialized = PETSC_FALSE;
49: QEPRegisterAllCalled = PETSC_FALSE;
50: return(0);
51: }
55: /*@C
56: QEPInitializePackage - This function initializes everything in the QEP package. It is called
57: from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to QEPCreate()
58: when using static libraries.
60: Level: developer
62: .seealso: SlepcInitialize()
63: @*/
64: PetscErrorCode QEPInitializePackage(void)
65: {
66: char logList[256];
67: char *className;
68: PetscBool opt;
72: if (QEPPackageInitialized) return(0);
73: QEPPackageInitialized = PETSC_TRUE;
74: /* Register Classes */
75: PetscClassIdRegister("Quadratic Eigenvalue Problem solver",&QEP_CLASSID);
76: /* Register Constructors */
77: QEPRegisterAll();
78: /* Register Events */
79: PetscLogEventRegister("QEPSetUp",QEP_CLASSID,&QEP_SetUp);
80: PetscLogEventRegister("QEPSolve",QEP_CLASSID,&QEP_Solve);
81: PetscLogEventRegister("QEPDense",QEP_CLASSID,&QEP_Dense);
82: /* Process info exclusions */
83: PetscOptionsGetString(NULL,"-info_exclude",logList,256,&opt);
84: if (opt) {
85: PetscStrstr(logList,"qep",&className);
86: if (className) {
87: PetscInfoDeactivateClass(QEP_CLASSID);
88: }
89: }
90: /* Process summary exclusions */
91: PetscOptionsGetString(NULL,"-log_summary_exclude",logList,256,&opt);
92: if (opt) {
93: PetscStrstr(logList,"qep",&className);
94: if (className) {
95: PetscLogEventDeactivateClass(QEP_CLASSID);
96: }
97: }
98: PetscRegisterFinalize(QEPFinalizePackage);
99: return(0);
100: }
104: /*@C
105: QEPView - Prints the QEP data structure.
107: Collective on QEP
109: Input Parameters:
110: + qep - the quadratic eigenproblem solver context
111: - viewer - optional visualization context
113: Options Database Key:
114: . -qep_view - Calls QEPView() at end of QEPSolve()
116: Note:
117: The available visualization contexts include
118: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
119: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
120: output where only the first processor opens
121: the file. All other processors send their
122: data to the first processor to print.
124: The user can open an alternative visualization context with
125: PetscViewerASCIIOpen() - output to a specified file.
127: Level: beginner
129: .seealso: PetscViewerASCIIOpen()
130: @*/
131: PetscErrorCode QEPView(QEP qep,PetscViewer viewer)
132: {
134: const char *type;
135: char str[50];
136: PetscBool isascii,islinear;
140: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)qep));
144: #if defined(PETSC_USE_COMPLEX)
145: #define HERM "hermitian"
146: #else
147: #define HERM "symmetric"
148: #endif
149: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
150: if (isascii) {
151: PetscObjectPrintClassNamePrefixType((PetscObject)qep,viewer,"QEP Object");
152: if (qep->ops->view) {
153: PetscViewerASCIIPushTab(viewer);
154: (*qep->ops->view)(qep,viewer);
155: PetscViewerASCIIPopTab(viewer);
156: }
157: if (qep->problem_type) {
158: switch (qep->problem_type) {
159: case QEP_GENERAL: type = "general quadratic eigenvalue problem"; break;
160: case QEP_HERMITIAN: type = HERM " quadratic eigenvalue problem"; break;
161: case QEP_GYROSCOPIC: type = "gyroscopic quadratic eigenvalue problem"; break;
162: default: SETERRQ(PetscObjectComm((PetscObject)qep),1,"Wrong value of qep->problem_type");
163: }
164: } else type = "not yet set";
165: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
166: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
167: SlepcSNPrintfScalar(str,50,qep->target,PETSC_FALSE);
168: if (!qep->which) {
169: PetscViewerASCIIPrintf(viewer,"not yet set\n");
170: } else switch (qep->which) {
171: case QEP_TARGET_MAGNITUDE:
172: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
173: break;
174: case QEP_TARGET_REAL:
175: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
176: break;
177: case QEP_TARGET_IMAGINARY:
178: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
179: break;
180: case QEP_LARGEST_MAGNITUDE:
181: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
182: break;
183: case QEP_SMALLEST_MAGNITUDE:
184: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
185: break;
186: case QEP_LARGEST_REAL:
187: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
188: break;
189: case QEP_SMALLEST_REAL:
190: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
191: break;
192: case QEP_LARGEST_IMAGINARY:
193: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
194: break;
195: case QEP_SMALLEST_IMAGINARY:
196: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
197: break;
198: default: SETERRQ(PetscObjectComm((PetscObject)qep),1,"Wrong value of qep->which");
199: }
200: if (qep->leftvecs) {
201: PetscViewerASCIIPrintf(viewer," computing left eigenvectors also\n");
202: }
203: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",qep->nev);
204: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",qep->ncv);
205: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",qep->mpd);
206: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",qep->max_it);
207: PetscViewerASCIIPrintf(viewer," tolerance: %G\n",qep->tol);
208: PetscViewerASCIIPrintf(viewer," scaling factor: %G\n",qep->sfactor);
209: if (qep->nini) {
210: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(qep->nini));
211: }
212: if (qep->ninil) {
213: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial left space: %D\n",PetscAbs(qep->ninil));
214: }
215: } else {
216: if (qep->ops->view) {
217: (*qep->ops->view)(qep,viewer);
218: }
219: }
220: PetscObjectTypeCompare((PetscObject)qep,QEPLINEAR,&islinear);
221: if (!islinear) {
222: if (!qep->ip) { QEPGetIP(qep,&qep->ip); }
223: IPView(qep->ip,viewer);
224: if (!qep->ds) { QEPGetDS(qep,&qep->ds); }
225: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
226: DSView(qep->ds,viewer);
227: PetscViewerPopFormat(viewer);
228: if (!qep->st) { QEPGetST(qep,&qep->st); }
229: STView(qep->st,viewer);
230: }
231: return(0);
232: }
236: /*@
237: QEPPrintSolution - Prints the computed eigenvalues.
239: Collective on QEP
241: Input Parameters:
242: + qep - the eigensolver context
243: - viewer - optional visualization context
245: Options Database Key:
246: . -qep_terse - print only minimal information
248: Note:
249: By default, this function prints a table with eigenvalues and associated
250: relative errors. With -qep_terse only the eigenvalues are printed.
252: Level: intermediate
254: .seealso: PetscViewerASCIIOpen()
255: @*/
256: PetscErrorCode QEPPrintSolution(QEP qep,PetscViewer viewer)
257: {
258: PetscBool terse,errok,isascii;
259: PetscReal error,re,im;
260: PetscScalar kr,ki;
261: PetscInt i,j;
266: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)qep));
269: if (!qep->eigr || !qep->eigi || !qep->V) SETERRQ(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_WRONGSTATE,"QEPSolve must be called first");
270: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
271: if (!isascii) return(0);
273: PetscOptionsHasName(NULL,"-qep_terse",&terse);
274: if (terse) {
275: if (qep->nconv<qep->nev) {
276: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",qep->nev);
277: } else {
278: errok = PETSC_TRUE;
279: for (i=0;i<qep->nev;i++) {
280: QEPComputeRelativeError(qep,i,&error);
281: errok = (errok && error<qep->tol)? PETSC_TRUE: PETSC_FALSE;
282: }
283: if (errok) {
284: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
285: for (i=0;i<=(qep->nev-1)/8;i++) {
286: PetscViewerASCIIPrintf(viewer,"\n ");
287: for (j=0;j<PetscMin(8,qep->nev-8*i);j++) {
288: QEPGetEigenpair(qep,8*i+j,&kr,&ki,NULL,NULL);
289: #if defined(PETSC_USE_COMPLEX)
290: re = PetscRealPart(kr);
291: im = PetscImaginaryPart(kr);
292: #else
293: re = kr;
294: im = ki;
295: #endif
296: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
297: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
298: if (im!=0.0) {
299: PetscViewerASCIIPrintf(viewer,"%.5F%+.5Fi",re,im);
300: } else {
301: PetscViewerASCIIPrintf(viewer,"%.5F",re);
302: }
303: if (8*i+j+1<qep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
304: }
305: }
306: PetscViewerASCIIPrintf(viewer,"\n\n");
307: } else {
308: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",qep->nev);
309: }
310: }
311: } else {
312: PetscViewerASCIIPrintf(viewer," Number of converged approximate eigenpairs: %D\n\n",qep->nconv);
313: if (qep->nconv>0) {
314: PetscViewerASCIIPrintf(viewer,
315: " k ||(k^2M+Ck+K)x||/||kx||\n"
316: " ----------------- -------------------------\n");
317: for (i=0;i<qep->nconv;i++) {
318: QEPGetEigenpair(qep,i,&kr,&ki,NULL,NULL);
319: QEPComputeRelativeError(qep,i,&error);
320: #if defined(PETSC_USE_COMPLEX)
321: re = PetscRealPart(kr);
322: im = PetscImaginaryPart(kr);
323: #else
324: re = kr;
325: im = ki;
326: #endif
327: if (im!=0.0) {
328: PetscViewerASCIIPrintf(viewer," % 9F%+9F i %12G\n",re,im,error);
329: } else {
330: PetscViewerASCIIPrintf(viewer," % 12F %12G\n",re,error);
331: }
332: }
333: PetscViewerASCIIPrintf(viewer,"\n");
334: }
335: }
336: return(0);
337: }
341: /*@C
342: QEPCreate - Creates the default QEP context.
344: Collective on MPI_Comm
346: Input Parameter:
347: . comm - MPI communicator
349: Output Parameter:
350: . qep - location to put the QEP context
352: Note:
353: The default QEP type is QEPLINEAR
355: Level: beginner
357: .seealso: QEPSetUp(), QEPSolve(), QEPDestroy(), QEP
358: @*/
359: PetscErrorCode QEPCreate(MPI_Comm comm,QEP *outqep)
360: {
362: QEP qep;
366: *outqep = 0;
367: #if !defined(PETSC_USE_DYNAMIC_LIBRARIES)
368: QEPInitializePackage();
369: #endif
371: SlepcHeaderCreate(qep,_p_QEP,struct _QEPOps,QEP_CLASSID,"QEP","Quadratic Eigenvalue Problem","QEP",comm,QEPDestroy,QEPView);
373: qep->M = 0;
374: qep->C = 0;
375: qep->K = 0;
376: qep->max_it = 0;
377: qep->nev = 1;
378: qep->ncv = 0;
379: qep->mpd = 0;
380: qep->nini = 0;
381: qep->ninil = 0;
382: qep->allocated_ncv = 0;
383: qep->ip = 0;
384: qep->ds = 0;
385: qep->tol = PETSC_DEFAULT;
386: qep->sfactor = 0.0;
387: qep->sfactor_set = PETSC_FALSE;
388: qep->converged = QEPConvergedDefault;
389: qep->convergedctx = NULL;
390: qep->which = (QEPWhich)0;
391: qep->comparison = NULL;
392: qep->comparisonctx = NULL;
393: qep->leftvecs = PETSC_FALSE;
394: qep->problem_type = (QEPProblemType)0;
395: qep->V = NULL;
396: qep->W = NULL;
397: qep->IS = NULL;
398: qep->ISL = NULL;
399: qep->eigr = NULL;
400: qep->eigi = NULL;
401: qep->errest = NULL;
402: qep->data = NULL;
403: qep->t = NULL;
404: qep->nconv = 0;
405: qep->its = 0;
406: qep->perm = NULL;
407: qep->matvecs = 0;
408: qep->linits = 0;
409: qep->nwork = 0;
410: qep->work = NULL;
411: qep->setupcalled = 0;
412: qep->reason = QEP_CONVERGED_ITERATING;
413: qep->numbermonitors = 0;
414: qep->trackall = PETSC_FALSE;
415: qep->rand = 0;
417: PetscRandomCreate(comm,&qep->rand);
418: PetscRandomSetSeed(qep->rand,0x12345678);
419: PetscLogObjectParent(qep,qep->rand);
420: *outqep = qep;
421: return(0);
422: }
426: /*@C
427: QEPSetType - Selects the particular solver to be used in the QEP object.
429: Logically Collective on QEP
431: Input Parameters:
432: + qep - the quadratic eigensolver context
433: - type - a known method
435: Options Database Key:
436: . -qep_type <method> - Sets the method; use -help for a list
437: of available methods
439: Notes:
440: See "slepc/include/slepcqep.h" for available methods. The default
441: is QEPLINEAR.
443: Normally, it is best to use the QEPSetFromOptions() command and
444: then set the QEP type from the options database rather than by using
445: this routine. Using the options database provides the user with
446: maximum flexibility in evaluating the different available methods.
447: The QEPSetType() routine is provided for those situations where it
448: is necessary to set the iterative solver independently of the command
449: line or options database.
451: Level: intermediate
453: .seealso: QEPType
454: @*/
455: PetscErrorCode QEPSetType(QEP qep,QEPType type)
456: {
457: PetscErrorCode ierr,(*r)(QEP);
458: PetscBool match;
464: PetscObjectTypeCompare((PetscObject)qep,type,&match);
465: if (match) return(0);
467: PetscFunctionListFind(QEPList,type,&r);
468: if (!r) SETERRQ1(PetscObjectComm((PetscObject)qep),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown QEP type given: %s",type);
470: if (qep->ops->destroy) { (*qep->ops->destroy)(qep); }
471: PetscMemzero(qep->ops,sizeof(struct _QEPOps));
473: qep->setupcalled = 0;
474: PetscObjectChangeTypeName((PetscObject)qep,type);
475: (*r)(qep);
476: return(0);
477: }
481: /*@C
482: QEPGetType - Gets the QEP type as a string from the QEP object.
484: Not Collective
486: Input Parameter:
487: . qep - the eigensolver context
489: Output Parameter:
490: . name - name of QEP method
492: Level: intermediate
494: .seealso: QEPSetType()
495: @*/
496: PetscErrorCode QEPGetType(QEP qep,QEPType *type)
497: {
501: *type = ((PetscObject)qep)->type_name;
502: return(0);
503: }
507: /*@C
508: QEPRegister - Adds a method to the quadratic eigenproblem solver package.
510: Not Collective
512: Input Parameters:
513: + name - name of a new user-defined solver
514: - function - routine to create the solver context
516: Notes:
517: QEPRegister() may be called multiple times to add several user-defined solvers.
519: Sample usage:
520: .vb
521: QEPRegister("my_solver",MySolverCreate);
522: .ve
524: Then, your solver can be chosen with the procedural interface via
525: $ QEPSetType(qep,"my_solver")
526: or at runtime via the option
527: $ -qep_type my_solver
529: Level: advanced
531: .seealso: QEPRegisterAll()
532: @*/
533: PetscErrorCode QEPRegister(const char *name,PetscErrorCode (*function)(QEP))
534: {
538: PetscFunctionListAdd(&QEPList,name,function);
539: return(0);
540: }
544: /*@
545: QEPReset - Resets the QEP context to the setupcalled=0 state and removes any
546: allocated objects.
548: Collective on QEP
550: Input Parameter:
551: . qep - eigensolver context obtained from QEPCreate()
553: Level: advanced
555: .seealso: QEPDestroy()
556: @*/
557: PetscErrorCode QEPReset(QEP qep)
558: {
563: if (qep->ops->reset) { (qep->ops->reset)(qep); }
564: if (qep->ip) { IPReset(qep->ip); }
565: if (qep->ds) { DSReset(qep->ds); }
566: MatDestroy(&qep->M);
567: MatDestroy(&qep->C);
568: MatDestroy(&qep->K);
569: VecDestroy(&qep->t);
570: QEPFreeSolution(qep);
571: qep->matvecs = 0;
572: qep->linits = 0;
573: qep->setupcalled = 0;
574: return(0);
575: }
579: /*@C
580: QEPDestroy - Destroys the QEP context.
582: Collective on QEP
584: Input Parameter:
585: . qep - eigensolver context obtained from QEPCreate()
587: Level: beginner
589: .seealso: QEPCreate(), QEPSetUp(), QEPSolve()
590: @*/
591: PetscErrorCode QEPDestroy(QEP *qep)
592: {
596: if (!*qep) return(0);
598: if (--((PetscObject)(*qep))->refct > 0) { *qep = 0; return(0); }
599: QEPReset(*qep);
600: if ((*qep)->ops->destroy) { (*(*qep)->ops->destroy)(*qep); }
601: STDestroy(&(*qep)->st);
602: IPDestroy(&(*qep)->ip);
603: DSDestroy(&(*qep)->ds);
604: PetscRandomDestroy(&(*qep)->rand);
605: /* just in case the initial vectors have not been used */
606: SlepcBasisDestroy_Private(&(*qep)->nini,&(*qep)->IS);
607: SlepcBasisDestroy_Private(&(*qep)->ninil,&(*qep)->ISL);
608: QEPMonitorCancel(*qep);
609: PetscHeaderDestroy(qep);
610: return(0);
611: }
615: /*@
616: QEPSetIP - Associates an inner product object to the quadratic eigensolver.
618: Collective on QEP
620: Input Parameters:
621: + qep - eigensolver context obtained from QEPCreate()
622: - ip - the inner product object
624: Note:
625: Use QEPGetIP() to retrieve the inner product context (for example,
626: to free it at the end of the computations).
628: Level: advanced
630: .seealso: QEPGetIP()
631: @*/
632: PetscErrorCode QEPSetIP(QEP qep,IP ip)
633: {
640: PetscObjectReference((PetscObject)ip);
641: IPDestroy(&qep->ip);
642: qep->ip = ip;
643: PetscLogObjectParent(qep,qep->ip);
644: return(0);
645: }
649: /*@C
650: QEPGetIP - Obtain the inner product object associated
651: to the quadratic eigensolver object.
653: Not Collective
655: Input Parameters:
656: . qep - eigensolver context obtained from QEPCreate()
658: Output Parameter:
659: . ip - inner product context
661: Level: advanced
663: .seealso: QEPSetIP()
664: @*/
665: PetscErrorCode QEPGetIP(QEP qep,IP *ip)
666: {
672: if (!qep->ip) {
673: IPCreate(PetscObjectComm((PetscObject)qep),&qep->ip);
674: PetscLogObjectParent(qep,qep->ip);
675: }
676: *ip = qep->ip;
677: return(0);
678: }
682: /*@
683: QEPSetDS - Associates a direct solver object to the quadratic eigensolver.
685: Collective on QEP
687: Input Parameters:
688: + qep - eigensolver context obtained from QEPCreate()
689: - ds - the direct solver object
691: Note:
692: Use QEPGetDS() to retrieve the direct solver context (for example,
693: to free it at the end of the computations).
695: Level: advanced
697: .seealso: QEPGetDS()
698: @*/
699: PetscErrorCode QEPSetDS(QEP qep,DS ds)
700: {
707: PetscObjectReference((PetscObject)ds);
708: DSDestroy(&qep->ds);
709: qep->ds = ds;
710: PetscLogObjectParent(qep,qep->ds);
711: return(0);
712: }
716: /*@C
717: QEPGetDS - Obtain the direct solver object associated to the
718: quadratic eigensolver object.
720: Not Collective
722: Input Parameters:
723: . qep - eigensolver context obtained from QEPCreate()
725: Output Parameter:
726: . ds - direct solver context
728: Level: advanced
730: .seealso: QEPSetDS()
731: @*/
732: PetscErrorCode QEPGetDS(QEP qep,DS *ds)
733: {
739: if (!qep->ds) {
740: DSCreate(PetscObjectComm((PetscObject)qep),&qep->ds);
741: PetscLogObjectParent(qep,qep->ds);
742: }
743: *ds = qep->ds;
744: return(0);
745: }
749: /*@
750: QEPSetST - Associates a spectral transformation object to the eigensolver.
752: Collective on QEP
754: Input Parameters:
755: + qep - eigensolver context obtained from QEPCreate()
756: - st - the spectral transformation object
758: Note:
759: Use QEPGetST() to retrieve the spectral transformation context (for example,
760: to free it at the end of the computations).
762: Level: developer
764: .seealso: QEPGetST()
765: @*/
766: PetscErrorCode QEPSetST(QEP qep,ST st)
767: {
774: PetscObjectReference((PetscObject)st);
775: STDestroy(&qep->st);
776: qep->st = st;
777: PetscLogObjectParent(qep,qep->st);
778: return(0);
779: }
783: /*@C
784: QEPGetST - Obtain the spectral transformation (ST) object associated
785: to the eigensolver object.
787: Not Collective
789: Input Parameters:
790: . qep - eigensolver context obtained from QEPCreate()
792: Output Parameter:
793: . st - spectral transformation context
795: Level: beginner
797: .seealso: QEPSetST()
798: @*/
799: PetscErrorCode QEPGetST(QEP qep,ST *st)
800: {
806: if (!qep->st) {
807: STCreate(PetscObjectComm((PetscObject)qep),&qep->st);
808: PetscLogObjectParent(qep,qep->st);
809: }
810: *st = qep->st;
811: return(0);
812: }
816: /*@
817: QEPSetTarget - Sets the value of the target.
819: Logically Collective on QEP
821: Input Parameters:
822: + qep - eigensolver context
823: - target - the value of the target
825: Notes:
826: The target is a scalar value used to determine the portion of the spectrum
827: of interest. It is used in combination with QEPSetWhichEigenpairs().
829: Level: beginner
831: .seealso: QEPGetTarget(), QEPSetWhichEigenpairs()
832: @*/
833: PetscErrorCode QEPSetTarget(QEP qep,PetscScalar target)
834: {
840: qep->target = target;
841: if (!qep->st) { QEPGetST(qep,&qep->st); }
842: STSetDefaultShift(qep->st,target);
843: return(0);
844: }
848: /*@
849: QEPGetTarget - Gets the value of the target.
851: Not Collective
853: Input Parameter:
854: . qep - eigensolver context
856: Output Parameter:
857: . target - the value of the target
859: Level: beginner
861: Note:
862: If the target was not set by the user, then zero is returned.
864: .seealso: QEPSetTarget()
865: @*/
866: PetscErrorCode QEPGetTarget(QEP qep,PetscScalar* target)
867: {
871: *target = qep->target;
872: return(0);
873: }