Actual source code: basic.c
1: /*
2: The basic EPS routines, Create, View, etc. are here.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2010, Universidad Politecnica de Valencia, Spain
8: This file is part of SLEPc.
9:
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 private/epsimpl.h
26: PetscFList EPSList = 0;
27: PetscCookie EPS_COOKIE = 0;
28: PetscLogEvent EPS_SetUp = 0, EPS_Solve = 0, EPS_Dense = 0;
29: static PetscTruth EPSPackageInitialized = PETSC_FALSE;
33: /*@C
34: EPSFinalizePackage - This function destroys everything in the Slepc interface to the EPS package. It is
35: called from SlepcFinalize().
37: Level: developer
39: .seealso: SlepcFinalize()
40: @*/
41: PetscErrorCode EPSFinalizePackage(void)
42: {
44: EPSPackageInitialized = PETSC_FALSE;
45: EPSList = 0;
46: return(0);
47: }
51: /*@C
52: EPSInitializePackage - This function initializes everything in the EPS package. It is called
53: from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to EPSCreate()
54: when using static libraries.
56: Input Parameter:
57: path - The dynamic library path, or PETSC_NULL
59: Level: developer
61: .seealso: SlepcInitialize()
62: @*/
63: PetscErrorCode EPSInitializePackage(char *path) {
64: char logList[256];
65: char *className;
66: PetscTruth opt;
70: if (EPSPackageInitialized) return(0);
71: EPSPackageInitialized = PETSC_TRUE;
72: /* Register Classes */
73: PetscCookieRegister("Eigenproblem Solver",&EPS_COOKIE);
74: /* Register Constructors */
75: EPSRegisterAll(path);
76: /* Register Events */
77: PetscLogEventRegister("EPSSetUp",EPS_COOKIE,&EPS_SetUp);
78: PetscLogEventRegister("EPSSolve",EPS_COOKIE,&EPS_Solve);
79: PetscLogEventRegister("EPSDense",EPS_COOKIE,&EPS_Dense);
80: /* Process info exclusions */
81: PetscOptionsGetString(PETSC_NULL, "-log_info_exclude", logList, 256, &opt);
82: if (opt) {
83: PetscStrstr(logList, "eps", &className);
84: if (className) {
85: PetscInfoDeactivateClass(EPS_COOKIE);
86: }
87: }
88: /* Process summary exclusions */
89: PetscOptionsGetString(PETSC_NULL, "-log_summary_exclude", logList, 256, &opt);
90: if (opt) {
91: PetscStrstr(logList, "eps", &className);
92: if (className) {
93: PetscLogEventDeactivateClass(EPS_COOKIE);
94: }
95: }
96: PetscRegisterFinalize(EPSFinalizePackage);
97: return(0);
98: }
102: /*@C
103: EPSView - Prints the EPS data structure.
105: Collective on EPS
107: Input Parameters:
108: + eps - the eigenproblem solver context
109: - viewer - optional visualization context
111: Options Database Key:
112: . -eps_view - Calls EPSView() at end of EPSSolve()
114: Note:
115: The available visualization contexts include
116: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
117: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
118: output where only the first processor opens
119: the file. All other processors send their
120: data to the first processor to print.
122: The user can open an alternative visualization context with
123: PetscViewerASCIIOpen() - output to a specified file.
125: Level: beginner
127: .seealso: STView(), PetscViewerASCIIOpen()
128: @*/
129: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
130: {
132: const char *type,*extr,*bal;
133: PetscTruth isascii;
137: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(((PetscObject)eps)->comm);
141: #if defined(PETSC_USE_COMPLEX)
142: #define HERM "hermitian"
143: #else
144: #define HERM "symmetric"
145: #endif
146: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
147: if (isascii) {
148: PetscViewerASCIIPrintf(viewer,"EPS Object:\n");
149: switch (eps->problem_type) {
150: case EPS_HEP: type = HERM " eigenvalue problem"; break;
151: case EPS_GHEP: type = "generalized " HERM " eigenvalue problem"; break;
152: case EPS_NHEP: type = "non-" HERM " eigenvalue problem"; break;
153: case EPS_GNHEP: type = "generalized non-" HERM " eigenvalue problem"; break;
154: case EPS_PGNHEP: type = "generalized non-" HERM " eigenvalue problem with " HERM " positive definite B"; break;
155: case 0: type = "not yet set"; break;
156: default: SETERRQ(1,"Wrong value of eps->problem_type");
157: }
158: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
159: EPSGetType(eps,&type);
160: if (type) {
161: PetscViewerASCIIPrintf(viewer," method: %s\n",type);
162: } else {
163: PetscViewerASCIIPrintf(viewer," method: not yet set\n");
164: }
165: if (eps->ops->view) {
166: PetscViewerASCIIPushTab(viewer);
167: (*eps->ops->view)(eps,viewer);
168: PetscViewerASCIIPopTab(viewer);
169: }
170: if (eps->extraction) {
171: switch (eps->extraction) {
172: case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
173: case EPS_HARMONIC: extr = "harmonic Ritz"; break;
174: case EPS_HARMONIC_RELATIVE:extr = "relative harmonic Ritz"; break;
175: case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
176: case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
177: case EPS_REFINED: extr = "refined Ritz"; break;
178: case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
179: default: SETERRQ(1,"Wrong value of eps->extraction");
180: }
181: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr);
182: }
183: if (eps->balance && !eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
184: switch (eps->balance) {
185: case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
186: case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
187: case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
188: default: SETERRQ(1,"Wrong value of eps->balance");
189: }
190: PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal);
191: if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) {
192: PetscViewerASCIIPrintf(viewer,", with its=%d",eps->balance_its);
193: }
194: if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) {
195: PetscViewerASCIIPrintf(viewer," and cutoff=%g",eps->balance_cutoff);
196: }
197: PetscViewerASCIIPrintf(viewer,"\n");
198: }
199: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
200: if (!eps->which) {
201: PetscViewerASCIIPrintf(viewer,"not yet set\n");
202: } else switch (eps->which) {
203: case EPS_WHICH_USER:
204: PetscViewerASCIIPrintf(viewer,"user defined\n");
205: break;
206: case EPS_TARGET_MAGNITUDE:
207: #if !defined(PETSC_USE_COMPLEX)
208: PetscViewerASCIIPrintf(viewer,"closest to target: %g (in magnitude)\n",eps->target);
209: #else
210: PetscViewerASCIIPrintf(viewer,"closest to target: %g+%g i (in magnitude)\n",PetscRealPart(eps->target),PetscImaginaryPart(eps->target));
211: #endif
212: break;
213: case EPS_TARGET_REAL:
214: #if !defined(PETSC_USE_COMPLEX)
215: PetscViewerASCIIPrintf(viewer,"closest to target: %g (along the real axis)\n",eps->target);
216: #else
217: PetscViewerASCIIPrintf(viewer,"closest to target: %g+%g i (along the real axis)\n",PetscRealPart(eps->target),PetscImaginaryPart(eps->target));
218: #endif
219: break;
220: #if defined(PETSC_USE_COMPLEX)
221: case EPS_TARGET_IMAGINARY:
222: PetscViewerASCIIPrintf(viewer,"closest to target: %g+%g i (along the imaginary axis)\n",PetscRealPart(eps->target),PetscImaginaryPart(eps->target));
223: break;
224: #endif
225: case EPS_LARGEST_MAGNITUDE:
226: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
227: break;
228: case EPS_SMALLEST_MAGNITUDE:
229: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
230: break;
231: case EPS_LARGEST_REAL:
232: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
233: break;
234: case EPS_SMALLEST_REAL:
235: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
236: break;
237: case EPS_LARGEST_IMAGINARY:
238: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
239: break;
240: case EPS_SMALLEST_IMAGINARY:
241: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
242: break;
243: default: SETERRQ(1,"Wrong value of eps->which");
244: }
245: if (eps->leftvecs) {
246: PetscViewerASCIIPrintf(viewer," computing left eigenvectors also\n");
247: }
248: if (eps->trueres) {
249: PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n");
250: }
251: if (eps->trackall) {
252: PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n");
253: }
254: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %d\n",eps->nev);
255: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %d\n",eps->ncv);
256: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %d\n",eps->mpd);
257: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %d\n", eps->max_it);
258: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",eps->tol);
259: PetscViewerASCIIPrintf(viewer," convergence test: ");
260: switch(eps->conv) {
261: case EPS_CONV_ABS:
262: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
263: case EPS_CONV_EIG:
264: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
265: case EPS_CONV_NORM:
266: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");break;
267: default:
268: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
269: }
270: if (eps->nini!=0) {
271: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %d\n",PetscAbs(eps->nini));
272: }
273: if (eps->ninil!=0) {
274: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial left space: %d\n",PetscAbs(eps->ninil));
275: }
276: if (eps->nds>0) {
277: PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %d\n",eps->nds);
278: }
279: PetscViewerASCIIPrintf(viewer," estimates of matrix norms (%s): norm(A)=%g",eps->adaptive?"adaptive":"constant",eps->nrma);
280: if (eps->isgeneralized) {
281: PetscViewerASCIIPrintf(viewer,", norm(B)=%g",eps->nrmb);
282: }
283: PetscViewerASCIIPrintf(viewer,"\n");
284: PetscViewerASCIIPushTab(viewer);
285: IPView(eps->ip,viewer);
286: STView(eps->OP,viewer);
287: PetscViewerASCIIPopTab(viewer);
288: } else {
289: if (eps->ops->view) {
290: (*eps->ops->view)(eps,viewer);
291: }
292: STView(eps->OP,viewer);
293: }
294: return(0);
295: }
299: /*@C
300: EPSCreate - Creates the default EPS context.
302: Collective on MPI_Comm
304: Input Parameter:
305: . comm - MPI communicator
307: Output Parameter:
308: . eps - location to put the EPS context
310: Note:
311: The default EPS type is EPSKRYLOVSCHUR
313: Level: beginner
315: .seealso: EPSSetUp(), EPSSolve(), EPSDestroy(), EPS
316: @*/
317: PetscErrorCode EPSCreate(MPI_Comm comm,EPS *outeps)
318: {
320: EPS eps;
324: *outeps = 0;
326: PetscHeaderCreate(eps,_p_EPS,struct _EPSOps,EPS_COOKIE,-1,"EPS",comm,EPSDestroy,EPSView);
327: *outeps = eps;
329: PetscMemzero(eps->ops,sizeof(struct _EPSOps));
331: eps->max_it = 0;
332: eps->nev = 1;
333: eps->ncv = 0;
334: eps->mpd = 0;
335: eps->allocated_ncv = 0;
336: eps->nini = 0;
337: eps->ninil = 0;
338: eps->nds = 0;
339: eps->tol = 1e-7;
340: eps->conv = EPS_CONV_EIG;
341: eps->conv_func = EPSEigRelativeConverged;
342: eps->conv_ctx = PETSC_NULL;
343: eps->which = (EPSWhich)0;
344: eps->which_func = PETSC_NULL;
345: eps->which_ctx = PETSC_NULL;
346: eps->leftvecs = PETSC_FALSE;
347: eps->trueres = PETSC_FALSE;
348: eps->trackall = PETSC_FALSE;
349: eps->target = 0.0;
350: eps->evecsavailable = PETSC_FALSE;
351: eps->problem_type = (EPSProblemType)0;
352: eps->extraction = (EPSExtraction)0;
353: eps->balance = (EPSBalance)0;
354: eps->balance_its = 5;
355: eps->balance_cutoff = 1e-8;
356: eps->nrma = 1.0;
357: eps->nrmb = 1.0;
358: eps->adaptive = PETSC_FALSE;
360: eps->V = 0;
361: eps->W = 0;
362: eps->T = 0;
363: eps->DS = 0;
364: eps->IS = 0;
365: eps->ISL = 0;
366: eps->ds_ortho = PETSC_FALSE;
367: eps->eigr = 0;
368: eps->eigi = 0;
369: eps->errest = 0;
370: eps->errest_left = 0;
371: eps->OP = 0;
372: eps->ip = 0;
373: eps->data = 0;
374: eps->nconv = 0;
375: eps->its = 0;
376: eps->perm = PETSC_NULL;
378: eps->nwork = 0;
379: eps->work = 0;
380: eps->isgeneralized = PETSC_FALSE;
381: eps->ishermitian = PETSC_FALSE;
382: eps->ispositive = PETSC_FALSE;
383: eps->setupcalled = 0;
384: eps->reason = EPS_CONVERGED_ITERATING;
386: eps->numbermonitors = 0;
388: STCreate(comm,&eps->OP);
389: PetscLogObjectParent(eps,eps->OP);
390: IPCreate(comm,&eps->ip);
391: IPSetOptionsPrefix(eps->ip,((PetscObject)eps)->prefix);
392: IPAppendOptionsPrefix(eps->ip,"eps_");
393: PetscLogObjectParent(eps,eps->ip);
394: PetscPublishAll(eps);
395: return(0);
396: }
397:
400: /*@C
401: EPSSetType - Selects the particular solver to be used in the EPS object.
403: Collective on EPS
405: Input Parameters:
406: + eps - the eigensolver context
407: - type - a known method
409: Options Database Key:
410: . -eps_type <method> - Sets the method; use -help for a list
411: of available methods
412:
413: Notes:
414: See "slepc/include/slepceps.h" for available methods. The default
415: is EPSKRYLOVSCHUR.
417: Normally, it is best to use the EPSSetFromOptions() command and
418: then set the EPS type from the options database rather than by using
419: this routine. Using the options database provides the user with
420: maximum flexibility in evaluating the different available methods.
421: The EPSSetType() routine is provided for those situations where it
422: is necessary to set the iterative solver independently of the command
423: line or options database.
425: Level: intermediate
427: .seealso: STSetType(), EPSType
428: @*/
429: PetscErrorCode EPSSetType(EPS eps,const EPSType type)
430: {
431: PetscErrorCode ierr,(*r)(EPS);
432: PetscTruth match;
438: PetscTypeCompare((PetscObject)eps,type,&match);
439: if (match) return(0);
441: if (eps->data) {
442: /* destroy the old private EPS context */
443: (*eps->ops->destroy)(eps);
444: eps->data = 0;
445: }
447: PetscFListFind(EPSList,((PetscObject)eps)->comm,type,(void (**)(void)) &r);
449: if (!r) SETERRQ1(1,"Unknown EPS type given: %s",type);
451: eps->setupcalled = 0;
452: PetscMemzero(eps->ops,sizeof(struct _EPSOps));
453: (*r)(eps);
455: PetscObjectChangeTypeName((PetscObject)eps,type);
456: return(0);
457: }
461: /*@C
462: EPSGetType - Gets the EPS type as a string from the EPS object.
464: Not Collective
466: Input Parameter:
467: . eps - the eigensolver context
469: Output Parameter:
470: . name - name of EPS method
472: Level: intermediate
474: .seealso: EPSSetType()
475: @*/
476: PetscErrorCode EPSGetType(EPS eps,const EPSType *type)
477: {
481: *type = ((PetscObject)eps)->type_name;
482: return(0);
483: }
485: /*MC
486: EPSRegisterDynamic - Adds a method to the eigenproblem solver package.
488: Synopsis:
489: EPSRegisterDynamic(char *name_solver,char *path,char *name_create,PetscErrorCode (*routine_create)(EPS))
491: Not Collective
493: Input Parameters:
494: + name_solver - name of a new user-defined solver
495: . path - path (either absolute or relative) the library containing this solver
496: . name_create - name of routine to create the solver context
497: - routine_create - routine to create the solver context
499: Notes:
500: EPSRegisterDynamic() may be called multiple times to add several user-defined solvers.
502: If dynamic libraries are used, then the fourth input argument (routine_create)
503: is ignored.
505: Sample usage:
506: .vb
507: EPSRegisterDynamic("my_solver",/home/username/my_lib/lib/libO/solaris/mylib.a,
508: "MySolverCreate",MySolverCreate);
509: .ve
511: Then, your solver can be chosen with the procedural interface via
512: $ EPSSetType(eps,"my_solver")
513: or at runtime via the option
514: $ -eps_type my_solver
516: Level: advanced
518: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
519: and others of the form ${any_environmental_variable} occuring in pathname will be
520: replaced with appropriate values.
522: .seealso: EPSRegisterDestroy(), EPSRegisterAll()
524: M*/
528: /*@C
529: EPSRegister - See EPSRegisterDynamic()
531: Level: advanced
532: @*/
533: PetscErrorCode EPSRegister(const char *sname,const char *path,const char *name,PetscErrorCode (*function)(EPS))
534: {
536: char fullname[256];
539: PetscFListConcat(path,name,fullname);
540: PetscFListAdd(&EPSList,sname,fullname,(void (*)(void))function);
541: return(0);
542: }
546: /*@
547: EPSRegisterDestroy - Frees the list of EPS methods that were
548: registered by EPSRegisterDynamic().
550: Not Collective
552: Level: advanced
554: .seealso: EPSRegisterDynamic(), EPSRegisterAll()
555: @*/
556: PetscErrorCode EPSRegisterDestroy(void)
557: {
561: PetscFListDestroy(&EPSList);
562: EPSRegisterAll(PETSC_NULL);
563: return(0);
564: }
568: /*@
569: EPSDestroy - Destroys the EPS context.
571: Collective on EPS
573: Input Parameter:
574: . eps - eigensolver context obtained from EPSCreate()
576: Level: beginner
578: .seealso: EPSCreate(), EPSSetUp(), EPSSolve()
579: @*/
580: PetscErrorCode EPSDestroy(EPS eps)
581: {
586: if (--((PetscObject)eps)->refct > 0) return(0);
588: /* if memory was published with AMS then destroy it */
589: PetscObjectDepublish(eps);
591: STDestroy(eps->OP);
592: IPDestroy(eps->ip);
594: if (eps->ops->destroy) {
595: (*eps->ops->destroy)(eps);
596: }
597:
598: PetscFree(eps->T);
599: PetscFree(eps->Tl);
600: PetscFree(eps->perm);
601: if (eps->rand) {
602: PetscRandomDestroy(eps->rand);
603: }
605: if (eps->D) {
606: VecDestroy(eps->D);
607: }
609: EPSRemoveDeflationSpace(eps);
610: EPSMonitorCancel(eps);
612: PetscHeaderDestroy(eps);
613: return(0);
614: }
618: /*@
619: EPSSetTarget - Sets the value of the target.
621: Collective on EPS
623: Input Parameters:
624: + eps - eigensolver context
625: - target - the value of the target
627: Notes:
628: The target is a scalar value used to determine the portion of the spectrum
629: of interest. It is used in combination with EPSSetWhichEigenpairs().
630:
631: Level: beginner
633: .seealso: EPSGetTarget(), EPSSetWhichEigenpairs()
634: @*/
635: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
636: {
641: eps->target = target;
642: STSetDefaultShift(eps->OP,target);
643: return(0);
644: }
648: /*@
649: EPSGetTarget - Gets the value of the target.
651: Not collective
653: Input Parameter:
654: . eps - eigensolver context
656: Output Parameter:
657: . target - the value of the target
659: Level: beginner
661: Note:
662: If the target was not set by the user, then zero is returned.
664: .seealso: EPSSetTarget()
665: @*/
666: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
667: {
670: if (target) {
671: *target = eps->target;
672: }
673: return(0);
674: }
678: /*@
679: EPSSetST - Associates a spectral transformation object to the
680: eigensolver.
682: Collective on EPS and ST
684: Input Parameters:
685: + eps - eigensolver context obtained from EPSCreate()
686: - st - the spectral transformation object
688: Note:
689: Use EPSGetST() to retrieve the spectral transformation context (for example,
690: to free it at the end of the computations).
692: Level: advanced
694: .seealso: EPSGetST()
695: @*/
696: PetscErrorCode EPSSetST(EPS eps,ST st)
697: {
704: PetscObjectReference((PetscObject)st);
705: STDestroy(eps->OP);
706: eps->OP = st;
707: return(0);
708: }
712: /*@C
713: EPSGetST - Obtain the spectral transformation (ST) object associated
714: to the eigensolver object.
716: Not Collective
718: Input Parameters:
719: . eps - eigensolver context obtained from EPSCreate()
721: Output Parameter:
722: . st - spectral transformation context
724: Level: beginner
726: .seealso: EPSSetST()
727: @*/
728: PetscErrorCode EPSGetST(EPS eps, ST *st)
729: {
733: *st = eps->OP;
734: return(0);
735: }
739: /*@
740: EPSSetIP - Associates an inner product object to the eigensolver.
742: Collective on EPS and IP
744: Input Parameters:
745: + eps - eigensolver context obtained from EPSCreate()
746: - ip - the inner product object
748: Note:
749: Use EPSGetIP() to retrieve the inner product context (for example,
750: to free it at the end of the computations).
752: Level: advanced
754: .seealso: EPSGetIP()
755: @*/
756: PetscErrorCode EPSSetIP(EPS eps,IP ip)
757: {
764: PetscObjectReference((PetscObject)ip);
765: IPDestroy(eps->ip);
766: eps->ip = ip;
767: return(0);
768: }
772: /*@C
773: EPSGetIP - Obtain the inner product object associated
774: to the eigensolver object.
776: Not Collective
778: Input Parameters:
779: . eps - eigensolver context obtained from EPSCreate()
781: Output Parameter:
782: . ip - inner product context
784: Level: advanced
786: .seealso: EPSSetIP()
787: @*/
788: PetscErrorCode EPSGetIP(EPS eps,IP *ip)
789: {
793: *ip = eps->ip;
794: return(0);
795: }
799: /*@
800: EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
801: eigenvalue problem.
803: Not collective
805: Input Parameter:
806: . eps - the eigenproblem solver context
808: Output Parameter:
809: . is - the answer
811: Level: intermediate
813: @*/
814: PetscErrorCode EPSIsGeneralized(EPS eps,PetscTruth* is)
815: {
817: Mat B;
821: STGetOperators(eps->OP,PETSC_NULL,&B);
822: if( B ) *is = PETSC_TRUE;
823: else *is = PETSC_FALSE;
824: if( eps->setupcalled ) {
825: if( eps->isgeneralized != *is ) {
826: SETERRQ(0,"Warning: Inconsistent EPS state");
827: }
828: }
829: return(0);
830: }
834: /*@
835: EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
836: eigenvalue problem.
838: Not collective
840: Input Parameter:
841: . eps - the eigenproblem solver context
843: Output Parameter:
844: . is - the answer
846: Level: intermediate
848: @*/
849: PetscErrorCode EPSIsHermitian(EPS eps,PetscTruth* is)
850: {
853: if( eps->ishermitian ) *is = PETSC_TRUE;
854: else *is = PETSC_FALSE;
855: return(0);
856: }