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-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/epsimpl.h> /*I "slepceps.h" I*/
26: PetscFunctionList 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
40: to the EPS package. It is called from SlepcFinalize().
42: Level: developer
44: .seealso: SlepcFinalize()
45: @*/
46: PetscErrorCode EPSFinalizePackage(void)
47: {
51: PetscFunctionListDestroy(&EPSList);
52: EPSPackageInitialized = PETSC_FALSE;
53: EPSRegisterAllCalled = PETSC_FALSE;
54: return(0);
55: }
59: /*@C
60: EPSInitializePackage - This function initializes everything in the EPS package.
61: It is called from PetscDLLibraryRegister() when using dynamic libraries, and
62: on the first call to EPSCreate() when using static libraries.
64: Level: developer
66: .seealso: SlepcInitialize()
67: @*/
68: PetscErrorCode EPSInitializePackage()
69: {
70: char logList[256];
71: char *className;
72: PetscBool opt;
76: if (EPSPackageInitialized) return(0);
77: EPSPackageInitialized = PETSC_TRUE;
78: /* Register Classes */
79: PetscClassIdRegister("Eigenvalue Problem Solver",&EPS_CLASSID);
80: /* Register Constructors */
81: EPSRegisterAll();
82: /* Register Events */
83: PetscLogEventRegister("EPSSetUp",EPS_CLASSID,&EPS_SetUp);
84: PetscLogEventRegister("EPSSolve",EPS_CLASSID,&EPS_Solve);
85: /* Process info exclusions */
86: PetscOptionsGetString(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(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: char str[50];
139: PetscBool isascii,ispower,isexternal;
143: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
147: #if defined(PETSC_USE_COMPLEX)
148: #define HERM "hermitian"
149: #else
150: #define HERM "symmetric"
151: #endif
152: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
153: if (isascii) {
154: PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer,"EPS Object");
155: if (eps->ops->view) {
156: PetscViewerASCIIPushTab(viewer);
157: (*eps->ops->view)(eps,viewer);
158: PetscViewerASCIIPopTab(viewer);
159: }
160: if (eps->problem_type) {
161: switch (eps->problem_type) {
162: case EPS_HEP: type = HERM " eigenvalue problem"; break;
163: case EPS_GHEP: type = "generalized " HERM " eigenvalue problem"; break;
164: case EPS_NHEP: type = "non-" HERM " eigenvalue problem"; break;
165: case EPS_GNHEP: type = "generalized non-" HERM " eigenvalue problem"; break;
166: case EPS_PGNHEP: type = "generalized non-" HERM " eigenvalue problem with " HERM " positive definite B"; break;
167: case EPS_GHIEP: type = "generalized " HERM "-indefinite eigenvalue problem"; break;
168: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->problem_type");
169: }
170: } else type = "not yet set";
171: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
172: if (eps->extraction) {
173: switch (eps->extraction) {
174: case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
175: case EPS_HARMONIC: extr = "harmonic Ritz"; break;
176: case EPS_HARMONIC_RELATIVE:extr = "relative harmonic Ritz"; break;
177: case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
178: case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
179: case EPS_REFINED: extr = "refined Ritz"; break;
180: case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
181: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->extraction");
182: }
183: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr);
184: }
185: if (eps->balance && !eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
186: switch (eps->balance) {
187: case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
188: case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
189: case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
190: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->balance");
191: }
192: PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal);
193: if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) {
194: PetscViewerASCIIPrintf(viewer,", with its=%D",eps->balance_its);
195: }
196: if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) {
197: PetscViewerASCIIPrintf(viewer," and cutoff=%G",eps->balance_cutoff);
198: }
199: PetscViewerASCIIPrintf(viewer,"\n");
200: }
201: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
202: SlepcSNPrintfScalar(str,50,eps->target,PETSC_FALSE);
203: if (!eps->which) {
204: PetscViewerASCIIPrintf(viewer,"not yet set\n");
205: } else switch (eps->which) {
206: case EPS_WHICH_USER:
207: PetscViewerASCIIPrintf(viewer,"user defined\n");
208: break;
209: case EPS_TARGET_MAGNITUDE:
210: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
211: break;
212: case EPS_TARGET_REAL:
213: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
214: break;
215: case EPS_TARGET_IMAGINARY:
216: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
217: break;
218: case EPS_LARGEST_MAGNITUDE:
219: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
220: break;
221: case EPS_SMALLEST_MAGNITUDE:
222: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
223: break;
224: case EPS_LARGEST_REAL:
225: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
226: break;
227: case EPS_SMALLEST_REAL:
228: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
229: break;
230: case EPS_LARGEST_IMAGINARY:
231: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
232: break;
233: case EPS_SMALLEST_IMAGINARY:
234: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
235: break;
236: case EPS_ALL:
237: PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%G,%G]\n",eps->inta,eps->intb);
238: break;
239: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->which");
240: }
241: if (eps->leftvecs) {
242: PetscViewerASCIIPrintf(viewer," computing left eigenvectors also\n");
243: }
244: if (eps->trueres) {
245: PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n");
246: }
247: if (eps->trackall) {
248: PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n");
249: }
250: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",eps->nev);
251: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",eps->ncv);
252: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",eps->mpd);
253: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",eps->max_it);
254: PetscViewerASCIIPrintf(viewer," tolerance: %G\n",eps->tol);
255: PetscViewerASCIIPrintf(viewer," convergence test: ");
256: switch (eps->conv) {
257: case EPS_CONV_ABS:
258: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
259: case EPS_CONV_EIG:
260: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
261: case EPS_CONV_NORM:
262: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");break;
263: default:
264: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
265: }
266: if (eps->nini) {
267: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(eps->nini));
268: }
269: if (eps->ninil) {
270: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial left space: %D\n",PetscAbs(eps->ninil));
271: }
272: if (eps->nds>0) {
273: PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %D\n",eps->nds);
274: }
275: PetscViewerASCIIPrintf(viewer," estimates of matrix norms (%s): norm(A)=%G",eps->adaptive?"adaptive":"constant",eps->nrma);
276: if (eps->isgeneralized) {
277: PetscViewerASCIIPrintf(viewer,", norm(B)=%G",eps->nrmb);
278: }
279: PetscViewerASCIIPrintf(viewer,"\n");
280: } else {
281: if (eps->ops->view) {
282: (*eps->ops->view)(eps,viewer);
283: }
284: }
285: PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLZPACK,EPSTRLAN,EPSBLOPEX,EPSPRIMME,"");
286: if (!isexternal) {
287: if (!eps->ip) { EPSGetIP(eps,&eps->ip); }
288: IPView(eps->ip,viewer);
289: PetscObjectTypeCompare((PetscObject)eps,EPSPOWER,&ispower);
290: if (!ispower) {
291: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
292: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
293: DSView(eps->ds,viewer);
294: PetscViewerPopFormat(viewer);
295: }
296: }
297: if (!eps->st) { EPSGetST(eps,&eps->st); }
298: STView(eps->st,viewer);
299: return(0);
300: }
304: /*@
305: EPSPrintSolution - Prints the computed eigenvalues.
307: Collective on EPS
309: Input Parameters:
310: + eps - the eigensolver context
311: - viewer - optional visualization context
313: Options Database Key:
314: . -eps_terse - print only minimal information
316: Note:
317: By default, this function prints a table with eigenvalues and associated
318: relative errors. With -eps_terse only the eigenvalues are printed.
320: Level: intermediate
322: .seealso: PetscViewerASCIIOpen()
323: @*/
324: PetscErrorCode EPSPrintSolution(EPS eps,PetscViewer viewer)
325: {
326: PetscBool terse,errok,isascii;
327: PetscReal error,re,im;
328: PetscScalar kr,ki;
329: PetscInt i,j;
334: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
337: if (!eps->eigr || !eps->eigi || !eps->V) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONGSTATE,"EPSSolve must be called first");
338: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
339: if (!isascii) return(0);
341: PetscOptionsHasName(NULL,"-eps_terse",&terse);
342: if (terse) {
343: if (eps->nconv<eps->nev) {
344: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
345: } else {
346: errok = PETSC_TRUE;
347: for (i=0;i<eps->nev;i++) {
348: EPSComputeRelativeError(eps,i,&error);
349: errok = (errok && error<eps->tol)? PETSC_TRUE: PETSC_FALSE;
350: }
351: if (errok) {
352: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
353: for (i=0;i<=(eps->nev-1)/8;i++) {
354: PetscViewerASCIIPrintf(viewer,"\n ");
355: for (j=0;j<PetscMin(8,eps->nev-8*i);j++) {
356: EPSGetEigenpair(eps,8*i+j,&kr,&ki,NULL,NULL);
357: #if defined(PETSC_USE_COMPLEX)
358: re = PetscRealPart(kr);
359: im = PetscImaginaryPart(kr);
360: #else
361: re = kr;
362: im = ki;
363: #endif
364: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
365: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
366: if (im!=0.0) {
367: PetscViewerASCIIPrintf(viewer,"%.5F%+.5Fi",re,im);
368: } else {
369: PetscViewerASCIIPrintf(viewer,"%.5F",re);
370: }
371: if (8*i+j+1<eps->nev) { PetscViewerASCIIPrintf(viewer,", "); }
372: }
373: }
374: PetscViewerASCIIPrintf(viewer,"\n\n");
375: } else {
376: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",eps->nev);
377: }
378: }
379: } else {
380: PetscViewerASCIIPrintf(viewer," Number of converged approximate eigenpairs: %D\n\n",eps->nconv);
381: if (eps->nconv>0) {
382: PetscViewerASCIIPrintf(viewer,
383: " k ||Ax-k%sx||/||kx||\n"
384: " ----------------- ------------------\n",eps->isgeneralized?"B":"");
385: for (i=0;i<eps->nconv;i++) {
386: EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
387: EPSComputeRelativeError(eps,i,&error);
388: #if defined(PETSC_USE_COMPLEX)
389: re = PetscRealPart(kr);
390: im = PetscImaginaryPart(kr);
391: #else
392: re = kr;
393: im = ki;
394: #endif
395: if (im!=0.0) {
396: PetscViewerASCIIPrintf(viewer," % 9F%+9F i %12G\n",re,im,error);
397: } else {
398: PetscViewerASCIIPrintf(viewer," % 12F %12G\n",re,error);
399: }
400: }
401: PetscViewerASCIIPrintf(viewer,"\n");
402: }
403: }
404: return(0);
405: }
409: /*@C
410: EPSCreate - Creates the default EPS context.
412: Collective on MPI_Comm
414: Input Parameter:
415: . comm - MPI communicator
417: Output Parameter:
418: . eps - location to put the EPS context
420: Note:
421: The default EPS type is EPSKRYLOVSCHUR
423: Level: beginner
425: .seealso: EPSSetUp(), EPSSolve(), EPSDestroy(), EPS
426: @*/
427: PetscErrorCode EPSCreate(MPI_Comm comm,EPS *outeps)
428: {
430: EPS eps;
434: *outeps = 0;
435: #if !defined(PETSC_USE_DYNAMIC_LIBRARIES)
436: EPSInitializePackage();
437: #endif
439: SlepcHeaderCreate(eps,_p_EPS,struct _EPSOps,EPS_CLASSID,"EPS","Eigenvalue Problem Solver","EPS",comm,EPSDestroy,EPSView);
441: eps->max_it = 0;
442: eps->nev = 1;
443: eps->ncv = 0;
444: eps->mpd = 0;
445: eps->allocated_ncv = 0;
446: eps->nini = 0;
447: eps->ninil = 0;
448: eps->nds = 0;
449: eps->tol = PETSC_DEFAULT;
450: eps->conv = EPS_CONV_EIG;
451: eps->converged = EPSConvergedEigRelative;
452: eps->convergedctx = NULL;
453: eps->which = (EPSWhich)0;
454: eps->comparison = NULL;
455: eps->comparisonctx = NULL;
456: eps->arbitrary = NULL;
457: eps->arbitraryctx = NULL;
458: eps->leftvecs = PETSC_FALSE;
459: eps->trueres = PETSC_FALSE;
460: eps->trackall = PETSC_FALSE;
461: eps->target = 0.0;
462: eps->inta = 0.0;
463: eps->intb = 0.0;
464: eps->evecsavailable = PETSC_FALSE;
465: eps->problem_type = (EPSProblemType)0;
466: eps->extraction = (EPSExtraction)0;
467: eps->balance = (EPSBalance)0;
468: eps->balance_its = 5;
469: eps->balance_cutoff = 1e-8;
470: eps->nrma = 1.0;
471: eps->nrmb = 1.0;
472: eps->adaptive = PETSC_FALSE;
474: eps->V = 0;
475: eps->W = 0;
476: eps->D = 0;
477: eps->defl = 0;
478: eps->IS = 0;
479: eps->ISL = 0;
480: eps->t = 0;
481: eps->ds_ortho = PETSC_FALSE;
482: eps->eigr = 0;
483: eps->eigi = 0;
484: eps->errest = 0;
485: eps->errest_left = 0;
486: eps->st = 0;
487: eps->ip = 0;
488: eps->ds = 0;
489: eps->rand = 0;
490: eps->data = 0;
491: eps->nconv = 0;
492: eps->its = 0;
493: eps->perm = NULL;
495: eps->nwork = 0;
496: eps->work = 0;
497: eps->isgeneralized = PETSC_FALSE;
498: eps->ishermitian = PETSC_FALSE;
499: eps->ispositive = PETSC_FALSE;
500: eps->setupcalled = 0;
501: eps->reason = EPS_CONVERGED_ITERATING;
502: eps->numbermonitors = 0;
504: PetscRandomCreate(comm,&eps->rand);
505: PetscRandomSetSeed(eps->rand,0x12345678);
506: PetscLogObjectParent(eps,eps->rand);
507: *outeps = eps;
508: return(0);
509: }
513: /*@C
514: EPSSetType - Selects the particular solver to be used in the EPS object.
516: Logically Collective on EPS
518: Input Parameters:
519: + eps - the eigensolver context
520: - type - a known method
522: Options Database Key:
523: . -eps_type <method> - Sets the method; use -help for a list
524: of available methods
526: Notes:
527: See "slepc/include/slepceps.h" for available methods. The default
528: is EPSKRYLOVSCHUR.
530: Normally, it is best to use the EPSSetFromOptions() command and
531: then set the EPS type from the options database rather than by using
532: this routine. Using the options database provides the user with
533: maximum flexibility in evaluating the different available methods.
534: The EPSSetType() routine is provided for those situations where it
535: is necessary to set the iterative solver independently of the command
536: line or options database.
538: Level: intermediate
540: .seealso: STSetType(), EPSType
541: @*/
542: PetscErrorCode EPSSetType(EPS eps,EPSType type)
543: {
544: PetscErrorCode ierr,(*r)(EPS);
545: PetscBool match;
551: PetscObjectTypeCompare((PetscObject)eps,type,&match);
552: if (match) return(0);
554: PetscFunctionListFind(EPSList,type,&r);
555: if (!r) SETERRQ1(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown EPS type given: %s",type);
557: if (eps->ops->destroy) { (*eps->ops->destroy)(eps); }
558: PetscMemzero(eps->ops,sizeof(struct _EPSOps));
560: eps->setupcalled = 0;
561: PetscObjectChangeTypeName((PetscObject)eps,type);
562: (*r)(eps);
563: return(0);
564: }
568: /*@C
569: EPSGetType - Gets the EPS type as a string from the EPS object.
571: Not Collective
573: Input Parameter:
574: . eps - the eigensolver context
576: Output Parameter:
577: . name - name of EPS method
579: Level: intermediate
581: .seealso: EPSSetType()
582: @*/
583: PetscErrorCode EPSGetType(EPS eps,EPSType *type)
584: {
588: *type = ((PetscObject)eps)->type_name;
589: return(0);
590: }
594: /*@C
595: EPSRegister - Adds a method to the eigenproblem solver package.
597: Not Collective
599: Input Parameters:
600: + name - name of a new user-defined solver
601: - function - routine to create the solver context
603: Notes:
604: EPSRegister() may be called multiple times to add several user-defined solvers.
606: Sample usage:
607: .vb
608: EPSRegister("my_solver",MySolverCreate);
609: .ve
611: Then, your solver can be chosen with the procedural interface via
612: $ EPSSetType(eps,"my_solver")
613: or at runtime via the option
614: $ -eps_type my_solver
616: Level: advanced
618: .seealso: EPSRegisterAll()
619: @*/
620: PetscErrorCode EPSRegister(const char *name,PetscErrorCode (*function)(EPS))
621: {
625: PetscFunctionListAdd(&EPSList,name,function);
626: return(0);
627: }
631: /*@
632: EPSReset - Resets the EPS context to the setupcalled=0 state and removes any
633: allocated objects.
635: Collective on EPS
637: Input Parameter:
638: . eps - eigensolver context obtained from EPSCreate()
640: Level: advanced
642: .seealso: EPSDestroy()
643: @*/
644: PetscErrorCode EPSReset(EPS eps)
645: {
650: if (eps->ops->reset) { (eps->ops->reset)(eps); }
651: if (eps->st) { STReset(eps->st); }
652: if (eps->ip) { IPReset(eps->ip); }
653: if (eps->ds) { DSReset(eps->ds); }
654: VecDestroy(&eps->t);
655: VecDestroy(&eps->D);
656: eps->setupcalled = 0;
657: return(0);
658: }
662: /*@C
663: EPSDestroy - Destroys the EPS context.
665: Collective on EPS
667: Input Parameter:
668: . eps - eigensolver context obtained from EPSCreate()
670: Level: beginner
672: .seealso: EPSCreate(), EPSSetUp(), EPSSolve()
673: @*/
674: PetscErrorCode EPSDestroy(EPS *eps)
675: {
679: if (!*eps) return(0);
681: if (--((PetscObject)(*eps))->refct > 0) { *eps = 0; return(0); }
682: EPSReset(*eps);
683: if ((*eps)->ops->destroy) { (*(*eps)->ops->destroy)(*eps); }
684: STDestroy(&(*eps)->st);
685: IPDestroy(&(*eps)->ip);
686: DSDestroy(&(*eps)->ds);
687: PetscRandomDestroy(&(*eps)->rand);
688: /* just in case the initial vectors have not been used */
689: SlepcBasisDestroy_Private(&(*eps)->nini,&(*eps)->IS);
690: SlepcBasisDestroy_Private(&(*eps)->ninil,&(*eps)->ISL);
691: EPSRemoveDeflationSpace(*eps);
692: EPSMonitorCancel(*eps);
693: PetscHeaderDestroy(eps);
694: return(0);
695: }
699: /*@
700: EPSSetTarget - Sets the value of the target.
702: Logically Collective on EPS
704: Input Parameters:
705: + eps - eigensolver context
706: - target - the value of the target
708: Notes:
709: The target is a scalar value used to determine the portion of the spectrum
710: of interest. It is used in combination with EPSSetWhichEigenpairs().
712: Level: beginner
714: .seealso: EPSGetTarget(), EPSSetWhichEigenpairs()
715: @*/
716: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
717: {
723: eps->target = target;
724: if (!eps->st) { EPSGetST(eps,&eps->st); }
725: STSetDefaultShift(eps->st,target);
726: return(0);
727: }
731: /*@
732: EPSGetTarget - Gets the value of the target.
734: Not Collective
736: Input Parameter:
737: . eps - eigensolver context
739: Output Parameter:
740: . target - the value of the target
742: Level: beginner
744: Note:
745: If the target was not set by the user, then zero is returned.
747: .seealso: EPSSetTarget()
748: @*/
749: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
750: {
754: *target = eps->target;
755: return(0);
756: }
760: /*@
761: EPSSetInterval - Defines the computational interval for spectrum slicing.
763: Logically Collective on EPS
765: Input Parameters:
766: + eps - eigensolver context
767: . inta - left end of the interval
768: - intb - right end of the interval
770: Options Database Key:
771: . -eps_interval <a,b> - set [a,b] as the interval of interest
773: Notes:
774: Spectrum slicing is a technique employed for computing all eigenvalues of
775: symmetric eigenproblems in a given interval. This function provides the
776: interval to be considered. It must be used in combination with EPS_ALL, see
777: EPSSetWhichEigenpairs().
779: In the command-line option, two values must be provided. For an open interval,
780: one can give an infinite, e.g., -eps_interval 1.0,inf or -eps_interval -inf,1.0.
781: An open interval in the programmatic interface can be specified with
782: PETSC_MAX_REAL and -PETSC_MAX_REAL.
784: Level: intermediate
786: .seealso: EPSGetInterval(), EPSSetWhichEigenpairs()
787: @*/
788: PetscErrorCode EPSSetInterval(EPS eps,PetscReal inta,PetscReal intb)
789: {
794: if (inta>=intb) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
795: eps->inta = inta;
796: eps->intb = intb;
797: return(0);
798: }
802: /*@
803: EPSGetInterval - Gets the computational interval for spectrum slicing.
805: Not Collective
807: Input Parameter:
808: . eps - eigensolver context
810: Output Parameters:
811: + inta - left end of the interval
812: - intb - right end of the interval
814: Level: intermediate
816: Note:
817: If the interval was not set by the user, then zeros are returned.
819: .seealso: EPSSetInterval()
820: @*/
821: PetscErrorCode EPSGetInterval(EPS eps,PetscReal* inta,PetscReal* intb)
822: {
827: if (inta) *inta = eps->inta;
828: if (intb) *intb = eps->intb;
829: return(0);
830: }
834: /*@
835: EPSSetST - Associates a spectral transformation object to the eigensolver.
837: Collective on EPS
839: Input Parameters:
840: + eps - eigensolver context obtained from EPSCreate()
841: - st - the spectral transformation object
843: Note:
844: Use EPSGetST() to retrieve the spectral transformation context (for example,
845: to free it at the end of the computations).
847: Level: developer
849: .seealso: EPSGetST()
850: @*/
851: PetscErrorCode EPSSetST(EPS eps,ST st)
852: {
859: PetscObjectReference((PetscObject)st);
860: STDestroy(&eps->st);
861: eps->st = st;
862: PetscLogObjectParent(eps,eps->st);
863: return(0);
864: }
868: /*@C
869: EPSGetST - Obtain the spectral transformation (ST) object associated
870: to the eigensolver object.
872: Not Collective
874: Input Parameters:
875: . eps - eigensolver context obtained from EPSCreate()
877: Output Parameter:
878: . st - spectral transformation context
880: Level: beginner
882: .seealso: EPSSetST()
883: @*/
884: PetscErrorCode EPSGetST(EPS eps,ST *st)
885: {
891: if (!eps->st) {
892: STCreate(PetscObjectComm((PetscObject)eps),&eps->st);
893: PetscLogObjectParent(eps,eps->st);
894: }
895: *st = eps->st;
896: return(0);
897: }
901: /*@
902: EPSSetIP - Associates an inner product object to the eigensolver.
904: Collective on EPS
906: Input Parameters:
907: + eps - eigensolver context obtained from EPSCreate()
908: - ip - the inner product object
910: Note:
911: Use EPSGetIP() to retrieve the inner product context (for example,
912: to free it at the end of the computations).
914: Level: advanced
916: .seealso: EPSGetIP()
917: @*/
918: PetscErrorCode EPSSetIP(EPS eps,IP ip)
919: {
926: PetscObjectReference((PetscObject)ip);
927: IPDestroy(&eps->ip);
928: eps->ip = ip;
929: PetscLogObjectParent(eps,eps->ip);
930: return(0);
931: }
935: /*@C
936: EPSGetIP - Obtain the inner product object associated to the eigensolver object.
938: Not Collective
940: Input Parameters:
941: . eps - eigensolver context obtained from EPSCreate()
943: Output Parameter:
944: . ip - inner product context
946: Level: advanced
948: .seealso: EPSSetIP()
949: @*/
950: PetscErrorCode EPSGetIP(EPS eps,IP *ip)
951: {
957: if (!eps->ip) {
958: IPCreate(PetscObjectComm((PetscObject)eps),&eps->ip);
959: PetscLogObjectParent(eps,eps->ip);
960: }
961: *ip = eps->ip;
962: return(0);
963: }
967: /*@
968: EPSSetDS - Associates a direct solver object to the eigensolver.
970: Collective on EPS
972: Input Parameters:
973: + eps - eigensolver context obtained from EPSCreate()
974: - ds - the direct solver object
976: Note:
977: Use EPSGetDS() to retrieve the direct solver context (for example,
978: to free it at the end of the computations).
980: Level: advanced
982: .seealso: EPSGetDS()
983: @*/
984: PetscErrorCode EPSSetDS(EPS eps,DS ds)
985: {
992: PetscObjectReference((PetscObject)ds);
993: DSDestroy(&eps->ds);
994: eps->ds = ds;
995: PetscLogObjectParent(eps,eps->ds);
996: return(0);
997: }
1001: /*@C
1002: EPSGetDS - Obtain the direct solver object associated to the eigensolver object.
1004: Not Collective
1006: Input Parameters:
1007: . eps - eigensolver context obtained from EPSCreate()
1009: Output Parameter:
1010: . ds - direct solver context
1012: Level: advanced
1014: .seealso: EPSSetDS()
1015: @*/
1016: PetscErrorCode EPSGetDS(EPS eps,DS *ds)
1017: {
1023: if (!eps->ds) {
1024: DSCreate(PetscObjectComm((PetscObject)eps),&eps->ds);
1025: PetscLogObjectParent(eps,eps->ds);
1026: }
1027: *ds = eps->ds;
1028: return(0);
1029: }
1033: /*@
1034: EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
1035: eigenvalue problem.
1037: Not collective
1039: Input Parameter:
1040: . eps - the eigenproblem solver context
1042: Output Parameter:
1043: . is - the answer
1045: Level: intermediate
1047: .seealso: EPSIsHermitian(), EPSIsPositive()
1048: @*/
1049: PetscErrorCode EPSIsGeneralized(EPS eps,PetscBool* is)
1050: {
1054: *is = eps->isgeneralized;
1055: return(0);
1056: }
1060: /*@
1061: EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
1062: eigenvalue problem.
1064: Not collective
1066: Input Parameter:
1067: . eps - the eigenproblem solver context
1069: Output Parameter:
1070: . is - the answer
1072: Level: intermediate
1074: .seealso: EPSIsGeneralized(), EPSIsPositive()
1075: @*/
1076: PetscErrorCode EPSIsHermitian(EPS eps,PetscBool* is)
1077: {
1081: *is = eps->ishermitian;
1082: return(0);
1083: }
1087: /*@
1088: EPSIsPositive - Ask if the EPS object corresponds to an eigenvalue
1089: problem type that requires a positive (semi-) definite matrix B.
1091: Not collective
1093: Input Parameter:
1094: . eps - the eigenproblem solver context
1096: Output Parameter:
1097: . is - the answer
1099: Level: intermediate
1101: .seealso: EPSIsGeneralized(), EPSIsHermitian()
1102: @*/
1103: PetscErrorCode EPSIsPositive(EPS eps,PetscBool* is)
1104: {
1108: *is = eps->ispositive;
1109: return(0);
1110: }