1: /*
2: EPS routines related to options that can be set via the command-line
3: or procedurally.
5: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
6: SLEPc - Scalable Library for Eigenvalue Problem Computations
7: Copyright (c) 2002-2014, Universitat Politecnica de Valencia, Spain
9: This file is part of SLEPc.
11: SLEPc is free software: you can redistribute it and/or modify it under the
12: terms of version 3 of the GNU Lesser General Public License as published by
13: the Free Software Foundation.
15: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
16: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
18: more details.
20: You should have received a copy of the GNU Lesser General Public License
21: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
22: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
23: */
25: #include <slepc-private/epsimpl.h> /*I "slepceps.h" I*/
29: /*@
30: EPSSetFromOptions - Sets EPS options from the options database.
31: This routine must be called before EPSSetUp() if the user is to be
32: allowed to set the solver type.
34: Collective on EPS 36: Input Parameters:
37: . eps - the eigensolver context
39: Notes:
40: To see all options, run your program with the -help option.
42: Level: beginner
43: @*/
44: PetscErrorCode EPSSetFromOptions(EPS eps) 45: {
46: PetscErrorCode ierr;
47: char type[256],monfilename[PETSC_MAX_PATH_LEN];
48: PetscBool flg,flg1,flg2,flg3;
49: PetscReal r,array[2]={0,0};
50: PetscScalar s;
51: PetscInt i,j,k;
52: PetscViewer monviewer;
53: SlepcConvMonitor ctx;
57: if (!EPSRegisterAllCalled) { EPSRegisterAll(); }
58: PetscObjectOptionsBegin((PetscObject)eps);
59: PetscOptionsFList("-eps_type","Eigenvalue Problem Solver method","EPSSetType",EPSList,(char*)(((PetscObject)eps)->type_name?((PetscObject)eps)->type_name:EPSKRYLOVSCHUR),type,256,&flg);
60: if (flg) {
61: EPSSetType(eps,type);
62: }
63: /*
64: Set the type if it was never set.
65: */
66: if (!((PetscObject)eps)->type_name) {
67: EPSSetType(eps,EPSKRYLOVSCHUR);
68: }
70: PetscOptionsBoolGroupBegin("-eps_hermitian","hermitian eigenvalue problem","EPSSetProblemType",&flg);
71: if (flg) { EPSSetProblemType(eps,EPS_HEP); }
72: PetscOptionsBoolGroup("-eps_gen_hermitian","generalized hermitian eigenvalue problem","EPSSetProblemType",&flg);
73: if (flg) { EPSSetProblemType(eps,EPS_GHEP); }
74: PetscOptionsBoolGroup("-eps_non_hermitian","non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
75: if (flg) { EPSSetProblemType(eps,EPS_NHEP); }
76: PetscOptionsBoolGroup("-eps_gen_non_hermitian","generalized non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
77: if (flg) { EPSSetProblemType(eps,EPS_GNHEP); }
78: PetscOptionsBoolGroup("-eps_pos_gen_non_hermitian","generalized non-hermitian eigenvalue problem with positive semi-definite B","EPSSetProblemType",&flg);
79: if (flg) { EPSSetProblemType(eps,EPS_PGNHEP); }
80: PetscOptionsBoolGroupEnd("-eps_gen_indefinite","generalized hermitian-indefinite eigenvalue problem","EPSSetProblemType",&flg);
81: if (flg) { EPSSetProblemType(eps,EPS_GHIEP); }
83: PetscOptionsBoolGroupBegin("-eps_ritz","Rayleigh-Ritz extraction","EPSSetExtraction",&flg);
84: if (flg) { EPSSetExtraction(eps,EPS_RITZ); }
85: PetscOptionsBoolGroup("-eps_harmonic","harmonic Ritz extraction","EPSSetExtraction",&flg);
86: if (flg) { EPSSetExtraction(eps,EPS_HARMONIC); }
87: PetscOptionsBoolGroup("-eps_harmonic_relative","relative harmonic Ritz extraction","EPSSetExtraction",&flg);
88: if (flg) { EPSSetExtraction(eps,EPS_HARMONIC_RELATIVE); }
89: PetscOptionsBoolGroup("-eps_harmonic_right","right harmonic Ritz extraction","EPSSetExtraction",&flg);
90: if (flg) { EPSSetExtraction(eps,EPS_HARMONIC_RIGHT); }
91: PetscOptionsBoolGroup("-eps_harmonic_largest","largest harmonic Ritz extraction","EPSSetExtraction",&flg);
92: if (flg) { EPSSetExtraction(eps,EPS_HARMONIC_LARGEST); }
93: PetscOptionsBoolGroup("-eps_refined","refined Ritz extraction","EPSSetExtraction",&flg);
94: if (flg) { EPSSetExtraction(eps,EPS_REFINED); }
95: PetscOptionsBoolGroupEnd("-eps_refined_harmonic","refined harmonic Ritz extraction","EPSSetExtraction",&flg);
96: if (flg) { EPSSetExtraction(eps,EPS_REFINED_HARMONIC); }
98: PetscOptionsEnum("-eps_balance","Balancing method","EPSSetBalance",EPSBalanceTypes,(PetscEnum)eps->balance,(PetscEnum*)&eps->balance,NULL);
100: j = eps->balance_its;
101: PetscOptionsInt("-eps_balance_its","Number of iterations in balancing","EPSSetBalance",eps->balance_its,&j,&flg1);
102: r = eps->balance_cutoff;
103: PetscOptionsReal("-eps_balance_cutoff","Cutoff value in balancing","EPSSetBalance",eps->balance_cutoff,&r,&flg2);
104: if (flg1 || flg2) {
105: EPSSetBalance(eps,eps->balance,j,r);
106: }
108: i = eps->max_it? eps->max_it: PETSC_DEFAULT;
109: PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,&flg1);
110: r = eps->tol;
111: PetscOptionsReal("-eps_tol","Tolerance","EPSSetTolerances",eps->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:eps->tol,&r,&flg2);
112: if (flg1 || flg2) {
113: EPSSetTolerances(eps,r,i);
114: }
116: PetscOptionsBoolGroupBegin("-eps_conv_eig","Relative error convergence test","EPSSetConvergenceTest",&flg);
117: if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_EIG); }
118: PetscOptionsBoolGroup("-eps_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","EPSSetConvergenceTest",&flg);
119: if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_NORM); }
120: PetscOptionsBoolGroup("-eps_conv_abs","Absolute error convergence test","EPSSetConvergenceTest",&flg);
121: if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_ABS); }
122: PetscOptionsBoolGroupEnd("-eps_conv_user","User-defined convergence test","EPSSetConvergenceTest",&flg);
123: if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_USER); }
125: i = eps->nev;
126: PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,&flg1);
127: j = eps->ncv? eps->ncv: PETSC_DEFAULT;
128: PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&j,&flg2);
129: k = eps->mpd? eps->mpd: PETSC_DEFAULT;
130: PetscOptionsInt("-eps_mpd","Maximum dimension of projected problem","EPSSetDimensions",eps->mpd,&k,&flg3);
131: if (flg1 || flg2 || flg3) {
132: EPSSetDimensions(eps,i,j,k);
133: }
135: /* -----------------------------------------------------------------------*/
136: /*
137: Cancels all monitors hardwired into code before call to EPSSetFromOptions()
138: */
139: flg = PETSC_FALSE;
140: PetscOptionsBool("-eps_monitor_cancel","Remove any hardwired monitor routines","EPSMonitorCancel",flg,&flg,NULL);
141: if (flg) {
142: EPSMonitorCancel(eps);
143: }
144: /*
145: Prints approximate eigenvalues and error estimates at each iteration
146: */
147: PetscOptionsString("-eps_monitor","Monitor first unconverged approximate eigenvalue and error estimate","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
148: if (flg) {
149: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)eps),monfilename,&monviewer);
150: EPSMonitorSet(eps,EPSMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
151: }
152: PetscOptionsString("-eps_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
153: if (flg) {
154: PetscNew(&ctx);
155: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)eps),monfilename,&ctx->viewer);
156: EPSMonitorSet(eps,EPSMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
157: }
158: PetscOptionsString("-eps_monitor_all","Monitor approximate eigenvalues and error estimates","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
159: if (flg) {
160: PetscViewerASCIIOpen(PetscObjectComm((PetscObject)eps),monfilename,&monviewer);
161: EPSMonitorSet(eps,EPSMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
162: EPSSetTrackAll(eps,PETSC_TRUE);
163: }
164: flg = PETSC_FALSE;
165: PetscOptionsBool("-eps_monitor_lg","Monitor first unconverged approximate eigenvalue and error estimate graphically","EPSMonitorSet",flg,&flg,NULL);
166: if (flg) {
167: EPSMonitorSet(eps,EPSMonitorLG,NULL,NULL);
168: }
169: flg = PETSC_FALSE;
170: PetscOptionsBool("-eps_monitor_lg_all","Monitor error estimates graphically","EPSMonitorSet",flg,&flg,NULL);
171: if (flg) {
172: EPSMonitorSet(eps,EPSMonitorLGAll,NULL,NULL);
173: EPSSetTrackAll(eps,PETSC_TRUE);
174: }
175: /* -----------------------------------------------------------------------*/
176: PetscOptionsBoolGroupBegin("-eps_largest_magnitude","compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
177: if (flg) { EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE); }
178: PetscOptionsBoolGroup("-eps_smallest_magnitude","compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
179: if (flg) { EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE); }
180: PetscOptionsBoolGroup("-eps_largest_real","compute largest real parts","EPSSetWhichEigenpairs",&flg);
181: if (flg) { EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL); }
182: PetscOptionsBoolGroup("-eps_smallest_real","compute smallest real parts","EPSSetWhichEigenpairs",&flg);
183: if (flg) { EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL); }
184: PetscOptionsBoolGroup("-eps_largest_imaginary","compute largest imaginary parts","EPSSetWhichEigenpairs",&flg);
185: if (flg) { EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY); }
186: PetscOptionsBoolGroup("-eps_smallest_imaginary","compute smallest imaginary parts","EPSSetWhichEigenpairs",&flg);
187: if (flg) { EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY); }
188: PetscOptionsBoolGroup("-eps_target_magnitude","compute nearest eigenvalues to target","EPSSetWhichEigenpairs",&flg);
189: if (flg) { EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE); }
190: PetscOptionsBoolGroup("-eps_target_real","compute eigenvalues with real parts close to target","EPSSetWhichEigenpairs",&flg);
191: if (flg) { EPSSetWhichEigenpairs(eps,EPS_TARGET_REAL); }
192: PetscOptionsBoolGroup("-eps_target_imaginary","compute eigenvalues with imaginary parts close to target","EPSSetWhichEigenpairs",&flg);
193: if (flg) { EPSSetWhichEigenpairs(eps,EPS_TARGET_IMAGINARY); }
194: PetscOptionsBoolGroupEnd("-eps_all","compute all eigenvalues in an interval","EPSSetWhichEigenpairs",&flg);
195: if (flg) { EPSSetWhichEigenpairs(eps,EPS_ALL); }
197: PetscOptionsScalar("-eps_target","Value of the target","EPSSetTarget",eps->target,&s,&flg);
198: if (flg) {
199: if (eps->which!=EPS_TARGET_REAL && eps->which!=EPS_TARGET_IMAGINARY) {
200: EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);
201: }
202: EPSSetTarget(eps,s);
203: }
204: k = 2;
205: PetscOptionsRealArray("-eps_interval","Computational interval (two real values separated with a comma without spaces)","EPSSetInterval",array,&k,&flg);
206: if (flg) {
207: if (k<2) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_SIZ,"Must pass two values in -eps_interval (comma-separated without spaces)");
208: EPSSetWhichEigenpairs(eps,EPS_ALL);
209: EPSSetInterval(eps,array[0],array[1]);
210: }
212: PetscOptionsBool("-eps_true_residual","Compute true residuals explicitly","EPSSetTrueResidual",eps->trueres,&eps->trueres,NULL);
214: PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",0);
215: PetscOptionsName("-eps_plot_eigs","Make a plot of the computed eigenvalues","EPSSolve",0);
217: if (eps->ops->setfromoptions) {
218: (*eps->ops->setfromoptions)(eps);
219: }
220: PetscObjectProcessOptionsHandlers((PetscObject)eps);
221: PetscOptionsEnd();
223: if (!eps->V) { EPSGetBV(eps,&eps->V); }
224: BVSetFromOptions(eps->V);
225: if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
226: RGSetFromOptions(eps->rg);
227: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
228: DSSetFromOptions(eps->ds);
229: STSetFromOptions(eps->st);
230: PetscRandomSetFromOptions(eps->rand);
231: return(0);
232: }
236: /*@
237: EPSGetTolerances - Gets the tolerance and maximum iteration count used
238: by the EPS convergence tests.
240: Not Collective
242: Input Parameter:
243: . eps - the eigensolver context
245: Output Parameters:
246: + tol - the convergence tolerance
247: - maxits - maximum number of iterations
249: Notes:
250: The user can specify NULL for any parameter that is not needed.
252: Level: intermediate
254: .seealso: EPSSetTolerances()
255: @*/
256: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,PetscInt *maxits)257: {
260: if (tol) *tol = eps->tol;
261: if (maxits) *maxits = eps->max_it;
262: return(0);
263: }
267: /*@
268: EPSSetTolerances - Sets the tolerance and maximum iteration count used
269: by the EPS convergence tests.
271: Logically Collective on EPS273: Input Parameters:
274: + eps - the eigensolver context
275: . tol - the convergence tolerance
276: - maxits - maximum number of iterations to use
278: Options Database Keys:
279: + -eps_tol <tol> - Sets the convergence tolerance
280: - -eps_max_it <maxits> - Sets the maximum number of iterations allowed
282: Notes:
283: Use PETSC_DEFAULT for either argument to assign a reasonably good value.
285: Level: intermediate
287: .seealso: EPSGetTolerances()
288: @*/
289: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,PetscInt maxits)290: {
295: if (tol == PETSC_DEFAULT) {
296: eps->tol = PETSC_DEFAULT;
297: eps->state = EPS_STATE_INITIAL;
298: } else {
299: if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
300: eps->tol = tol;
301: }
302: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
303: eps->max_it = 0;
304: eps->state = EPS_STATE_INITIAL;
305: } else {
306: if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
307: eps->max_it = maxits;
308: }
309: return(0);
310: }
314: /*@
315: EPSGetDimensions - Gets the number of eigenvalues to compute
316: and the dimension of the subspace.
318: Not Collective
320: Input Parameter:
321: . eps - the eigensolver context
323: Output Parameters:
324: + nev - number of eigenvalues to compute
325: . ncv - the maximum dimension of the subspace to be used by the solver
326: - mpd - the maximum dimension allowed for the projected problem
328: Level: intermediate
330: .seealso: EPSSetDimensions()
331: @*/
332: PetscErrorCode EPSGetDimensions(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)333: {
336: if (nev) *nev = eps->nev;
337: if (ncv) *ncv = eps->ncv;
338: if (mpd) *mpd = eps->mpd;
339: return(0);
340: }
344: /*@
345: EPSSetDimensions - Sets the number of eigenvalues to compute
346: and the dimension of the subspace.
348: Logically Collective on EPS350: Input Parameters:
351: + eps - the eigensolver context
352: . nev - number of eigenvalues to compute
353: . ncv - the maximum dimension of the subspace to be used by the solver
354: - mpd - the maximum dimension allowed for the projected problem
356: Options Database Keys:
357: + -eps_nev <nev> - Sets the number of eigenvalues
358: . -eps_ncv <ncv> - Sets the dimension of the subspace
359: - -eps_mpd <mpd> - Sets the maximum projected dimension
361: Notes:
362: Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
363: dependent on the solution method.
365: The parameters ncv and mpd are intimately related, so that the user is advised
366: to set one of them at most. Normal usage is that
367: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
368: (b) in cases where nev is large, the user sets mpd.
370: The value of ncv should always be between nev and (nev+mpd), typically
371: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
372: a smaller value should be used.
374: When computing all eigenvalues in an interval, see EPSSetInterval(), these
375: parameters lose relevance, and tuning must be done with
376: EPSKrylovSchurSetDimensions().
378: Level: intermediate
380: .seealso: EPSGetDimensions(), EPSSetInterval(), EPSKrylovSchurSetDimensions()
381: @*/
382: PetscErrorCode EPSSetDimensions(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)383: {
389: if (nev<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
390: eps->nev = nev;
391: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
392: eps->ncv = 0;
393: } else {
394: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
395: eps->ncv = ncv;
396: }
397: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
398: eps->mpd = 0;
399: } else {
400: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
401: eps->mpd = mpd;
402: }
403: eps->state = EPS_STATE_INITIAL;
404: return(0);
405: }
409: /*@
410: EPSSetWhichEigenpairs - Specifies which portion of the spectrum is
411: to be sought.
413: Logically Collective on EPS415: Input Parameters:
416: + eps - eigensolver context obtained from EPSCreate()
417: - which - the portion of the spectrum to be sought
419: Possible values:
420: The parameter 'which' can have one of these values
422: + EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
423: . EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
424: . EPS_LARGEST_REAL - largest real parts
425: . EPS_SMALLEST_REAL - smallest real parts
426: . EPS_LARGEST_IMAGINARY - largest imaginary parts
427: . EPS_SMALLEST_IMAGINARY - smallest imaginary parts
428: . EPS_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
429: . EPS_TARGET_REAL - eigenvalues with real part closest to target
430: . EPS_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
431: . EPS_ALL - all eigenvalues contained in a given interval
432: - EPS_WHICH_USER - user defined ordering set with EPSSetEigenvalueComparison()
434: Options Database Keys:
435: + -eps_largest_magnitude - Sets largest eigenvalues in magnitude
436: . -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
437: . -eps_largest_real - Sets largest real parts
438: . -eps_smallest_real - Sets smallest real parts
439: . -eps_largest_imaginary - Sets largest imaginary parts
440: . -eps_smallest_imaginary - Sets smallest imaginary parts
441: . -eps_target_magnitude - Sets eigenvalues closest to target
442: . -eps_target_real - Sets real parts closest to target
443: . -eps_target_imaginary - Sets imaginary parts closest to target
444: - -eps_all - Sets all eigenvalues in an interval
446: Notes:
447: Not all eigensolvers implemented in EPS account for all the possible values
448: stated above. Also, some values make sense only for certain types of
449: problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
450: and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part
451: for eigenvalue selection.
453: The target is a scalar value provided with EPSSetTarget().
455: The criterion EPS_TARGET_IMAGINARY is available only in case PETSc and
456: SLEPc have been built with complex scalars.
458: EPS_ALL is intended for use in combination with an interval (see
459: EPSSetInterval()), when all eigenvalues within the interval are requested.
460: In that case, the number of eigenvalues is unknown, so the nev parameter
461: has a different sense, see EPSSetDimensions().
463: Level: intermediate
465: .seealso: EPSGetWhichEigenpairs(), EPSSetTarget(), EPSSetInterval(),
466: EPSSetDimensions(), EPSSetEigenvalueComparison(), EPSWhich467: @*/
468: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)469: {
473: switch (which) {
474: case EPS_LARGEST_MAGNITUDE:
475: case EPS_SMALLEST_MAGNITUDE:
476: case EPS_LARGEST_REAL:
477: case EPS_SMALLEST_REAL:
478: case EPS_LARGEST_IMAGINARY:
479: case EPS_SMALLEST_IMAGINARY:
480: case EPS_TARGET_MAGNITUDE:
481: case EPS_TARGET_REAL:
482: #if defined(PETSC_USE_COMPLEX)
483: case EPS_TARGET_IMAGINARY:
484: #endif
485: case EPS_ALL:
486: case EPS_WHICH_USER:
487: if (eps->which != which) {
488: eps->state = EPS_STATE_INITIAL;
489: eps->which = which;
490: }
491: break;
492: default:493: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
494: }
495: return(0);
496: }
500: /*@C
501: EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be
502: sought.
504: Not Collective
506: Input Parameter:
507: . eps - eigensolver context obtained from EPSCreate()
509: Output Parameter:
510: . which - the portion of the spectrum to be sought
512: Notes:
513: See EPSSetWhichEigenpairs() for possible values of 'which'.
515: Level: intermediate
517: .seealso: EPSSetWhichEigenpairs(), EPSWhich518: @*/
519: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)520: {
524: *which = eps->which;
525: return(0);
526: }
530: /*@C
531: EPSSetEigenvalueComparison - Specifies the eigenvalue comparison function
532: when EPSSetWhichEigenpairs() is set to EPS_WHICH_USER.
534: Logically Collective on EPS536: Input Parameters:
537: + eps - eigensolver context obtained from EPSCreate()
538: . func - a pointer to the comparison function
539: - ctx - a context pointer (the last parameter to the comparison function)
541: Calling Sequence of func:
542: $ func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
544: + ar - real part of the 1st eigenvalue
545: . ai - imaginary part of the 1st eigenvalue
546: . br - real part of the 2nd eigenvalue
547: . bi - imaginary part of the 2nd eigenvalue
548: . res - result of comparison
549: - ctx - optional context, as set by EPSSetEigenvalueComparison()
551: Note:
552: The returning parameter 'res' can be:
553: + negative - if the 1st eigenvalue is preferred to the 2st one
554: . zero - if both eigenvalues are equally preferred
555: - positive - if the 2st eigenvalue is preferred to the 1st one
557: Level: advanced
559: .seealso: EPSSetWhichEigenpairs(), EPSWhich560: @*/
561: PetscErrorCode EPSSetEigenvalueComparison(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)562: {
565: eps->sc->comparison = func;
566: eps->sc->comparisonctx = ctx;
567: eps->which = EPS_WHICH_USER;
568: return(0);
569: }
573: /*@C
574: EPSSetArbitrarySelection - Specifies a function intended to look for
575: eigenvalues according to an arbitrary selection criterion. This criterion
576: can be based on a computation involving the current eigenvector approximation.
578: Logically Collective on EPS580: Input Parameters:
581: + eps - eigensolver context obtained from EPSCreate()
582: . func - a pointer to the evaluation function
583: - ctx - a context pointer (the last parameter to the evaluation function)
585: Calling Sequence of func:
586: $ func(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx)
588: + er - real part of the current eigenvalue approximation
589: . ei - imaginary part of the current eigenvalue approximation
590: . xr - real part of the current eigenvector approximation
591: . xi - imaginary part of the current eigenvector approximation
592: . rr - result of evaluation (real part)
593: . ri - result of evaluation (imaginary part)
594: - ctx - optional context, as set by EPSSetArbitrarySelection()
596: Notes:
597: This provides a mechanism to select eigenpairs by evaluating a user-defined
598: function. When a function has been provided, the default selection based on
599: sorting the eigenvalues is replaced by the sorting of the results of this
600: function (with the same sorting criterion given in EPSSetWhichEigenpairs()).
602: For instance, suppose you want to compute those eigenvectors that maximize
603: a certain computable expression. Then implement the computation using
604: the arguments xr and xi, and return the result in rr. Then set the standard
605: sorting by magnitude so that the eigenpair with largest value of rr is
606: selected.
608: This evaluation function is collective, that is, all processes call it and
609: it can use collective operations; furthermore, the computed result must
610: be the same in all processes.
612: The result of func is expressed as a complex number so that it is possible to
613: use the standard eigenvalue sorting functions, but normally only rr is used.
614: Set ri to zero unless it is meaningful in your application.
616: Level: advanced
618: .seealso: EPSSetWhichEigenpairs()
619: @*/
620: PetscErrorCode EPSSetArbitrarySelection(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,Vec,Vec,PetscScalar*,PetscScalar*,void*),void* ctx)621: {
624: eps->arbitrary = func;
625: eps->arbitraryctx = ctx;
626: eps->state = EPS_STATE_INITIAL;
627: return(0);
628: }
632: /*@C
633: EPSSetConvergenceTestFunction - Sets a function to compute the error estimate
634: used in the convergence test.
636: Logically Collective on EPS638: Input Parameters:
639: + eps - eigensolver context obtained from EPSCreate()
640: . func - a pointer to the convergence test function
641: . ctx - [optional] context for private data for the convergence routine
642: - destroy - [optional] destructor for the context (may be NULL;
643: PETSC_NULL_FUNCTION in Fortran)
645: Calling Sequence of func:
646: $ func(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
648: + eps - eigensolver context obtained from EPSCreate()
649: . eigr - real part of the eigenvalue
650: . eigi - imaginary part of the eigenvalue
651: . res - residual norm associated to the eigenpair
652: . errest - (output) computed error estimate
653: - ctx - optional context, as set by EPSSetConvergenceTest()
655: Note:
656: If the error estimate returned by the convergence test function is less than
657: the tolerance, then the eigenvalue is accepted as converged.
659: Level: advanced
661: .seealso: EPSSetConvergenceTest(), EPSSetTolerances()
662: @*/
663: PetscErrorCode EPSSetConvergenceTestFunction(EPS eps,PetscErrorCode (*func)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))664: {
669: if (eps->convergeddestroy) {
670: (*eps->convergeddestroy)(eps->convergedctx);
671: }
672: eps->converged = func;
673: eps->convergeddestroy = destroy;
674: eps->convergedctx = ctx;
675: if (func == EPSConvergedEigRelative) eps->conv = EPS_CONV_EIG;
676: else if (func == EPSConvergedNormRelative) eps->conv = EPS_CONV_NORM;
677: else if (func == EPSConvergedAbsolute) eps->conv = EPS_CONV_ABS;
678: else eps->conv = EPS_CONV_USER;
679: return(0);
680: }
684: /*@
685: EPSSetConvergenceTest - Specifies how to compute the error estimate
686: used in the convergence test.
688: Logically Collective on EPS690: Input Parameters:
691: + eps - eigensolver context obtained from EPSCreate()
692: - conv - the type of convergence test
694: Options Database Keys:
695: + -eps_conv_abs - Sets the absolute convergence test
696: . -eps_conv_eig - Sets the convergence test relative to the eigenvalue
697: . -eps_conv_norm - Sets the convergence test relative to the matrix norms
698: - -eps_conv_user - Selects the user-defined convergence test
700: Note:
701: The parameter 'conv' can have one of these values
702: + EPS_CONV_ABS - absolute error ||r||
703: . EPS_CONV_EIG - error relative to the eigenvalue l, ||r||/|l|
704: . EPS_CONV_NORM - error relative to the matrix norms, ||r||/(||A||+|l|*||B||)
705: - EPS_CONV_USER - function set by EPSSetConvergenceTestFunction()
707: Level: intermediate
709: .seealso: EPSGetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSConv710: @*/
711: PetscErrorCode EPSSetConvergenceTest(EPS eps,EPSConv conv)712: {
716: switch (conv) {
717: case EPS_CONV_ABS: eps->converged = EPSConvergedAbsolute; break;
718: case EPS_CONV_EIG: eps->converged = EPSConvergedEigRelative; break;
719: case EPS_CONV_NORM: eps->converged = EPSConvergedNormRelative; break;
720: case EPS_CONV_USER: break;
721: default:722: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
723: }
724: eps->conv = conv;
725: return(0);
726: }
730: /*@C
731: EPSGetConvergenceTest - Gets the method used to compute the error estimate
732: used in the convergence test.
734: Not Collective
736: Input Parameters:
737: . eps - eigensolver context obtained from EPSCreate()
739: Output Parameters:
740: . conv - the type of convergence test
742: Level: intermediate
744: .seealso: EPSSetConvergenceTest(), EPSConv745: @*/
746: PetscErrorCode EPSGetConvergenceTest(EPS eps,EPSConv *conv)747: {
751: *conv = eps->conv;
752: return(0);
753: }
757: /*@
758: EPSSetProblemType - Specifies the type of the eigenvalue problem.
760: Logically Collective on EPS762: Input Parameters:
763: + eps - the eigensolver context
764: - type - a known type of eigenvalue problem
766: Options Database Keys:
767: + -eps_hermitian - Hermitian eigenvalue problem
768: . -eps_gen_hermitian - generalized Hermitian eigenvalue problem
769: . -eps_non_hermitian - non-Hermitian eigenvalue problem
770: . -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
771: - -eps_pos_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
772: with positive semi-definite B
774: Notes:
775: Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
776: (EPS_NHEP), generalized Hermitian (EPS_GHEP), generalized non-Hermitian
777: (EPS_GNHEP), generalized non-Hermitian with positive semi-definite B
778: (EPS_PGNHEP), and generalized Hermitian-indefinite (EPS_GHIEP).
780: This function must be used to instruct SLEPc to exploit symmetry. If no
781: problem type is specified, by default a non-Hermitian problem is assumed
782: (either standard or generalized). If the user knows that the problem is
783: Hermitian (i.e. A=A^H) or generalized Hermitian (i.e. A=A^H, B=B^H, and
784: B positive definite) then it is recommended to set the problem type so
785: that eigensolver can exploit these properties.
787: Level: beginner
789: .seealso: EPSSetOperators(), EPSSetType(), EPSGetProblemType(), EPSProblemType790: @*/
791: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)792: {
796: switch (type) {
797: case EPS_HEP:
798: eps->isgeneralized = PETSC_FALSE;
799: eps->ishermitian = PETSC_TRUE;
800: eps->ispositive = PETSC_FALSE;
801: break;
802: case EPS_NHEP:
803: eps->isgeneralized = PETSC_FALSE;
804: eps->ishermitian = PETSC_FALSE;
805: eps->ispositive = PETSC_FALSE;
806: break;
807: case EPS_GHEP:
808: eps->isgeneralized = PETSC_TRUE;
809: eps->ishermitian = PETSC_TRUE;
810: eps->ispositive = PETSC_TRUE;
811: break;
812: case EPS_GNHEP:
813: eps->isgeneralized = PETSC_TRUE;
814: eps->ishermitian = PETSC_FALSE;
815: eps->ispositive = PETSC_FALSE;
816: break;
817: case EPS_PGNHEP:
818: eps->isgeneralized = PETSC_TRUE;
819: eps->ishermitian = PETSC_FALSE;
820: eps->ispositive = PETSC_TRUE;
821: break;
822: case EPS_GHIEP:
823: eps->isgeneralized = PETSC_TRUE;
824: eps->ishermitian = PETSC_TRUE;
825: eps->ispositive = PETSC_FALSE;
826: break;
827: default:828: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
829: }
830: eps->problem_type = type;
831: return(0);
832: }
836: /*@C
837: EPSGetProblemType - Gets the problem type from the EPS object.
839: Not Collective
841: Input Parameter:
842: . eps - the eigensolver context
844: Output Parameter:
845: . type - name of EPS problem type
847: Level: intermediate
849: .seealso: EPSSetProblemType(), EPSProblemType850: @*/
851: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)852: {
856: *type = eps->problem_type;
857: return(0);
858: }
862: /*@
863: EPSSetExtraction - Specifies the type of extraction technique to be employed
864: by the eigensolver.
866: Logically Collective on EPS868: Input Parameters:
869: + eps - the eigensolver context
870: - extr - a known type of extraction
872: Options Database Keys:
873: + -eps_ritz - Rayleigh-Ritz extraction
874: . -eps_harmonic - harmonic Ritz extraction
875: . -eps_harmonic_relative - harmonic Ritz extraction relative to the eigenvalue
876: . -eps_harmonic_right - harmonic Ritz extraction for rightmost eigenvalues
877: . -eps_harmonic_largest - harmonic Ritz extraction for largest magnitude
878: (without target)
879: . -eps_refined - refined Ritz extraction
880: - -eps_refined_harmonic - refined harmonic Ritz extraction
882: Notes:
883: Not all eigensolvers support all types of extraction. See the SLEPc
884: Users Manual for details.
886: By default, a standard Rayleigh-Ritz extraction is used. Other extractions
887: may be useful when computing interior eigenvalues.
889: Harmonic-type extractions are used in combination with a 'target'.
891: Level: beginner
893: .seealso: EPSSetTarget(), EPSGetExtraction(), EPSExtraction894: @*/
895: PetscErrorCode EPSSetExtraction(EPS eps,EPSExtraction extr)896: {
900: eps->extraction = extr;
901: return(0);
902: }
906: /*@C
907: EPSGetExtraction - Gets the extraction type used by the EPS object.
909: Not Collective
911: Input Parameter:
912: . eps - the eigensolver context
914: Output Parameter:
915: . extr - name of extraction type
917: Level: intermediate
919: .seealso: EPSSetExtraction(), EPSExtraction920: @*/
921: PetscErrorCode EPSGetExtraction(EPS eps,EPSExtraction *extr)922: {
926: *extr = eps->extraction;
927: return(0);
928: }
932: /*@
933: EPSSetBalance - Specifies the balancing technique to be employed by the
934: eigensolver, and some parameters associated to it.
936: Logically Collective on EPS938: Input Parameters:
939: + eps - the eigensolver context
940: . bal - the balancing method, one of EPS_BALANCE_NONE, EPS_BALANCE_ONESIDE,
941: EPS_BALANCE_TWOSIDE, or EPS_BALANCE_USER
942: . its - number of iterations of the balancing algorithm
943: - cutoff - cutoff value
945: Options Database Keys:
946: + -eps_balance <method> - the balancing method, where <method> is one of
947: 'none', 'oneside', 'twoside', or 'user'
948: . -eps_balance_its <its> - number of iterations
949: - -eps_balance_cutoff <cutoff> - cutoff value
951: Notes:
952: When balancing is enabled, the solver works implicitly with matrix DAD^-1,
953: where D is an appropriate diagonal matrix. This improves the accuracy of
954: the computed results in some cases. See the SLEPc Users Manual for details.
956: Balancing makes sense only for non-Hermitian problems when the required
957: precision is high (i.e. a small tolerance such as 1e-15).
959: By default, balancing is disabled. The two-sided method is much more
960: effective than the one-sided counterpart, but it requires the system
961: matrices to have the MatMultTranspose operation defined.
963: The parameter 'its' is the number of iterations performed by the method. The
964: cutoff value is used only in the two-side variant. Use PETSC_DEFAULT to assign
965: a reasonably good value.
967: User-defined balancing is allowed provided that the corresponding matrix
968: is set via STSetBalanceMatrix.
970: Level: intermediate
972: .seealso: EPSGetBalance(), EPSBalance, STSetBalanceMatrix()
973: @*/
974: PetscErrorCode EPSSetBalance(EPS eps,EPSBalance bal,PetscInt its,PetscReal cutoff)975: {
981: switch (bal) {
982: case EPS_BALANCE_NONE:
983: case EPS_BALANCE_ONESIDE:
984: case EPS_BALANCE_TWOSIDE:
985: case EPS_BALANCE_USER:
986: eps->balance = bal;
987: break;
988: default:989: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid value of argument 'bal'");
990: }
991: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) eps->balance_its = 5;
992: else {
993: if (its <= 0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be > 0");
994: eps->balance_its = its;
995: }
996: if (cutoff==PETSC_DECIDE || cutoff==PETSC_DEFAULT) eps->balance_cutoff = 1e-8;
997: else {
998: if (cutoff <= 0.0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of cutoff. Must be > 0");
999: eps->balance_cutoff = cutoff;
1000: }
1001: return(0);
1002: }
1006: /*@C
1007: EPSGetBalance - Gets the balancing type used by the EPS object, and the
1008: associated parameters.
1010: Not Collective
1012: Input Parameter:
1013: . eps - the eigensolver context
1015: Output Parameters:
1016: + bal - the balancing method
1017: . its - number of iterations of the balancing algorithm
1018: - cutoff - cutoff value
1020: Level: intermediate
1022: Note:
1023: The user can specify NULL for any parameter that is not needed.
1025: .seealso: EPSSetBalance(), EPSBalance1026: @*/
1027: PetscErrorCode EPSGetBalance(EPS eps,EPSBalance *bal,PetscInt *its,PetscReal *cutoff)1028: {
1031: if (bal) *bal = eps->balance;
1032: if (its) *its = eps->balance_its;
1033: if (cutoff) *cutoff = eps->balance_cutoff;
1034: return(0);
1035: }
1039: /*@
1040: EPSSetTrueResidual - Specifies if the solver must compute the true residual
1041: explicitly or not.
1043: Logically Collective on EPS1045: Input Parameters:
1046: + eps - the eigensolver context
1047: - trueres - whether true residuals are required or not
1049: Options Database Keys:
1050: . -eps_true_residual <boolean> - Sets/resets the boolean flag 'trueres'
1052: Notes:
1053: If the user sets trueres=PETSC_TRUE then the solver explicitly computes
1054: the true residual for each eigenpair approximation, and uses it for
1055: convergence testing. Computing the residual is usually an expensive
1056: operation. Some solvers (e.g., Krylov solvers) can avoid this computation
1057: by using a cheap estimate of the residual norm, but this may sometimes
1058: give inaccurate results (especially if a spectral transform is being
1059: used). On the contrary, preconditioned eigensolvers (e.g., Davidson solvers)
1060: do rely on computing the true residual, so this option is irrelevant for them.
1062: Level: intermediate
1064: .seealso: EPSGetTrueResidual()
1065: @*/
1066: PetscErrorCode EPSSetTrueResidual(EPS eps,PetscBool trueres)1067: {
1071: eps->trueres = trueres;
1072: return(0);
1073: }
1077: /*@C
1078: EPSGetTrueResidual - Returns the flag indicating whether true
1079: residuals must be computed explicitly or not.
1081: Not Collective
1083: Input Parameter:
1084: . eps - the eigensolver context
1086: Output Parameter:
1087: . trueres - the returned flag
1089: Level: intermediate
1091: .seealso: EPSSetTrueResidual()
1092: @*/
1093: PetscErrorCode EPSGetTrueResidual(EPS eps,PetscBool *trueres)1094: {
1098: *trueres = eps->trueres;
1099: return(0);
1100: }
1104: /*@
1105: EPSSetTrackAll - Specifies if the solver must compute the residual norm of all
1106: approximate eigenpairs or not.
1108: Logically Collective on EPS1110: Input Parameters:
1111: + eps - the eigensolver context
1112: - trackall - whether to compute all residuals or not
1114: Notes:
1115: If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
1116: the residual norm for each eigenpair approximation. Computing the residual is
1117: usually an expensive operation and solvers commonly compute only the residual
1118: associated to the first unconverged eigenpair.
1120: The options '-eps_monitor_all' and '-eps_monitor_lg_all' automatically
1121: activate this option.
1123: Level: intermediate
1125: .seealso: EPSGetTrackAll()
1126: @*/
1127: PetscErrorCode EPSSetTrackAll(EPS eps,PetscBool trackall)1128: {
1132: eps->trackall = trackall;
1133: return(0);
1134: }
1138: /*@
1139: EPSGetTrackAll - Returns the flag indicating whether all residual norms must
1140: be computed or not.
1142: Not Collective
1144: Input Parameter:
1145: . eps - the eigensolver context
1147: Output Parameter:
1148: . trackall - the returned flag
1150: Level: intermediate
1152: .seealso: EPSSetTrackAll()
1153: @*/
1154: PetscErrorCode EPSGetTrackAll(EPS eps,PetscBool *trackall)1155: {
1159: *trackall = eps->trackall;
1160: return(0);
1161: }
1165: /*@C
1166: EPSSetOptionsPrefix - Sets the prefix used for searching for all
1167: EPS options in the database.
1169: Logically Collective on EPS1171: Input Parameters:
1172: + eps - the eigensolver context
1173: - prefix - the prefix string to prepend to all EPS option requests
1175: Notes:
1176: A hyphen (-) must NOT be given at the beginning of the prefix name.
1177: The first character of all runtime options is AUTOMATICALLY the
1178: hyphen.
1180: For example, to distinguish between the runtime options for two
1181: different EPS contexts, one could call
1182: .vb
1183: EPSSetOptionsPrefix(eps1,"eig1_")
1184: EPSSetOptionsPrefix(eps2,"eig2_")
1185: .ve
1187: Level: advanced
1189: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
1190: @*/
1191: PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix)1192: {
1197: if (!eps->st) { EPSGetST(eps,&eps->st); }
1198: STSetOptionsPrefix(eps->st,prefix);
1199: if (!eps->V) { EPSGetBV(eps,&eps->V); }
1200: BVSetOptionsPrefix(eps->V,prefix);
1201: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1202: DSSetOptionsPrefix(eps->ds,prefix);
1203: if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
1204: RGSetOptionsPrefix(eps->rg,prefix);
1205: PetscObjectSetOptionsPrefix((PetscObject)eps,prefix);
1206: return(0);
1207: }
1211: /*@C
1212: EPSAppendOptionsPrefix - Appends to the prefix used for searching for all
1213: EPS options in the database.
1215: Logically Collective on EPS1217: Input Parameters:
1218: + eps - the eigensolver context
1219: - prefix - the prefix string to prepend to all EPS option requests
1221: Notes:
1222: A hyphen (-) must NOT be given at the beginning of the prefix name.
1223: The first character of all runtime options is AUTOMATICALLY the hyphen.
1225: Level: advanced
1227: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
1228: @*/
1229: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix)1230: {
1235: if (!eps->st) { EPSGetST(eps,&eps->st); }
1236: STAppendOptionsPrefix(eps->st,prefix);
1237: if (!eps->V) { EPSGetBV(eps,&eps->V); }
1238: BVSetOptionsPrefix(eps->V,prefix);
1239: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1240: DSSetOptionsPrefix(eps->ds,prefix);
1241: if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
1242: RGSetOptionsPrefix(eps->rg,prefix);
1243: PetscObjectAppendOptionsPrefix((PetscObject)eps,prefix);
1244: return(0);
1245: }
1249: /*@C
1250: EPSGetOptionsPrefix - Gets the prefix used for searching for all
1251: EPS options in the database.
1253: Not Collective
1255: Input Parameters:
1256: . eps - the eigensolver context
1258: Output Parameters:
1259: . prefix - pointer to the prefix string used is returned
1261: Notes: On the fortran side, the user should pass in a string 'prefix' of
1262: sufficient length to hold the prefix.
1264: Level: advanced
1266: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
1267: @*/
1268: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])1269: {
1275: PetscObjectGetOptionsPrefix((PetscObject)eps,prefix);
1276: return(0);
1277: }