Actual source code: opts.c
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-2012, Universitat Politecnica de Valencia, Spain
9: This file is part of SLEPc.
10:
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,val;
49: PetscReal r,nrma,nrmb,array[2]={0,0};
50: PetscScalar s;
51: PetscInt i,j,k;
52: const char *bal_list[4] = {"none","oneside","twoside","user"};
53: PetscViewer monviewer;
54: SlepcConvMonitor ctx;
58: if (!EPSRegisterAllCalled) { EPSRegisterAll(PETSC_NULL); }
59: PetscOptionsBegin(((PetscObject)eps)->comm,((PetscObject)eps)->prefix,"Eigenproblem Solver (EPS) Options","EPS");
60: PetscOptionsList("-eps_type","Eigenproblem Solver method","EPSSetType",EPSList,(char*)(((PetscObject)eps)->type_name?((PetscObject)eps)->type_name:EPSKRYLOVSCHUR),type,256,&flg);
61: if (flg) {
62: EPSSetType(eps,type);
63: }
64: /*
65: Set the type if it was never set.
66: */
67: if (!((PetscObject)eps)->type_name) {
68: EPSSetType(eps,EPSKRYLOVSCHUR);
69: }
71: PetscOptionsBoolGroupBegin("-eps_hermitian","hermitian eigenvalue problem","EPSSetProblemType",&flg);
72: if (flg) {EPSSetProblemType(eps,EPS_HEP);}
73: PetscOptionsBoolGroup("-eps_gen_hermitian","generalized hermitian eigenvalue problem","EPSSetProblemType",&flg);
74: if (flg) {EPSSetProblemType(eps,EPS_GHEP);}
75: PetscOptionsBoolGroup("-eps_non_hermitian","non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
76: if (flg) {EPSSetProblemType(eps,EPS_NHEP);}
77: PetscOptionsBoolGroup("-eps_gen_non_hermitian","generalized non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
78: if (flg) {EPSSetProblemType(eps,EPS_GNHEP);}
79: PetscOptionsBoolGroup("-eps_pos_gen_non_hermitian","generalized non-hermitian eigenvalue problem with positive semi-definite B","EPSSetProblemType",&flg);
80: if (flg) {EPSSetProblemType(eps,EPS_PGNHEP);}
81: PetscOptionsBoolGroupEnd("-eps_gen_indefinite","generalized hermitian-indefinite eigenvalue problem","EPSSetProblemType",&flg);
82: if (flg) {EPSSetProblemType(eps,EPS_GHIEP);}
84: PetscOptionsBoolGroupBegin("-eps_ritz","Rayleigh-Ritz extraction","EPSSetExtraction",&flg);
85: if (flg) {EPSSetExtraction(eps,EPS_RITZ);}
86: PetscOptionsBoolGroup("-eps_harmonic","harmonic Ritz extraction","EPSSetExtraction",&flg);
87: if (flg) {EPSSetExtraction(eps,EPS_HARMONIC);}
88: PetscOptionsBoolGroup("-eps_harmonic_relative","relative harmonic Ritz extraction","EPSSetExtraction",&flg);
89: if (flg) {EPSSetExtraction(eps,EPS_HARMONIC_RELATIVE);}
90: PetscOptionsBoolGroup("-eps_harmonic_right","right harmonic Ritz extraction","EPSSetExtraction",&flg);
91: if (flg) {EPSSetExtraction(eps,EPS_HARMONIC_RIGHT);}
92: PetscOptionsBoolGroup("-eps_harmonic_largest","largest harmonic Ritz extraction","EPSSetExtraction",&flg);
93: if (flg) {EPSSetExtraction(eps,EPS_HARMONIC_LARGEST);}
94: PetscOptionsBoolGroup("-eps_refined","refined Ritz extraction","EPSSetExtraction",&flg);
95: if (flg) {EPSSetExtraction(eps,EPS_REFINED);}
96: PetscOptionsBoolGroupEnd("-eps_refined_harmonic","refined harmonic Ritz extraction","EPSSetExtraction",&flg);
97: if (flg) {EPSSetExtraction(eps,EPS_REFINED_HARMONIC);}
99: if (!eps->balance) eps->balance = EPS_BALANCE_NONE;
100: PetscOptionsEList("-eps_balance","Balancing method","EPSSetBalance",bal_list,4,bal_list[eps->balance-EPS_BALANCE_NONE],&i,&flg);
101: if (flg) { eps->balance = (EPSBalance)(i+EPS_BALANCE_NONE); }
102: r = j = PETSC_IGNORE;
103: PetscOptionsInt("-eps_balance_its","Number of iterations in balancing","EPSSetBalance",eps->balance_its,&j,PETSC_NULL);
104: PetscOptionsReal("-eps_balance_cutoff","Cutoff value in balancing","EPSSetBalance",eps->balance_cutoff,&r,PETSC_NULL);
105: EPSSetBalance(eps,(EPSBalance)PETSC_IGNORE,j,r);
107: r = i = PETSC_IGNORE;
108: PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,PETSC_NULL);
109: PetscOptionsReal("-eps_tol","Tolerance","EPSSetTolerances",eps->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:eps->tol,&r,PETSC_NULL);
110: EPSSetTolerances(eps,r,i);
111: PetscOptionsBoolGroupBegin("-eps_conv_eig","relative error convergence test","EPSSetConvergenceTest",&flg);
112: if (flg) {EPSSetConvergenceTest(eps,EPS_CONV_EIG);}
113: PetscOptionsBoolGroup("-eps_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","EPSSetConvergenceTest",&flg);
114: if (flg) {EPSSetConvergenceTest(eps,EPS_CONV_NORM);}
115: PetscOptionsBoolGroupEnd("-eps_conv_abs","Absolute error convergence test","EPSSetConvergenceTest",&flg);
116: if (flg) {EPSSetConvergenceTest(eps,EPS_CONV_ABS);}
118: i = j = k = PETSC_IGNORE;
119: PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,PETSC_NULL);
120: PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&j,PETSC_NULL);
121: PetscOptionsInt("-eps_mpd","Maximum dimension of projected problem","EPSSetDimensions",eps->mpd,&k,PETSC_NULL);
122: EPSSetDimensions(eps,i,j,k);
123:
124: /* -----------------------------------------------------------------------*/
125: /*
126: Cancels all monitors hardwired into code before call to EPSSetFromOptions()
127: */
128: flg = PETSC_FALSE;
129: PetscOptionsBool("-eps_monitor_cancel","Remove any hardwired monitor routines","EPSMonitorCancel",flg,&flg,PETSC_NULL);
130: if (flg) {
131: EPSMonitorCancel(eps);
132: }
133: /*
134: Prints approximate eigenvalues and error estimates at each iteration
135: */
136: PetscOptionsString("-eps_monitor","Monitor first unconverged approximate eigenvalue and error estimate","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
137: if (flg) {
138: PetscViewerASCIIOpen(((PetscObject)eps)->comm,monfilename,&monviewer);
139: EPSMonitorSet(eps,EPSMonitorFirst,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
140: }
141: PetscOptionsString("-eps_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
142: if (flg) {
143: PetscNew(struct _n_SlepcConvMonitor,&ctx);
144: PetscViewerASCIIOpen(((PetscObject)eps)->comm,monfilename,&ctx->viewer);
145: EPSMonitorSet(eps,EPSMonitorConverged,ctx,(PetscErrorCode (*)(void**))SlepcConvMonitorDestroy);
146: }
147: PetscOptionsString("-eps_monitor_all","Monitor approximate eigenvalues and error estimates","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
148: if (flg) {
149: PetscViewerASCIIOpen(((PetscObject)eps)->comm,monfilename,&monviewer);
150: EPSMonitorSet(eps,EPSMonitorAll,monviewer,(PetscErrorCode (*)(void**))PetscViewerDestroy);
151: EPSSetTrackAll(eps,PETSC_TRUE);
152: }
153: flg = PETSC_FALSE;
154: PetscOptionsBool("-eps_monitor_draw","Monitor first unconverged approximate eigenvalue and error estimate graphically","EPSMonitorSet",flg,&flg,PETSC_NULL);
155: if (flg) {
156: EPSMonitorSet(eps,EPSMonitorLG,PETSC_NULL,PETSC_NULL);
157: }
158: flg = PETSC_FALSE;
159: PetscOptionsBool("-eps_monitor_draw_all","Monitor error estimates graphically","EPSMonitorSet",flg,&flg,PETSC_NULL);
160: if (flg) {
161: EPSMonitorSet(eps,EPSMonitorLGAll,PETSC_NULL,PETSC_NULL);
162: EPSSetTrackAll(eps,PETSC_TRUE);
163: }
164: /* -----------------------------------------------------------------------*/
165: PetscOptionsScalar("-eps_target","Value of the target","EPSSetTarget",eps->target,&s,&flg);
166: if (flg) {
167: EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);
168: EPSSetTarget(eps,s);
169: }
170: k = 2;
171: PetscOptionsRealArray("-eps_interval","Computational interval (two real values separated with a comma without spaces)","EPSSetInterval",array,&k,&flg);
172: if (flg) {
173: if (k<2) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_SIZ,"Must pass two values in -eps_interval (comma-separated without spaces)");
174: EPSSetWhichEigenpairs(eps,EPS_ALL);
175: EPSSetInterval(eps,array[0],array[1]);
176: }
178: PetscOptionsBoolGroupBegin("-eps_largest_magnitude","compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
179: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE);}
180: PetscOptionsBoolGroup("-eps_smallest_magnitude","compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
181: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE);}
182: PetscOptionsBoolGroup("-eps_largest_real","compute largest real parts","EPSSetWhichEigenpairs",&flg);
183: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL);}
184: PetscOptionsBoolGroup("-eps_smallest_real","compute smallest real parts","EPSSetWhichEigenpairs",&flg);
185: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);}
186: PetscOptionsBoolGroup("-eps_largest_imaginary","compute largest imaginary parts","EPSSetWhichEigenpairs",&flg);
187: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY);}
188: PetscOptionsBoolGroup("-eps_smallest_imaginary","compute smallest imaginary parts","EPSSetWhichEigenpairs",&flg);
189: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY);}
190: PetscOptionsBoolGroup("-eps_target_magnitude","compute nearest eigenvalues to target","EPSSetWhichEigenpairs",&flg);
191: if (flg) {EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);}
192: PetscOptionsBoolGroup("-eps_target_real","compute eigenvalues with real parts close to target","EPSSetWhichEigenpairs",&flg);
193: if (flg) {EPSSetWhichEigenpairs(eps,EPS_TARGET_REAL);}
194: PetscOptionsBoolGroup("-eps_target_imaginary","compute eigenvalues with imaginary parts close to target","EPSSetWhichEigenpairs",&flg);
195: if (flg) {EPSSetWhichEigenpairs(eps,EPS_TARGET_IMAGINARY);}
196: PetscOptionsBoolGroupEnd("-eps_all","compute all eigenvalues in an interval","EPSSetWhichEigenpairs",&flg);
197: if (flg) {EPSSetWhichEigenpairs(eps,EPS_ALL);}
199: PetscOptionsBool("-eps_left_vectors","Compute left eigenvectors also","EPSSetLeftVectorsWanted",eps->leftvecs,&val,&flg);
200: if (flg) {
201: EPSSetLeftVectorsWanted(eps,val);
202: }
203: PetscOptionsBool("-eps_true_residual","Compute true residuals explicitly","EPSSetTrueResidual",eps->trueres,&val,&flg);
204: if (flg) {
205: EPSSetTrueResidual(eps,val);
206: }
208: nrma = nrmb = PETSC_IGNORE;
209: PetscOptionsReal("-eps_norm_a","Norm of matrix A","EPSSetMatrixNorms",eps->nrma,&nrma,PETSC_NULL);
210: PetscOptionsReal("-eps_norm_b","Norm of matrix B","EPSSetMatrixNorms",eps->nrmb,&nrmb,PETSC_NULL);
211: EPSSetMatrixNorms(eps,nrma,nrmb,eps->adaptive);
212: PetscOptionsBool("-eps_norms_adaptive","Update the value of matrix norms adaptively","EPSSetMatrixNorms",eps->adaptive,&val,&flg);
213: if (flg) {
214: EPSSetMatrixNorms(eps,PETSC_IGNORE,PETSC_IGNORE,val);
215: }
217: PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",0);
218: PetscOptionsName("-eps_view_binary","Save the matrices associated to the eigenproblem","EPSSetFromOptions",0);
219: PetscOptionsName("-eps_plot_eigs","Make a plot of the computed eigenvalues","EPSSolve",0);
220:
221: if (eps->ops->setfromoptions) {
222: (*eps->ops->setfromoptions)(eps);
223: }
224: PetscObjectProcessOptionsHandlers((PetscObject)eps);
225: PetscOptionsEnd();
227: if (!eps->ip) { EPSGetIP(eps,&eps->ip); }
228: IPSetFromOptions(eps->ip);
229: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
230: DSSetFromOptions(eps->ds);
231: STSetFromOptions(eps->OP);
232: PetscRandomSetFromOptions(eps->rand);
233: return(0);
234: }
238: /*@
239: EPSGetTolerances - Gets the tolerance and maximum iteration count used
240: by the EPS convergence tests.
242: Not Collective
244: Input Parameter:
245: . eps - the eigensolver context
246:
247: Output Parameters:
248: + tol - the convergence tolerance
249: - maxits - maximum number of iterations
251: Notes:
252: The user can specify PETSC_NULL for any parameter that is not needed.
254: Level: intermediate
256: .seealso: EPSSetTolerances()
257: @*/
258: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,PetscInt *maxits)
259: {
262: if (tol) *tol = eps->tol;
263: if (maxits) *maxits = eps->max_it;
264: return(0);
265: }
269: /*@
270: EPSSetTolerances - Sets the tolerance and maximum iteration count used
271: by the EPS convergence tests.
273: Logically Collective on EPS
275: Input Parameters:
276: + eps - the eigensolver context
277: . tol - the convergence tolerance
278: - maxits - maximum number of iterations to use
280: Options Database Keys:
281: + -eps_tol <tol> - Sets the convergence tolerance
282: - -eps_max_it <maxits> - Sets the maximum number of iterations allowed
284: Notes:
285: Use PETSC_IGNORE for an argument that need not be changed.
287: Use PETSC_DECIDE for maxits to assign a reasonably good value, which is
288: dependent on the solution method.
290: Level: intermediate
292: .seealso: EPSGetTolerances()
293: @*/
294: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,PetscInt maxits)
295: {
300: if (tol != PETSC_IGNORE) {
301: if (tol == PETSC_DEFAULT) {
302: eps->tol = PETSC_DEFAULT;
303: } else {
304: if (tol < 0.0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
305: eps->tol = tol;
306: }
307: }
308: if (maxits != PETSC_IGNORE) {
309: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
310: eps->max_it = 0;
311: eps->setupcalled = 0;
312: } else {
313: if (maxits < 0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
314: eps->max_it = maxits;
315: }
316: }
317: return(0);
318: }
322: /*@
323: EPSGetDimensions - Gets the number of eigenvalues to compute
324: and the dimension of the subspace.
326: Not Collective
328: Input Parameter:
329: . eps - the eigensolver context
330:
331: Output Parameters:
332: + nev - number of eigenvalues to compute
333: . ncv - the maximum dimension of the subspace to be used by the solver
334: - mpd - the maximum dimension allowed for the projected problem
336: Notes:
337: The user can specify PETSC_NULL for any parameter that is not needed.
339: Level: intermediate
341: .seealso: EPSSetDimensions()
342: @*/
343: PetscErrorCode EPSGetDimensions(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
344: {
347: if (nev) *nev = eps->nev;
348: if (ncv) *ncv = eps->ncv;
349: if (mpd) *mpd = eps->mpd;
350: return(0);
351: }
355: /*@
356: EPSSetDimensions - Sets the number of eigenvalues to compute
357: and the dimension of the subspace.
359: Logically Collective on EPS
361: Input Parameters:
362: + eps - the eigensolver context
363: . nev - number of eigenvalues to compute
364: . ncv - the maximum dimension of the subspace to be used by the solver
365: - mpd - the maximum dimension allowed for the projected problem
367: Options Database Keys:
368: + -eps_nev <nev> - Sets the number of eigenvalues
369: . -eps_ncv <ncv> - Sets the dimension of the subspace
370: - -eps_mpd <mpd> - Sets the maximum projected dimension
372: Notes:
373: Use PETSC_IGNORE to retain the previous value of any parameter.
375: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
376: dependent on the solution method.
378: The parameters ncv and mpd are intimately related, so that the user is advised
379: to set one of them at most. Normal usage is that
380: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
381: (b) in cases where nev is large, the user sets mpd.
383: The value of ncv should always be between nev and (nev+mpd), typically
384: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
385: a smaller value should be used.
387: When computing all eigenvalues in an interval, see EPSSetInterval(), the
388: meaning of nev changes. In that case, the number of eigenvalues in the
389: interval is not known a priori; the meaning of nev is then the number of
390: eigenvalues that are computed at a time when sweeping the interval from one
391: end to the other. The value of nev in this case may have an impact on overall
392: performance. The value of ncv should not be assigned in this case.
394: Level: intermediate
396: .seealso: EPSGetDimensions(), EPSSetInterval()
397: @*/
398: PetscErrorCode EPSSetDimensions(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)
399: {
405: if (nev != PETSC_IGNORE) {
406: if (nev<1) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
407: eps->nev = nev;
408: eps->setupcalled = 0;
409: }
410: if (ncv != PETSC_IGNORE) {
411: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
412: eps->ncv = 0;
413: } else {
414: if (ncv<1) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
415: eps->ncv = ncv;
416: }
417: eps->setupcalled = 0;
418: }
419: if (mpd != PETSC_IGNORE) {
420: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
421: eps->mpd = 0;
422: } else {
423: if (mpd<1) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
424: eps->mpd = mpd;
425: }
426: }
427: return(0);
428: }
432: /*@
433: EPSSetWhichEigenpairs - Specifies which portion of the spectrum is
434: to be sought.
436: Logically Collective on EPS
438: Input Parameters:
439: + eps - eigensolver context obtained from EPSCreate()
440: - which - the portion of the spectrum to be sought
442: Possible values:
443: The parameter 'which' can have one of these values
444:
445: + EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
446: . EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
447: . EPS_LARGEST_REAL - largest real parts
448: . EPS_SMALLEST_REAL - smallest real parts
449: . EPS_LARGEST_IMAGINARY - largest imaginary parts
450: . EPS_SMALLEST_IMAGINARY - smallest imaginary parts
451: . EPS_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
452: . EPS_TARGET_REAL - eigenvalues with real part closest to target
453: . EPS_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
454: . EPS_ALL - all eigenvalues contained in a given interval
455: - EPS_WHICH_USER - user defined ordering set with EPSSetEigenvalueComparison()
457: Options Database Keys:
458: + -eps_largest_magnitude - Sets largest eigenvalues in magnitude
459: . -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
460: . -eps_largest_real - Sets largest real parts
461: . -eps_smallest_real - Sets smallest real parts
462: . -eps_largest_imaginary - Sets largest imaginary parts
463: . -eps_smallest_imaginary - Sets smallest imaginary parts
464: . -eps_target_magnitude - Sets eigenvalues closest to target
465: . -eps_target_real - Sets real parts closest to target
466: . -eps_target_imaginary - Sets imaginary parts closest to target
467: - -eps_all - Sets all eigenvalues in an interval
469: Notes:
470: Not all eigensolvers implemented in EPS account for all the possible values
471: stated above. Also, some values make sense only for certain types of
472: problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
473: and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part
474: for eigenvalue selection.
475:
476: The target is a scalar value provided with EPSSetTarget().
478: The criterion EPS_TARGET_IMAGINARY is available only in case PETSc and
479: SLEPc have been built with complex scalars.
481: EPS_ALL is intended for use in combination with an interval (see
482: EPSSetInterval()), when all eigenvalues within the interval are requested.
483: In that case, the number of eigenvalues is unknown, so the nev parameter
484: has a different sense, see EPSSetDimensions().
486: Level: intermediate
488: .seealso: EPSGetWhichEigenpairs(), EPSSetTarget(), EPSSetInterval(),
489: EPSSetDimensions(), EPSSetEigenvalueComparison(),
490: EPSSortEigenvalues(), EPSWhich
491: @*/
492: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)
493: {
497: if (which!=PETSC_IGNORE) {
498: if (which==PETSC_DECIDE || which==PETSC_DEFAULT) eps->which = (EPSWhich)0;
499: else switch (which) {
500: case EPS_LARGEST_MAGNITUDE:
501: case EPS_SMALLEST_MAGNITUDE:
502: case EPS_LARGEST_REAL:
503: case EPS_SMALLEST_REAL:
504: case EPS_LARGEST_IMAGINARY:
505: case EPS_SMALLEST_IMAGINARY:
506: case EPS_TARGET_MAGNITUDE:
507: case EPS_TARGET_REAL:
508: #if defined(PETSC_USE_COMPLEX)
509: case EPS_TARGET_IMAGINARY:
510: #endif
511: case EPS_ALL:
512: case EPS_WHICH_USER:
513: if (eps->which != which) {
514: eps->setupcalled = 0;
515: eps->which = which;
516: }
517: break;
518: default:
519: SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
520: }
521: }
522: return(0);
523: }
527: /*@C
528: EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be
529: sought.
531: Not Collective
533: Input Parameter:
534: . eps - eigensolver context obtained from EPSCreate()
536: Output Parameter:
537: . which - the portion of the spectrum to be sought
539: Notes:
540: See EPSSetWhichEigenpairs() for possible values of 'which'.
542: Level: intermediate
544: .seealso: EPSSetWhichEigenpairs(), EPSWhich
545: @*/
546: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
547: {
551: *which = eps->which;
552: return(0);
553: }
557: /*@
558: EPSSetLeftVectorsWanted - Specifies which eigenvectors are required.
560: Logically Collective on EPS
562: Input Parameters:
563: + eps - the eigensolver context
564: - leftvecs - whether left eigenvectors are required or not
566: Options Database Keys:
567: . -eps_left_vectors <boolean> - Sets/resets the boolean flag 'leftvecs'
569: Notes:
570: If the user sets leftvecs=PETSC_TRUE then the solver uses a variant of
571: the algorithm that computes both right and left eigenvectors. This is
572: usually much more costly. This option is not available in all solvers.
574: Level: intermediate
576: .seealso: EPSGetLeftVectorsWanted(), EPSGetEigenvectorLeft()
577: @*/
578: PetscErrorCode EPSSetLeftVectorsWanted(EPS eps,PetscBool leftvecs)
579: {
583: if (eps->leftvecs != leftvecs) {
584: eps->leftvecs = leftvecs;
585: eps->setupcalled = 0;
586: }
587: return(0);
588: }
592: /*@
593: EPSGetLeftVectorsWanted - Returns the flag indicating whether left
594: eigenvectors are required or not.
596: Not Collective
598: Input Parameter:
599: . eps - the eigensolver context
601: Output Parameter:
602: . leftvecs - the returned flag
604: Level: intermediate
606: .seealso: EPSSetLeftVectorsWanted(), EPSWhich
607: @*/
608: PetscErrorCode EPSGetLeftVectorsWanted(EPS eps,PetscBool *leftvecs)
609: {
613: *leftvecs = eps->leftvecs;
614: return(0);
615: }
619: /*@
620: EPSSetMatrixNorms - Gives the reference values of the matrix norms
621: and specifies whether these values should be improved adaptively.
623: Logically Collective on EPS
625: Input Parameters:
626: + eps - the eigensolver context
627: . nrma - a reference value for the norm of matrix A
628: . nrmb - a reference value for the norm of matrix B
629: - adaptive - whether matrix norms are improved adaptively
631: Options Database Keys:
632: + -eps_norm_a <nrma> - norm of A
633: . -eps_norm_b <nrma> - norm of B
634: - -eps_norms_adaptive <boolean> - Sets/resets the boolean flag 'adaptive'
636: Notes:
637: If the user sets adaptive=PETSC_FALSE then the solver uses the values
638: of nrma and nrmb for the matrix norms, and these values do not change
639: throughout the iteration.
641: If the user sets adaptive=PETSC_TRUE then the solver tries to adaptively
642: improve the supplied values, with the numerical information generated
643: during the iteration. This option is not available in all solvers.
645: If a passed value is PETSC_DEFAULT, the corresponding norm will be set to 1.
646: If a passed value is PETSC_DETERMINE, the corresponding norm will be computed
647: as the NORM_INFINITY with MatNorm().
649: Level: intermediate
651: .seealso: EPSGetMatrixNorms()
652: @*/
653: PetscErrorCode EPSSetMatrixNorms(EPS eps,PetscReal nrma,PetscReal nrmb,PetscBool adaptive)
654: {
660: if (nrma != PETSC_IGNORE) {
661: if (nrma == PETSC_DEFAULT) eps->nrma = 1.0;
662: else if (nrma == PETSC_DETERMINE) {
663: eps->nrma = nrma;
664: eps->setupcalled = 0;
665: } else {
666: if (nrma < 0.0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nrma. Must be > 0");
667: eps->nrma = nrma;
668: }
669: }
670: if (nrmb != PETSC_IGNORE) {
671: if (!eps->isgeneralized) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_WRONG,"Norm of B only allowed in generalized problems");
672: if (nrmb == PETSC_DEFAULT) eps->nrmb = 1.0;
673: else if (nrmb == PETSC_DETERMINE) {
674: eps->nrmb = nrmb;
675: eps->setupcalled = 0;
676: } else {
677: if (nrmb < 0.0) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nrmb. Must be > 0");
678: eps->nrmb = nrmb;
679: }
680: }
681: if (eps->adaptive != adaptive) {
682: eps->adaptive = adaptive;
683: eps->setupcalled = 0;
684: if (adaptive) SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_SUP,"Sorry, adaptive norms are not implemented in this release");
685: }
686: return(0);
687: }
691: /*@
692: EPSGetMatrixNorms - Returns the value of the matrix norms (either set
693: by the user or estimated by the solver) and the flag indicating whether
694: the norms are being adaptively improved.
696: Not Collective
698: Input Parameter:
699: . eps - the eigensolver context
701: Output Parameters:
702: + nrma - the norm of matrix A
703: . nrmb - the norm of matrix B
704: - adaptive - whether matrix norms are improved adaptively
706: Level: intermediate
708: .seealso: EPSSetMatrixNorms()
709: @*/
710: PetscErrorCode EPSGetMatrixNorms(EPS eps,PetscReal *nrma,PetscReal *nrmb,PetscBool *adaptive)
711: {
714: if (nrma) *nrma = eps->nrma;
715: if (nrmb) *nrmb = eps->nrmb;
716: if (adaptive) *adaptive = eps->adaptive;
717: return(0);
718: }
722: /*@C
723: EPSSetEigenvalueComparison - Specifies the eigenvalue comparison function
724: when EPSSetWhichEigenpairs() is set to EPS_WHICH_USER.
726: Logically Collective on EPS
728: Input Parameters:
729: + eps - eigensolver context obtained from EPSCreate()
730: . func - a pointer to the comparison function
731: - ctx - a context pointer (the last parameter to the comparison function)
733: Calling Sequence of func:
734: $ func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
736: + ar - real part of the 1st eigenvalue
737: . ai - imaginary part of the 1st eigenvalue
738: . br - real part of the 2nd eigenvalue
739: . bi - imaginary part of the 2nd eigenvalue
740: . res - result of comparison
741: - ctx - optional context, as set by EPSSetEigenvalueComparison()
743: Note:
744: The returning parameter 'res' can be:
745: + negative - if the 1st eigenvalue is preferred to the 2st one
746: . zero - if both eigenvalues are equally preferred
747: - positive - if the 2st eigenvalue is preferred to the 1st one
749: Level: advanced
751: .seealso: EPSSetWhichEigenpairs(), EPSSortEigenvalues(), EPSWhich
752: @*/
753: PetscErrorCode EPSSetEigenvalueComparison(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)
754: {
757: eps->which_func = func;
758: eps->which_ctx = ctx;
759: eps->which = EPS_WHICH_USER;
760: return(0);
761: }
765: /*@C
766: EPSSetArbitrarySelection - Specifies a function intended to look for
767: eigenvalues according to an arbitrary selection criterion. This criterion
768: can be based on a computation involving the current eigenvector approximation.
770: Logically Collective on EPS
772: Input Parameters:
773: + eps - eigensolver context obtained from EPSCreate()
774: . func - a pointer to the evaluation function
775: - ctx - a context pointer (the last parameter to the evaluation function)
777: Calling Sequence of func:
778: $ func(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx)
780: + er - real part of the current eigenvalue approximation
781: . ei - imaginary part of the current eigenvalue approximation
782: . xr - real part of the current eigenvector approximation
783: . xi - imaginary part of the current eigenvector approximation
784: . rr - result of evaluation (real part)
785: . ri - result of evaluation (imaginary part)
786: - ctx - optional context, as set by EPSSetArbitrarySelection()
788: Note:
789: This provides a mechanism to select eigenpairs by evaluating a user-defined
790: function. When a function has been provided, the default selection based on
791: sorting the eigenvalues is replaced by the sorting of the results of this
792: function (with the same sorting criterion given in EPSSortEigenvalues()).
794: For instance, suppose you want to compute those eigenvectors that maximize
795: a certain computable expression. Then implement the computation using
796: the arguments xr and xi, and return the result in rr. Then set the standard
797: sorting by magnitude so that the eigenpair with largest value of rr is
798: selected.
800: The result of func is expressed as a complex number so that it is possible to
801: use the standard eigenvalue sorting functions, but normally only rr is used.
802: Set ri to zero unless it is meaningful in your application.
804: Level: advanced
806: .seealso: EPSSetWhichEigenpairs(), EPSSortEigenvalues()
807: @*/
808: PetscErrorCode EPSSetArbitrarySelection(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,Vec,Vec,PetscScalar*,PetscScalar*,void*),void* ctx)
809: {
812: eps->arbit_func = func;
813: eps->arbit_ctx = ctx;
814: eps->setupcalled = 0;
815: return(0);
816: }
820: /*@C
821: EPSSetConvergenceTestFunction - Sets a function to compute the error estimate
822: used in the convergence test.
824: Logically Collective on EPS
826: Input Parameters:
827: + eps - eigensolver context obtained from EPSCreate()
828: . func - a pointer to the convergence test function
829: - ctx - a context pointer (the last parameter to the convergence test function)
831: Calling Sequence of func:
832: $ func(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
834: + eps - eigensolver context obtained from EPSCreate()
835: . eigr - real part of the eigenvalue
836: . eigi - imaginary part of the eigenvalue
837: . res - residual norm associated to the eigenpair
838: . errest - (output) computed error estimate
839: - ctx - optional context, as set by EPSSetConvergenceTest()
841: Note:
842: If the error estimate returned by the convergence test function is less than
843: the tolerance, then the eigenvalue is accepted as converged.
845: Level: advanced
847: .seealso: EPSSetConvergenceTest(),EPSSetTolerances()
848: @*/
849: extern PetscErrorCode EPSSetConvergenceTestFunction(EPS eps,PetscErrorCode (*func)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx)
850: {
853: eps->conv_func = func;
854: eps->conv_ctx = ctx;
855: if (func == EPSConvergedEigRelative) eps->conv = EPS_CONV_EIG;
856: else if (func == EPSConvergedNormRelative) eps->conv = EPS_CONV_NORM;
857: else if (func == EPSConvergedAbsolute) eps->conv = EPS_CONV_ABS;
858: else eps->conv = EPS_CONV_USER;
859: return(0);
860: }
864: /*@
865: EPSSetConvergenceTest - Specifies how to compute the error estimate
866: used in the convergence test.
868: Logically Collective on EPS
870: Input Parameters:
871: + eps - eigensolver context obtained from EPSCreate()
872: - conv - the type of convergence test
874: Options Database Keys:
875: + -eps_conv_abs - Sets the absolute convergence test
876: . -eps_conv_eig - Sets the convergence test relative to the eigenvalue
877: - -eps_conv_norm - Sets the convergence test relative to the matrix norms
879: Note:
880: The parameter 'conv' can have one of these values
881: + EPS_CONV_ABS - absolute error ||r||
882: . EPS_CONV_EIG - error relative to the eigenvalue l, ||r||/|l|
883: . EPS_CONV_NORM - error relative to the matrix norms, ||r||/(||A||+|l|*||B||)
884: - EPS_CONV_USER - function set by EPSSetConvergenceTestFunction()
886: Level: intermediate
888: .seealso: EPSGetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSConv
889: @*/
890: PetscErrorCode EPSSetConvergenceTest(EPS eps,EPSConv conv)
891: {
895: switch(conv) {
896: case EPS_CONV_EIG: eps->conv_func = EPSConvergedEigRelative; break;
897: case EPS_CONV_NORM: eps->conv_func = EPSConvergedNormRelative; break;
898: case EPS_CONV_ABS: eps->conv_func = EPSConvergedAbsolute; break;
899: case EPS_CONV_USER: break;
900: default:
901: SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
902: }
903: eps->conv = conv;
904: return(0);
905: }
909: /*@
910: EPSGetConvergenceTest - Gets the method used to compute the error estimate
911: used in the convergence test.
913: Not Collective
915: Input Parameters:
916: . eps - eigensolver context obtained from EPSCreate()
918: Output Parameters:
919: . conv - the type of convergence test
921: Level: intermediate
923: .seealso: EPSSetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSConv
924: @*/
925: PetscErrorCode EPSGetConvergenceTest(EPS eps,EPSConv *conv)
926: {
930: *conv = eps->conv;
931: return(0);
932: }
936: /*@
937: EPSSetProblemType - Specifies the type of the eigenvalue problem.
939: Logically Collective on EPS
941: Input Parameters:
942: + eps - the eigensolver context
943: - type - a known type of eigenvalue problem
945: Options Database Keys:
946: + -eps_hermitian - Hermitian eigenvalue problem
947: . -eps_gen_hermitian - generalized Hermitian eigenvalue problem
948: . -eps_non_hermitian - non-Hermitian eigenvalue problem
949: . -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
950: - -eps_pos_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
951: with positive semi-definite B
952:
953: Notes:
954: Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
955: (EPS_NHEP), generalized Hermitian (EPS_GHEP), generalized non-Hermitian
956: (EPS_GNHEP), generalized non-Hermitian with positive semi-definite B
957: (EPS_PGNHEP), and generalized Hermitian-indefinite (EPS_GHIEP).
959: This function must be used to instruct SLEPc to exploit symmetry. If no
960: problem type is specified, by default a non-Hermitian problem is assumed
961: (either standard or generalized). If the user knows that the problem is
962: Hermitian (i.e. A=A^H) or generalized Hermitian (i.e. A=A^H, B=B^H, and
963: B positive definite) then it is recommended to set the problem type so
964: that eigensolver can exploit these properties.
966: Level: beginner
968: .seealso: EPSSetOperators(), EPSSetType(), EPSGetProblemType(), EPSProblemType
969: @*/
970: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
971: {
975: switch (type) {
976: case EPS_HEP:
977: eps->isgeneralized = PETSC_FALSE;
978: eps->ishermitian = PETSC_TRUE;
979: eps->ispositive = PETSC_FALSE;
980: break;
981: case EPS_NHEP:
982: eps->isgeneralized = PETSC_FALSE;
983: eps->ishermitian = PETSC_FALSE;
984: eps->ispositive = PETSC_FALSE;
985: break;
986: case EPS_GHEP:
987: eps->isgeneralized = PETSC_TRUE;
988: eps->ishermitian = PETSC_TRUE;
989: eps->ispositive = PETSC_TRUE;
990: break;
991: case EPS_GNHEP:
992: eps->isgeneralized = PETSC_TRUE;
993: eps->ishermitian = PETSC_FALSE;
994: eps->ispositive = PETSC_FALSE;
995: break;
996: case EPS_PGNHEP:
997: eps->isgeneralized = PETSC_TRUE;
998: eps->ishermitian = PETSC_FALSE;
999: eps->ispositive = PETSC_TRUE;
1000: break;
1001: case EPS_GHIEP:
1002: eps->isgeneralized = PETSC_TRUE;
1003: eps->ishermitian = PETSC_TRUE;
1004: eps->ispositive = PETSC_FALSE;
1005: break;
1006: default:
1007: SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
1008: }
1009: eps->problem_type = type;
1010: return(0);
1011: }
1015: /*@C
1016: EPSGetProblemType - Gets the problem type from the EPS object.
1018: Not Collective
1020: Input Parameter:
1021: . eps - the eigensolver context
1023: Output Parameter:
1024: . type - name of EPS problem type
1026: Level: intermediate
1028: .seealso: EPSSetProblemType(), EPSProblemType
1029: @*/
1030: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
1031: {
1035: *type = eps->problem_type;
1036: return(0);
1037: }
1041: /*@
1042: EPSSetExtraction - Specifies the type of extraction technique to be employed
1043: by the eigensolver.
1045: Logically Collective on EPS
1047: Input Parameters:
1048: + eps - the eigensolver context
1049: - extr - a known type of extraction
1051: Options Database Keys:
1052: + -eps_ritz - Rayleigh-Ritz extraction
1053: . -eps_harmonic - harmonic Ritz extraction
1054: . -eps_harmonic_relative - harmonic Ritz extraction relative to the eigenvalue
1055: . -eps_harmonic_right - harmonic Ritz extraction for rightmost eigenvalues
1056: . -eps_harmonic_largest - harmonic Ritz extraction for largest magnitude
1057: (without target)
1058: . -eps_refined - refined Ritz extraction
1059: - -eps_refined_harmonic - refined harmonic Ritz extraction
1060:
1061: Notes:
1062: Not all eigensolvers support all types of extraction. See the SLEPc
1063: Users Manual for details.
1065: By default, a standard Rayleigh-Ritz extraction is used. Other extractions
1066: may be useful when computing interior eigenvalues.
1068: Harmonic-type extractions are used in combination with a 'target'.
1070: Level: beginner
1072: .seealso: EPSSetTarget(), EPSGetExtraction(), EPSExtraction
1073: @*/
1074: PetscErrorCode EPSSetExtraction(EPS eps,EPSExtraction extr)
1075: {
1079: eps->extraction = extr;
1080: return(0);
1081: }
1085: /*@C
1086: EPSGetExtraction - Gets the extraction type used by the EPS object.
1088: Not Collective
1090: Input Parameter:
1091: . eps - the eigensolver context
1093: Output Parameter:
1094: . extr - name of extraction type
1096: Level: intermediate
1098: .seealso: EPSSetExtraction(), EPSExtraction
1099: @*/
1100: PetscErrorCode EPSGetExtraction(EPS eps,EPSExtraction *extr)
1101: {
1105: *extr = eps->extraction;
1106: return(0);
1107: }
1111: /*@
1112: EPSSetBalance - Specifies the balancing technique to be employed by the
1113: eigensolver, and some parameters associated to it.
1115: Logically Collective on EPS
1117: Input Parameters:
1118: + eps - the eigensolver context
1119: . bal - the balancing method, one of EPS_BALANCE_NONE, EPS_BALANCE_ONESIDE,
1120: EPS_BALANCE_TWOSIDE, or EPS_BALANCE_USER
1121: . its - number of iterations of the balancing algorithm
1122: - cutoff - cutoff value
1124: Options Database Keys:
1125: + -eps_balance <method> - the balancing method, where <method> is one of
1126: 'none', 'oneside', 'twoside', or 'user'
1127: . -eps_balance_its <its> - number of iterations
1128: - -eps_balance_cutoff <cutoff> - cutoff value
1129:
1130: Notes:
1131: When balancing is enabled, the solver works implicitly with matrix DAD^-1,
1132: where D is an appropriate diagonal matrix. This improves the accuracy of
1133: the computed results in some cases. See the SLEPc Users Manual for details.
1135: Balancing makes sense only for non-Hermitian problems when the required
1136: precision is high (i.e. a small tolerance such as 1e-15).
1138: By default, balancing is disabled. The two-sided method is much more
1139: effective than the one-sided counterpart, but it requires the system
1140: matrices to have the MatMultTranspose operation defined.
1142: The parameter 'its' is the number of iterations performed by the method. The
1143: cutoff value is used only in the two-side variant. Use PETSC_IGNORE for an
1144: argument that need not be changed. Use PETSC_DECIDE to assign a reasonably
1145: good value.
1147: User-defined balancing is allowed provided that the corresponding matrix
1148: is set via STSetBalanceMatrix.
1150: Level: intermediate
1152: .seealso: EPSGetBalance(), EPSBalance, STSetBalanceMatrix()
1153: @*/
1154: PetscErrorCode EPSSetBalance(EPS eps,EPSBalance bal,PetscInt its,PetscReal cutoff)
1155: {
1161: if (bal!=PETSC_IGNORE) {
1162: if (bal==PETSC_DECIDE || bal==PETSC_DEFAULT) eps->balance = EPS_BALANCE_TWOSIDE;
1163: else switch (bal) {
1164: case EPS_BALANCE_NONE:
1165: case EPS_BALANCE_ONESIDE:
1166: case EPS_BALANCE_TWOSIDE:
1167: case EPS_BALANCE_USER:
1168: eps->balance = bal;
1169: break;
1170: default:
1171: SETERRQ(((PetscObject)eps)->comm,PETSC_ERR_ARG_OUTOFRANGE,"Invalid value of argument 'bal'");
1172: }
1173: }
1174: if (its!=PETSC_IGNORE) {
1175: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) eps->balance_its = 5;
1176: else eps->balance_its = its;
1177: }
1178: if (cutoff!=PETSC_IGNORE) {
1179: if (cutoff==PETSC_DECIDE || cutoff==PETSC_DEFAULT) eps->balance_cutoff = 1e-8;
1180: else eps->balance_cutoff = cutoff;
1181: }
1182: return(0);
1183: }
1187: /*@
1188: EPSGetBalance - Gets the balancing type used by the EPS object, and the associated
1189: parameters.
1191: Not Collective
1193: Input Parameter:
1194: . eps - the eigensolver context
1196: Output Parameters:
1197: + bal - the balancing method
1198: . its - number of iterations of the balancing algorithm
1199: - cutoff - cutoff value
1201: Level: intermediate
1203: Note:
1204: The user can specify PETSC_NULL for any parameter that is not needed.
1206: .seealso: EPSSetBalance(), EPSBalance
1207: @*/
1208: PetscErrorCode EPSGetBalance(EPS eps,EPSBalance *bal,PetscInt *its,PetscReal *cutoff)
1209: {
1212: if (bal) *bal = eps->balance;
1213: if (its) *its = eps->balance_its;
1214: if (cutoff) *cutoff = eps->balance_cutoff;
1215: return(0);
1216: }
1220: /*@
1221: EPSSetTrueResidual - Specifies if the solver must compute the true residual
1222: explicitly or not.
1224: Logically Collective on EPS
1226: Input Parameters:
1227: + eps - the eigensolver context
1228: - trueres - whether true residuals are required or not
1230: Options Database Keys:
1231: . -eps_true_residual <boolean> - Sets/resets the boolean flag 'trueres'
1233: Notes:
1234: If the user sets trueres=PETSC_TRUE then the solver explicitly computes
1235: the true residual for each eigenpair approximation, and uses it for
1236: convergence testing. Computing the residual is usually an expensive
1237: operation. Some solvers (e.g., Krylov solvers) can avoid this computation
1238: by using a cheap estimate of the residual norm, but this may sometimes
1239: give inaccurate results (especially if a spectral transform is being
1240: used). On the contrary, preconditioned eigensolvers (e.g., Davidson solvers)
1241: do rely on computing the true residual, so this option is irrelevant for them.
1243: Level: intermediate
1245: .seealso: EPSGetTrueResidual()
1246: @*/
1247: PetscErrorCode EPSSetTrueResidual(EPS eps,PetscBool trueres)
1248: {
1252: eps->trueres = trueres;
1253: return(0);
1254: }
1258: /*@C
1259: EPSGetTrueResidual - Returns the flag indicating whether true
1260: residuals must be computed explicitly or not.
1262: Not Collective
1264: Input Parameter:
1265: . eps - the eigensolver context
1267: Output Parameter:
1268: . trueres - the returned flag
1270: Level: intermediate
1272: .seealso: EPSSetTrueResidual()
1273: @*/
1274: PetscErrorCode EPSGetTrueResidual(EPS eps,PetscBool *trueres)
1275: {
1279: *trueres = eps->trueres;
1280: return(0);
1281: }
1285: /*@
1286: EPSSetTrackAll - Specifies if the solver must compute the residual norm of all
1287: approximate eigenpairs or not.
1289: Logically Collective on EPS
1291: Input Parameters:
1292: + eps - the eigensolver context
1293: - trackall - whether to compute all residuals or not
1295: Notes:
1296: If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
1297: the residual norm for each eigenpair approximation. Computing the residual is
1298: usually an expensive operation and solvers commonly compute only the residual
1299: associated to the first unconverged eigenpair.
1301: The options '-eps_monitor_all' and '-eps_monitor_draw_all' automatically
1302: activate this option.
1304: Level: intermediate
1306: .seealso: EPSGetTrackAll()
1307: @*/
1308: PetscErrorCode EPSSetTrackAll(EPS eps,PetscBool trackall)
1309: {
1313: eps->trackall = trackall;
1314: return(0);
1315: }
1319: /*@
1320: EPSGetTrackAll - Returns the flag indicating whether all residual norms must
1321: be computed or not.
1323: Not Collective
1325: Input Parameter:
1326: . eps - the eigensolver context
1328: Output Parameter:
1329: . trackall - the returned flag
1331: Level: intermediate
1333: .seealso: EPSSetTrackAll()
1334: @*/
1335: PetscErrorCode EPSGetTrackAll(EPS eps,PetscBool *trackall)
1336: {
1340: *trackall = eps->trackall;
1341: return(0);
1342: }
1346: /*@C
1347: EPSSetOptionsPrefix - Sets the prefix used for searching for all
1348: EPS options in the database.
1350: Logically Collective on EPS
1352: Input Parameters:
1353: + eps - the eigensolver context
1354: - prefix - the prefix string to prepend to all EPS option requests
1356: Notes:
1357: A hyphen (-) must NOT be given at the beginning of the prefix name.
1358: The first character of all runtime options is AUTOMATICALLY the
1359: hyphen.
1361: For example, to distinguish between the runtime options for two
1362: different EPS contexts, one could call
1363: .vb
1364: EPSSetOptionsPrefix(eps1,"eig1_")
1365: EPSSetOptionsPrefix(eps2,"eig2_")
1366: .ve
1368: Level: advanced
1370: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
1371: @*/
1372: PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix)
1373: {
1378: if (!eps->OP) { EPSGetST(eps,&eps->OP); }
1379: STSetOptionsPrefix(eps->OP,prefix);
1380: if (!eps->ip) { EPSGetIP(eps,&eps->ip); }
1381: IPSetOptionsPrefix(eps->ip,prefix);
1382: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1383: DSSetOptionsPrefix(eps->ds,prefix);
1384: PetscObjectSetOptionsPrefix((PetscObject)eps,prefix);
1385: return(0);
1386: }
1387:
1390: /*@C
1391: EPSAppendOptionsPrefix - Appends to the prefix used for searching for all
1392: EPS options in the database.
1394: Logically Collective on EPS
1396: Input Parameters:
1397: + eps - the eigensolver context
1398: - prefix - the prefix string to prepend to all EPS option requests
1400: Notes:
1401: A hyphen (-) must NOT be given at the beginning of the prefix name.
1402: The first character of all runtime options is AUTOMATICALLY the hyphen.
1404: Level: advanced
1406: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
1407: @*/
1408: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix)
1409: {
1414: if (!eps->OP) { EPSGetST(eps,&eps->OP); }
1415: STAppendOptionsPrefix(eps->OP,prefix);
1416: if (!eps->ip) { EPSGetIP(eps,&eps->ip); }
1417: IPSetOptionsPrefix(eps->ip,prefix);
1418: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1419: DSSetOptionsPrefix(eps->ds,prefix);
1420: PetscObjectAppendOptionsPrefix((PetscObject)eps,prefix);
1421: return(0);
1422: }
1426: /*@C
1427: EPSGetOptionsPrefix - Gets the prefix used for searching for all
1428: EPS options in the database.
1430: Not Collective
1432: Input Parameters:
1433: . eps - the eigensolver context
1435: Output Parameters:
1436: . prefix - pointer to the prefix string used is returned
1438: Notes: On the fortran side, the user should pass in a string 'prefix' of
1439: sufficient length to hold the prefix.
1441: Level: advanced
1443: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
1444: @*/
1445: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])
1446: {
1452: PetscObjectGetOptionsPrefix((PetscObject)eps,prefix);
1453: return(0);
1454: }