1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2021, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: EPS routines related to options that can be set via the command-line
12: or procedurally.
13: */
15: #include <slepc/private/epsimpl.h> 16: #include <petscdraw.h>
18: /*@C
19: EPSMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
20: indicated by the user.
22: Collective on eps
24: Input Parameters:
25: + eps - the eigensolver context
26: . opt - the command line option for this monitor
27: . name - the monitor type one is seeking
28: . ctx - an optional user context for the monitor, or NULL
29: - trackall - whether this monitor tracks all eigenvalues or not
31: Level: developer
33: .seealso: EPSMonitorSet(), EPSSetTrackAll()
34: @*/
35: PetscErrorCode EPSMonitorSetFromOptions(EPS eps,const char opt[],const char name[],void *ctx,PetscBool trackall) 36: {
37: PetscErrorCode (*mfunc)(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
38: PetscErrorCode (*cfunc)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**);
39: PetscErrorCode (*dfunc)(PetscViewerAndFormat**);
40: PetscViewerAndFormat *vf;
41: PetscViewer viewer;
42: PetscViewerFormat format;
43: PetscViewerType vtype;
44: char key[PETSC_MAX_PATH_LEN];
45: PetscBool flg;
46: PetscErrorCode ierr;
49: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,opt,&viewer,&format,&flg);
50: if (!flg) return(0);
52: PetscViewerGetType(viewer,&vtype);
53: SlepcMonitorMakeKey_Internal(name,vtype,format,key);
54: PetscFunctionListFind(EPSMonitorList,key,&mfunc);
55: if (!mfunc) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Specified viewer and format not supported");
56: PetscFunctionListFind(EPSMonitorCreateList,key,&cfunc);
57: PetscFunctionListFind(EPSMonitorDestroyList,key,&dfunc);
58: if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
59: if (!dfunc) dfunc = PetscViewerAndFormatDestroy;
61: (*cfunc)(viewer,format,ctx,&vf);
62: PetscObjectDereference((PetscObject)viewer);
63: EPSMonitorSet(eps,mfunc,vf,(PetscErrorCode(*)(void **))dfunc);
64: if (trackall) {
65: EPSSetTrackAll(eps,PETSC_TRUE);
66: }
67: return(0);
68: }
70: /*@
71: EPSSetFromOptions - Sets EPS options from the options database.
72: This routine must be called before EPSSetUp() if the user is to be
73: allowed to set the solver type.
75: Collective on eps
77: Input Parameters:
78: . eps - the eigensolver context
80: Notes:
81: To see all options, run your program with the -help option.
83: Level: beginner
84: @*/
85: PetscErrorCode EPSSetFromOptions(EPS eps) 86: {
88: char type[256];
89: PetscBool set,flg,flg1,flg2,flg3,bval;
90: PetscReal r,array[2]={0,0};
91: PetscScalar s;
92: PetscInt i,j,k;
93: EPSBalance bal;
97: EPSRegisterAll();
98: PetscObjectOptionsBegin((PetscObject)eps);
99: PetscOptionsFList("-eps_type","Eigensolver method","EPSSetType",EPSList,(char*)(((PetscObject)eps)->type_name?((PetscObject)eps)->type_name:EPSKRYLOVSCHUR),type,sizeof(type),&flg);
100: if (flg) {
101: EPSSetType(eps,type);
102: } else if (!((PetscObject)eps)->type_name) {
103: EPSSetType(eps,EPSKRYLOVSCHUR);
104: }
106: PetscOptionsBoolGroupBegin("-eps_hermitian","Hermitian eigenvalue problem","EPSSetProblemType",&flg);
107: if (flg) { EPSSetProblemType(eps,EPS_HEP); }
108: PetscOptionsBoolGroup("-eps_gen_hermitian","Generalized Hermitian eigenvalue problem","EPSSetProblemType",&flg);
109: if (flg) { EPSSetProblemType(eps,EPS_GHEP); }
110: PetscOptionsBoolGroup("-eps_non_hermitian","Non-Hermitian eigenvalue problem","EPSSetProblemType",&flg);
111: if (flg) { EPSSetProblemType(eps,EPS_NHEP); }
112: PetscOptionsBoolGroup("-eps_gen_non_hermitian","Generalized non-Hermitian eigenvalue problem","EPSSetProblemType",&flg);
113: if (flg) { EPSSetProblemType(eps,EPS_GNHEP); }
114: PetscOptionsBoolGroup("-eps_pos_gen_non_hermitian","Generalized non-Hermitian eigenvalue problem with positive semi-definite B","EPSSetProblemType",&flg);
115: if (flg) { EPSSetProblemType(eps,EPS_PGNHEP); }
116: PetscOptionsBoolGroupEnd("-eps_gen_indefinite","Generalized Hermitian-indefinite eigenvalue problem","EPSSetProblemType",&flg);
117: if (flg) { EPSSetProblemType(eps,EPS_GHIEP); }
119: PetscOptionsBoolGroupBegin("-eps_ritz","Rayleigh-Ritz extraction","EPSSetExtraction",&flg);
120: if (flg) { EPSSetExtraction(eps,EPS_RITZ); }
121: PetscOptionsBoolGroup("-eps_harmonic","Harmonic Ritz extraction","EPSSetExtraction",&flg);
122: if (flg) { EPSSetExtraction(eps,EPS_HARMONIC); }
123: PetscOptionsBoolGroup("-eps_harmonic_relative","Relative harmonic Ritz extraction","EPSSetExtraction",&flg);
124: if (flg) { EPSSetExtraction(eps,EPS_HARMONIC_RELATIVE); }
125: PetscOptionsBoolGroup("-eps_harmonic_right","Right harmonic Ritz extraction","EPSSetExtraction",&flg);
126: if (flg) { EPSSetExtraction(eps,EPS_HARMONIC_RIGHT); }
127: PetscOptionsBoolGroup("-eps_harmonic_largest","Largest harmonic Ritz extraction","EPSSetExtraction",&flg);
128: if (flg) { EPSSetExtraction(eps,EPS_HARMONIC_LARGEST); }
129: PetscOptionsBoolGroup("-eps_refined","Refined Ritz extraction","EPSSetExtraction",&flg);
130: if (flg) { EPSSetExtraction(eps,EPS_REFINED); }
131: PetscOptionsBoolGroupEnd("-eps_refined_harmonic","Refined harmonic Ritz extraction","EPSSetExtraction",&flg);
132: if (flg) { EPSSetExtraction(eps,EPS_REFINED_HARMONIC); }
134: bal = eps->balance;
135: PetscOptionsEnum("-eps_balance","Balancing method","EPSSetBalance",EPSBalanceTypes,(PetscEnum)bal,(PetscEnum*)&bal,&flg1);
136: j = eps->balance_its;
137: PetscOptionsInt("-eps_balance_its","Number of iterations in balancing","EPSSetBalance",eps->balance_its,&j,&flg2);
138: r = eps->balance_cutoff;
139: PetscOptionsReal("-eps_balance_cutoff","Cutoff value in balancing","EPSSetBalance",eps->balance_cutoff,&r,&flg3);
140: if (flg1 || flg2 || flg3) { EPSSetBalance(eps,bal,j,r); }
142: i = eps->max_it;
143: PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,&flg1);
144: r = eps->tol;
145: PetscOptionsReal("-eps_tol","Tolerance","EPSSetTolerances",SlepcDefaultTol(eps->tol),&r,&flg2);
146: if (flg1 || flg2) { EPSSetTolerances(eps,r,i); }
148: PetscOptionsBoolGroupBegin("-eps_conv_rel","Relative error convergence test","EPSSetConvergenceTest",&flg);
149: if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_REL); }
150: PetscOptionsBoolGroup("-eps_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","EPSSetConvergenceTest",&flg);
151: if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_NORM); }
152: PetscOptionsBoolGroup("-eps_conv_abs","Absolute error convergence test","EPSSetConvergenceTest",&flg);
153: if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_ABS); }
154: PetscOptionsBoolGroupEnd("-eps_conv_user","User-defined convergence test","EPSSetConvergenceTest",&flg);
155: if (flg) { EPSSetConvergenceTest(eps,EPS_CONV_USER); }
157: PetscOptionsBoolGroupBegin("-eps_stop_basic","Stop iteration if all eigenvalues converged or max_it reached","EPSSetStoppingTest",&flg);
158: if (flg) { EPSSetStoppingTest(eps,EPS_STOP_BASIC); }
159: PetscOptionsBoolGroupEnd("-eps_stop_user","User-defined stopping test","EPSSetStoppingTest",&flg);
160: if (flg) { EPSSetStoppingTest(eps,EPS_STOP_USER); }
162: i = eps->nev;
163: PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,&flg1);
164: j = eps->ncv;
165: PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&j,&flg2);
166: k = eps->mpd;
167: PetscOptionsInt("-eps_mpd","Maximum dimension of projected problem","EPSSetDimensions",eps->mpd,&k,&flg3);
168: if (flg1 || flg2 || flg3) { EPSSetDimensions(eps,i,j,k); }
170: PetscOptionsBoolGroupBegin("-eps_largest_magnitude","Compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
171: if (flg) { EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE); }
172: PetscOptionsBoolGroup("-eps_smallest_magnitude","Compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
173: if (flg) { EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE); }
174: PetscOptionsBoolGroup("-eps_largest_real","Compute eigenvalues with largest real parts","EPSSetWhichEigenpairs",&flg);
175: if (flg) { EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL); }
176: PetscOptionsBoolGroup("-eps_smallest_real","Compute eigenvalues with smallest real parts","EPSSetWhichEigenpairs",&flg);
177: if (flg) { EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL); }
178: PetscOptionsBoolGroup("-eps_largest_imaginary","Compute eigenvalues with largest imaginary parts","EPSSetWhichEigenpairs",&flg);
179: if (flg) { EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY); }
180: PetscOptionsBoolGroup("-eps_smallest_imaginary","Compute eigenvalues with smallest imaginary parts","EPSSetWhichEigenpairs",&flg);
181: if (flg) { EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY); }
182: PetscOptionsBoolGroup("-eps_target_magnitude","Compute eigenvalues closest to target","EPSSetWhichEigenpairs",&flg);
183: if (flg) { EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE); }
184: PetscOptionsBoolGroup("-eps_target_real","Compute eigenvalues with real parts closest to target","EPSSetWhichEigenpairs",&flg);
185: if (flg) { EPSSetWhichEigenpairs(eps,EPS_TARGET_REAL); }
186: PetscOptionsBoolGroup("-eps_target_imaginary","Compute eigenvalues with imaginary parts closest to target","EPSSetWhichEigenpairs",&flg);
187: if (flg) { EPSSetWhichEigenpairs(eps,EPS_TARGET_IMAGINARY); }
188: PetscOptionsBoolGroupEnd("-eps_all","Compute all eigenvalues in an interval or a region","EPSSetWhichEigenpairs",&flg);
189: if (flg) { EPSSetWhichEigenpairs(eps,EPS_ALL); }
191: PetscOptionsScalar("-eps_target","Value of the target","EPSSetTarget",eps->target,&s,&flg);
192: if (flg) {
193: if (eps->which!=EPS_TARGET_REAL && eps->which!=EPS_TARGET_IMAGINARY) {
194: EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);
195: }
196: EPSSetTarget(eps,s);
197: }
199: k = 2;
200: PetscOptionsRealArray("-eps_interval","Computational interval (two real values separated with a comma without spaces)","EPSSetInterval",array,&k,&flg);
201: if (flg) {
202: if (k<2) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_SIZ,"Must pass two values in -eps_interval (comma-separated without spaces)");
203: EPSSetWhichEigenpairs(eps,EPS_ALL);
204: EPSSetInterval(eps,array[0],array[1]);
205: }
207: PetscOptionsBool("-eps_true_residual","Compute true residuals explicitly","EPSSetTrueResidual",eps->trueres,&eps->trueres,NULL);
208: PetscOptionsBool("-eps_purify","Postprocess eigenvectors for purification","EPSSetPurify",eps->purify,&bval,&flg);
209: if (flg) { EPSSetPurify(eps,bval); }
210: PetscOptionsBool("-eps_two_sided","Use two-sided variant (to compute left eigenvectors)","EPSSetTwoSided",eps->twosided,&bval,&flg);
211: if (flg) { EPSSetTwoSided(eps,bval); }
213: /* -----------------------------------------------------------------------*/
214: /*
215: Cancels all monitors hardwired into code before call to EPSSetFromOptions()
216: */
217: PetscOptionsBool("-eps_monitor_cancel","Remove any hardwired monitor routines","EPSMonitorCancel",PETSC_FALSE,&flg,&set);
218: if (set && flg) { EPSMonitorCancel(eps); }
219: EPSMonitorSetFromOptions(eps,"-eps_monitor","first_approximation",NULL,PETSC_FALSE);
220: EPSMonitorSetFromOptions(eps,"-eps_monitor_all","all_approximations",NULL,PETSC_TRUE);
221: EPSMonitorSetFromOptions(eps,"-eps_monitor_conv","convergence_history",NULL,PETSC_FALSE);
223: /* -----------------------------------------------------------------------*/
224: PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",NULL);
225: PetscOptionsName("-eps_view_vectors","View computed eigenvectors","EPSVectorsView",NULL);
226: PetscOptionsName("-eps_view_values","View computed eigenvalues","EPSValuesView",NULL);
227: PetscOptionsName("-eps_converged_reason","Print reason for convergence, and number of iterations","EPSConvergedReasonView",NULL);
228: PetscOptionsName("-eps_error_absolute","Print absolute errors of each eigenpair","EPSErrorView",NULL);
229: PetscOptionsName("-eps_error_relative","Print relative errors of each eigenpair","EPSErrorView",NULL);
230: PetscOptionsName("-eps_error_backward","Print backward errors of each eigenpair","EPSErrorView",NULL);
232: if (eps->ops->setfromoptions) {
233: (*eps->ops->setfromoptions)(PetscOptionsObject,eps);
234: }
235: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)eps);
236: PetscOptionsEnd();
238: if (!eps->V) { EPSGetBV(eps,&eps->V); }
239: BVSetFromOptions(eps->V);
240: if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
241: RGSetFromOptions(eps->rg);
242: if (eps->useds) {
243: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
244: DSSetFromOptions(eps->ds);
245: }
246: if (!eps->st) { EPSGetST(eps,&eps->st); }
247: EPSSetDefaultST(eps);
248: STSetFromOptions(eps->st);
249: return(0);
250: }
252: /*@C
253: EPSGetTolerances - Gets the tolerance and maximum iteration count used
254: by the EPS convergence tests.
256: Not Collective
258: Input Parameter:
259: . eps - the eigensolver context
261: Output Parameters:
262: + tol - the convergence tolerance
263: - maxits - maximum number of iterations
265: Notes:
266: The user can specify NULL for any parameter that is not needed.
268: Level: intermediate
270: .seealso: EPSSetTolerances()
271: @*/
272: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,PetscInt *maxits)273: {
276: if (tol) *tol = eps->tol;
277: if (maxits) *maxits = eps->max_it;
278: return(0);
279: }
281: /*@
282: EPSSetTolerances - Sets the tolerance and maximum iteration count used
283: by the EPS convergence tests.
285: Logically Collective on eps
287: Input Parameters:
288: + eps - the eigensolver context
289: . tol - the convergence tolerance
290: - maxits - maximum number of iterations to use
292: Options Database Keys:
293: + -eps_tol <tol> - Sets the convergence tolerance
294: - -eps_max_it <maxits> - Sets the maximum number of iterations allowed
296: Notes:
297: Use PETSC_DEFAULT for either argument to assign a reasonably good value.
299: Level: intermediate
301: .seealso: EPSGetTolerances()
302: @*/
303: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,PetscInt maxits)304: {
309: if (tol == PETSC_DEFAULT) {
310: eps->tol = PETSC_DEFAULT;
311: eps->state = EPS_STATE_INITIAL;
312: } else {
313: if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
314: eps->tol = tol;
315: }
316: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
317: eps->max_it = PETSC_DEFAULT;
318: eps->state = EPS_STATE_INITIAL;
319: } else {
320: if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
321: eps->max_it = maxits;
322: }
323: return(0);
324: }
326: /*@C
327: EPSGetDimensions - Gets the number of eigenvalues to compute
328: and the dimension of the subspace.
330: Not Collective
332: Input Parameter:
333: . eps - the eigensolver context
335: Output Parameters:
336: + nev - number of eigenvalues to compute
337: . ncv - the maximum dimension of the subspace to be used by the solver
338: - mpd - the maximum dimension allowed for the projected problem
340: Level: intermediate
342: .seealso: EPSSetDimensions()
343: @*/
344: PetscErrorCode EPSGetDimensions(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)345: {
348: if (nev) *nev = eps->nev;
349: if (ncv) *ncv = eps->ncv;
350: if (mpd) *mpd = eps->mpd;
351: return(0);
352: }
354: /*@
355: EPSSetDimensions - Sets the number of eigenvalues to compute
356: and the dimension of the subspace.
358: Logically Collective on eps
360: Input Parameters:
361: + eps - the eigensolver context
362: . nev - number of eigenvalues to compute
363: . ncv - the maximum dimension of the subspace to be used by the solver
364: - mpd - the maximum dimension allowed for the projected problem
366: Options Database Keys:
367: + -eps_nev <nev> - Sets the number of eigenvalues
368: . -eps_ncv <ncv> - Sets the dimension of the subspace
369: - -eps_mpd <mpd> - Sets the maximum projected dimension
371: Notes:
372: Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
373: dependent on the solution method.
375: The parameters ncv and mpd are intimately related, so that the user is advised
376: to set one of them at most. Normal usage is that
377: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
378: (b) in cases where nev is large, the user sets mpd.
380: The value of ncv should always be between nev and (nev+mpd), typically
381: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
382: a smaller value should be used.
384: When computing all eigenvalues in an interval, see EPSSetInterval(), these
385: parameters lose relevance, and tuning must be done with
386: EPSKrylovSchurSetDimensions().
388: Level: intermediate
390: .seealso: EPSGetDimensions(), EPSSetInterval(), EPSKrylovSchurSetDimensions()
391: @*/
392: PetscErrorCode EPSSetDimensions(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)393: {
399: if (nev<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
400: eps->nev = nev;
401: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
402: eps->ncv = PETSC_DEFAULT;
403: } else {
404: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
405: eps->ncv = ncv;
406: }
407: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
408: eps->mpd = PETSC_DEFAULT;
409: } else {
410: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
411: eps->mpd = mpd;
412: }
413: eps->state = EPS_STATE_INITIAL;
414: return(0);
415: }
417: /*@
418: EPSSetWhichEigenpairs - Specifies which portion of the spectrum is
419: to be sought.
421: Logically Collective on eps
423: Input Parameters:
424: + eps - eigensolver context obtained from EPSCreate()
425: - which - the portion of the spectrum to be sought
427: Possible values:
428: The parameter 'which' can have one of these values
430: + EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
431: . EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
432: . EPS_LARGEST_REAL - largest real parts
433: . EPS_SMALLEST_REAL - smallest real parts
434: . EPS_LARGEST_IMAGINARY - largest imaginary parts
435: . EPS_SMALLEST_IMAGINARY - smallest imaginary parts
436: . EPS_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
437: . EPS_TARGET_REAL - eigenvalues with real part closest to target
438: . EPS_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
439: . EPS_ALL - all eigenvalues contained in a given interval or region
440: - EPS_WHICH_USER - user defined ordering set with EPSSetEigenvalueComparison()
442: Options Database Keys:
443: + -eps_largest_magnitude - Sets largest eigenvalues in magnitude
444: . -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
445: . -eps_largest_real - Sets largest real parts
446: . -eps_smallest_real - Sets smallest real parts
447: . -eps_largest_imaginary - Sets largest imaginary parts
448: . -eps_smallest_imaginary - Sets smallest imaginary parts
449: . -eps_target_magnitude - Sets eigenvalues closest to target
450: . -eps_target_real - Sets real parts closest to target
451: . -eps_target_imaginary - Sets imaginary parts closest to target
452: - -eps_all - Sets all eigenvalues in an interval or region
454: Notes:
455: Not all eigensolvers implemented in EPS account for all the possible values
456: stated above. Also, some values make sense only for certain types of
457: problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY458: and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part
459: for eigenvalue selection.
461: The target is a scalar value provided with EPSSetTarget().
463: The criterion EPS_TARGET_IMAGINARY is available only in case PETSc and
464: SLEPc have been built with complex scalars.
466: EPS_ALL is intended for use in combination with an interval (see
467: EPSSetInterval()), when all eigenvalues within the interval are requested,
468: or in the context of the CISS solver for computing all eigenvalues in a region.
469: In those cases, the number of eigenvalues is unknown, so the nev parameter
470: has a different sense, see EPSSetDimensions().
472: Level: intermediate
474: .seealso: EPSGetWhichEigenpairs(), EPSSetTarget(), EPSSetInterval(),
475: EPSSetDimensions(), EPSSetEigenvalueComparison(), EPSWhich476: @*/
477: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)478: {
482: switch (which) {
483: case EPS_LARGEST_MAGNITUDE:
484: case EPS_SMALLEST_MAGNITUDE:
485: case EPS_LARGEST_REAL:
486: case EPS_SMALLEST_REAL:
487: case EPS_LARGEST_IMAGINARY:
488: case EPS_SMALLEST_IMAGINARY:
489: case EPS_TARGET_MAGNITUDE:
490: case EPS_TARGET_REAL:
491: #if defined(PETSC_USE_COMPLEX)
492: case EPS_TARGET_IMAGINARY:
493: #endif
494: case EPS_ALL:
495: case EPS_WHICH_USER:
496: if (eps->which != which) {
497: eps->state = EPS_STATE_INITIAL;
498: eps->which = which;
499: }
500: break;
501: #if !defined(PETSC_USE_COMPLEX)
502: case EPS_TARGET_IMAGINARY:
503: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"EPS_TARGET_IMAGINARY can be used only with complex scalars");
504: #endif
505: default:506: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
507: }
508: return(0);
509: }
511: /*@
512: EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be
513: sought.
515: Not Collective
517: Input Parameter:
518: . eps - eigensolver context obtained from EPSCreate()
520: Output Parameter:
521: . which - the portion of the spectrum to be sought
523: Notes:
524: See EPSSetWhichEigenpairs() for possible values of 'which'.
526: Level: intermediate
528: .seealso: EPSSetWhichEigenpairs(), EPSWhich529: @*/
530: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)531: {
535: *which = eps->which;
536: return(0);
537: }
539: /*@C
540: EPSSetEigenvalueComparison - Specifies the eigenvalue comparison function
541: when EPSSetWhichEigenpairs() is set to EPS_WHICH_USER.
543: Logically Collective on eps
545: Input Parameters:
546: + eps - eigensolver context obtained from EPSCreate()
547: . func - a pointer to the comparison function
548: - ctx - a context pointer (the last parameter to the comparison function)
550: Calling Sequence of func:
551: $ func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
553: + ar - real part of the 1st eigenvalue
554: . ai - imaginary part of the 1st eigenvalue
555: . br - real part of the 2nd eigenvalue
556: . bi - imaginary part of the 2nd eigenvalue
557: . res - result of comparison
558: - ctx - optional context, as set by EPSSetEigenvalueComparison()
560: Note:
561: The returning parameter 'res' can be
562: + negative - if the 1st eigenvalue is preferred to the 2st one
563: . zero - if both eigenvalues are equally preferred
564: - positive - if the 2st eigenvalue is preferred to the 1st one
566: Level: advanced
568: .seealso: EPSSetWhichEigenpairs(), EPSWhich569: @*/
570: PetscErrorCode EPSSetEigenvalueComparison(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)571: {
574: eps->sc->comparison = func;
575: eps->sc->comparisonctx = ctx;
576: eps->which = EPS_WHICH_USER;
577: return(0);
578: }
580: /*@C
581: EPSSetArbitrarySelection - Specifies a function intended to look for
582: eigenvalues according to an arbitrary selection criterion. This criterion
583: can be based on a computation involving the current eigenvector approximation.
585: Logically Collective on eps
587: Input Parameters:
588: + eps - eigensolver context obtained from EPSCreate()
589: . func - a pointer to the evaluation function
590: - ctx - a context pointer (the last parameter to the evaluation function)
592: Calling Sequence of func:
593: $ func(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx)
595: + er - real part of the current eigenvalue approximation
596: . ei - imaginary part of the current eigenvalue approximation
597: . xr - real part of the current eigenvector approximation
598: . xi - imaginary part of the current eigenvector approximation
599: . rr - result of evaluation (real part)
600: . ri - result of evaluation (imaginary part)
601: - ctx - optional context, as set by EPSSetArbitrarySelection()
603: Notes:
604: This provides a mechanism to select eigenpairs by evaluating a user-defined
605: function. When a function has been provided, the default selection based on
606: sorting the eigenvalues is replaced by the sorting of the results of this
607: function (with the same sorting criterion given in EPSSetWhichEigenpairs()).
609: For instance, suppose you want to compute those eigenvectors that maximize
610: a certain computable expression. Then implement the computation using
611: the arguments xr and xi, and return the result in rr. Then set the standard
612: sorting by magnitude so that the eigenpair with largest value of rr is
613: selected.
615: This evaluation function is collective, that is, all processes call it and
616: it can use collective operations; furthermore, the computed result must
617: be the same in all processes.
619: The result of func is expressed as a complex number so that it is possible to
620: use the standard eigenvalue sorting functions, but normally only rr is used.
621: Set ri to zero unless it is meaningful in your application.
623: Level: advanced
625: .seealso: EPSSetWhichEigenpairs()
626: @*/
627: PetscErrorCode EPSSetArbitrarySelection(EPS eps,PetscErrorCode (*func)(PetscScalar,PetscScalar,Vec,Vec,PetscScalar*,PetscScalar*,void*),void* ctx)628: {
631: eps->arbitrary = func;
632: eps->arbitraryctx = ctx;
633: eps->state = EPS_STATE_INITIAL;
634: return(0);
635: }
637: /*@C
638: EPSSetConvergenceTestFunction - Sets a function to compute the error estimate
639: used in the convergence test.
641: Logically Collective on eps
643: Input Parameters:
644: + eps - eigensolver context obtained from EPSCreate()
645: . func - a pointer to the convergence test function
646: . ctx - context for private data for the convergence routine (may be null)
647: - destroy - a routine for destroying the context (may be null)
649: Calling Sequence of func:
650: $ func(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
652: + eps - eigensolver context obtained from EPSCreate()
653: . eigr - real part of the eigenvalue
654: . eigi - imaginary part of the eigenvalue
655: . res - residual norm associated to the eigenpair
656: . errest - (output) computed error estimate
657: - ctx - optional context, as set by EPSSetConvergenceTestFunction()
659: Note:
660: If the error estimate returned by the convergence test function is less than
661: the tolerance, then the eigenvalue is accepted as converged.
663: Level: advanced
665: .seealso: EPSSetConvergenceTest(), EPSSetTolerances()
666: @*/
667: PetscErrorCode EPSSetConvergenceTestFunction(EPS eps,PetscErrorCode (*func)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))668: {
673: if (eps->convergeddestroy) {
674: (*eps->convergeddestroy)(eps->convergedctx);
675: }
676: eps->convergeduser = func;
677: eps->convergeddestroy = destroy;
678: eps->convergedctx = ctx;
679: if (func == EPSConvergedRelative) eps->conv = EPS_CONV_REL;
680: else if (func == EPSConvergedNorm) eps->conv = EPS_CONV_NORM;
681: else if (func == EPSConvergedAbsolute) eps->conv = EPS_CONV_ABS;
682: else {
683: eps->conv = EPS_CONV_USER;
684: eps->converged = eps->convergeduser;
685: }
686: return(0);
687: }
689: /*@
690: EPSSetConvergenceTest - Specifies how to compute the error estimate
691: used in the convergence test.
693: Logically Collective on eps
695: Input Parameters:
696: + eps - eigensolver context obtained from EPSCreate()
697: - conv - the type of convergence test
699: Options Database Keys:
700: + -eps_conv_abs - Sets the absolute convergence test
701: . -eps_conv_rel - Sets the convergence test relative to the eigenvalue
702: . -eps_conv_norm - Sets the convergence test relative to the matrix norms
703: - -eps_conv_user - Selects the user-defined convergence test
705: Note:
706: The parameter 'conv' can have one of these values
707: + EPS_CONV_ABS - absolute error ||r||
708: . EPS_CONV_REL - error relative to the eigenvalue l, ||r||/|l|
709: . EPS_CONV_NORM - error relative to the matrix norms, ||r||/(||A||+|l|*||B||)
710: - EPS_CONV_USER - function set by EPSSetConvergenceTestFunction()
712: Level: intermediate
714: .seealso: EPSGetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSSetStoppingTest(), EPSConv715: @*/
716: PetscErrorCode EPSSetConvergenceTest(EPS eps,EPSConv conv)717: {
721: switch (conv) {
722: case EPS_CONV_ABS: eps->converged = EPSConvergedAbsolute; break;
723: case EPS_CONV_REL: eps->converged = EPSConvergedRelative; break;
724: case EPS_CONV_NORM: eps->converged = EPSConvergedNorm; break;
725: case EPS_CONV_USER:
726: if (!eps->convergeduser) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ORDER,"Must call EPSSetConvergenceTestFunction() first");
727: eps->converged = eps->convergeduser;
728: break;
729: default:730: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
731: }
732: eps->conv = conv;
733: return(0);
734: }
736: /*@
737: EPSGetConvergenceTest - Gets the method used to compute the error estimate
738: used in the convergence test.
740: Not Collective
742: Input Parameters:
743: . eps - eigensolver context obtained from EPSCreate()
745: Output Parameters:
746: . conv - the type of convergence test
748: Level: intermediate
750: .seealso: EPSSetConvergenceTest(), EPSConv751: @*/
752: PetscErrorCode EPSGetConvergenceTest(EPS eps,EPSConv *conv)753: {
757: *conv = eps->conv;
758: return(0);
759: }
761: /*@C
762: EPSSetStoppingTestFunction - Sets a function to decide when to stop the outer
763: iteration of the eigensolver.
765: Logically Collective on eps
767: Input Parameters:
768: + eps - eigensolver context obtained from EPSCreate()
769: . func - pointer to the stopping test function
770: . ctx - context for private data for the stopping routine (may be null)
771: - destroy - a routine for destroying the context (may be null)
773: Calling Sequence of func:
774: $ func(EPS eps,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,EPSConvergedReason *reason,void *ctx)
776: + eps - eigensolver context obtained from EPSCreate()
777: . its - current number of iterations
778: . max_it - maximum number of iterations
779: . nconv - number of currently converged eigenpairs
780: . nev - number of requested eigenpairs
781: . reason - (output) result of the stopping test
782: - ctx - optional context, as set by EPSSetStoppingTestFunction()
784: Note:
785: Normal usage is to first call the default routine EPSStoppingBasic() and then
786: set reason to EPS_CONVERGED_USER if some user-defined conditions have been
787: met. To let the eigensolver continue iterating, the result must be left as
788: EPS_CONVERGED_ITERATING.
790: Level: advanced
792: .seealso: EPSSetStoppingTest(), EPSStoppingBasic()
793: @*/
794: PetscErrorCode EPSSetStoppingTestFunction(EPS eps,PetscErrorCode (*func)(EPS,PetscInt,PetscInt,PetscInt,PetscInt,EPSConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))795: {
800: if (eps->stoppingdestroy) {
801: (*eps->stoppingdestroy)(eps->stoppingctx);
802: }
803: eps->stoppinguser = func;
804: eps->stoppingdestroy = destroy;
805: eps->stoppingctx = ctx;
806: if (func == EPSStoppingBasic) eps->stop = EPS_STOP_BASIC;
807: else {
808: eps->stop = EPS_STOP_USER;
809: eps->stopping = eps->stoppinguser;
810: }
811: return(0);
812: }
814: /*@
815: EPSSetStoppingTest - Specifies how to decide the termination of the outer
816: loop of the eigensolver.
818: Logically Collective on eps
820: Input Parameters:
821: + eps - eigensolver context obtained from EPSCreate()
822: - stop - the type of stopping test
824: Options Database Keys:
825: + -eps_stop_basic - Sets the default stopping test
826: - -eps_stop_user - Selects the user-defined stopping test
828: Note:
829: The parameter 'stop' can have one of these values
830: + EPS_STOP_BASIC - default stopping test
831: - EPS_STOP_USER - function set by EPSSetStoppingTestFunction()
833: Level: advanced
835: .seealso: EPSGetStoppingTest(), EPSSetStoppingTestFunction(), EPSSetConvergenceTest(), EPSStop836: @*/
837: PetscErrorCode EPSSetStoppingTest(EPS eps,EPSStop stop)838: {
842: switch (stop) {
843: case EPS_STOP_BASIC: eps->stopping = EPSStoppingBasic; break;
844: case EPS_STOP_USER:
845: if (!eps->stoppinguser) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ORDER,"Must call EPSSetStoppingTestFunction() first");
846: eps->stopping = eps->stoppinguser;
847: break;
848: default:849: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
850: }
851: eps->stop = stop;
852: return(0);
853: }
855: /*@
856: EPSGetStoppingTest - Gets the method used to decide the termination of the outer
857: loop of the eigensolver.
859: Not Collective
861: Input Parameters:
862: . eps - eigensolver context obtained from EPSCreate()
864: Output Parameters:
865: . stop - the type of stopping test
867: Level: advanced
869: .seealso: EPSSetStoppingTest(), EPSStop870: @*/
871: PetscErrorCode EPSGetStoppingTest(EPS eps,EPSStop *stop)872: {
876: *stop = eps->stop;
877: return(0);
878: }
880: /*@
881: EPSSetProblemType - Specifies the type of the eigenvalue problem.
883: Logically Collective on eps
885: Input Parameters:
886: + eps - the eigensolver context
887: - type - a known type of eigenvalue problem
889: Options Database Keys:
890: + -eps_hermitian - Hermitian eigenvalue problem
891: . -eps_gen_hermitian - generalized Hermitian eigenvalue problem
892: . -eps_non_hermitian - non-Hermitian eigenvalue problem
893: . -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
894: . -eps_pos_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
895: with positive semi-definite B
896: - -eps_gen_indefinite - generalized Hermitian-indefinite eigenvalue problem
898: Notes:
899: Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
900: (EPS_NHEP), generalized Hermitian (EPS_GHEP), generalized non-Hermitian
901: (EPS_GNHEP), generalized non-Hermitian with positive semi-definite B
902: (EPS_PGNHEP), and generalized Hermitian-indefinite (EPS_GHIEP).
904: This function must be used to instruct SLEPc to exploit symmetry. If no
905: problem type is specified, by default a non-Hermitian problem is assumed
906: (either standard or generalized). If the user knows that the problem is
907: Hermitian (i.e. A=A^H) or generalized Hermitian (i.e. A=A^H, B=B^H, and
908: B positive definite) then it is recommended to set the problem type so
909: that eigensolver can exploit these properties.
911: Level: intermediate
913: .seealso: EPSSetOperators(), EPSSetType(), EPSGetProblemType(), EPSProblemType914: @*/
915: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)916: {
920: if (type == eps->problem_type) return(0);
921: switch (type) {
922: case EPS_HEP:
923: eps->isgeneralized = PETSC_FALSE;
924: eps->ishermitian = PETSC_TRUE;
925: eps->ispositive = PETSC_FALSE;
926: break;
927: case EPS_NHEP:
928: eps->isgeneralized = PETSC_FALSE;
929: eps->ishermitian = PETSC_FALSE;
930: eps->ispositive = PETSC_FALSE;
931: break;
932: case EPS_GHEP:
933: eps->isgeneralized = PETSC_TRUE;
934: eps->ishermitian = PETSC_TRUE;
935: eps->ispositive = PETSC_TRUE;
936: break;
937: case EPS_GNHEP:
938: eps->isgeneralized = PETSC_TRUE;
939: eps->ishermitian = PETSC_FALSE;
940: eps->ispositive = PETSC_FALSE;
941: break;
942: case EPS_PGNHEP:
943: eps->isgeneralized = PETSC_TRUE;
944: eps->ishermitian = PETSC_FALSE;
945: eps->ispositive = PETSC_TRUE;
946: break;
947: case EPS_GHIEP:
948: eps->isgeneralized = PETSC_TRUE;
949: eps->ishermitian = PETSC_TRUE;
950: eps->ispositive = PETSC_FALSE;
951: break;
952: default:953: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
954: }
955: eps->problem_type = type;
956: eps->state = EPS_STATE_INITIAL;
957: return(0);
958: }
960: /*@
961: EPSGetProblemType - Gets the problem type from the EPS object.
963: Not Collective
965: Input Parameter:
966: . eps - the eigensolver context
968: Output Parameter:
969: . type - the problem type
971: Level: intermediate
973: .seealso: EPSSetProblemType(), EPSProblemType974: @*/
975: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)976: {
980: *type = eps->problem_type;
981: return(0);
982: }
984: /*@
985: EPSSetExtraction - Specifies the type of extraction technique to be employed
986: by the eigensolver.
988: Logically Collective on eps
990: Input Parameters:
991: + eps - the eigensolver context
992: - extr - a known type of extraction
994: Options Database Keys:
995: + -eps_ritz - Rayleigh-Ritz extraction
996: . -eps_harmonic - harmonic Ritz extraction
997: . -eps_harmonic_relative - harmonic Ritz extraction relative to the eigenvalue
998: . -eps_harmonic_right - harmonic Ritz extraction for rightmost eigenvalues
999: . -eps_harmonic_largest - harmonic Ritz extraction for largest magnitude
1000: (without target)
1001: . -eps_refined - refined Ritz extraction
1002: - -eps_refined_harmonic - refined harmonic Ritz extraction
1004: Notes:
1005: Not all eigensolvers support all types of extraction. See the SLEPc
1006: Users Manual for details.
1008: By default, a standard Rayleigh-Ritz extraction is used. Other extractions
1009: may be useful when computing interior eigenvalues.
1011: Harmonic-type extractions are used in combination with a 'target'.
1013: Level: advanced
1015: .seealso: EPSSetTarget(), EPSGetExtraction(), EPSExtraction1016: @*/
1017: PetscErrorCode EPSSetExtraction(EPS eps,EPSExtraction extr)1018: {
1022: eps->extraction = extr;
1023: return(0);
1024: }
1026: /*@
1027: EPSGetExtraction - Gets the extraction type used by the EPS object.
1029: Not Collective
1031: Input Parameter:
1032: . eps - the eigensolver context
1034: Output Parameter:
1035: . extr - name of extraction type
1037: Level: advanced
1039: .seealso: EPSSetExtraction(), EPSExtraction1040: @*/
1041: PetscErrorCode EPSGetExtraction(EPS eps,EPSExtraction *extr)1042: {
1046: *extr = eps->extraction;
1047: return(0);
1048: }
1050: /*@
1051: EPSSetBalance - Specifies the balancing technique to be employed by the
1052: eigensolver, and some parameters associated to it.
1054: Logically Collective on eps
1056: Input Parameters:
1057: + eps - the eigensolver context
1058: . bal - the balancing method, one of EPS_BALANCE_NONE, EPS_BALANCE_ONESIDE,
1059: EPS_BALANCE_TWOSIDE, or EPS_BALANCE_USER1060: . its - number of iterations of the balancing algorithm
1061: - cutoff - cutoff value
1063: Options Database Keys:
1064: + -eps_balance <method> - the balancing method, where <method> is one of
1065: 'none', 'oneside', 'twoside', or 'user'
1066: . -eps_balance_its <its> - number of iterations
1067: - -eps_balance_cutoff <cutoff> - cutoff value
1069: Notes:
1070: When balancing is enabled, the solver works implicitly with matrix DAD^-1,
1071: where D is an appropriate diagonal matrix. This improves the accuracy of
1072: the computed results in some cases. See the SLEPc Users Manual for details.
1074: Balancing makes sense only for non-Hermitian problems when the required
1075: precision is high (i.e. a small tolerance such as 1e-15).
1077: By default, balancing is disabled. The two-sided method is much more
1078: effective than the one-sided counterpart, but it requires the system
1079: matrices to have the MatMultTranspose operation defined.
1081: The parameter 'its' is the number of iterations performed by the method. The
1082: cutoff value is used only in the two-side variant. Use PETSC_DEFAULT to assign
1083: a reasonably good value.
1085: User-defined balancing is allowed provided that the corresponding matrix
1086: is set via STSetBalanceMatrix.
1088: Level: intermediate
1090: .seealso: EPSGetBalance(), EPSBalance, STSetBalanceMatrix()
1091: @*/
1092: PetscErrorCode EPSSetBalance(EPS eps,EPSBalance bal,PetscInt its,PetscReal cutoff)1093: {
1099: switch (bal) {
1100: case EPS_BALANCE_NONE:
1101: case EPS_BALANCE_ONESIDE:
1102: case EPS_BALANCE_TWOSIDE:
1103: case EPS_BALANCE_USER:
1104: if (eps->balance != bal) {
1105: eps->state = EPS_STATE_INITIAL;
1106: eps->balance = bal;
1107: }
1108: break;
1109: default:1110: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid value of argument 'bal'");
1111: }
1112: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) eps->balance_its = 5;
1113: else {
1114: if (its <= 0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be > 0");
1115: eps->balance_its = its;
1116: }
1117: if (cutoff==PETSC_DECIDE || cutoff==PETSC_DEFAULT) eps->balance_cutoff = 1e-8;
1118: else {
1119: if (cutoff <= 0.0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of cutoff. Must be > 0");
1120: eps->balance_cutoff = cutoff;
1121: }
1122: return(0);
1123: }
1125: /*@C
1126: EPSGetBalance - Gets the balancing type used by the EPS object, and the
1127: associated parameters.
1129: Not Collective
1131: Input Parameter:
1132: . eps - the eigensolver context
1134: Output Parameters:
1135: + bal - the balancing method
1136: . its - number of iterations of the balancing algorithm
1137: - cutoff - cutoff value
1139: Level: intermediate
1141: Note:
1142: The user can specify NULL for any parameter that is not needed.
1144: .seealso: EPSSetBalance(), EPSBalance1145: @*/
1146: PetscErrorCode EPSGetBalance(EPS eps,EPSBalance *bal,PetscInt *its,PetscReal *cutoff)1147: {
1150: if (bal) *bal = eps->balance;
1151: if (its) *its = eps->balance_its;
1152: if (cutoff) *cutoff = eps->balance_cutoff;
1153: return(0);
1154: }
1156: /*@
1157: EPSSetTwoSided - Sets the solver to use a two-sided variant so that left
1158: eigenvectors are also computed.
1160: Logically Collective on eps
1162: Input Parameters:
1163: + eps - the eigensolver context
1164: - twosided - whether the two-sided variant is to be used or not
1166: Options Database Keys:
1167: . -eps_two_sided <boolean> - Sets/resets the twosided flag
1169: Notes:
1170: If the user sets twosided=PETSC_TRUE then the solver uses a variant of
1171: the algorithm that computes both right and left eigenvectors. This is
1172: usually much more costly. This option is not available in all solvers.
1174: When using two-sided solvers, the problem matrices must have both the
1175: MatMult and MatMultTranspose operations defined.
1177: Level: advanced
1179: .seealso: EPSGetTwoSided(), EPSGetLeftEigenvector()
1180: @*/
1181: PetscErrorCode EPSSetTwoSided(EPS eps,PetscBool twosided)1182: {
1186: if (twosided!=eps->twosided) {
1187: eps->twosided = twosided;
1188: eps->state = EPS_STATE_INITIAL;
1189: }
1190: return(0);
1191: }
1193: /*@
1194: EPSGetTwoSided - Returns the flag indicating whether a two-sided variant
1195: of the algorithm is being used or not.
1197: Not Collective
1199: Input Parameter:
1200: . eps - the eigensolver context
1202: Output Parameter:
1203: . twosided - the returned flag
1205: Level: advanced
1207: .seealso: EPSSetTwoSided()
1208: @*/
1209: PetscErrorCode EPSGetTwoSided(EPS eps,PetscBool *twosided)1210: {
1214: *twosided = eps->twosided;
1215: return(0);
1216: }
1218: /*@
1219: EPSSetTrueResidual - Specifies if the solver must compute the true residual
1220: explicitly or not.
1222: Logically Collective on eps
1224: Input Parameters:
1225: + eps - the eigensolver context
1226: - trueres - whether true residuals are required or not
1228: Options Database Keys:
1229: . -eps_true_residual <boolean> - Sets/resets the boolean flag 'trueres'
1231: Notes:
1232: If the user sets trueres=PETSC_TRUE then the solver explicitly computes
1233: the true residual for each eigenpair approximation, and uses it for
1234: convergence testing. Computing the residual is usually an expensive
1235: operation. Some solvers (e.g., Krylov solvers) can avoid this computation
1236: by using a cheap estimate of the residual norm, but this may sometimes
1237: give inaccurate results (especially if a spectral transform is being
1238: used). On the contrary, preconditioned eigensolvers (e.g., Davidson solvers)
1239: do rely on computing the true residual, so this option is irrelevant for them.
1241: Level: advanced
1243: .seealso: EPSGetTrueResidual()
1244: @*/
1245: PetscErrorCode EPSSetTrueResidual(EPS eps,PetscBool trueres)1246: {
1250: eps->trueres = trueres;
1251: return(0);
1252: }
1254: /*@
1255: EPSGetTrueResidual - Returns the flag indicating whether true
1256: residuals must be computed explicitly or not.
1258: Not Collective
1260: Input Parameter:
1261: . eps - the eigensolver context
1263: Output Parameter:
1264: . trueres - the returned flag
1266: Level: advanced
1268: .seealso: EPSSetTrueResidual()
1269: @*/
1270: PetscErrorCode EPSGetTrueResidual(EPS eps,PetscBool *trueres)1271: {
1275: *trueres = eps->trueres;
1276: return(0);
1277: }
1279: /*@
1280: EPSSetTrackAll - Specifies if the solver must compute the residual norm of all
1281: approximate eigenpairs or not.
1283: Logically Collective on eps
1285: Input Parameters:
1286: + eps - the eigensolver context
1287: - trackall - whether to compute all residuals or not
1289: Notes:
1290: If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
1291: the residual norm for each eigenpair approximation. Computing the residual is
1292: usually an expensive operation and solvers commonly compute only the residual
1293: associated to the first unconverged eigenpair.
1295: The option '-eps_monitor_all' automatically activates this option.
1297: Level: developer
1299: .seealso: EPSGetTrackAll()
1300: @*/
1301: PetscErrorCode EPSSetTrackAll(EPS eps,PetscBool trackall)1302: {
1306: eps->trackall = trackall;
1307: return(0);
1308: }
1310: /*@
1311: EPSGetTrackAll - Returns the flag indicating whether all residual norms must
1312: be computed or not.
1314: Not Collective
1316: Input Parameter:
1317: . eps - the eigensolver context
1319: Output Parameter:
1320: . trackall - the returned flag
1322: Level: developer
1324: .seealso: EPSSetTrackAll()
1325: @*/
1326: PetscErrorCode EPSGetTrackAll(EPS eps,PetscBool *trackall)1327: {
1331: *trackall = eps->trackall;
1332: return(0);
1333: }
1335: /*@
1336: EPSSetPurify - Deactivate eigenvector purification (which is activated by default).
1338: Logically Collective on eps
1340: Input Parameters:
1341: + eps - the eigensolver context
1342: - purify - whether purification is required or not
1344: Options Database Keys:
1345: . -eps_purify <boolean> - Sets/resets the boolean flag 'purify'
1347: Notes:
1348: By default, eigenvectors of generalized symmetric eigenproblems are purified
1349: in order to purge directions in the nullspace of matrix B. If the user knows
1350: that B is non-singular, then purification can be safely deactivated and some
1351: computational cost is avoided (this is particularly important in interval computations).
1353: Level: intermediate
1355: .seealso: EPSGetPurify(), EPSSetInterval()
1356: @*/
1357: PetscErrorCode EPSSetPurify(EPS eps,PetscBool purify)1358: {
1362: if (purify!=eps->purify) {
1363: eps->purify = purify;
1364: eps->state = EPS_STATE_INITIAL;
1365: }
1366: return(0);
1367: }
1369: /*@
1370: EPSGetPurify - Returns the flag indicating whether purification is activated
1371: or not.
1373: Not Collective
1375: Input Parameter:
1376: . eps - the eigensolver context
1378: Output Parameter:
1379: . purify - the returned flag
1381: Level: intermediate
1383: .seealso: EPSSetPurify()
1384: @*/
1385: PetscErrorCode EPSGetPurify(EPS eps,PetscBool *purify)1386: {
1390: *purify = eps->purify;
1391: return(0);
1392: }
1394: /*@C
1395: EPSSetOptionsPrefix - Sets the prefix used for searching for all
1396: EPS options in the database.
1398: Logically Collective on eps
1400: Input Parameters:
1401: + eps - the eigensolver context
1402: - prefix - the prefix string to prepend to all EPS option requests
1404: Notes:
1405: A hyphen (-) must NOT be given at the beginning of the prefix name.
1406: The first character of all runtime options is AUTOMATICALLY the
1407: hyphen.
1409: For example, to distinguish between the runtime options for two
1410: different EPS contexts, one could call
1411: .vb
1412: EPSSetOptionsPrefix(eps1,"eig1_")
1413: EPSSetOptionsPrefix(eps2,"eig2_")
1414: .ve
1416: Level: advanced
1418: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
1419: @*/
1420: PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix)1421: {
1426: if (!eps->st) { EPSGetST(eps,&eps->st); }
1427: STSetOptionsPrefix(eps->st,prefix);
1428: if (!eps->V) { EPSGetBV(eps,&eps->V); }
1429: BVSetOptionsPrefix(eps->V,prefix);
1430: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1431: DSSetOptionsPrefix(eps->ds,prefix);
1432: if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
1433: RGSetOptionsPrefix(eps->rg,prefix);
1434: PetscObjectSetOptionsPrefix((PetscObject)eps,prefix);
1435: return(0);
1436: }
1438: /*@C
1439: EPSAppendOptionsPrefix - Appends to the prefix used for searching for all
1440: EPS options in the database.
1442: Logically Collective on eps
1444: Input Parameters:
1445: + eps - the eigensolver context
1446: - prefix - the prefix string to prepend to all EPS option requests
1448: Notes:
1449: A hyphen (-) must NOT be given at the beginning of the prefix name.
1450: The first character of all runtime options is AUTOMATICALLY the hyphen.
1452: Level: advanced
1454: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
1455: @*/
1456: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix)1457: {
1462: if (!eps->st) { EPSGetST(eps,&eps->st); }
1463: STAppendOptionsPrefix(eps->st,prefix);
1464: if (!eps->V) { EPSGetBV(eps,&eps->V); }
1465: BVAppendOptionsPrefix(eps->V,prefix);
1466: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
1467: DSAppendOptionsPrefix(eps->ds,prefix);
1468: if (!eps->rg) { EPSGetRG(eps,&eps->rg); }
1469: RGAppendOptionsPrefix(eps->rg,prefix);
1470: PetscObjectAppendOptionsPrefix((PetscObject)eps,prefix);
1471: return(0);
1472: }
1474: /*@C
1475: EPSGetOptionsPrefix - Gets the prefix used for searching for all
1476: EPS options in the database.
1478: Not Collective
1480: Input Parameters:
1481: . eps - the eigensolver context
1483: Output Parameters:
1484: . prefix - pointer to the prefix string used is returned
1486: Note:
1487: On the Fortran side, the user should pass in a string 'prefix' of
1488: sufficient length to hold the prefix.
1490: Level: advanced
1492: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
1493: @*/
1494: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])1495: {
1501: PetscObjectGetOptionsPrefix((PetscObject)eps,prefix);
1502: return(0);
1503: }