Actual source code: epsbasic.c
slepc-3.5.2 2014-10-10
1: /*
2: The basic EPS routines, Create, View, etc. are here.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2014, 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;
33: /*@C
34: EPSView - Prints the EPS data structure.
36: Collective on EPS
38: Input Parameters:
39: + eps - the eigenproblem solver context
40: - viewer - optional visualization context
42: Options Database Key:
43: . -eps_view - Calls EPSView() at end of EPSSolve()
45: Note:
46: The available visualization contexts include
47: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
48: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
49: output where only the first processor opens
50: the file. All other processors send their
51: data to the first processor to print.
53: The user can open an alternative visualization context with
54: PetscViewerASCIIOpen() - output to a specified file.
56: Level: beginner
58: .seealso: STView(), PetscViewerASCIIOpen()
59: @*/
60: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
61: {
63: const char *type,*extr,*bal;
64: char str[50];
65: PetscBool isascii,ispower,isexternal,istrivial;
69: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
73: #if defined(PETSC_USE_COMPLEX)
74: #define HERM "hermitian"
75: #else
76: #define HERM "symmetric"
77: #endif
78: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
79: if (isascii) {
80: PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer);
81: if (eps->ops->view) {
82: PetscViewerASCIIPushTab(viewer);
83: (*eps->ops->view)(eps,viewer);
84: PetscViewerASCIIPopTab(viewer);
85: }
86: if (eps->problem_type) {
87: switch (eps->problem_type) {
88: case EPS_HEP: type = HERM " eigenvalue problem"; break;
89: case EPS_GHEP: type = "generalized " HERM " eigenvalue problem"; break;
90: case EPS_NHEP: type = "non-" HERM " eigenvalue problem"; break;
91: case EPS_GNHEP: type = "generalized non-" HERM " eigenvalue problem"; break;
92: case EPS_PGNHEP: type = "generalized non-" HERM " eigenvalue problem with " HERM " positive definite B"; break;
93: case EPS_GHIEP: type = "generalized " HERM "-indefinite eigenvalue problem"; break;
94: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->problem_type");
95: }
96: } else type = "not yet set";
97: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
98: if (eps->extraction) {
99: switch (eps->extraction) {
100: case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
101: case EPS_HARMONIC: extr = "harmonic Ritz"; break;
102: case EPS_HARMONIC_RELATIVE:extr = "relative harmonic Ritz"; break;
103: case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
104: case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
105: case EPS_REFINED: extr = "refined Ritz"; break;
106: case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
107: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->extraction");
108: }
109: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr);
110: }
111: if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
112: switch (eps->balance) {
113: case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
114: case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
115: case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
116: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->balance");
117: }
118: PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal);
119: if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) {
120: PetscViewerASCIIPrintf(viewer,", with its=%D",eps->balance_its);
121: }
122: if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) {
123: PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff);
124: }
125: PetscViewerASCIIPrintf(viewer,"\n");
126: }
127: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
128: SlepcSNPrintfScalar(str,50,eps->target,PETSC_FALSE);
129: if (!eps->which) {
130: PetscViewerASCIIPrintf(viewer,"not yet set\n");
131: } else switch (eps->which) {
132: case EPS_WHICH_USER:
133: PetscViewerASCIIPrintf(viewer,"user defined\n");
134: break;
135: case EPS_TARGET_MAGNITUDE:
136: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
137: break;
138: case EPS_TARGET_REAL:
139: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
140: break;
141: case EPS_TARGET_IMAGINARY:
142: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
143: break;
144: case EPS_LARGEST_MAGNITUDE:
145: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
146: break;
147: case EPS_SMALLEST_MAGNITUDE:
148: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
149: break;
150: case EPS_LARGEST_REAL:
151: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
152: break;
153: case EPS_SMALLEST_REAL:
154: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
155: break;
156: case EPS_LARGEST_IMAGINARY:
157: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
158: break;
159: case EPS_SMALLEST_IMAGINARY:
160: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
161: break;
162: case EPS_ALL:
163: PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb);
164: break;
165: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->which");
166: }
167: if (eps->trueres) {
168: PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n");
169: }
170: if (eps->trackall) {
171: PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n");
172: }
173: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",eps->nev);
174: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",eps->ncv);
175: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",eps->mpd);
176: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",eps->max_it);
177: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)eps->tol);
178: PetscViewerASCIIPrintf(viewer," convergence test: ");
179: switch (eps->conv) {
180: case EPS_CONV_ABS:
181: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
182: case EPS_CONV_EIG:
183: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
184: case EPS_CONV_NORM:
185: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");
186: PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)eps->nrma);
187: if (eps->isgeneralized) {
188: PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb);
189: }
190: PetscViewerASCIIPrintf(viewer,"\n");
191: break;
192: case EPS_CONV_USER:
193: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
194: }
195: if (eps->nini) {
196: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(eps->nini));
197: }
198: if (eps->nds) {
199: PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %D\n",PetscAbs(eps->nds));
200: }
201: } else {
202: if (eps->ops->view) {
203: (*eps->ops->view)(eps,viewer);
204: }
205: }
206: PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLZPACK,EPSTRLAN,EPSBLOPEX,EPSPRIMME,"");
207: if (!isexternal) {
208: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
209: if (!eps->V) { EPSGetBV(eps,&eps->V); }
210: BVView(eps->V,viewer);
211: if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
212: RGIsTrivial(eps->rg,&istrivial);
213: if (!istrivial) { RGView(eps->rg,viewer); }
214: PetscObjectTypeCompare((PetscObject)eps,EPSPOWER,&ispower);
215: if (!ispower) {
216: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
217: DSView(eps->ds,viewer);
218: }
219: PetscViewerPopFormat(viewer);
220: }
221: if (!eps->st) { EPSGetST(eps,&eps->st); }
222: STView(eps->st,viewer);
223: return(0);
224: }
228: /*@
229: EPSPrintSolution - Prints the computed eigenvalues.
231: Collective on EPS
233: Input Parameters:
234: + eps - the eigensolver context
235: - viewer - optional visualization context
237: Options Database Key:
238: . -eps_terse - print only minimal information
240: Note:
241: By default, this function prints a table with eigenvalues and associated
242: relative errors. With -eps_terse only the eigenvalues are printed.
244: Level: intermediate
246: .seealso: PetscViewerASCIIOpen()
247: @*/
248: PetscErrorCode EPSPrintSolution(EPS eps,PetscViewer viewer)
249: {
250: PetscBool terse,errok,isascii;
251: PetscReal error,re,im;
252: PetscScalar kr,ki;
253: PetscInt i,j;
258: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
261: EPSCheckSolved(eps,1);
262: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
263: if (!isascii) return(0);
265: PetscOptionsHasName(NULL,"-eps_terse",&terse);
266: if (terse) {
267: if (eps->nconv<eps->nev) {
268: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
269: } else {
270: errok = PETSC_TRUE;
271: for (i=0;i<eps->nev;i++) {
272: EPSComputeRelativeError(eps,i,&error);
273: errok = (errok && error<5.0*eps->tol)? PETSC_TRUE: PETSC_FALSE;
274: }
275: if (errok) {
276: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
277: for (i=0;i<=(eps->nev-1)/8;i++) {
278: PetscViewerASCIIPrintf(viewer,"\n ");
279: for (j=0;j<PetscMin(8,eps->nev-8*i);j++) {
280: EPSGetEigenpair(eps,8*i+j,&kr,&ki,NULL,NULL);
281: #if defined(PETSC_USE_COMPLEX)
282: re = PetscRealPart(kr);
283: im = PetscImaginaryPart(kr);
284: #else
285: re = kr;
286: im = ki;
287: #endif
288: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
289: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
290: if (im!=0.0) {
291: PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
292: } else {
293: PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
294: }
295: if (8*i+j+1<eps->nev) { PetscViewerASCIIPrintf(viewer,", "); }
296: }
297: }
298: PetscViewerASCIIPrintf(viewer,"\n\n");
299: } else {
300: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",eps->nev);
301: }
302: }
303: } else {
304: PetscViewerASCIIPrintf(viewer," Number of converged approximate eigenpairs: %D\n\n",eps->nconv);
305: if (eps->nconv>0) {
306: PetscViewerASCIIPrintf(viewer,
307: " k ||Ax-k%sx||/||kx||\n"
308: " ----------------- ------------------\n",eps->isgeneralized?"B":"");
309: for (i=0;i<eps->nconv;i++) {
310: EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
311: EPSComputeRelativeError(eps,i,&error);
312: #if defined(PETSC_USE_COMPLEX)
313: re = PetscRealPart(kr);
314: im = PetscImaginaryPart(kr);
315: #else
316: re = kr;
317: im = ki;
318: #endif
319: if (im!=0.0) {
320: PetscViewerASCIIPrintf(viewer," % 9f%+9f i %12g\n",(double)re,(double)im,(double)error);
321: } else {
322: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
323: }
324: }
325: PetscViewerASCIIPrintf(viewer,"\n");
326: }
327: }
328: return(0);
329: }
333: /*@C
334: EPSCreate - Creates the default EPS context.
336: Collective on MPI_Comm
338: Input Parameter:
339: . comm - MPI communicator
341: Output Parameter:
342: . eps - location to put the EPS context
344: Note:
345: The default EPS type is EPSKRYLOVSCHUR
347: Level: beginner
349: .seealso: EPSSetUp(), EPSSolve(), EPSDestroy(), EPS
350: @*/
351: PetscErrorCode EPSCreate(MPI_Comm comm,EPS *outeps)
352: {
354: EPS eps;
358: *outeps = 0;
359: EPSInitializePackage();
360: SlepcHeaderCreate(eps,_p_EPS,struct _EPSOps,EPS_CLASSID,"EPS","Eigenvalue Problem Solver","EPS",comm,EPSDestroy,EPSView);
362: eps->max_it = 0;
363: eps->nev = 1;
364: eps->ncv = 0;
365: eps->mpd = 0;
366: eps->nini = 0;
367: eps->nds = 0;
368: eps->target = 0.0;
369: eps->tol = PETSC_DEFAULT;
370: eps->conv = EPS_CONV_EIG;
371: eps->which = (EPSWhich)0;
372: eps->inta = 0.0;
373: eps->intb = 0.0;
374: eps->problem_type = (EPSProblemType)0;
375: eps->extraction = EPS_RITZ;
376: eps->balance = EPS_BALANCE_NONE;
377: eps->balance_its = 5;
378: eps->balance_cutoff = 1e-8;
379: eps->trueres = PETSC_FALSE;
380: eps->trackall = PETSC_FALSE;
382: eps->converged = EPSConvergedEigRelative;
383: eps->convergeddestroy= NULL;
384: eps->arbitrary = NULL;
385: eps->convergedctx = NULL;
386: eps->arbitraryctx = NULL;
387: eps->numbermonitors = 0;
389: eps->st = NULL;
390: eps->ds = NULL;
391: eps->V = NULL;
392: eps->rg = NULL;
393: eps->rand = NULL;
394: eps->D = NULL;
395: eps->IS = NULL;
396: eps->defl = NULL;
397: eps->eigr = NULL;
398: eps->eigi = NULL;
399: eps->errest = NULL;
400: eps->rr = NULL;
401: eps->ri = NULL;
402: eps->perm = NULL;
403: eps->nwork = 0;
404: eps->work = NULL;
405: eps->data = NULL;
407: eps->state = EPS_STATE_INITIAL;
408: eps->nconv = 0;
409: eps->its = 0;
410: eps->nloc = 0;
411: eps->nrma = 0.0;
412: eps->nrmb = 0.0;
413: eps->isgeneralized = PETSC_FALSE;
414: eps->ispositive = PETSC_FALSE;
415: eps->ishermitian = PETSC_FALSE;
416: eps->reason = EPS_CONVERGED_ITERATING;
418: PetscNewLog(eps,&eps->sc);
419: PetscRandomCreate(comm,&eps->rand);
420: PetscRandomSetSeed(eps->rand,0x12345678);
421: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rand);
422: *outeps = eps;
423: return(0);
424: }
428: /*@C
429: EPSSetType - Selects the particular solver to be used in the EPS object.
431: Logically Collective on EPS
433: Input Parameters:
434: + eps - the eigensolver context
435: - type - a known method
437: Options Database Key:
438: . -eps_type <method> - Sets the method; use -help for a list
439: of available methods
441: Notes:
442: See "slepc/include/slepceps.h" for available methods. The default
443: is EPSKRYLOVSCHUR.
445: Normally, it is best to use the EPSSetFromOptions() command and
446: then set the EPS type from the options database rather than by using
447: this routine. Using the options database provides the user with
448: maximum flexibility in evaluating the different available methods.
449: The EPSSetType() routine is provided for those situations where it
450: is necessary to set the iterative solver independently of the command
451: line or options database.
453: Level: intermediate
455: .seealso: STSetType(), EPSType
456: @*/
457: PetscErrorCode EPSSetType(EPS eps,EPSType type)
458: {
459: PetscErrorCode ierr,(*r)(EPS);
460: PetscBool match;
466: PetscObjectTypeCompare((PetscObject)eps,type,&match);
467: if (match) return(0);
469: PetscFunctionListFind(EPSList,type,&r);
470: if (!r) SETERRQ1(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown EPS type given: %s",type);
472: if (eps->ops->destroy) { (*eps->ops->destroy)(eps); }
473: PetscMemzero(eps->ops,sizeof(struct _EPSOps));
475: eps->state = EPS_STATE_INITIAL;
476: PetscObjectChangeTypeName((PetscObject)eps,type);
477: (*r)(eps);
478: return(0);
479: }
483: /*@C
484: EPSGetType - Gets the EPS type as a string from the EPS object.
486: Not Collective
488: Input Parameter:
489: . eps - the eigensolver context
491: Output Parameter:
492: . name - name of EPS method
494: Level: intermediate
496: .seealso: EPSSetType()
497: @*/
498: PetscErrorCode EPSGetType(EPS eps,EPSType *type)
499: {
503: *type = ((PetscObject)eps)->type_name;
504: return(0);
505: }
509: /*@C
510: EPSRegister - Adds a method to the eigenproblem solver package.
512: Not Collective
514: Input Parameters:
515: + name - name of a new user-defined solver
516: - function - routine to create the solver context
518: Notes:
519: EPSRegister() may be called multiple times to add several user-defined solvers.
521: Sample usage:
522: .vb
523: EPSRegister("my_solver",MySolverCreate);
524: .ve
526: Then, your solver can be chosen with the procedural interface via
527: $ EPSSetType(eps,"my_solver")
528: or at runtime via the option
529: $ -eps_type my_solver
531: Level: advanced
533: .seealso: EPSRegisterAll()
534: @*/
535: PetscErrorCode EPSRegister(const char *name,PetscErrorCode (*function)(EPS))
536: {
540: PetscFunctionListAdd(&EPSList,name,function);
541: return(0);
542: }
546: /*@
547: EPSReset - Resets the EPS context to the initial state and removes any
548: allocated objects.
550: Collective on EPS
552: Input Parameter:
553: . eps - eigensolver context obtained from EPSCreate()
555: Level: advanced
557: .seealso: EPSDestroy()
558: @*/
559: PetscErrorCode EPSReset(EPS eps)
560: {
562: PetscInt ncols;
566: if (eps->ops->reset) { (eps->ops->reset)(eps); }
567: if (eps->st) { STReset(eps->st); }
568: if (eps->ds) { DSReset(eps->ds); }
569: VecDestroy(&eps->D);
570: BVGetSizes(eps->V,NULL,NULL,&ncols);
571: if (ncols) {
572: PetscFree4(eps->eigr,eps->eigi,eps->errest,eps->perm);
573: PetscFree2(eps->rr,eps->ri);
574: }
575: BVDestroy(&eps->V);
576: VecDestroyVecs(eps->nwork,&eps->work);
577: eps->nwork = 0;
578: eps->state = EPS_STATE_INITIAL;
579: return(0);
580: }
584: /*@C
585: EPSDestroy - Destroys the EPS context.
587: Collective on EPS
589: Input Parameter:
590: . eps - eigensolver context obtained from EPSCreate()
592: Level: beginner
594: .seealso: EPSCreate(), EPSSetUp(), EPSSolve()
595: @*/
596: PetscErrorCode EPSDestroy(EPS *eps)
597: {
601: if (!*eps) return(0);
603: if (--((PetscObject)(*eps))->refct > 0) { *eps = 0; return(0); }
604: EPSReset(*eps);
605: if ((*eps)->ops->destroy) { (*(*eps)->ops->destroy)(*eps); }
606: STDestroy(&(*eps)->st);
607: RGDestroy(&(*eps)->rg);
608: DSDestroy(&(*eps)->ds);
609: PetscRandomDestroy(&(*eps)->rand);
610: PetscFree((*eps)->sc);
611: /* just in case the initial vectors have not been used */
612: SlepcBasisDestroy_Private(&(*eps)->nds,&(*eps)->defl);
613: SlepcBasisDestroy_Private(&(*eps)->nini,&(*eps)->IS);
614: if ((*eps)->convergeddestroy) {
615: (*(*eps)->convergeddestroy)((*eps)->convergedctx);
616: }
617: EPSMonitorCancel(*eps);
618: PetscHeaderDestroy(eps);
619: return(0);
620: }
624: /*@
625: EPSSetTarget - Sets the value of the target.
627: Logically Collective on EPS
629: Input Parameters:
630: + eps - eigensolver context
631: - target - the value of the target
633: Options Database Key:
634: . -eps_target <scalar> - the value of the target
636: Notes:
637: The target is a scalar value used to determine the portion of the spectrum
638: of interest. It is used in combination with EPSSetWhichEigenpairs().
640: In the case of complex scalars, a complex value can be provided in the
641: command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
642: -eps_target 1.0+2.0i
644: Level: beginner
646: .seealso: EPSGetTarget(), EPSSetWhichEigenpairs()
647: @*/
648: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
649: {
655: eps->target = target;
656: if (!eps->st) { EPSGetST(eps,&eps->st); }
657: STSetDefaultShift(eps->st,target);
658: return(0);
659: }
663: /*@
664: EPSGetTarget - Gets the value of the target.
666: Not Collective
668: Input Parameter:
669: . eps - eigensolver context
671: Output Parameter:
672: . target - the value of the target
674: Note:
675: If the target was not set by the user, then zero is returned.
677: Level: beginner
679: .seealso: EPSSetTarget()
680: @*/
681: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
682: {
686: *target = eps->target;
687: return(0);
688: }
692: /*@
693: EPSSetInterval - Defines the computational interval for spectrum slicing.
695: Logically Collective on EPS
697: Input Parameters:
698: + eps - eigensolver context
699: . inta - left end of the interval
700: - intb - right end of the interval
702: Options Database Key:
703: . -eps_interval <a,b> - set [a,b] as the interval of interest
705: Notes:
706: Spectrum slicing is a technique employed for computing all eigenvalues of
707: symmetric eigenproblems in a given interval. This function provides the
708: interval to be considered. It must be used in combination with EPS_ALL, see
709: EPSSetWhichEigenpairs().
711: In the command-line option, two values must be provided. For an open interval,
712: one can give an infinite, e.g., -eps_interval 1.0,inf or -eps_interval -inf,1.0.
713: An open interval in the programmatic interface can be specified with
714: PETSC_MAX_REAL and -PETSC_MAX_REAL.
716: Level: intermediate
718: .seealso: EPSGetInterval(), EPSSetWhichEigenpairs()
719: @*/
720: PetscErrorCode EPSSetInterval(EPS eps,PetscReal inta,PetscReal intb)
721: {
726: if (inta>=intb) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
727: eps->inta = inta;
728: eps->intb = intb;
729: eps->state = EPS_STATE_INITIAL;
730: return(0);
731: }
735: /*@
736: EPSGetInterval - Gets the computational interval for spectrum slicing.
738: Not Collective
740: Input Parameter:
741: . eps - eigensolver context
743: Output Parameters:
744: + inta - left end of the interval
745: - intb - right end of the interval
747: Level: intermediate
749: Note:
750: If the interval was not set by the user, then zeros are returned.
752: .seealso: EPSSetInterval()
753: @*/
754: PetscErrorCode EPSGetInterval(EPS eps,PetscReal* inta,PetscReal* intb)
755: {
760: if (inta) *inta = eps->inta;
761: if (intb) *intb = eps->intb;
762: return(0);
763: }
767: /*@
768: EPSSetST - Associates a spectral transformation object to the eigensolver.
770: Collective on EPS
772: Input Parameters:
773: + eps - eigensolver context obtained from EPSCreate()
774: - st - the spectral transformation object
776: Note:
777: Use EPSGetST() to retrieve the spectral transformation context (for example,
778: to free it at the end of the computations).
780: Level: developer
782: .seealso: EPSGetST()
783: @*/
784: PetscErrorCode EPSSetST(EPS eps,ST st)
785: {
792: PetscObjectReference((PetscObject)st);
793: STDestroy(&eps->st);
794: eps->st = st;
795: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
796: return(0);
797: }
801: /*@C
802: EPSGetST - Obtain the spectral transformation (ST) object associated
803: to the eigensolver object.
805: Not Collective
807: Input Parameters:
808: . eps - eigensolver context obtained from EPSCreate()
810: Output Parameter:
811: . st - spectral transformation context
813: Level: beginner
815: .seealso: EPSSetST()
816: @*/
817: PetscErrorCode EPSGetST(EPS eps,ST *st)
818: {
824: if (!eps->st) {
825: STCreate(PetscObjectComm((PetscObject)eps),&eps->st);
826: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
827: }
828: *st = eps->st;
829: return(0);
830: }
834: /*@
835: EPSSetBV - Associates a basis vectors object to the eigensolver.
837: Collective on EPS
839: Input Parameters:
840: + eps - eigensolver context obtained from EPSCreate()
841: - V - the basis vectors object
843: Note:
844: Use EPSGetBV() to retrieve the basis vectors context (for example,
845: to free them at the end of the computations).
847: Level: advanced
849: .seealso: EPSGetBV()
850: @*/
851: PetscErrorCode EPSSetBV(EPS eps,BV V)
852: {
859: PetscObjectReference((PetscObject)V);
860: BVDestroy(&eps->V);
861: eps->V = V;
862: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
863: return(0);
864: }
868: /*@C
869: EPSGetBV - Obtain the basis vectors object associated to the eigensolver object.
871: Not Collective
873: Input Parameters:
874: . eps - eigensolver context obtained from EPSCreate()
876: Output Parameter:
877: . V - basis vectors context
879: Level: advanced
881: .seealso: EPSSetBV()
882: @*/
883: PetscErrorCode EPSGetBV(EPS eps,BV *V)
884: {
890: if (!eps->V) {
891: BVCreate(PetscObjectComm((PetscObject)eps),&eps->V);
892: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
893: }
894: *V = eps->V;
895: return(0);
896: }
900: /*@
901: EPSSetRG - Associates a region object to the eigensolver.
903: Collective on EPS
905: Input Parameters:
906: + eps - eigensolver context obtained from EPSCreate()
907: - rg - the region object
909: Note:
910: Use EPSGetRG() to retrieve the region context (for example,
911: to free it at the end of the computations).
913: Level: advanced
915: .seealso: EPSGetRG()
916: @*/
917: PetscErrorCode EPSSetRG(EPS eps,RG rg)
918: {
925: PetscObjectReference((PetscObject)rg);
926: RGDestroy(&eps->rg);
927: eps->rg = rg;
928: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
929: return(0);
930: }
934: /*@C
935: EPSGetRG - Obtain the region object associated to the eigensolver.
937: Not Collective
939: Input Parameters:
940: . eps - eigensolver context obtained from EPSCreate()
942: Output Parameter:
943: . rg - region context
945: Level: advanced
947: .seealso: EPSSetRG()
948: @*/
949: PetscErrorCode EPSGetRG(EPS eps,RG *rg)
950: {
956: if (!eps->rg) {
957: RGCreate(PetscObjectComm((PetscObject)eps),&eps->rg);
958: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
959: }
960: *rg = eps->rg;
961: return(0);
962: }
966: /*@
967: EPSSetDS - Associates a direct solver object to the eigensolver.
969: Collective on EPS
971: Input Parameters:
972: + eps - eigensolver context obtained from EPSCreate()
973: - ds - the direct solver object
975: Note:
976: Use EPSGetDS() to retrieve the direct solver context (for example,
977: to free it at the end of the computations).
979: Level: advanced
981: .seealso: EPSGetDS()
982: @*/
983: PetscErrorCode EPSSetDS(EPS eps,DS ds)
984: {
991: PetscObjectReference((PetscObject)ds);
992: DSDestroy(&eps->ds);
993: eps->ds = ds;
994: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
995: return(0);
996: }
1000: /*@C
1001: EPSGetDS - Obtain the direct solver object associated to the eigensolver object.
1003: Not Collective
1005: Input Parameters:
1006: . eps - eigensolver context obtained from EPSCreate()
1008: Output Parameter:
1009: . ds - direct solver context
1011: Level: advanced
1013: .seealso: EPSSetDS()
1014: @*/
1015: PetscErrorCode EPSGetDS(EPS eps,DS *ds)
1016: {
1022: if (!eps->ds) {
1023: DSCreate(PetscObjectComm((PetscObject)eps),&eps->ds);
1024: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
1025: }
1026: *ds = eps->ds;
1027: return(0);
1028: }
1032: /*@
1033: EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
1034: eigenvalue problem.
1036: Not collective
1038: Input Parameter:
1039: . eps - the eigenproblem solver context
1041: Output Parameter:
1042: . is - the answer
1044: Level: intermediate
1046: .seealso: EPSIsHermitian(), EPSIsPositive()
1047: @*/
1048: PetscErrorCode EPSIsGeneralized(EPS eps,PetscBool* is)
1049: {
1053: *is = eps->isgeneralized;
1054: return(0);
1055: }
1059: /*@
1060: EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
1061: eigenvalue problem.
1063: Not collective
1065: Input Parameter:
1066: . eps - the eigenproblem solver context
1068: Output Parameter:
1069: . is - the answer
1071: Level: intermediate
1073: .seealso: EPSIsGeneralized(), EPSIsPositive()
1074: @*/
1075: PetscErrorCode EPSIsHermitian(EPS eps,PetscBool* is)
1076: {
1080: *is = eps->ishermitian;
1081: return(0);
1082: }
1086: /*@
1087: EPSIsPositive - Ask if the EPS object corresponds to an eigenvalue
1088: problem type that requires a positive (semi-) definite matrix B.
1090: Not collective
1092: Input Parameter:
1093: . eps - the eigenproblem solver context
1095: Output Parameter:
1096: . is - the answer
1098: Level: intermediate
1100: .seealso: EPSIsGeneralized(), EPSIsHermitian()
1101: @*/
1102: PetscErrorCode EPSIsPositive(EPS eps,PetscBool* is)
1103: {
1107: *is = eps->ispositive;
1108: return(0);
1109: }