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-2012, Universitat 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 <slepc-private/epsimpl.h> /*I "slepceps.h" I*/
26: PetscFList EPSList = 0;
27: PetscBool EPSRegisterAllCalled = PETSC_FALSE;
28: PetscClassId EPS_CLASSID = 0;
29: PetscLogEvent EPS_SetUp = 0,EPS_Solve = 0;
30: static PetscBool EPSPackageInitialized = PETSC_FALSE;
32: const char *EPSPowerShiftTypes[] = {"CONSTANT","RAYLEIGH","WILKINSON","EPSPowerShiftType","EPS_POWER_SHIFT_",0};
33: const char *EPSLanczosReorthogTypes[] = {"LOCAL","FULL","SELECTIVE","PERIODIC","PARTIAL","DELAYED","EPSLanczosReorthogType","EPS_LANCZOS_REORTHOG_",0};
34: const char *EPSPRIMMEMethods[] = {"DYNAMIC","DEFAULT_MIN_TIME","DEFAULT_MIN_MATVECS","ARNOLDI","GD","GD_PLUSK","GD_OLSEN_PLUSK","JD_OLSEN_PLUSK","RQI","JDQR","JDQMR","JDQMR_ETOL","SUBSPACE_ITERATION","LOBPCG_ORTHOBASIS","LOBPCG_ORTHOBASISW","EPSPRIMMEMethod","EPS_PRIMME_",0};
38: /*@C
39: EPSFinalizePackage - This function destroys everything in the Slepc interface to the EPS package. It is
40: called from SlepcFinalize().
42: Level: developer
44: .seealso: SlepcFinalize()
45: @*/
46: PetscErrorCode EPSFinalizePackage(void)
47: {
49: EPSPackageInitialized = PETSC_FALSE;
50: EPSList = 0;
51: EPSRegisterAllCalled = PETSC_FALSE;
52: return(0);
53: }
57: /*@C
58: EPSInitializePackage - This function initializes everything in the EPS package. It is called
59: from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to EPSCreate()
60: when using static libraries.
62: Input Parameter:
63: path - The dynamic library path, or PETSC_NULL
65: Level: developer
67: .seealso: SlepcInitialize()
68: @*/
69: PetscErrorCode EPSInitializePackage(const char *path) {
70: char logList[256];
71: char *className;
72: PetscBool opt;
76: if (EPSPackageInitialized) return(0);
77: EPSPackageInitialized = PETSC_TRUE;
78: /* Register Classes */
79: PetscClassIdRegister("Eigenproblem Solver",&EPS_CLASSID);
80: /* Register Constructors */
81: EPSRegisterAll(path);
82: /* Register Events */
83: PetscLogEventRegister("EPSSetUp",EPS_CLASSID,&EPS_SetUp);
84: PetscLogEventRegister("EPSSolve",EPS_CLASSID,&EPS_Solve);
85: /* Process info exclusions */
86: PetscOptionsGetString(PETSC_NULL,"-info_exclude",logList,256,&opt);
87: if (opt) {
88: PetscStrstr(logList,"eps",&className);
89: if (className) {
90: PetscInfoDeactivateClass(EPS_CLASSID);
91: }
92: }
93: /* Process summary exclusions */
94: PetscOptionsGetString(PETSC_NULL,"-log_summary_exclude",logList,256,&opt);
95: if (opt) {
96: PetscStrstr(logList,"eps",&className);
97: if (className) {
98: PetscLogEventDeactivateClass(EPS_CLASSID);
99: }
100: }
101: PetscRegisterFinalize(EPSFinalizePackage);
102: return(0);
103: }
107: /*@C
108: EPSView - Prints the EPS data structure.
110: Collective on EPS
112: Input Parameters:
113: + eps - the eigenproblem solver context
114: - viewer - optional visualization context
116: Options Database Key:
117: . -eps_view - Calls EPSView() at end of EPSSolve()
119: Note:
120: The available visualization contexts include
121: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
122: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
123: output where only the first processor opens
124: the file. All other processors send their
125: data to the first processor to print.
127: The user can open an alternative visualization context with
128: PetscViewerASCIIOpen() - output to a specified file.
130: Level: beginner
132: .seealso: STView(), PetscViewerASCIIOpen()
133: @*/
134: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
135: {
137: const char *type,*extr,*bal;
138: PetscBool isascii,ispower,isexternal;
142: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(((PetscObject)eps)->comm);
146: #if defined(PETSC_USE_COMPLEX)
147: #define HERM "hermitian"
148: #else
149: #define HERM "symmetric"
150: #endif
151: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
152: if (isascii) {
153: PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer,"EPS Object");
154: if (eps->ops->view) {
155: PetscViewerASCIIPushTab(viewer);
156: (*eps->ops->view)(eps,viewer);
157: PetscViewerASCIIPopTab(viewer);
158: }
159: if (eps->problem_type) {
160: switch (eps->problem_type) {
161: case EPS_HEP: type = HERM " eigenvalue problem"; break;
162: case EPS_GHEP: type = "generalized " HERM " eigenvalue problem"; break;
163: case EPS_NHEP: type = "non-" HERM " eigenvalue problem"; break;
164: case EPS_GNHEP: type = "generalized non-" HERM " eigenvalue problem"; break;
165: case EPS_PGNHEP: type = "generalized non-" HERM " eigenvalue problem with " HERM " positive definite B"; break;
166: case EPS_GHIEP: type = "generalized " HERM "-indefinite eigenvalue problem"; break;
167: default: SETERRQ(((PetscObject)eps)->comm,1,"Wrong value of eps->problem_type");
168: }
169: } else type = "not yet set";
170: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
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(((PetscObject)eps)->comm,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(((PetscObject)eps)->comm,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: case EPS_ALL:
245: PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%G,%G]\n",eps->inta,eps->intb);
246: break;
247: default: SETERRQ(((PetscObject)eps)->comm,1,"Wrong value of eps->which");
248: }
249: if (eps->leftvecs) {
250: PetscViewerASCIIPrintf(viewer," computing left eigenvectors also\n");
251: }
252: if (eps->trueres) {
253: PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n");
254: }
255: if (eps->trackall) {
256: PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n");
257: }
258: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",eps->nev);
259: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",eps->ncv);
260: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",eps->mpd);
261: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",eps->max_it);
262: PetscViewerASCIIPrintf(viewer," tolerance: %G\n",eps->tol);
263: PetscViewerASCIIPrintf(viewer," convergence test: ");
264: switch(eps->conv) {
265: case EPS_CONV_ABS:
266: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
267: case EPS_CONV_EIG:
268: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
269: case EPS_CONV_NORM:
270: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");break;
271: default:
272: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
273: }
274: if (eps->nini!=0) {
275: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(eps->nini));
276: }
277: if (eps->ninil!=0) {
278: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial left space: %D\n",PetscAbs(eps->ninil));
279: }
280: if (eps->nds>0) {
281: PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %D\n",eps->nds);
282: }
283: PetscViewerASCIIPrintf(viewer," estimates of matrix norms (%s): norm(A)=%G",eps->adaptive?"adaptive":"constant",eps->nrma);
284: if (eps->isgeneralized) {
285: PetscViewerASCIIPrintf(viewer,", norm(B)=%G",eps->nrmb);
286: }
287: PetscViewerASCIIPrintf(viewer,"\n");
288: } else {
289: if (eps->ops->view) {
290: (*eps->ops->view)(eps,viewer);
291: }
292: }
293: PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLZPACK,EPSTRLAN,EPSBLOPEX,EPSPRIMME,"");
294: if (!isexternal) {
295: if (!eps->ip) { EPSGetIP(eps,&eps->ip); }
296: IPView(eps->ip,viewer);
297: PetscObjectTypeCompare((PetscObject)eps,EPSPOWER,&ispower);
298: if (!ispower) {
299: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
300: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
301: DSView(eps->ds,viewer);
302: PetscViewerPopFormat(viewer);
303: }
304: }
305: if (!eps->OP) { EPSGetST(eps,&eps->OP); }
306: STView(eps->OP,viewer);
307: return(0);
308: }
312: /*@
313: EPSPrintSolution - Prints the computed eigenvalues.
315: Collective on EPS
317: Input Parameters:
318: + eps - the eigensolver context
319: - viewer - optional visualization context
321: Options Database:
322: . -eps_terse - print only minimal information
324: Note:
325: By default, this function prints a table with eigenvalues and associated
326: relative errors. With -eps_terse only the eigenvalues are printed.
328: Level: intermediate
330: .seealso: PetscViewerASCIIOpen()
331: @*/
332: PetscErrorCode EPSPrintSolution(EPS eps,PetscViewer viewer)
333: {
334: PetscBool terse,errok,isascii;
335: PetscReal error,re,im;
336: PetscScalar kr,ki;
337: PetscInt i,j;
342: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(((PetscObject)eps)->comm);
345: if (!eps->eigr || !eps->eigi || !eps->V) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_WRONGSTATE,"EPSSolve must be called first");
346: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
347: if (!isascii) return(0);
349: PetscOptionsHasName(PETSC_NULL,"-eps_terse",&terse);
350: if (terse) {
351: if (eps->nconv<eps->nev) {
352: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
353: } else {
354: errok = PETSC_TRUE;
355: for (i=0;i<eps->nev;i++) {
356: EPSComputeRelativeError(eps,i,&error);
357: errok = (errok && error<eps->tol)? PETSC_TRUE: PETSC_FALSE;
358: }
359: if (errok) {
360: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
361: for (i=0;i<=(eps->nev-1)/8;i++) {
362: PetscViewerASCIIPrintf(viewer,"\n ");
363: for (j=0;j<PetscMin(8,eps->nev-8*i);j++) {
364: EPSGetEigenpair(eps,8*i+j,&kr,&ki,PETSC_NULL,PETSC_NULL);
365: #if defined(PETSC_USE_COMPLEX)
366: re = PetscRealPart(kr);
367: im = PetscImaginaryPart(kr);
368: #else
369: re = kr;
370: im = ki;
371: #endif
372: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
373: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
374: if (im!=0.0) {
375: PetscViewerASCIIPrintf(viewer,"%.5F%+.5Fi",re,im);
376: } else {
377: PetscViewerASCIIPrintf(viewer,"%.5F",re);
378: }
379: if (8*i+j+1<eps->nev) { PetscViewerASCIIPrintf(viewer,", "); }
380: }
381: }
382: PetscViewerASCIIPrintf(viewer,"\n\n");
383: } else {
384: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",eps->nev);
385: }
386: }
387: } else {
388: PetscViewerASCIIPrintf(viewer," Number of converged approximate eigenpairs: %D\n\n",eps->nconv);
389: if (eps->nconv>0) {
390: PetscViewerASCIIPrintf(viewer,
391: " k ||Ax-k%sx||/||kx||\n"
392: " ----------------- ------------------\n",eps->isgeneralized?"B":"");
393: for (i=0;i<eps->nconv;i++) {
394: EPSGetEigenpair(eps,i,&kr,&ki,PETSC_NULL,PETSC_NULL);
395: EPSComputeRelativeError(eps,i,&error);
396: #if defined(PETSC_USE_COMPLEX)
397: re = PetscRealPart(kr);
398: im = PetscImaginaryPart(kr);
399: #else
400: re = kr;
401: im = ki;
402: #endif
403: if (im!=0.0) {
404: PetscViewerASCIIPrintf(viewer," % 9F%+9F i %12G\n",re,im,error);
405: } else {
406: PetscViewerASCIIPrintf(viewer," % 12F %12G\n",re,error);
407: }
408: }
409: PetscViewerASCIIPrintf(viewer,"\n");
410: }
411: }
412: return(0);
413: }
417: /*@C
418: EPSCreate - Creates the default EPS context.
420: Collective on MPI_Comm
422: Input Parameter:
423: . comm - MPI communicator
425: Output Parameter:
426: . eps - location to put the EPS context
428: Note:
429: The default EPS type is EPSKRYLOVSCHUR
431: Level: beginner
433: .seealso: EPSSetUp(), EPSSolve(), EPSDestroy(), EPS
434: @*/
435: PetscErrorCode EPSCreate(MPI_Comm comm,EPS *outeps)
436: {
438: EPS eps;
442: *outeps = 0;
444: SlepcHeaderCreate(eps,_p_EPS,struct _EPSOps,EPS_CLASSID,-1,"EPS","Eigenvalue Problem Solver","EPS",comm,EPSDestroy,EPSView);
446: eps->max_it = 0;
447: eps->nev = 1;
448: eps->ncv = 0;
449: eps->mpd = 0;
450: eps->allocated_ncv = 0;
451: eps->nini = 0;
452: eps->ninil = 0;
453: eps->nds = 0;
454: eps->tol = PETSC_DEFAULT;
455: eps->conv = EPS_CONV_EIG;
456: eps->conv_func = EPSConvergedEigRelative;
457: eps->conv_ctx = PETSC_NULL;
458: eps->which = (EPSWhich)0;
459: eps->which_func = PETSC_NULL;
460: eps->which_ctx = PETSC_NULL;
461: eps->arbit_func = PETSC_NULL;
462: eps->arbit_ctx = PETSC_NULL;
463: eps->leftvecs = PETSC_FALSE;
464: eps->trueres = PETSC_FALSE;
465: eps->trackall = PETSC_FALSE;
466: eps->target = 0.0;
467: eps->inta = 0.0;
468: eps->intb = 0.0;
469: eps->evecsavailable = PETSC_FALSE;
470: eps->problem_type = (EPSProblemType)0;
471: eps->extraction = (EPSExtraction)0;
472: eps->balance = (EPSBalance)0;
473: eps->balance_its = 5;
474: eps->balance_cutoff = 1e-8;
475: eps->nrma = 1.0;
476: eps->nrmb = 1.0;
477: eps->adaptive = PETSC_FALSE;
479: eps->V = 0;
480: eps->W = 0;
481: eps->D = 0;
482: eps->defl = 0;
483: eps->IS = 0;
484: eps->ISL = 0;
485: eps->t = 0;
486: eps->ds_ortho = PETSC_FALSE;
487: eps->eigr = 0;
488: eps->eigi = 0;
489: eps->errest = 0;
490: eps->errest_left = 0;
491: eps->OP = 0;
492: eps->ip = 0;
493: eps->ds = 0;
494: eps->rand = 0;
495: eps->data = 0;
496: eps->nconv = 0;
497: eps->its = 0;
498: eps->perm = PETSC_NULL;
500: eps->nwork = 0;
501: eps->work = 0;
502: eps->isgeneralized = PETSC_FALSE;
503: eps->ishermitian = PETSC_FALSE;
504: eps->ispositive = PETSC_FALSE;
505: eps->setupcalled = 0;
506: eps->reason = EPS_CONVERGED_ITERATING;
507: eps->numbermonitors = 0;
509: PetscRandomCreate(comm,&eps->rand);
510: PetscRandomSetSeed(eps->rand,0x12345678);
511: PetscLogObjectParent(eps,eps->rand);
512: *outeps = eps;
513: return(0);
514: }
515:
518: /*@C
519: EPSSetType - Selects the particular solver to be used in the EPS object.
521: Logically Collective on EPS
523: Input Parameters:
524: + eps - the eigensolver context
525: - type - a known method
527: Options Database Key:
528: . -eps_type <method> - Sets the method; use -help for a list
529: of available methods
530:
531: Notes:
532: See "slepc/include/slepceps.h" for available methods. The default
533: is EPSKRYLOVSCHUR.
535: Normally, it is best to use the EPSSetFromOptions() command and
536: then set the EPS type from the options database rather than by using
537: this routine. Using the options database provides the user with
538: maximum flexibility in evaluating the different available methods.
539: The EPSSetType() routine is provided for those situations where it
540: is necessary to set the iterative solver independently of the command
541: line or options database.
543: Level: intermediate
545: .seealso: STSetType(), EPSType
546: @*/
547: PetscErrorCode EPSSetType(EPS eps,const EPSType type)
548: {
549: PetscErrorCode ierr,(*r)(EPS);
550: PetscBool match;
556: PetscObjectTypeCompare((PetscObject)eps,type,&match);
557: if (match) return(0);
559: PetscFListFind(EPSList,((PetscObject)eps)->comm,type,PETSC_TRUE,(void (**)(void)) &r);
560: if (!r) SETERRQ1(((PetscObject)eps)->comm,PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown EPS type given: %s",type);
562: if (eps->ops->destroy) { (*eps->ops->destroy)(eps); }
563: PetscMemzero(eps->ops,sizeof(struct _EPSOps));
565: eps->setupcalled = 0;
566: PetscObjectChangeTypeName((PetscObject)eps,type);
567: (*r)(eps);
568: return(0);
569: }
573: /*@C
574: EPSGetType - Gets the EPS type as a string from the EPS object.
576: Not Collective
578: Input Parameter:
579: . eps - the eigensolver context
581: Output Parameter:
582: . name - name of EPS method
584: Level: intermediate
586: .seealso: EPSSetType()
587: @*/
588: PetscErrorCode EPSGetType(EPS eps,const EPSType *type)
589: {
593: *type = ((PetscObject)eps)->type_name;
594: return(0);
595: }
599: /*@C
600: EPSRegister - See EPSRegisterDynamic()
602: Level: advanced
603: @*/
604: PetscErrorCode EPSRegister(const char *sname,const char *path,const char *name,PetscErrorCode (*function)(EPS))
605: {
607: char fullname[PETSC_MAX_PATH_LEN];
610: PetscFListConcat(path,name,fullname);
611: PetscFListAdd(&EPSList,sname,fullname,(void (*)(void))function);
612: return(0);
613: }
617: /*@
618: EPSRegisterDestroy - Frees the list of EPS methods that were
619: registered by EPSRegisterDynamic().
621: Not Collective
623: Level: advanced
625: .seealso: EPSRegisterDynamic(), EPSRegisterAll()
626: @*/
627: PetscErrorCode EPSRegisterDestroy(void)
628: {
632: PetscFListDestroy(&EPSList);
633: EPSRegisterAllCalled = PETSC_FALSE;
634: return(0);
635: }
639: /*@
640: EPSReset - Resets the EPS context to the setupcalled=0 state and removes any
641: allocated objects.
643: Collective on EPS
645: Input Parameter:
646: . eps - eigensolver context obtained from EPSCreate()
648: Level: advanced
650: .seealso: EPSDestroy()
651: @*/
652: PetscErrorCode EPSReset(EPS eps)
653: {
658: if (eps->ops->reset) { (eps->ops->reset)(eps); }
659: if (eps->OP) { STReset(eps->OP); }
660: if (eps->ip) { IPReset(eps->ip); }
661: if (eps->ds) { DSReset(eps->ds); }
662: VecDestroy(&eps->t);
663: VecDestroy(&eps->D);
664: eps->setupcalled = 0;
665: return(0);
666: }
670: /*@C
671: EPSDestroy - Destroys the EPS context.
673: Collective on EPS
675: Input Parameter:
676: . eps - eigensolver context obtained from EPSCreate()
678: Level: beginner
680: .seealso: EPSCreate(), EPSSetUp(), EPSSolve()
681: @*/
682: PetscErrorCode EPSDestroy(EPS *eps)
683: {
687: if (!*eps) return(0);
689: if (--((PetscObject)(*eps))->refct > 0) { *eps = 0; return(0); }
690: EPSReset(*eps);
691: PetscObjectDepublish(*eps);
692: if ((*eps)->ops->destroy) { (*(*eps)->ops->destroy)(*eps); }
693: STDestroy(&(*eps)->OP);
694: IPDestroy(&(*eps)->ip);
695: DSDestroy(&(*eps)->ds);
696: PetscRandomDestroy(&(*eps)->rand);
697: EPSRemoveDeflationSpace(*eps);
698: EPSMonitorCancel(*eps);
699: PetscHeaderDestroy(eps);
700: return(0);
701: }
705: /*@
706: EPSSetTarget - Sets the value of the target.
708: Logically Collective on EPS
710: Input Parameters:
711: + eps - eigensolver context
712: - target - the value of the target
714: Notes:
715: The target is a scalar value used to determine the portion of the spectrum
716: of interest. It is used in combination with EPSSetWhichEigenpairs().
717:
718: Level: beginner
720: .seealso: EPSGetTarget(), EPSSetWhichEigenpairs()
721: @*/
722: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
723: {
729: eps->target = target;
730: if (!eps->OP) { EPSGetST(eps,&eps->OP); }
731: STSetDefaultShift(eps->OP,target);
732: return(0);
733: }
737: /*@
738: EPSGetTarget - Gets the value of the target.
740: Not Collective
742: Input Parameter:
743: . eps - eigensolver context
745: Output Parameter:
746: . target - the value of the target
748: Level: beginner
750: Note:
751: If the target was not set by the user, then zero is returned.
753: .seealso: EPSSetTarget()
754: @*/
755: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
756: {
760: *target = eps->target;
761: return(0);
762: }
766: /*@
767: EPSSetInterval - Defines the computational interval for spectrum slicing.
769: Logically Collective on EPS
771: Input Parameters:
772: + eps - eigensolver context
773: . inta - left end of the interval
774: - intb - right end of the interval
776: Options Database Key:
777: . -eps_interval <a,b> - set [a,b] as the interval of interest
779: Notes:
780: Spectrum slicing is a technique employed for computing all eigenvalues of
781: symmetric eigenproblems in a given interval. This function provides the
782: interval to be considered. It must be used in combination with EPS_ALL, see
783: EPSSetWhichEigenpairs().
785: In the command-line option, two values must be provided. For an open interval,
786: one can give an infinite, e.g., -eps_interval 1.0,inf or -eps_interval -inf,1.0.
787: An open interval in the programmatic interface can be specified with
788: PETSC_MAX_REAL and -PETSC_MAX_REAL.
789:
790: Level: intermediate
792: .seealso: EPSGetInterval(), EPSSetWhichEigenpairs()
793: @*/
794: PetscErrorCode EPSSetInterval(EPS eps,PetscReal inta,PetscReal intb)
795: {
800: if (inta>=intb) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
801: eps->inta = inta;
802: eps->intb = intb;
803: return(0);
804: }
808: /*@
809: EPSGetInterval - Gets the computational interval for spectrum slicing.
811: Not Collective
813: Input Parameter:
814: . eps - eigensolver context
816: Output Parameters:
817: + inta - left end of the interval
818: - intb - right end of the interval
820: Level: intermediate
822: Note:
823: If the interval was not set by the user, then zeros are returned.
825: .seealso: EPSSetInterval()
826: @*/
827: PetscErrorCode EPSGetInterval(EPS eps,PetscReal* inta,PetscReal* intb)
828: {
833: if (inta) *inta = eps->inta;
834: if (intb) *intb = eps->intb;
835: return(0);
836: }
840: /*@
841: EPSSetST - Associates a spectral transformation object to the eigensolver.
843: Collective on EPS
845: Input Parameters:
846: + eps - eigensolver context obtained from EPSCreate()
847: - st - the spectral transformation object
849: Note:
850: Use EPSGetST() to retrieve the spectral transformation context (for example,
851: to free it at the end of the computations).
853: Level: developer
855: .seealso: EPSGetST()
856: @*/
857: PetscErrorCode EPSSetST(EPS eps,ST st)
858: {
865: PetscObjectReference((PetscObject)st);
866: STDestroy(&eps->OP);
867: eps->OP = st;
868: PetscLogObjectParent(eps,eps->OP);
869: return(0);
870: }
874: /*@C
875: EPSGetST - Obtain the spectral transformation (ST) object associated
876: to the eigensolver object.
878: Not Collective
880: Input Parameters:
881: . eps - eigensolver context obtained from EPSCreate()
883: Output Parameter:
884: . st - spectral transformation context
886: Level: beginner
888: .seealso: EPSSetST()
889: @*/
890: PetscErrorCode EPSGetST(EPS eps,ST *st)
891: {
897: if (!eps->OP) {
898: STCreate(((PetscObject)eps)->comm,&eps->OP);
899: PetscLogObjectParent(eps,eps->OP);
900: }
901: *st = eps->OP;
902: return(0);
903: }
907: /*@
908: EPSSetIP - Associates an inner product object to the eigensolver.
910: Collective on EPS
912: Input Parameters:
913: + eps - eigensolver context obtained from EPSCreate()
914: - ip - the inner product object
916: Note:
917: Use EPSGetIP() to retrieve the inner product context (for example,
918: to free it at the end of the computations).
920: Level: advanced
922: .seealso: EPSGetIP()
923: @*/
924: PetscErrorCode EPSSetIP(EPS eps,IP ip)
925: {
932: PetscObjectReference((PetscObject)ip);
933: IPDestroy(&eps->ip);
934: eps->ip = ip;
935: PetscLogObjectParent(eps,eps->ip);
936: return(0);
937: }
941: /*@C
942: EPSGetIP - Obtain the inner product object associated to the eigensolver object.
944: Not Collective
946: Input Parameters:
947: . eps - eigensolver context obtained from EPSCreate()
949: Output Parameter:
950: . ip - inner product context
952: Level: advanced
954: .seealso: EPSSetIP()
955: @*/
956: PetscErrorCode EPSGetIP(EPS eps,IP *ip)
957: {
963: if (!eps->ip) {
964: IPCreate(((PetscObject)eps)->comm,&eps->ip);
965: PetscLogObjectParent(eps,eps->ip);
966: }
967: *ip = eps->ip;
968: return(0);
969: }
973: /*@
974: EPSSetDS - Associates a direct solver object to the eigensolver.
976: Collective on EPS
978: Input Parameters:
979: + eps - eigensolver context obtained from EPSCreate()
980: - ds - the direct solver object
982: Note:
983: Use EPSGetDS() to retrieve the direct solver context (for example,
984: to free it at the end of the computations).
986: Level: advanced
988: .seealso: EPSGetDS()
989: @*/
990: PetscErrorCode EPSSetDS(EPS eps,DS ds)
991: {
998: PetscObjectReference((PetscObject)ds);
999: DSDestroy(&eps->ds);
1000: eps->ds = ds;
1001: PetscLogObjectParent(eps,eps->ds);
1002: return(0);
1003: }
1007: /*@C
1008: EPSGetDS - Obtain the direct solver object associated to the eigensolver object.
1010: Not Collective
1012: Input Parameters:
1013: . eps - eigensolver context obtained from EPSCreate()
1015: Output Parameter:
1016: . ds - direct solver context
1018: Level: advanced
1020: .seealso: EPSSetDS()
1021: @*/
1022: PetscErrorCode EPSGetDS(EPS eps,DS *ds)
1023: {
1029: if (!eps->ds) {
1030: DSCreate(((PetscObject)eps)->comm,&eps->ds);
1031: PetscLogObjectParent(eps,eps->ds);
1032: }
1033: *ds = eps->ds;
1034: return(0);
1035: }
1039: /*@
1040: EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
1041: eigenvalue problem.
1043: Not collective
1045: Input Parameter:
1046: . eps - the eigenproblem solver context
1048: Output Parameter:
1049: . is - the answer
1051: Level: intermediate
1053: .seealso: EPSIsHermitian(), EPSIsPositive()
1054: @*/
1055: PetscErrorCode EPSIsGeneralized(EPS eps,PetscBool* is)
1056: {
1060: *is = eps->isgeneralized;
1061: return(0);
1062: }
1066: /*@
1067: EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
1068: eigenvalue problem.
1070: Not collective
1072: Input Parameter:
1073: . eps - the eigenproblem solver context
1075: Output Parameter:
1076: . is - the answer
1078: Level: intermediate
1080: .seealso: EPSIsGeneralized(), EPSIsPositive()
1081: @*/
1082: PetscErrorCode EPSIsHermitian(EPS eps,PetscBool* is)
1083: {
1087: *is = eps->ishermitian;
1088: return(0);
1089: }
1093: /*@
1094: EPSIsPositive - Ask if the EPS object corresponds to an eigenvalue
1095: problem type that requires a positive (semi-) definite matrix B.
1097: Not collective
1099: Input Parameter:
1100: . eps - the eigenproblem solver context
1102: Output Parameter:
1103: . is - the answer
1105: Level: intermediate
1107: .seealso: EPSIsGeneralized(), EPSIsHermitian()
1108: @*/
1109: PetscErrorCode EPSIsPositive(EPS eps,PetscBool* is)
1110: {
1114: *is = eps->ispositive;
1115: return(0);
1116: }