Actual source code: epsopts.c
slepc-3.20.2 2024-03-15
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, 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
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;
47: PetscFunctionBegin;
48: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,opt,&viewer,&format,&flg));
49: if (!flg) PetscFunctionReturn(PETSC_SUCCESS);
51: PetscCall(PetscViewerGetType(viewer,&vtype));
52: PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
53: PetscCall(PetscFunctionListFind(EPSMonitorList,key,&mfunc));
54: PetscCheck(mfunc,PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Specified viewer and format not supported");
55: PetscCall(PetscFunctionListFind(EPSMonitorCreateList,key,&cfunc));
56: PetscCall(PetscFunctionListFind(EPSMonitorDestroyList,key,&dfunc));
57: if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
58: if (!dfunc) dfunc = PetscViewerAndFormatDestroy;
60: PetscCall((*cfunc)(viewer,format,ctx,&vf));
61: PetscCall(PetscObjectDereference((PetscObject)viewer));
62: PetscCall(EPSMonitorSet(eps,mfunc,vf,(PetscErrorCode(*)(void **))dfunc));
63: if (trackall) PetscCall(EPSSetTrackAll(eps,PETSC_TRUE));
64: PetscFunctionReturn(PETSC_SUCCESS);
65: }
67: /*@
68: EPSSetFromOptions - Sets EPS options from the options database.
69: This routine must be called before EPSSetUp() if the user is to be
70: allowed to set the solver type.
72: Collective
74: Input Parameters:
75: . eps - the eigensolver context
77: Notes:
78: To see all options, run your program with the -help option.
80: Level: beginner
82: .seealso: EPSSetOptionsPrefix()
83: @*/
84: PetscErrorCode EPSSetFromOptions(EPS eps)
85: {
86: char type[256];
87: PetscBool set,flg,flg1,flg2,flg3,bval;
88: PetscReal r,array[2]={0,0};
89: PetscScalar s;
90: PetscInt i,j,k;
91: EPSBalance bal;
93: PetscFunctionBegin;
95: PetscCall(EPSRegisterAll());
96: PetscObjectOptionsBegin((PetscObject)eps);
97: PetscCall(PetscOptionsFList("-eps_type","Eigensolver method","EPSSetType",EPSList,(char*)(((PetscObject)eps)->type_name?((PetscObject)eps)->type_name:EPSKRYLOVSCHUR),type,sizeof(type),&flg));
98: if (flg) PetscCall(EPSSetType(eps,type));
99: else if (!((PetscObject)eps)->type_name) PetscCall(EPSSetType(eps,EPSKRYLOVSCHUR));
101: PetscCall(PetscOptionsBoolGroupBegin("-eps_hermitian","Hermitian eigenvalue problem","EPSSetProblemType",&flg));
102: if (flg) PetscCall(EPSSetProblemType(eps,EPS_HEP));
103: PetscCall(PetscOptionsBoolGroup("-eps_gen_hermitian","Generalized Hermitian eigenvalue problem","EPSSetProblemType",&flg));
104: if (flg) PetscCall(EPSSetProblemType(eps,EPS_GHEP));
105: PetscCall(PetscOptionsBoolGroup("-eps_non_hermitian","Non-Hermitian eigenvalue problem","EPSSetProblemType",&flg));
106: if (flg) PetscCall(EPSSetProblemType(eps,EPS_NHEP));
107: PetscCall(PetscOptionsBoolGroup("-eps_gen_non_hermitian","Generalized non-Hermitian eigenvalue problem","EPSSetProblemType",&flg));
108: if (flg) PetscCall(EPSSetProblemType(eps,EPS_GNHEP));
109: PetscCall(PetscOptionsBoolGroup("-eps_pos_gen_non_hermitian","Generalized non-Hermitian eigenvalue problem with positive semi-definite B","EPSSetProblemType",&flg));
110: if (flg) PetscCall(EPSSetProblemType(eps,EPS_PGNHEP));
111: PetscCall(PetscOptionsBoolGroupEnd("-eps_gen_indefinite","Generalized Hermitian-indefinite eigenvalue problem","EPSSetProblemType",&flg));
112: if (flg) PetscCall(EPSSetProblemType(eps,EPS_GHIEP));
114: PetscCall(PetscOptionsBoolGroupBegin("-eps_ritz","Rayleigh-Ritz extraction","EPSSetExtraction",&flg));
115: if (flg) PetscCall(EPSSetExtraction(eps,EPS_RITZ));
116: PetscCall(PetscOptionsBoolGroup("-eps_harmonic","Harmonic Ritz extraction","EPSSetExtraction",&flg));
117: if (flg) PetscCall(EPSSetExtraction(eps,EPS_HARMONIC));
118: PetscCall(PetscOptionsBoolGroup("-eps_harmonic_relative","Relative harmonic Ritz extraction","EPSSetExtraction",&flg));
119: if (flg) PetscCall(EPSSetExtraction(eps,EPS_HARMONIC_RELATIVE));
120: PetscCall(PetscOptionsBoolGroup("-eps_harmonic_right","Right harmonic Ritz extraction","EPSSetExtraction",&flg));
121: if (flg) PetscCall(EPSSetExtraction(eps,EPS_HARMONIC_RIGHT));
122: PetscCall(PetscOptionsBoolGroup("-eps_harmonic_largest","Largest harmonic Ritz extraction","EPSSetExtraction",&flg));
123: if (flg) PetscCall(EPSSetExtraction(eps,EPS_HARMONIC_LARGEST));
124: PetscCall(PetscOptionsBoolGroup("-eps_refined","Refined Ritz extraction","EPSSetExtraction",&flg));
125: if (flg) PetscCall(EPSSetExtraction(eps,EPS_REFINED));
126: PetscCall(PetscOptionsBoolGroupEnd("-eps_refined_harmonic","Refined harmonic Ritz extraction","EPSSetExtraction",&flg));
127: if (flg) PetscCall(EPSSetExtraction(eps,EPS_REFINED_HARMONIC));
129: bal = eps->balance;
130: PetscCall(PetscOptionsEnum("-eps_balance","Balancing method","EPSSetBalance",EPSBalanceTypes,(PetscEnum)bal,(PetscEnum*)&bal,&flg1));
131: j = eps->balance_its;
132: PetscCall(PetscOptionsInt("-eps_balance_its","Number of iterations in balancing","EPSSetBalance",eps->balance_its,&j,&flg2));
133: r = eps->balance_cutoff;
134: PetscCall(PetscOptionsReal("-eps_balance_cutoff","Cutoff value in balancing","EPSSetBalance",eps->balance_cutoff,&r,&flg3));
135: if (flg1 || flg2 || flg3) PetscCall(EPSSetBalance(eps,bal,j,r));
137: i = eps->max_it;
138: PetscCall(PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,&flg1));
139: r = eps->tol;
140: PetscCall(PetscOptionsReal("-eps_tol","Tolerance","EPSSetTolerances",SlepcDefaultTol(eps->tol),&r,&flg2));
141: if (flg1 || flg2) PetscCall(EPSSetTolerances(eps,r,i));
143: PetscCall(PetscOptionsBoolGroupBegin("-eps_conv_rel","Relative error convergence test","EPSSetConvergenceTest",&flg));
144: if (flg) PetscCall(EPSSetConvergenceTest(eps,EPS_CONV_REL));
145: PetscCall(PetscOptionsBoolGroup("-eps_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","EPSSetConvergenceTest",&flg));
146: if (flg) PetscCall(EPSSetConvergenceTest(eps,EPS_CONV_NORM));
147: PetscCall(PetscOptionsBoolGroup("-eps_conv_abs","Absolute error convergence test","EPSSetConvergenceTest",&flg));
148: if (flg) PetscCall(EPSSetConvergenceTest(eps,EPS_CONV_ABS));
149: PetscCall(PetscOptionsBoolGroupEnd("-eps_conv_user","User-defined convergence test","EPSSetConvergenceTest",&flg));
150: if (flg) PetscCall(EPSSetConvergenceTest(eps,EPS_CONV_USER));
152: PetscCall(PetscOptionsBoolGroupBegin("-eps_stop_basic","Stop iteration if all eigenvalues converged or max_it reached","EPSSetStoppingTest",&flg));
153: if (flg) PetscCall(EPSSetStoppingTest(eps,EPS_STOP_BASIC));
154: PetscCall(PetscOptionsBoolGroupEnd("-eps_stop_user","User-defined stopping test","EPSSetStoppingTest",&flg));
155: if (flg) PetscCall(EPSSetStoppingTest(eps,EPS_STOP_USER));
157: i = eps->nev;
158: PetscCall(PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,&flg1));
159: j = eps->ncv;
160: PetscCall(PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&j,&flg2));
161: k = eps->mpd;
162: PetscCall(PetscOptionsInt("-eps_mpd","Maximum dimension of projected problem","EPSSetDimensions",eps->mpd,&k,&flg3));
163: if (flg1 || flg2 || flg3) PetscCall(EPSSetDimensions(eps,i,j,k));
165: PetscCall(PetscOptionsBoolGroupBegin("-eps_largest_magnitude","Compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg));
166: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE));
167: PetscCall(PetscOptionsBoolGroup("-eps_smallest_magnitude","Compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg));
168: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE));
169: PetscCall(PetscOptionsBoolGroup("-eps_largest_real","Compute eigenvalues with largest real parts","EPSSetWhichEigenpairs",&flg));
170: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL));
171: PetscCall(PetscOptionsBoolGroup("-eps_smallest_real","Compute eigenvalues with smallest real parts","EPSSetWhichEigenpairs",&flg));
172: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL));
173: PetscCall(PetscOptionsBoolGroup("-eps_largest_imaginary","Compute eigenvalues with largest imaginary parts","EPSSetWhichEigenpairs",&flg));
174: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY));
175: PetscCall(PetscOptionsBoolGroup("-eps_smallest_imaginary","Compute eigenvalues with smallest imaginary parts","EPSSetWhichEigenpairs",&flg));
176: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY));
177: PetscCall(PetscOptionsBoolGroup("-eps_target_magnitude","Compute eigenvalues closest to target","EPSSetWhichEigenpairs",&flg));
178: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE));
179: PetscCall(PetscOptionsBoolGroup("-eps_target_real","Compute eigenvalues with real parts closest to target","EPSSetWhichEigenpairs",&flg));
180: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_REAL));
181: PetscCall(PetscOptionsBoolGroup("-eps_target_imaginary","Compute eigenvalues with imaginary parts closest to target","EPSSetWhichEigenpairs",&flg));
182: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_IMAGINARY));
183: PetscCall(PetscOptionsBoolGroupEnd("-eps_all","Compute all eigenvalues in an interval or a region","EPSSetWhichEigenpairs",&flg));
184: if (flg) PetscCall(EPSSetWhichEigenpairs(eps,EPS_ALL));
186: PetscCall(PetscOptionsScalar("-eps_target","Value of the target","EPSSetTarget",eps->target,&s,&flg));
187: if (flg) {
188: if (eps->which!=EPS_TARGET_REAL && eps->which!=EPS_TARGET_IMAGINARY) PetscCall(EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE));
189: PetscCall(EPSSetTarget(eps,s));
190: }
192: k = 2;
193: PetscCall(PetscOptionsRealArray("-eps_interval","Computational interval (two real values separated with a comma without spaces)","EPSSetInterval",array,&k,&flg));
194: if (flg) {
195: PetscCheck(k>1,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_SIZ,"Must pass two values in -eps_interval (comma-separated without spaces)");
196: PetscCall(EPSSetWhichEigenpairs(eps,EPS_ALL));
197: PetscCall(EPSSetInterval(eps,array[0],array[1]));
198: }
200: PetscCall(PetscOptionsBool("-eps_true_residual","Compute true residuals explicitly","EPSSetTrueResidual",eps->trueres,&eps->trueres,NULL));
201: PetscCall(PetscOptionsBool("-eps_purify","Postprocess eigenvectors for purification","EPSSetPurify",eps->purify,&bval,&flg));
202: if (flg) PetscCall(EPSSetPurify(eps,bval));
203: PetscCall(PetscOptionsBool("-eps_two_sided","Use two-sided variant (to compute left eigenvectors)","EPSSetTwoSided",eps->twosided,&bval,&flg));
204: if (flg) PetscCall(EPSSetTwoSided(eps,bval));
206: /* -----------------------------------------------------------------------*/
207: /*
208: Cancels all monitors hardwired into code before call to EPSSetFromOptions()
209: */
210: PetscCall(PetscOptionsBool("-eps_monitor_cancel","Remove any hardwired monitor routines","EPSMonitorCancel",PETSC_FALSE,&flg,&set));
211: if (set && flg) PetscCall(EPSMonitorCancel(eps));
212: PetscCall(EPSMonitorSetFromOptions(eps,"-eps_monitor","first_approximation",NULL,PETSC_FALSE));
213: PetscCall(EPSMonitorSetFromOptions(eps,"-eps_monitor_all","all_approximations",NULL,PETSC_TRUE));
214: PetscCall(EPSMonitorSetFromOptions(eps,"-eps_monitor_conv","convergence_history",NULL,PETSC_FALSE));
216: /* -----------------------------------------------------------------------*/
217: PetscCall(PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",&set));
218: PetscCall(PetscOptionsName("-eps_view_vectors","View computed eigenvectors","EPSVectorsView",&set));
219: PetscCall(PetscOptionsName("-eps_view_values","View computed eigenvalues","EPSValuesView",&set));
220: PetscCall(PetscOptionsName("-eps_converged_reason","Print reason for convergence, and number of iterations","EPSConvergedReasonView",&set));
221: PetscCall(PetscOptionsName("-eps_error_absolute","Print absolute errors of each eigenpair","EPSErrorView",&set));
222: PetscCall(PetscOptionsName("-eps_error_relative","Print relative errors of each eigenpair","EPSErrorView",&set));
223: PetscCall(PetscOptionsName("-eps_error_backward","Print backward errors of each eigenpair","EPSErrorView",&set));
225: PetscTryTypeMethod(eps,setfromoptions,PetscOptionsObject);
226: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)eps,PetscOptionsObject));
227: PetscOptionsEnd();
229: if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
230: PetscCall(BVSetFromOptions(eps->V));
231: if (!eps->rg) PetscCall(EPSGetRG(eps,&eps->rg));
232: PetscCall(RGSetFromOptions(eps->rg));
233: if (eps->useds) {
234: if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
235: PetscCall(EPSSetDSType(eps));
236: PetscCall(DSSetFromOptions(eps->ds));
237: }
238: if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
239: PetscCall(EPSSetDefaultST(eps));
240: PetscCall(STSetFromOptions(eps->st));
241: PetscFunctionReturn(PETSC_SUCCESS);
242: }
244: /*@C
245: EPSGetTolerances - Gets the tolerance and maximum iteration count used
246: by the EPS convergence tests.
248: Not Collective
250: Input Parameter:
251: . eps - the eigensolver context
253: Output Parameters:
254: + tol - the convergence tolerance
255: - maxits - maximum number of iterations
257: Notes:
258: The user can specify NULL for any parameter that is not needed.
260: Level: intermediate
262: .seealso: EPSSetTolerances()
263: @*/
264: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,PetscInt *maxits)
265: {
266: PetscFunctionBegin;
268: if (tol) *tol = eps->tol;
269: if (maxits) *maxits = eps->max_it;
270: PetscFunctionReturn(PETSC_SUCCESS);
271: }
273: /*@
274: EPSSetTolerances - Sets the tolerance and maximum iteration count used
275: by the EPS convergence tests.
277: Logically Collective
279: Input Parameters:
280: + eps - the eigensolver context
281: . tol - the convergence tolerance
282: - maxits - maximum number of iterations to use
284: Options Database Keys:
285: + -eps_tol <tol> - Sets the convergence tolerance
286: - -eps_max_it <maxits> - Sets the maximum number of iterations allowed
288: Notes:
289: Use PETSC_DEFAULT for either argument to assign a reasonably good value.
291: Level: intermediate
293: .seealso: EPSGetTolerances()
294: @*/
295: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,PetscInt maxits)
296: {
297: PetscFunctionBegin;
301: if (tol == (PetscReal)PETSC_DEFAULT) {
302: eps->tol = PETSC_DEFAULT;
303: eps->state = EPS_STATE_INITIAL;
304: } else {
305: PetscCheck(tol>0.0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
306: eps->tol = tol;
307: }
308: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
309: eps->max_it = PETSC_DEFAULT;
310: eps->state = EPS_STATE_INITIAL;
311: } else {
312: PetscCheck(maxits>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
313: eps->max_it = maxits;
314: }
315: PetscFunctionReturn(PETSC_SUCCESS);
316: }
318: /*@C
319: EPSGetDimensions - Gets the number of eigenvalues to compute
320: and the dimension of the subspace.
322: Not Collective
324: Input Parameter:
325: . eps - the eigensolver context
327: Output Parameters:
328: + nev - number of eigenvalues to compute
329: . ncv - the maximum dimension of the subspace to be used by the solver
330: - mpd - the maximum dimension allowed for the projected problem
332: Level: intermediate
334: .seealso: EPSSetDimensions()
335: @*/
336: PetscErrorCode EPSGetDimensions(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
337: {
338: PetscFunctionBegin;
340: if (nev) *nev = eps->nev;
341: if (ncv) *ncv = eps->ncv;
342: if (mpd) *mpd = eps->mpd;
343: PetscFunctionReturn(PETSC_SUCCESS);
344: }
346: /*@
347: EPSSetDimensions - Sets the number of eigenvalues to compute
348: and the dimension of the subspace.
350: Logically Collective
352: Input Parameters:
353: + eps - the eigensolver context
354: . nev - number of eigenvalues to compute
355: . ncv - the maximum dimension of the subspace to be used by the solver
356: - mpd - the maximum dimension allowed for the projected problem
358: Options Database Keys:
359: + -eps_nev <nev> - Sets the number of eigenvalues
360: . -eps_ncv <ncv> - Sets the dimension of the subspace
361: - -eps_mpd <mpd> - Sets the maximum projected dimension
363: Notes:
364: Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
365: dependent on the solution method.
367: The parameters ncv and mpd are intimately related, so that the user is advised
368: to set one of them at most. Normal usage is that
369: (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
370: (b) in cases where nev is large, the user sets mpd.
372: The value of ncv should always be between nev and (nev+mpd), typically
373: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
374: a smaller value should be used.
376: When computing all eigenvalues in an interval, see EPSSetInterval(), these
377: parameters lose relevance, and tuning must be done with
378: EPSKrylovSchurSetDimensions().
380: Level: intermediate
382: .seealso: EPSGetDimensions(), EPSSetInterval(), EPSKrylovSchurSetDimensions()
383: @*/
384: PetscErrorCode EPSSetDimensions(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)
385: {
386: PetscFunctionBegin;
391: PetscCheck(nev>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
392: eps->nev = nev;
393: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
394: eps->ncv = PETSC_DEFAULT;
395: } else {
396: PetscCheck(ncv>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
397: eps->ncv = ncv;
398: }
399: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
400: eps->mpd = PETSC_DEFAULT;
401: } else {
402: PetscCheck(mpd>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
403: eps->mpd = mpd;
404: }
405: eps->state = EPS_STATE_INITIAL;
406: PetscFunctionReturn(PETSC_SUCCESS);
407: }
409: /*@
410: EPSSetWhichEigenpairs - Specifies which portion of the spectrum is
411: to be sought.
413: Logically Collective
415: Input Parameters:
416: + eps - eigensolver context obtained from EPSCreate()
417: - which - the portion of the spectrum to be sought
419: Options Database Keys:
420: + -eps_largest_magnitude - Sets largest eigenvalues in magnitude
421: . -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
422: . -eps_largest_real - Sets largest real parts
423: . -eps_smallest_real - Sets smallest real parts
424: . -eps_largest_imaginary - Sets largest imaginary parts
425: . -eps_smallest_imaginary - Sets smallest imaginary parts
426: . -eps_target_magnitude - Sets eigenvalues closest to target
427: . -eps_target_real - Sets real parts closest to target
428: . -eps_target_imaginary - Sets imaginary parts closest to target
429: - -eps_all - Sets all eigenvalues in an interval or region
431: Notes:
432: The parameter 'which' can have one of these values
434: + EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
435: . EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
436: . EPS_LARGEST_REAL - largest real parts
437: . EPS_SMALLEST_REAL - smallest real parts
438: . EPS_LARGEST_IMAGINARY - largest imaginary parts
439: . EPS_SMALLEST_IMAGINARY - smallest imaginary parts
440: . EPS_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
441: . EPS_TARGET_REAL - eigenvalues with real part closest to target
442: . EPS_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
443: . EPS_ALL - all eigenvalues contained in a given interval or region
444: - EPS_WHICH_USER - user defined ordering set with EPSSetEigenvalueComparison()
446: Not all eigensolvers implemented in EPS account for all the possible values
447: stated above. Also, some values make sense only for certain types of
448: problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
449: and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part
450: for eigenvalue selection.
452: The target is a scalar value provided with EPSSetTarget().
454: The criterion EPS_TARGET_IMAGINARY is available only in case PETSc and
455: SLEPc have been built with complex scalars.
457: EPS_ALL is intended for use in combination with an interval (see
458: EPSSetInterval()), when all eigenvalues within the interval are requested,
459: or in the context of the CISS solver for computing all eigenvalues in a region.
460: In those cases, the number of eigenvalues is unknown, so the nev parameter
461: has a different sense, see EPSSetDimensions().
463: Level: intermediate
465: .seealso: EPSGetWhichEigenpairs(), EPSSetTarget(), EPSSetInterval(),
466: EPSSetDimensions(), EPSSetEigenvalueComparison(), EPSWhich
467: @*/
468: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)
469: {
470: PetscFunctionBegin;
473: switch (which) {
474: case EPS_LARGEST_MAGNITUDE:
475: case EPS_SMALLEST_MAGNITUDE:
476: case EPS_LARGEST_REAL:
477: case EPS_SMALLEST_REAL:
478: case EPS_LARGEST_IMAGINARY:
479: case EPS_SMALLEST_IMAGINARY:
480: case EPS_TARGET_MAGNITUDE:
481: case EPS_TARGET_REAL:
482: #if defined(PETSC_USE_COMPLEX)
483: case EPS_TARGET_IMAGINARY:
484: #endif
485: case EPS_ALL:
486: case EPS_WHICH_USER:
487: if (eps->which != which) {
488: eps->state = EPS_STATE_INITIAL;
489: eps->which = which;
490: }
491: break;
492: #if !defined(PETSC_USE_COMPLEX)
493: case EPS_TARGET_IMAGINARY:
494: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"EPS_TARGET_IMAGINARY can be used only with complex scalars");
495: #endif
496: default:
497: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
498: }
499: PetscFunctionReturn(PETSC_SUCCESS);
500: }
502: /*@
503: EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be
504: sought.
506: Not Collective
508: Input Parameter:
509: . eps - eigensolver context obtained from EPSCreate()
511: Output Parameter:
512: . which - the portion of the spectrum to be sought
514: Notes:
515: See EPSSetWhichEigenpairs() for possible values of 'which'.
517: Level: intermediate
519: .seealso: EPSSetWhichEigenpairs(), EPSWhich
520: @*/
521: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
522: {
523: PetscFunctionBegin;
525: PetscAssertPointer(which,2);
526: *which = eps->which;
527: PetscFunctionReturn(PETSC_SUCCESS);
528: }
530: /*@C
531: EPSSetEigenvalueComparison - Specifies the eigenvalue comparison function
532: when EPSSetWhichEigenpairs() is set to EPS_WHICH_USER.
534: Logically Collective
536: Input Parameters:
537: + eps - eigensolver context obtained from EPSCreate()
538: . func - a pointer to the comparison function
539: - ctx - a context pointer (the last parameter to the comparison function)
541: Calling sequence of func:
542: $ PetscErrorCode func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
543: + ar - real part of the 1st eigenvalue
544: . ai - imaginary part of the 1st eigenvalue
545: . br - real part of the 2nd eigenvalue
546: . bi - imaginary part of the 2nd eigenvalue
547: . res - result of comparison
548: - ctx - optional context, as set by EPSSetEigenvalueComparison()
550: Note:
551: The returning parameter 'res' can be
552: + negative - if the 1st eigenvalue is preferred to the 2st one
553: . zero - if both eigenvalues are equally preferred
554: - positive - if the 2st eigenvalue is preferred to the 1st one
556: Level: advanced
558: .seealso: EPSSetWhichEigenpairs(), EPSWhich
559: @*/
560: PetscErrorCode EPSSetEigenvalueComparison(EPS eps,PetscErrorCode (*func)(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx),void* ctx)
561: {
562: PetscFunctionBegin;
564: eps->sc->comparison = func;
565: eps->sc->comparisonctx = ctx;
566: eps->which = EPS_WHICH_USER;
567: PetscFunctionReturn(PETSC_SUCCESS);
568: }
570: /*@C
571: EPSSetArbitrarySelection - Specifies a function intended to look for
572: eigenvalues according to an arbitrary selection criterion. This criterion
573: can be based on a computation involving the current eigenvector approximation.
575: Logically Collective
577: Input Parameters:
578: + eps - eigensolver context obtained from EPSCreate()
579: . func - a pointer to the evaluation function
580: - ctx - a context pointer (the last parameter to the evaluation function)
582: Calling sequence of func:
583: $ PetscErrorCode func(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx)
584: + er - real part of the current eigenvalue approximation
585: . ei - imaginary part of the current eigenvalue approximation
586: . xr - real part of the current eigenvector approximation
587: . xi - imaginary part of the current eigenvector approximation
588: . rr - result of evaluation (real part)
589: . ri - result of evaluation (imaginary part)
590: - ctx - optional context, as set by EPSSetArbitrarySelection()
592: Notes:
593: This provides a mechanism to select eigenpairs by evaluating a user-defined
594: function. When a function has been provided, the default selection based on
595: sorting the eigenvalues is replaced by the sorting of the results of this
596: function (with the same sorting criterion given in EPSSetWhichEigenpairs()).
598: For instance, suppose you want to compute those eigenvectors that maximize
599: a certain computable expression. Then implement the computation using
600: the arguments xr and xi, and return the result in rr. Then set the standard
601: sorting by magnitude so that the eigenpair with largest value of rr is
602: selected.
604: This evaluation function is collective, that is, all processes call it and
605: it can use collective operations; furthermore, the computed result must
606: be the same in all processes.
608: The result of func is expressed as a complex number so that it is possible to
609: use the standard eigenvalue sorting functions, but normally only rr is used.
610: Set ri to zero unless it is meaningful in your application.
612: Level: advanced
614: .seealso: EPSSetWhichEigenpairs()
615: @*/
616: PetscErrorCode EPSSetArbitrarySelection(EPS eps,PetscErrorCode (*func)(PetscScalar er,PetscScalar ei,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx),void* ctx)
617: {
618: PetscFunctionBegin;
620: eps->arbitrary = func;
621: eps->arbitraryctx = ctx;
622: eps->state = EPS_STATE_INITIAL;
623: PetscFunctionReturn(PETSC_SUCCESS);
624: }
626: /*@C
627: EPSSetConvergenceTestFunction - Sets a function to compute the error estimate
628: used in the convergence test.
630: Logically Collective
632: Input Parameters:
633: + eps - eigensolver context obtained from EPSCreate()
634: . func - a pointer to the convergence test function
635: . ctx - context for private data for the convergence routine (may be null)
636: - destroy - a routine for destroying the context (may be null)
638: Calling sequence of func:
639: $ PetscErrorCode func(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
640: + eps - eigensolver context obtained from EPSCreate()
641: . eigr - real part of the eigenvalue
642: . eigi - imaginary part of the eigenvalue
643: . res - residual norm associated to the eigenpair
644: . errest - (output) computed error estimate
645: - ctx - optional context, as set by EPSSetConvergenceTestFunction()
647: Note:
648: If the error estimate returned by the convergence test function is less than
649: the tolerance, then the eigenvalue is accepted as converged.
651: Level: advanced
653: .seealso: EPSSetConvergenceTest(), EPSSetTolerances()
654: @*/
655: PetscErrorCode EPSSetConvergenceTestFunction(EPS eps,PetscErrorCode (*func)(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx),void* ctx,PetscErrorCode (*destroy)(void*))
656: {
657: PetscFunctionBegin;
659: if (eps->convergeddestroy) PetscCall((*eps->convergeddestroy)(eps->convergedctx));
660: eps->convergeduser = func;
661: eps->convergeddestroy = destroy;
662: eps->convergedctx = ctx;
663: if (func == EPSConvergedRelative) eps->conv = EPS_CONV_REL;
664: else if (func == EPSConvergedNorm) eps->conv = EPS_CONV_NORM;
665: else if (func == EPSConvergedAbsolute) eps->conv = EPS_CONV_ABS;
666: else {
667: eps->conv = EPS_CONV_USER;
668: eps->converged = eps->convergeduser;
669: }
670: PetscFunctionReturn(PETSC_SUCCESS);
671: }
673: /*@
674: EPSSetConvergenceTest - Specifies how to compute the error estimate
675: used in the convergence test.
677: Logically Collective
679: Input Parameters:
680: + eps - eigensolver context obtained from EPSCreate()
681: - conv - the type of convergence test
683: Options Database Keys:
684: + -eps_conv_abs - Sets the absolute convergence test
685: . -eps_conv_rel - Sets the convergence test relative to the eigenvalue
686: . -eps_conv_norm - Sets the convergence test relative to the matrix norms
687: - -eps_conv_user - Selects the user-defined convergence test
689: Note:
690: The parameter 'conv' can have one of these values
691: + EPS_CONV_ABS - absolute error ||r||
692: . EPS_CONV_REL - error relative to the eigenvalue l, ||r||/|l|
693: . EPS_CONV_NORM - error relative to the matrix norms, ||r||/(||A||+|l|*||B||)
694: - EPS_CONV_USER - function set by EPSSetConvergenceTestFunction()
696: Level: intermediate
698: .seealso: EPSGetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSSetStoppingTest(), EPSConv
699: @*/
700: PetscErrorCode EPSSetConvergenceTest(EPS eps,EPSConv conv)
701: {
702: PetscFunctionBegin;
705: switch (conv) {
706: case EPS_CONV_ABS: eps->converged = EPSConvergedAbsolute; break;
707: case EPS_CONV_REL: eps->converged = EPSConvergedRelative; break;
708: case EPS_CONV_NORM: eps->converged = EPSConvergedNorm; break;
709: case EPS_CONV_USER:
710: PetscCheck(eps->convergeduser,PetscObjectComm((PetscObject)eps),PETSC_ERR_ORDER,"Must call EPSSetConvergenceTestFunction() first");
711: eps->converged = eps->convergeduser;
712: break;
713: default:
714: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
715: }
716: eps->conv = conv;
717: PetscFunctionReturn(PETSC_SUCCESS);
718: }
720: /*@
721: EPSGetConvergenceTest - Gets the method used to compute the error estimate
722: used in the convergence test.
724: Not Collective
726: Input Parameters:
727: . eps - eigensolver context obtained from EPSCreate()
729: Output Parameters:
730: . conv - the type of convergence test
732: Level: intermediate
734: .seealso: EPSSetConvergenceTest(), EPSConv
735: @*/
736: PetscErrorCode EPSGetConvergenceTest(EPS eps,EPSConv *conv)
737: {
738: PetscFunctionBegin;
740: PetscAssertPointer(conv,2);
741: *conv = eps->conv;
742: PetscFunctionReturn(PETSC_SUCCESS);
743: }
745: /*@C
746: EPSSetStoppingTestFunction - Sets a function to decide when to stop the outer
747: iteration of the eigensolver.
749: Logically Collective
751: Input Parameters:
752: + eps - eigensolver context obtained from EPSCreate()
753: . func - pointer to the stopping test function
754: . ctx - context for private data for the stopping routine (may be null)
755: - destroy - a routine for destroying the context (may be null)
757: Calling sequence of func:
758: $ PetscErrorCode func(EPS eps,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,EPSConvergedReason *reason,void *ctx)
759: + eps - eigensolver context obtained from EPSCreate()
760: . its - current number of iterations
761: . max_it - maximum number of iterations
762: . nconv - number of currently converged eigenpairs
763: . nev - number of requested eigenpairs
764: . reason - (output) result of the stopping test
765: - ctx - optional context, as set by EPSSetStoppingTestFunction()
767: Note:
768: Normal usage is to first call the default routine EPSStoppingBasic() and then
769: set reason to EPS_CONVERGED_USER if some user-defined conditions have been
770: met. To let the eigensolver continue iterating, the result must be left as
771: EPS_CONVERGED_ITERATING.
773: Level: advanced
775: .seealso: EPSSetStoppingTest(), EPSStoppingBasic()
776: @*/
777: PetscErrorCode EPSSetStoppingTestFunction(EPS eps,PetscErrorCode (*func)(EPS eps,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,EPSConvergedReason *reason,void *ctx),void* ctx,PetscErrorCode (*destroy)(void*))
778: {
779: PetscFunctionBegin;
781: if (eps->stoppingdestroy) PetscCall((*eps->stoppingdestroy)(eps->stoppingctx));
782: eps->stoppinguser = func;
783: eps->stoppingdestroy = destroy;
784: eps->stoppingctx = ctx;
785: if (func == EPSStoppingBasic) eps->stop = EPS_STOP_BASIC;
786: else {
787: eps->stop = EPS_STOP_USER;
788: eps->stopping = eps->stoppinguser;
789: }
790: PetscFunctionReturn(PETSC_SUCCESS);
791: }
793: /*@
794: EPSSetStoppingTest - Specifies how to decide the termination of the outer
795: loop of the eigensolver.
797: Logically Collective
799: Input Parameters:
800: + eps - eigensolver context obtained from EPSCreate()
801: - stop - the type of stopping test
803: Options Database Keys:
804: + -eps_stop_basic - Sets the default stopping test
805: - -eps_stop_user - Selects the user-defined stopping test
807: Note:
808: The parameter 'stop' can have one of these values
809: + EPS_STOP_BASIC - default stopping test
810: - EPS_STOP_USER - function set by EPSSetStoppingTestFunction()
812: Level: advanced
814: .seealso: EPSGetStoppingTest(), EPSSetStoppingTestFunction(), EPSSetConvergenceTest(), EPSStop
815: @*/
816: PetscErrorCode EPSSetStoppingTest(EPS eps,EPSStop stop)
817: {
818: PetscFunctionBegin;
821: switch (stop) {
822: case EPS_STOP_BASIC: eps->stopping = EPSStoppingBasic; break;
823: case EPS_STOP_USER:
824: PetscCheck(eps->stoppinguser,PetscObjectComm((PetscObject)eps),PETSC_ERR_ORDER,"Must call EPSSetStoppingTestFunction() first");
825: eps->stopping = eps->stoppinguser;
826: break;
827: default:
828: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
829: }
830: eps->stop = stop;
831: PetscFunctionReturn(PETSC_SUCCESS);
832: }
834: /*@
835: EPSGetStoppingTest - Gets the method used to decide the termination of the outer
836: loop of the eigensolver.
838: Not Collective
840: Input Parameters:
841: . eps - eigensolver context obtained from EPSCreate()
843: Output Parameters:
844: . stop - the type of stopping test
846: Level: advanced
848: .seealso: EPSSetStoppingTest(), EPSStop
849: @*/
850: PetscErrorCode EPSGetStoppingTest(EPS eps,EPSStop *stop)
851: {
852: PetscFunctionBegin;
854: PetscAssertPointer(stop,2);
855: *stop = eps->stop;
856: PetscFunctionReturn(PETSC_SUCCESS);
857: }
859: /*@
860: EPSSetProblemType - Specifies the type of the eigenvalue problem.
862: Logically Collective
864: Input Parameters:
865: + eps - the eigensolver context
866: - type - a known type of eigenvalue problem
868: Options Database Keys:
869: + -eps_hermitian - Hermitian eigenvalue problem
870: . -eps_gen_hermitian - generalized Hermitian eigenvalue problem
871: . -eps_non_hermitian - non-Hermitian eigenvalue problem
872: . -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
873: . -eps_pos_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
874: with positive semi-definite B
875: - -eps_gen_indefinite - generalized Hermitian-indefinite eigenvalue problem
877: Notes:
878: Allowed values for the problem type are Hermitian (EPS_HEP), non-Hermitian
879: (EPS_NHEP), generalized Hermitian (EPS_GHEP), generalized non-Hermitian
880: (EPS_GNHEP), generalized non-Hermitian with positive semi-definite B
881: (EPS_PGNHEP), and generalized Hermitian-indefinite (EPS_GHIEP).
883: This function must be used to instruct SLEPc to exploit symmetry. If no
884: problem type is specified, by default a non-Hermitian problem is assumed
885: (either standard or generalized). If the user knows that the problem is
886: Hermitian (i.e. A=A^H) or generalized Hermitian (i.e. A=A^H, B=B^H, and
887: B positive definite) then it is recommended to set the problem type so
888: that eigensolver can exploit these properties.
890: Level: intermediate
892: .seealso: EPSSetOperators(), EPSSetType(), EPSGetProblemType(), EPSProblemType
893: @*/
894: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
895: {
896: PetscFunctionBegin;
899: if (type == eps->problem_type) PetscFunctionReturn(PETSC_SUCCESS);
900: switch (type) {
901: case EPS_HEP:
902: eps->isgeneralized = PETSC_FALSE;
903: eps->ishermitian = PETSC_TRUE;
904: eps->ispositive = PETSC_FALSE;
905: break;
906: case EPS_NHEP:
907: eps->isgeneralized = PETSC_FALSE;
908: eps->ishermitian = PETSC_FALSE;
909: eps->ispositive = PETSC_FALSE;
910: break;
911: case EPS_GHEP:
912: eps->isgeneralized = PETSC_TRUE;
913: eps->ishermitian = PETSC_TRUE;
914: eps->ispositive = PETSC_TRUE;
915: break;
916: case EPS_GNHEP:
917: eps->isgeneralized = PETSC_TRUE;
918: eps->ishermitian = PETSC_FALSE;
919: eps->ispositive = PETSC_FALSE;
920: break;
921: case EPS_PGNHEP:
922: eps->isgeneralized = PETSC_TRUE;
923: eps->ishermitian = PETSC_FALSE;
924: eps->ispositive = PETSC_TRUE;
925: break;
926: case EPS_GHIEP:
927: eps->isgeneralized = PETSC_TRUE;
928: eps->ishermitian = PETSC_TRUE;
929: eps->ispositive = PETSC_FALSE;
930: break;
931: default:
932: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
933: }
934: eps->problem_type = type;
935: eps->state = EPS_STATE_INITIAL;
936: PetscFunctionReturn(PETSC_SUCCESS);
937: }
939: /*@
940: EPSGetProblemType - Gets the problem type from the EPS object.
942: Not Collective
944: Input Parameter:
945: . eps - the eigensolver context
947: Output Parameter:
948: . type - the problem type
950: Level: intermediate
952: .seealso: EPSSetProblemType(), EPSProblemType
953: @*/
954: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
955: {
956: PetscFunctionBegin;
958: PetscAssertPointer(type,2);
959: *type = eps->problem_type;
960: PetscFunctionReturn(PETSC_SUCCESS);
961: }
963: /*@
964: EPSSetExtraction - Specifies the type of extraction technique to be employed
965: by the eigensolver.
967: Logically Collective
969: Input Parameters:
970: + eps - the eigensolver context
971: - extr - a known type of extraction
973: Options Database Keys:
974: + -eps_ritz - Rayleigh-Ritz extraction
975: . -eps_harmonic - harmonic Ritz extraction
976: . -eps_harmonic_relative - harmonic Ritz extraction relative to the eigenvalue
977: . -eps_harmonic_right - harmonic Ritz extraction for rightmost eigenvalues
978: . -eps_harmonic_largest - harmonic Ritz extraction for largest magnitude
979: (without target)
980: . -eps_refined - refined Ritz extraction
981: - -eps_refined_harmonic - refined harmonic Ritz extraction
983: Notes:
984: Not all eigensolvers support all types of extraction. See the SLEPc
985: Users Manual for details.
987: By default, a standard Rayleigh-Ritz extraction is used. Other extractions
988: may be useful when computing interior eigenvalues.
990: Harmonic-type extractions are used in combination with a 'target'.
992: Level: advanced
994: .seealso: EPSSetTarget(), EPSGetExtraction(), EPSExtraction
995: @*/
996: PetscErrorCode EPSSetExtraction(EPS eps,EPSExtraction extr)
997: {
998: PetscFunctionBegin;
1001: eps->extraction = extr;
1002: PetscFunctionReturn(PETSC_SUCCESS);
1003: }
1005: /*@
1006: EPSGetExtraction - Gets the extraction type used by the EPS object.
1008: Not Collective
1010: Input Parameter:
1011: . eps - the eigensolver context
1013: Output Parameter:
1014: . extr - name of extraction type
1016: Level: advanced
1018: .seealso: EPSSetExtraction(), EPSExtraction
1019: @*/
1020: PetscErrorCode EPSGetExtraction(EPS eps,EPSExtraction *extr)
1021: {
1022: PetscFunctionBegin;
1024: PetscAssertPointer(extr,2);
1025: *extr = eps->extraction;
1026: PetscFunctionReturn(PETSC_SUCCESS);
1027: }
1029: /*@
1030: EPSSetBalance - Specifies the balancing technique to be employed by the
1031: eigensolver, and some parameters associated to it.
1033: Logically Collective
1035: Input Parameters:
1036: + eps - the eigensolver context
1037: . bal - the balancing method, one of EPS_BALANCE_NONE, EPS_BALANCE_ONESIDE,
1038: EPS_BALANCE_TWOSIDE, or EPS_BALANCE_USER
1039: . its - number of iterations of the balancing algorithm
1040: - cutoff - cutoff value
1042: Options Database Keys:
1043: + -eps_balance <method> - the balancing method, where <method> is one of
1044: 'none', 'oneside', 'twoside', or 'user'
1045: . -eps_balance_its <its> - number of iterations
1046: - -eps_balance_cutoff <cutoff> - cutoff value
1048: Notes:
1049: When balancing is enabled, the solver works implicitly with matrix DAD^-1,
1050: where D is an appropriate diagonal matrix. This improves the accuracy of
1051: the computed results in some cases. See the SLEPc Users Manual for details.
1053: Balancing makes sense only for non-Hermitian problems when the required
1054: precision is high (i.e. a small tolerance such as 1e-15).
1056: By default, balancing is disabled. The two-sided method is much more
1057: effective than the one-sided counterpart, but it requires the system
1058: matrices to have the MatMultTranspose operation defined.
1060: The parameter 'its' is the number of iterations performed by the method. The
1061: cutoff value is used only in the two-side variant. Use PETSC_DEFAULT to assign
1062: a reasonably good value.
1064: User-defined balancing is allowed provided that the corresponding matrix
1065: is set via STSetBalanceMatrix.
1067: Level: intermediate
1069: .seealso: EPSGetBalance(), EPSBalance, STSetBalanceMatrix()
1070: @*/
1071: PetscErrorCode EPSSetBalance(EPS eps,EPSBalance bal,PetscInt its,PetscReal cutoff)
1072: {
1073: PetscFunctionBegin;
1078: switch (bal) {
1079: case EPS_BALANCE_NONE:
1080: case EPS_BALANCE_ONESIDE:
1081: case EPS_BALANCE_TWOSIDE:
1082: case EPS_BALANCE_USER:
1083: if (eps->balance != bal) {
1084: eps->state = EPS_STATE_INITIAL;
1085: eps->balance = bal;
1086: }
1087: break;
1088: default:
1089: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid value of argument 'bal'");
1090: }
1091: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) eps->balance_its = 5;
1092: else {
1093: PetscCheck(its>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be > 0");
1094: eps->balance_its = its;
1095: }
1096: if (cutoff==(PetscReal)PETSC_DECIDE || cutoff==(PetscReal)PETSC_DEFAULT) eps->balance_cutoff = 1e-8;
1097: else {
1098: PetscCheck(cutoff>0.0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of cutoff. Must be > 0");
1099: eps->balance_cutoff = cutoff;
1100: }
1101: PetscFunctionReturn(PETSC_SUCCESS);
1102: }
1104: /*@C
1105: EPSGetBalance - Gets the balancing type used by the EPS object, and the
1106: associated parameters.
1108: Not Collective
1110: Input Parameter:
1111: . eps - the eigensolver context
1113: Output Parameters:
1114: + bal - the balancing method
1115: . its - number of iterations of the balancing algorithm
1116: - cutoff - cutoff value
1118: Level: intermediate
1120: Note:
1121: The user can specify NULL for any parameter that is not needed.
1123: .seealso: EPSSetBalance(), EPSBalance
1124: @*/
1125: PetscErrorCode EPSGetBalance(EPS eps,EPSBalance *bal,PetscInt *its,PetscReal *cutoff)
1126: {
1127: PetscFunctionBegin;
1129: if (bal) *bal = eps->balance;
1130: if (its) *its = eps->balance_its;
1131: if (cutoff) *cutoff = eps->balance_cutoff;
1132: PetscFunctionReturn(PETSC_SUCCESS);
1133: }
1135: /*@
1136: EPSSetTwoSided - Sets the solver to use a two-sided variant so that left
1137: eigenvectors are also computed.
1139: Logically Collective
1141: Input Parameters:
1142: + eps - the eigensolver context
1143: - twosided - whether the two-sided variant is to be used or not
1145: Options Database Keys:
1146: . -eps_two_sided <boolean> - Sets/resets the twosided flag
1148: Notes:
1149: If the user sets twosided=PETSC_TRUE then the solver uses a variant of
1150: the algorithm that computes both right and left eigenvectors. This is
1151: usually much more costly. This option is not available in all solvers.
1153: When using two-sided solvers, the problem matrices must have both the
1154: MatMult and MatMultTranspose operations defined.
1156: Level: advanced
1158: .seealso: EPSGetTwoSided(), EPSGetLeftEigenvector()
1159: @*/
1160: PetscErrorCode EPSSetTwoSided(EPS eps,PetscBool twosided)
1161: {
1162: PetscFunctionBegin;
1165: if (twosided!=eps->twosided) {
1166: eps->twosided = twosided;
1167: eps->state = EPS_STATE_INITIAL;
1168: }
1169: PetscFunctionReturn(PETSC_SUCCESS);
1170: }
1172: /*@
1173: EPSGetTwoSided - Returns the flag indicating whether a two-sided variant
1174: of the algorithm is being used or not.
1176: Not Collective
1178: Input Parameter:
1179: . eps - the eigensolver context
1181: Output Parameter:
1182: . twosided - the returned flag
1184: Level: advanced
1186: .seealso: EPSSetTwoSided()
1187: @*/
1188: PetscErrorCode EPSGetTwoSided(EPS eps,PetscBool *twosided)
1189: {
1190: PetscFunctionBegin;
1192: PetscAssertPointer(twosided,2);
1193: *twosided = eps->twosided;
1194: PetscFunctionReturn(PETSC_SUCCESS);
1195: }
1197: /*@
1198: EPSSetTrueResidual - Specifies if the solver must compute the true residual
1199: explicitly or not.
1201: Logically Collective
1203: Input Parameters:
1204: + eps - the eigensolver context
1205: - trueres - whether true residuals are required or not
1207: Options Database Keys:
1208: . -eps_true_residual <boolean> - Sets/resets the boolean flag 'trueres'
1210: Notes:
1211: If the user sets trueres=PETSC_TRUE then the solver explicitly computes
1212: the true residual for each eigenpair approximation, and uses it for
1213: convergence testing. Computing the residual is usually an expensive
1214: operation. Some solvers (e.g., Krylov solvers) can avoid this computation
1215: by using a cheap estimate of the residual norm, but this may sometimes
1216: give inaccurate results (especially if a spectral transform is being
1217: used). On the contrary, preconditioned eigensolvers (e.g., Davidson solvers)
1218: do rely on computing the true residual, so this option is irrelevant for them.
1220: Level: advanced
1222: .seealso: EPSGetTrueResidual()
1223: @*/
1224: PetscErrorCode EPSSetTrueResidual(EPS eps,PetscBool trueres)
1225: {
1226: PetscFunctionBegin;
1229: eps->trueres = trueres;
1230: PetscFunctionReturn(PETSC_SUCCESS);
1231: }
1233: /*@
1234: EPSGetTrueResidual - Returns the flag indicating whether true
1235: residuals must be computed explicitly or not.
1237: Not Collective
1239: Input Parameter:
1240: . eps - the eigensolver context
1242: Output Parameter:
1243: . trueres - the returned flag
1245: Level: advanced
1247: .seealso: EPSSetTrueResidual()
1248: @*/
1249: PetscErrorCode EPSGetTrueResidual(EPS eps,PetscBool *trueres)
1250: {
1251: PetscFunctionBegin;
1253: PetscAssertPointer(trueres,2);
1254: *trueres = eps->trueres;
1255: PetscFunctionReturn(PETSC_SUCCESS);
1256: }
1258: /*@
1259: EPSSetTrackAll - Specifies if the solver must compute the residual norm of all
1260: approximate eigenpairs or not.
1262: Logically Collective
1264: Input Parameters:
1265: + eps - the eigensolver context
1266: - trackall - whether to compute all residuals or not
1268: Notes:
1269: If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
1270: the residual norm for each eigenpair approximation. Computing the residual is
1271: usually an expensive operation and solvers commonly compute only the residual
1272: associated to the first unconverged eigenpair.
1274: The option '-eps_monitor_all' automatically activates this option.
1276: Level: developer
1278: .seealso: EPSGetTrackAll()
1279: @*/
1280: PetscErrorCode EPSSetTrackAll(EPS eps,PetscBool trackall)
1281: {
1282: PetscFunctionBegin;
1285: eps->trackall = trackall;
1286: PetscFunctionReturn(PETSC_SUCCESS);
1287: }
1289: /*@
1290: EPSGetTrackAll - Returns the flag indicating whether all residual norms must
1291: be computed or not.
1293: Not Collective
1295: Input Parameter:
1296: . eps - the eigensolver context
1298: Output Parameter:
1299: . trackall - the returned flag
1301: Level: developer
1303: .seealso: EPSSetTrackAll()
1304: @*/
1305: PetscErrorCode EPSGetTrackAll(EPS eps,PetscBool *trackall)
1306: {
1307: PetscFunctionBegin;
1309: PetscAssertPointer(trackall,2);
1310: *trackall = eps->trackall;
1311: PetscFunctionReturn(PETSC_SUCCESS);
1312: }
1314: /*@
1315: EPSSetPurify - Deactivate eigenvector purification (which is activated by default).
1317: Logically Collective
1319: Input Parameters:
1320: + eps - the eigensolver context
1321: - purify - whether purification is required or not
1323: Options Database Keys:
1324: . -eps_purify <boolean> - Sets/resets the boolean flag 'purify'
1326: Notes:
1327: By default, eigenvectors of generalized symmetric eigenproblems are purified
1328: in order to purge directions in the nullspace of matrix B. If the user knows
1329: that B is non-singular, then purification can be safely deactivated and some
1330: computational cost is avoided (this is particularly important in interval computations).
1332: Level: intermediate
1334: .seealso: EPSGetPurify(), EPSSetInterval()
1335: @*/
1336: PetscErrorCode EPSSetPurify(EPS eps,PetscBool purify)
1337: {
1338: PetscFunctionBegin;
1341: if (purify!=eps->purify) {
1342: eps->purify = purify;
1343: eps->state = EPS_STATE_INITIAL;
1344: }
1345: PetscFunctionReturn(PETSC_SUCCESS);
1346: }
1348: /*@
1349: EPSGetPurify - Returns the flag indicating whether purification is activated
1350: or not.
1352: Not Collective
1354: Input Parameter:
1355: . eps - the eigensolver context
1357: Output Parameter:
1358: . purify - the returned flag
1360: Level: intermediate
1362: .seealso: EPSSetPurify()
1363: @*/
1364: PetscErrorCode EPSGetPurify(EPS eps,PetscBool *purify)
1365: {
1366: PetscFunctionBegin;
1368: PetscAssertPointer(purify,2);
1369: *purify = eps->purify;
1370: PetscFunctionReturn(PETSC_SUCCESS);
1371: }
1373: /*@C
1374: EPSSetOptionsPrefix - Sets the prefix used for searching for all
1375: EPS options in the database.
1377: Logically Collective
1379: Input Parameters:
1380: + eps - the eigensolver context
1381: - prefix - the prefix string to prepend to all EPS option requests
1383: Notes:
1384: A hyphen (-) must NOT be given at the beginning of the prefix name.
1385: The first character of all runtime options is AUTOMATICALLY the
1386: hyphen.
1388: For example, to distinguish between the runtime options for two
1389: different EPS contexts, one could call
1390: .vb
1391: EPSSetOptionsPrefix(eps1,"eig1_")
1392: EPSSetOptionsPrefix(eps2,"eig2_")
1393: .ve
1395: Level: advanced
1397: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
1398: @*/
1399: PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix)
1400: {
1401: PetscFunctionBegin;
1403: if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
1404: PetscCall(STSetOptionsPrefix(eps->st,prefix));
1405: if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
1406: PetscCall(BVSetOptionsPrefix(eps->V,prefix));
1407: if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
1408: PetscCall(DSSetOptionsPrefix(eps->ds,prefix));
1409: if (!eps->rg) PetscCall(EPSGetRG(eps,&eps->rg));
1410: PetscCall(RGSetOptionsPrefix(eps->rg,prefix));
1411: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)eps,prefix));
1412: PetscFunctionReturn(PETSC_SUCCESS);
1413: }
1415: /*@C
1416: EPSAppendOptionsPrefix - Appends to the prefix used for searching for all
1417: EPS options in the database.
1419: Logically Collective
1421: Input Parameters:
1422: + eps - the eigensolver context
1423: - prefix - the prefix string to prepend to all EPS option requests
1425: Notes:
1426: A hyphen (-) must NOT be given at the beginning of the prefix name.
1427: The first character of all runtime options is AUTOMATICALLY the hyphen.
1429: Level: advanced
1431: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
1432: @*/
1433: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix)
1434: {
1435: PetscFunctionBegin;
1437: if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
1438: PetscCall(STAppendOptionsPrefix(eps->st,prefix));
1439: if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
1440: PetscCall(BVAppendOptionsPrefix(eps->V,prefix));
1441: if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
1442: PetscCall(DSAppendOptionsPrefix(eps->ds,prefix));
1443: if (!eps->rg) PetscCall(EPSGetRG(eps,&eps->rg));
1444: PetscCall(RGAppendOptionsPrefix(eps->rg,prefix));
1445: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)eps,prefix));
1446: PetscFunctionReturn(PETSC_SUCCESS);
1447: }
1449: /*@C
1450: EPSGetOptionsPrefix - Gets the prefix used for searching for all
1451: EPS options in the database.
1453: Not Collective
1455: Input Parameters:
1456: . eps - the eigensolver context
1458: Output Parameters:
1459: . prefix - pointer to the prefix string used is returned
1461: Note:
1462: On the Fortran side, the user should pass in a string 'prefix' of
1463: sufficient length to hold the prefix.
1465: Level: advanced
1467: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
1468: @*/
1469: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])
1470: {
1471: PetscFunctionBegin;
1473: PetscAssertPointer(prefix,2);
1474: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)eps,prefix));
1475: PetscFunctionReturn(PETSC_SUCCESS);
1476: }