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: if (eps->problem_type) {
150: switch (eps->problem_type) {
151: case EPS_HEP: type = HERM " eigenvalue problem"; break;
152: case EPS_GHEP: type = "generalized " HERM " eigenvalue problem"; break;
153: case EPS_NHEP: type = "non-" HERM " eigenvalue problem"; break;
154: case EPS_GNHEP: type = "generalized non-" HERM " eigenvalue problem"; break;
155: case EPS_PGNHEP: type = "generalized non-" HERM " eigenvalue problem with " HERM " positive definite B"; break;
156: default: SETERRQ(1,"Wrong value of eps->problem_type");
157: }
158: } else type = "not yet set";
159: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
160: EPSGetType(eps,&type);
161: if (type) {
162: PetscViewerASCIIPrintf(viewer," method: %s\n",type);
163: } else {
164: PetscViewerASCIIPrintf(viewer," method: not yet set\n");
165: }
166: if (eps->ops->view) {
167: PetscViewerASCIIPushTab(viewer);
168: (*eps->ops->view)(eps,viewer);
169: PetscViewerASCIIPopTab(viewer);
170: }
171: if (eps->extraction) {
172: switch (eps->extraction) {
173: case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
174: case EPS_HARMONIC: extr = "harmonic Ritz"; break;
175: case EPS_HARMONIC_RELATIVE:extr = "relative harmonic Ritz"; break;
176: case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
177: case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
178: case EPS_REFINED: extr = "refined Ritz"; break;
179: case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
180: default: SETERRQ(1,"Wrong value of eps->extraction");
181: }
182: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr);
183: }
184: if (eps->balance && !eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
185: switch (eps->balance) {
186: case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
187: case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
188: case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
189: default: SETERRQ(1,"Wrong value of eps->balance");
190: }
191: PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal);
192: if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) {
193: PetscViewerASCIIPrintf(viewer,", with its=%d",eps->balance_its);
194: }
195: if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) {
196: PetscViewerASCIIPrintf(viewer," and cutoff=%g",eps->balance_cutoff);
197: }
198: PetscViewerASCIIPrintf(viewer,"\n");
199: }
200: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
201: if (!eps->which) {
202: PetscViewerASCIIPrintf(viewer,"not yet set\n");
203: } else switch (eps->which) {
204: case EPS_WHICH_USER:
205: PetscViewerASCIIPrintf(viewer,"user defined\n");
206: break;
207: case EPS_TARGET_MAGNITUDE:
208: #if !defined(PETSC_USE_COMPLEX)
209: PetscViewerASCIIPrintf(viewer,"closest to target: %g (in magnitude)\n",eps->target);
210: #else
211: PetscViewerASCIIPrintf(viewer,"closest to target: %g+%g i (in magnitude)\n",PetscRealPart(eps->target),PetscImaginaryPart(eps->target));
212: #endif
213: break;
214: case EPS_TARGET_REAL:
215: #if !defined(PETSC_USE_COMPLEX)
216: PetscViewerASCIIPrintf(viewer,"closest to target: %g (along the real axis)\n",eps->target);
217: #else
218: PetscViewerASCIIPrintf(viewer,"closest to target: %g+%g i (along the real axis)\n",PetscRealPart(eps->target),PetscImaginaryPart(eps->target));
219: #endif
220: break;
221: #if defined(PETSC_USE_COMPLEX)
222: case EPS_TARGET_IMAGINARY:
223: PetscViewerASCIIPrintf(viewer,"closest to target: %g+%g i (along the imaginary axis)\n",PetscRealPart(eps->target),PetscImaginaryPart(eps->target));
224: break;
225: #endif
226: case EPS_LARGEST_MAGNITUDE:
227: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
228: break;
229: case EPS_SMALLEST_MAGNITUDE:
230: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
231: break;
232: case EPS_LARGEST_REAL:
233: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
234: break;
235: case EPS_SMALLEST_REAL:
236: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
237: break;
238: case EPS_LARGEST_IMAGINARY:
239: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
240: break;
241: case EPS_SMALLEST_IMAGINARY:
242: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
243: break;
244: default: SETERRQ(1,"Wrong value of eps->which");
245: }
246: if (eps->leftvecs) {
247: PetscViewerASCIIPrintf(viewer," computing left eigenvectors also\n");
248: }
249: if (eps->trueres) {
250: PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n");
251: }
252: if (eps->trackall) {
253: PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n");
254: }
255: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %d\n",eps->nev);
256: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %d\n",eps->ncv);
257: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %d\n",eps->mpd);
258: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %d\n", eps->max_it);
259: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",eps->tol);
260: PetscViewerASCIIPrintf(viewer," convergence test: ");
261: switch(eps->conv) {
262: case EPS_CONV_ABS:
263: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
264: case EPS_CONV_EIG:
265: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
266: case EPS_CONV_NORM:
267: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");break;
268: default:
269: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
270: }
271: if (eps->nini!=0) {
272: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %d\n",PetscAbs(eps->nini));
273: }
274: if (eps->ninil!=0) {
275: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial left space: %d\n",PetscAbs(eps->ninil));
276: }
277: if (eps->nds>0) {
278: PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %d\n",eps->nds);
279: }
280: PetscViewerASCIIPrintf(viewer," estimates of matrix norms (%s): norm(A)=%g",eps->adaptive?"adaptive":"constant",eps->nrma);
281: if (eps->isgeneralized) {
282: PetscViewerASCIIPrintf(viewer,", norm(B)=%g",eps->nrmb);
283: }
284: PetscViewerASCIIPrintf(viewer,"\n");
285: PetscViewerASCIIPushTab(viewer);
286: IPView(eps->ip,viewer);
287: STView(eps->OP,viewer);
288: PetscViewerASCIIPopTab(viewer);
289: } else {
290: if (eps->ops->view) {
291: (*eps->ops->view)(eps,viewer);
292: }
293: STView(eps->OP,viewer);
294: }
295: return(0);
296: }
300: /*@C
301: EPSCreate - Creates the default EPS context.
303: Collective on MPI_Comm
305: Input Parameter:
306: . comm - MPI communicator
308: Output Parameter:
309: . eps - location to put the EPS context
311: Note:
312: The default EPS type is EPSKRYLOVSCHUR
314: Level: beginner
316: .seealso: EPSSetUp(), EPSSolve(), EPSDestroy(), EPS
317: @*/
318: PetscErrorCode EPSCreate(MPI_Comm comm,EPS *outeps)
319: {
321: EPS eps;
325: *outeps = 0;
327: PetscHeaderCreate(eps,_p_EPS,struct _EPSOps,EPS_COOKIE,-1,"EPS",comm,EPSDestroy,EPSView);
328: *outeps = eps;
330: PetscMemzero(eps->ops,sizeof(struct _EPSOps));
332: eps->max_it = 0;
333: eps->nev = 1;
334: eps->ncv = 0;
335: eps->mpd = 0;
336: eps->allocated_ncv = 0;
337: eps->nini = 0;
338: eps->ninil = 0;
339: eps->nds = 0;
340: eps->tol = 1e-7;
341: eps->conv = EPS_CONV_EIG;
342: eps->conv_func = EPSEigRelativeConverged;
343: eps->conv_ctx = PETSC_NULL;
344: eps->which = (EPSWhich)0;
345: eps->which_func = PETSC_NULL;
346: eps->which_ctx = PETSC_NULL;
347: eps->leftvecs = PETSC_FALSE;
348: eps->trueres = PETSC_FALSE;
349: eps->trackall = PETSC_FALSE;
350: eps->target = 0.0;
351: eps->evecsavailable = PETSC_FALSE;
352: eps->problem_type = (EPSProblemType)0;
353: eps->extraction = (EPSExtraction)0;
354: eps->balance = (EPSBalance)0;
355: eps->balance_its = 5;
356: eps->balance_cutoff = 1e-8;
357: eps->nrma = 1.0;
358: eps->nrmb = 1.0;
359: eps->adaptive = PETSC_FALSE;
361: eps->V = 0;
362: eps->W = 0;
363: eps->T = 0;
364: eps->DS = 0;
365: eps->IS = 0;
366: eps->ISL = 0;
367: eps->ds_ortho = PETSC_FALSE;
368: eps->eigr = 0;
369: eps->eigi = 0;
370: eps->errest = 0;
371: eps->errest_left = 0;
372: eps->OP = 0;
373: eps->ip = 0;
374: eps->data = 0;
375: eps->nconv = 0;
376: eps->its = 0;
377: eps->perm = PETSC_NULL;
379: eps->nwork = 0;
380: eps->work = 0;
381: eps->isgeneralized = PETSC_FALSE;
382: eps->ishermitian = PETSC_FALSE;
383: eps->ispositive = PETSC_FALSE;
384: eps->setupcalled = 0;
385: eps->reason = EPS_CONVERGED_ITERATING;
387: eps->numbermonitors = 0;
389: STCreate(comm,&eps->OP);
390: PetscLogObjectParent(eps,eps->OP);
391: IPCreate(comm,&eps->ip);
392: IPSetOptionsPrefix(eps->ip,((PetscObject)eps)->prefix);
393: IPAppendOptionsPrefix(eps->ip,"eps_");
394: PetscLogObjectParent(eps,eps->ip);
395: PetscPublishAll(eps);
396: return(0);
397: }
398:
401: /*@C
402: EPSSetType - Selects the particular solver to be used in the EPS object.
404: Collective on EPS
406: Input Parameters:
407: + eps - the eigensolver context
408: - type - a known method
410: Options Database Key:
411: . -eps_type <method> - Sets the method; use -help for a list
412: of available methods
413:
414: Notes:
415: See "slepc/include/slepceps.h" for available methods. The default
416: is EPSKRYLOVSCHUR.
418: Normally, it is best to use the EPSSetFromOptions() command and
419: then set the EPS type from the options database rather than by using
420: this routine. Using the options database provides the user with
421: maximum flexibility in evaluating the different available methods.
422: The EPSSetType() routine is provided for those situations where it
423: is necessary to set the iterative solver independently of the command
424: line or options database.
426: Level: intermediate
428: .seealso: STSetType(), EPSType
429: @*/
430: PetscErrorCode EPSSetType(EPS eps,const EPSType type)
431: {
432: PetscErrorCode ierr,(*r)(EPS);
433: PetscTruth match;
439: PetscTypeCompare((PetscObject)eps,type,&match);
440: if (match) return(0);
442: if (eps->data) {
443: /* destroy the old private EPS context */
444: (*eps->ops->destroy)(eps);
445: eps->data = 0;
446: }
448: PetscFListFind(EPSList,((PetscObject)eps)->comm,type,(void (**)(void)) &r);
450: if (!r) SETERRQ1(1,"Unknown EPS type given: %s",type);
452: eps->setupcalled = 0;
453: PetscMemzero(eps->ops,sizeof(struct _EPSOps));
454: (*r)(eps);
456: PetscObjectChangeTypeName((PetscObject)eps,type);
457: return(0);
458: }
462: /*@C
463: EPSGetType - Gets the EPS type as a string from the EPS object.
465: Not Collective
467: Input Parameter:
468: . eps - the eigensolver context
470: Output Parameter:
471: . name - name of EPS method
473: Level: intermediate
475: .seealso: EPSSetType()
476: @*/
477: PetscErrorCode EPSGetType(EPS eps,const EPSType *type)
478: {
482: *type = ((PetscObject)eps)->type_name;
483: return(0);
484: }
486: /*MC
487: EPSRegisterDynamic - Adds a method to the eigenproblem solver package.
489: Synopsis:
490: EPSRegisterDynamic(char *name_solver,char *path,char *name_create,PetscErrorCode (*routine_create)(EPS))
492: Not Collective
494: Input Parameters:
495: + name_solver - name of a new user-defined solver
496: . path - path (either absolute or relative) the library containing this solver
497: . name_create - name of routine to create the solver context
498: - routine_create - routine to create the solver context
500: Notes:
501: EPSRegisterDynamic() may be called multiple times to add several user-defined solvers.
503: If dynamic libraries are used, then the fourth input argument (routine_create)
504: is ignored.
506: Sample usage:
507: .vb
508: EPSRegisterDynamic("my_solver",/home/username/my_lib/lib/libO/solaris/mylib.a,
509: "MySolverCreate",MySolverCreate);
510: .ve
512: Then, your solver can be chosen with the procedural interface via
513: $ EPSSetType(eps,"my_solver")
514: or at runtime via the option
515: $ -eps_type my_solver
517: Level: advanced
519: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
520: and others of the form ${any_environmental_variable} occuring in pathname will be
521: replaced with appropriate values.
523: .seealso: EPSRegisterDestroy(), EPSRegisterAll()
525: M*/
529: /*@C
530: EPSRegister - See EPSRegisterDynamic()
532: Level: advanced
533: @*/
534: PetscErrorCode EPSRegister(const char *sname,const char *path,const char *name,PetscErrorCode (*function)(EPS))
535: {
537: char fullname[256];
540: PetscFListConcat(path,name,fullname);
541: PetscFListAdd(&EPSList,sname,fullname,(void (*)(void))function);
542: return(0);
543: }
547: /*@
548: EPSRegisterDestroy - Frees the list of EPS methods that were
549: registered by EPSRegisterDynamic().
551: Not Collective
553: Level: advanced
555: .seealso: EPSRegisterDynamic(), EPSRegisterAll()
556: @*/
557: PetscErrorCode EPSRegisterDestroy(void)
558: {
562: PetscFListDestroy(&EPSList);
563: EPSRegisterAll(PETSC_NULL);
564: return(0);
565: }
569: /*@
570: EPSDestroy - Destroys the EPS context.
572: Collective on EPS
574: Input Parameter:
575: . eps - eigensolver context obtained from EPSCreate()
577: Level: beginner
579: .seealso: EPSCreate(), EPSSetUp(), EPSSolve()
580: @*/
581: PetscErrorCode EPSDestroy(EPS eps)
582: {
587: if (--((PetscObject)eps)->refct > 0) return(0);
589: /* if memory was published with AMS then destroy it */
590: PetscObjectDepublish(eps);
592: STDestroy(eps->OP);
593: IPDestroy(eps->ip);
595: if (eps->ops->destroy) {
596: (*eps->ops->destroy)(eps);
597: }
598:
599: PetscFree(eps->T);
600: PetscFree(eps->Tl);
601: PetscFree(eps->perm);
602: if (eps->rand) {
603: PetscRandomDestroy(eps->rand);
604: }
606: if (eps->D) {
607: VecDestroy(eps->D);
608: }
610: EPSRemoveDeflationSpace(eps);
611: EPSMonitorCancel(eps);
613: PetscHeaderDestroy(eps);
614: return(0);
615: }
619: /*@
620: EPSSetTarget - Sets the value of the target.
622: Collective on EPS
624: Input Parameters:
625: + eps - eigensolver context
626: - target - the value of the target
628: Notes:
629: The target is a scalar value used to determine the portion of the spectrum
630: of interest. It is used in combination with EPSSetWhichEigenpairs().
631:
632: Level: beginner
634: .seealso: EPSGetTarget(), EPSSetWhichEigenpairs()
635: @*/
636: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
637: {
642: eps->target = target;
643: STSetDefaultShift(eps->OP,target);
644: return(0);
645: }
649: /*@
650: EPSGetTarget - Gets the value of the target.
652: Not collective
654: Input Parameter:
655: . eps - eigensolver context
657: Output Parameter:
658: . target - the value of the target
660: Level: beginner
662: Note:
663: If the target was not set by the user, then zero is returned.
665: .seealso: EPSSetTarget()
666: @*/
667: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
668: {
671: if (target) {
672: *target = eps->target;
673: }
674: return(0);
675: }
679: /*@
680: EPSSetST - Associates a spectral transformation object to the
681: eigensolver.
683: Collective on EPS and ST
685: Input Parameters:
686: + eps - eigensolver context obtained from EPSCreate()
687: - st - the spectral transformation object
689: Note:
690: Use EPSGetST() to retrieve the spectral transformation context (for example,
691: to free it at the end of the computations).
693: Level: advanced
695: .seealso: EPSGetST()
696: @*/
697: PetscErrorCode EPSSetST(EPS eps,ST st)
698: {
705: PetscObjectReference((PetscObject)st);
706: STDestroy(eps->OP);
707: eps->OP = st;
708: return(0);
709: }
713: /*@C
714: EPSGetST - Obtain the spectral transformation (ST) object associated
715: to the eigensolver object.
717: Not Collective
719: Input Parameters:
720: . eps - eigensolver context obtained from EPSCreate()
722: Output Parameter:
723: . st - spectral transformation context
725: Level: beginner
727: .seealso: EPSSetST()
728: @*/
729: PetscErrorCode EPSGetST(EPS eps, ST *st)
730: {
734: *st = eps->OP;
735: return(0);
736: }
740: /*@
741: EPSSetIP - Associates an inner product object to the eigensolver.
743: Collective on EPS and IP
745: Input Parameters:
746: + eps - eigensolver context obtained from EPSCreate()
747: - ip - the inner product object
749: Note:
750: Use EPSGetIP() to retrieve the inner product context (for example,
751: to free it at the end of the computations).
753: Level: advanced
755: .seealso: EPSGetIP()
756: @*/
757: PetscErrorCode EPSSetIP(EPS eps,IP ip)
758: {
765: PetscObjectReference((PetscObject)ip);
766: IPDestroy(eps->ip);
767: eps->ip = ip;
768: return(0);
769: }
773: /*@C
774: EPSGetIP - Obtain the inner product object associated
775: to the eigensolver object.
777: Not Collective
779: Input Parameters:
780: . eps - eigensolver context obtained from EPSCreate()
782: Output Parameter:
783: . ip - inner product context
785: Level: advanced
787: .seealso: EPSSetIP()
788: @*/
789: PetscErrorCode EPSGetIP(EPS eps,IP *ip)
790: {
794: *ip = eps->ip;
795: return(0);
796: }
800: /*@
801: EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
802: eigenvalue problem.
804: Not collective
806: Input Parameter:
807: . eps - the eigenproblem solver context
809: Output Parameter:
810: . is - the answer
812: Level: intermediate
814: @*/
815: PetscErrorCode EPSIsGeneralized(EPS eps,PetscTruth* is)
816: {
818: Mat B;
822: STGetOperators(eps->OP,PETSC_NULL,&B);
823: if( B ) *is = PETSC_TRUE;
824: else *is = PETSC_FALSE;
825: if( eps->setupcalled ) {
826: if( eps->isgeneralized != *is ) {
827: SETERRQ(0,"Warning: Inconsistent EPS state");
828: }
829: }
830: return(0);
831: }
835: /*@
836: EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
837: eigenvalue problem.
839: Not collective
841: Input Parameter:
842: . eps - the eigenproblem solver context
844: Output Parameter:
845: . is - the answer
847: Level: intermediate
849: @*/
850: PetscErrorCode EPSIsHermitian(EPS eps,PetscTruth* is)
851: {
854: if( eps->ishermitian ) *is = PETSC_TRUE;
855: else *is = PETSC_FALSE;
856: return(0);
857: }