Actual source code: opts.c
1: /*
2: EPS routines related to options that can be set via the command-line
3: or procedurally.
5: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
6: SLEPc - Scalable Library for Eigenvalue Problem Computations
7: Copyright (c) 2002-2010, Universidad Politecnica de Valencia, Spain
9: This file is part of SLEPc.
10:
11: SLEPc is free software: you can redistribute it and/or modify it under the
12: terms of version 3 of the GNU Lesser General Public License as published by
13: the Free Software Foundation.
15: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
16: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
18: more details.
20: You should have received a copy of the GNU Lesser General Public License
21: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
22: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
23: */
25: #include private/epsimpl.h
29: /*@
30: EPSSetFromOptions - Sets EPS options from the options database.
31: This routine must be called before EPSSetUp() if the user is to be
32: allowed to set the solver type.
34: Collective on EPS
36: Input Parameters:
37: . eps - the eigensolver context
39: Notes:
40: To see all options, run your program with the -help option.
42: Level: beginner
43: @*/
44: PetscErrorCode EPSSetFromOptions(EPS eps)
45: {
47: char type[256],monfilename[PETSC_MAX_PATH_LEN];
48: PetscTruth flg,val;
49: PetscReal r,nrma,nrmb;
50: PetscScalar s;
51: PetscInt i,j,k;
52: const char *bal_list[4] = { "none", "oneside", "twoside","user" };
53: PetscViewerASCIIMonitor monviewer;
54: EPSMONITOR_CONV *ctx;
58: PetscOptionsBegin(((PetscObject)eps)->comm,((PetscObject)eps)->prefix,"Eigenproblem Solver (EPS) Options","EPS");
59: PetscOptionsList("-eps_type","Eigenproblem Solver method","EPSSetType",EPSList,(char*)(((PetscObject)eps)->type_name?((PetscObject)eps)->type_name:EPSKRYLOVSCHUR),type,256,&flg);
60: if (flg) {
61: EPSSetType(eps,type);
62: }
64: PetscOptionsTruthGroupBegin("-eps_hermitian","hermitian eigenvalue problem","EPSSetProblemType",&flg);
65: if (flg) {EPSSetProblemType(eps,EPS_HEP);}
66: PetscOptionsTruthGroup("-eps_gen_hermitian","generalized hermitian eigenvalue problem","EPSSetProblemType",&flg);
67: if (flg) {EPSSetProblemType(eps,EPS_GHEP);}
68: PetscOptionsTruthGroup("-eps_non_hermitian","non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
69: if (flg) {EPSSetProblemType(eps,EPS_NHEP);}
70: PetscOptionsTruthGroup("-eps_gen_non_hermitian","generalized non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
71: if (flg) {EPSSetProblemType(eps,EPS_GNHEP);}
72: PetscOptionsTruthGroupEnd("-eps_pos_gen_non_hermitian","generalized non-hermitian eigenvalue problem with positive semi-definite B","EPSSetProblemType",&flg);
73: if (flg) {EPSSetProblemType(eps,EPS_PGNHEP);}
75: /*
76: Set the type if it was never set.
77: */
78: if (!((PetscObject)eps)->type_name) {
79: EPSSetType(eps,EPSKRYLOVSCHUR);
80: }
82: PetscOptionsTruthGroupBegin("-eps_ritz","Rayleigh-Ritz extraction","EPSSetExtraction",&flg);
83: if (flg) {EPSSetExtraction(eps,EPS_RITZ);}
84: PetscOptionsTruthGroup("-eps_harmonic","harmonic Ritz extraction","EPSSetExtraction",&flg);
85: if (flg) {EPSSetExtraction(eps,EPS_HARMONIC);}
86: PetscOptionsTruthGroup("-eps_harmonic_relative","relative harmonic Ritz extraction","EPSSetExtraction",&flg);
87: if (flg) {EPSSetExtraction(eps,EPS_HARMONIC_RELATIVE);}
88: PetscOptionsTruthGroup("-eps_harmonic_right","right harmonic Ritz extraction","EPSSetExtraction",&flg);
89: if (flg) {EPSSetExtraction(eps,EPS_HARMONIC_RIGHT);}
90: PetscOptionsTruthGroup("-eps_harmonic_largest","largest harmonic Ritz extraction","EPSSetExtraction",&flg);
91: if (flg) {EPSSetExtraction(eps,EPS_HARMONIC_LARGEST);}
92: PetscOptionsTruthGroup("-eps_refined","refined Ritz extraction","EPSSetExtraction",&flg);
93: if (flg) {EPSSetExtraction(eps,EPS_REFINED);}
94: PetscOptionsTruthGroupEnd("-eps_refined_harmonic","refined harmonic Ritz extraction","EPSSetExtraction",&flg);
95: if (flg) {EPSSetExtraction(eps,EPS_REFINED_HARMONIC);}
97: if (!eps->balance) eps->balance = EPS_BALANCE_NONE;
98: PetscOptionsEList("-eps_balance", "Balancing method","EPSSetBalance",bal_list,4,bal_list[eps->balance-EPS_BALANCE_NONE],&i,&flg);
99: if (flg) { eps->balance = (EPSBalance)(i+EPS_BALANCE_NONE); }
100: r = j = PETSC_IGNORE;
101: PetscOptionsInt("-eps_balance_its","Number of iterations in balancing","EPSSetBalance",eps->balance_its,&j,PETSC_NULL);
102: PetscOptionsReal("-eps_balance_cutoff","Cutoff value in balancing","EPSSetBalance",eps->balance_cutoff,&r,PETSC_NULL);
103: EPSSetBalance(eps,(EPSBalance)PETSC_IGNORE,j,r);
105: r = i = PETSC_IGNORE;
106: PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,PETSC_NULL);
107: PetscOptionsReal("-eps_tol","Tolerance","EPSSetTolerances",eps->tol,&r,PETSC_NULL);
108: EPSSetTolerances(eps,r,i);
109: PetscOptionsTruthGroupBegin("-eps_conv_eig","relative error convergence test","EPSSetConvergenceTest",&flg);
110: if (flg) {EPSSetConvergenceTest(eps,EPS_CONV_EIG);}
111: PetscOptionsTruthGroup("-eps_conv_norm","Convergence test relative to the eigenvalue and the matrix norms","EPSSetConvergenceTest",&flg);
112: if (flg) {EPSSetConvergenceTest(eps,EPS_CONV_NORM);}
113: PetscOptionsTruthGroupEnd("-eps_conv_abs","Absolute error convergence test","EPSSetConvergenceTest",&flg);
114: if (flg) {EPSSetConvergenceTest(eps,EPS_CONV_ABS);}
116: i = j = k = PETSC_IGNORE;
117: PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,PETSC_NULL);
118: PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&j,PETSC_NULL);
119: PetscOptionsInt("-eps_mpd","Maximum dimension of projected problem","EPSSetDimensions",eps->mpd,&k,PETSC_NULL);
120: EPSSetDimensions(eps,i,j,k);
121:
122: /* -----------------------------------------------------------------------*/
123: /*
124: Cancels all monitors hardwired into code before call to EPSSetFromOptions()
125: */
126: flg = PETSC_FALSE;
127: PetscOptionsTruth("-eps_monitor_cancel","Remove any hardwired monitor routines","EPSMonitorCancel",flg,&flg,PETSC_NULL);
128: if (flg) {
129: EPSMonitorCancel(eps);
130: }
131: /*
132: Prints approximate eigenvalues and error estimates at each iteration
133: */
134: PetscOptionsString("-eps_monitor","Monitor first unconverged approximate eigenvalue and error estimate","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
135: if (flg) {
136: PetscViewerASCIIMonitorCreate(((PetscObject)eps)->comm,monfilename,((PetscObject)eps)->tablevel,&monviewer);
137: EPSMonitorSet(eps,EPSMonitorFirst,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
138: }
139: PetscOptionsString("-eps_monitor_conv","Monitor approximate eigenvalues and error estimates as they converge","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
140: if (flg) {
141: PetscNew(EPSMONITOR_CONV,&ctx);
142: PetscViewerASCIIMonitorCreate(((PetscObject)eps)->comm,monfilename,((PetscObject)eps)->tablevel,&ctx->viewer);
143: EPSMonitorSet(eps,EPSMonitorConverged,ctx,(PetscErrorCode (*)(void*))EPSMonitorDestroy_Converged);
144: }
145: PetscOptionsString("-eps_monitor_all","Monitor approximate eigenvalues and error estimates","EPSMonitorSet","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
146: if (flg) {
147: PetscViewerASCIIMonitorCreate(((PetscObject)eps)->comm,monfilename,((PetscObject)eps)->tablevel,&monviewer);
148: EPSMonitorSet(eps,EPSMonitorAll,monviewer,(PetscErrorCode (*)(void*))PetscViewerASCIIMonitorDestroy);
149: EPSSetTrackAll(eps,PETSC_TRUE);
150: }
151: flg = PETSC_FALSE;
152: PetscOptionsTruth("-eps_monitor_draw","Monitor first unconverged approximate eigenvalue and error estimate graphically","EPSMonitorSet",flg,&flg,PETSC_NULL);
153: if (flg) {
154: EPSMonitorSet(eps,EPSMonitorLG,PETSC_NULL,PETSC_NULL);
155: }
156: flg = PETSC_FALSE;
157: PetscOptionsTruth("-eps_monitor_draw_all","Monitor error estimates graphically","EPSMonitorSet",flg,&flg,PETSC_NULL);
158: if (flg) {
159: EPSMonitorSet(eps,EPSMonitorLGAll,PETSC_NULL,PETSC_NULL);
160: EPSSetTrackAll(eps,PETSC_TRUE);
161: }
162: /* -----------------------------------------------------------------------*/
163: PetscOptionsScalar("-eps_target","Value of the target","EPSSetTarget",eps->target,&s,&flg);
164: if (flg) {
165: EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);
166: EPSSetTarget(eps,s);
167: }
169: PetscOptionsTruthGroupBegin("-eps_largest_magnitude","compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
170: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE);}
171: PetscOptionsTruthGroup("-eps_smallest_magnitude","compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
172: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE);}
173: PetscOptionsTruthGroup("-eps_largest_real","compute largest real parts","EPSSetWhichEigenpairs",&flg);
174: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL);}
175: PetscOptionsTruthGroup("-eps_smallest_real","compute smallest real parts","EPSSetWhichEigenpairs",&flg);
176: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);}
177: PetscOptionsTruthGroup("-eps_largest_imaginary","compute largest imaginary parts","EPSSetWhichEigenpairs",&flg);
178: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY);}
179: PetscOptionsTruthGroup("-eps_smallest_imaginary","compute smallest imaginary parts","EPSSetWhichEigenpairs",&flg);
180: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY);}
181: PetscOptionsTruthGroup("-eps_target_magnitude","compute nearest eigenvalues to target","EPSSetWhichEigenpairs",&flg);
182: if (flg) {EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);}
183: PetscOptionsTruthGroup("-eps_target_real","compute eigenvalues with real parts close to target","EPSSetWhichEigenpairs",&flg);
184: if (flg) {EPSSetWhichEigenpairs(eps,EPS_TARGET_REAL);}
185: PetscOptionsTruthGroupEnd("-eps_target_imaginary","compute eigenvalues with imaginary parts close to target","EPSSetWhichEigenpairs",&flg);
186: if (flg) {EPSSetWhichEigenpairs(eps,EPS_TARGET_IMAGINARY);}
188: PetscOptionsTruth("-eps_left_vectors","Compute left eigenvectors also","EPSSetLeftVectorsWanted",eps->leftvecs,&val,&flg);
189: if (flg) {
190: EPSSetLeftVectorsWanted(eps,val);
191: }
192: PetscOptionsTruth("-eps_true_residual","Compute true residuals explicitly","EPSSetTrueResidual",eps->trueres,&val,&flg);
193: if (flg) {
194: EPSSetTrueResidual(eps,val);
195: }
197: nrma = nrmb = PETSC_IGNORE;
198: PetscOptionsReal("-eps_norm_a","Norm of matrix A","EPSSetMatrixNorms",eps->nrma,&nrma,PETSC_NULL);
199: PetscOptionsReal("-eps_norm_b","Norm of matrix B","EPSSetMatrixNorms",eps->nrmb,&nrmb,PETSC_NULL);
200: EPSSetMatrixNorms(eps,nrma,nrmb,eps->adaptive);
201: PetscOptionsTruth("-eps_norms_adaptive","Update the value of matrix norms adaptively","EPSSetMatrixNorms",eps->adaptive,&val,&flg);
202: if (flg) {
203: EPSSetMatrixNorms(eps,PETSC_IGNORE,PETSC_IGNORE,val);
204: }
206: PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",0);
207: PetscOptionsName("-eps_view_binary","Save the matrices associated to the eigenproblem","EPSSetFromOptions",0);
208: PetscOptionsName("-eps_plot_eigs","Make a plot of the computed eigenvalues","EPSSolve",0);
209:
210: if (eps->ops->setfromoptions) {
211: (*eps->ops->setfromoptions)(eps);
212: }
213: PetscOptionsEnd();
215: IPSetFromOptions(eps->ip);
216: STSetFromOptions(eps->OP);
217: return(0);
218: }
222: /*@
223: EPSGetTolerances - Gets the tolerance and maximum iteration count used
224: by the EPS convergence tests.
226: Not Collective
228: Input Parameter:
229: . eps - the eigensolver context
230:
231: Output Parameters:
232: + tol - the convergence tolerance
233: - maxits - maximum number of iterations
235: Notes:
236: The user can specify PETSC_NULL for any parameter that is not needed.
238: Level: intermediate
240: .seealso: EPSSetTolerances()
241: @*/
242: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,PetscInt *maxits)
243: {
246: if (tol) *tol = eps->tol;
247: if (maxits) *maxits = eps->max_it;
248: return(0);
249: }
253: /*@
254: EPSSetTolerances - Sets the tolerance and maximum iteration count used
255: by the EPS convergence tests.
257: Collective on EPS
259: Input Parameters:
260: + eps - the eigensolver context
261: . tol - the convergence tolerance
262: - maxits - maximum number of iterations to use
264: Options Database Keys:
265: + -eps_tol <tol> - Sets the convergence tolerance
266: - -eps_max_it <maxits> - Sets the maximum number of iterations allowed
268: Notes:
269: Use PETSC_IGNORE for an argument that need not be changed.
271: Use PETSC_DECIDE for maxits to assign a reasonably good value, which is
272: dependent on the solution method.
274: Level: intermediate
276: .seealso: EPSGetTolerances()
277: @*/
278: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,PetscInt maxits)
279: {
282: if (tol != PETSC_IGNORE) {
283: if (tol == PETSC_DEFAULT) {
284: eps->tol = 1e-7;
285: } else {
286: if (tol < 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
287: eps->tol = tol;
288: }
289: }
290: if (maxits != PETSC_IGNORE) {
291: if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
292: eps->max_it = 0;
293: eps->setupcalled = 0;
294: } else {
295: if (maxits < 0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
296: eps->max_it = maxits;
297: }
298: }
299: return(0);
300: }
304: /*@
305: EPSGetDimensions - Gets the number of eigenvalues to compute
306: and the dimension of the subspace.
308: Not Collective
310: Input Parameter:
311: . eps - the eigensolver context
312:
313: Output Parameters:
314: + nev - number of eigenvalues to compute
315: . ncv - the maximum dimension of the subspace to be used by the solver
316: - mpd - the maximum dimension allowed for the projected problem
318: Notes:
319: The user can specify PETSC_NULL for any parameter that is not needed.
321: Level: intermediate
323: .seealso: EPSSetDimensions()
324: @*/
325: PetscErrorCode EPSGetDimensions(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
326: {
329: if (nev) *nev = eps->nev;
330: if (ncv) *ncv = eps->ncv;
331: if (mpd) *mpd = eps->mpd;
332: return(0);
333: }
337: /*@
338: EPSSetDimensions - Sets the number of eigenvalues to compute
339: and the dimension of the subspace.
341: Collective on EPS
343: Input Parameters:
344: + eps - the eigensolver context
345: . nev - number of eigenvalues to compute
346: . ncv - the maximum dimension of the subspace to be used by the solver
347: - mpd - the maximum dimension allowed for the projected problem
349: Options Database Keys:
350: + -eps_nev <nev> - Sets the number of eigenvalues
351: . -eps_ncv <ncv> - Sets the dimension of the subspace
352: - -eps_mpd <mpd> - Sets the maximum projected dimension
354: Notes:
355: Use PETSC_IGNORE to retain the previous value of any parameter.
357: Use PETSC_DECIDE for ncv and mpd to assign a reasonably good value, which is
358: dependent on the solution method.
360: The parameters ncv and mpd are intimately related, so that the user is advised
361: to set one of them at most. Normal usage is the following
362: + - In cases where nev is small, the user sets ncv (a reasonable default is 2*nev).
363: - - In cases where nev is large, the user sets mpd.
365: The value of ncv should always be between nev and (nev+mpd), typically
366: ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
367: a smaller value should be used.
369: Level: intermediate
371: .seealso: EPSGetDimensions()
372: @*/
373: PetscErrorCode EPSSetDimensions(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)
374: {
378: if( nev != PETSC_IGNORE ) {
379: if (nev<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
380: eps->nev = nev;
381: eps->setupcalled = 0;
382: }
383: if( ncv != PETSC_IGNORE ) {
384: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
385: eps->ncv = 0;
386: } else {
387: if (ncv<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
388: eps->ncv = ncv;
389: }
390: eps->setupcalled = 0;
391: }
392: if( mpd != PETSC_IGNORE ) {
393: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
394: eps->mpd = 0;
395: } else {
396: if (mpd<1) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
397: eps->mpd = mpd;
398: }
399: }
400: return(0);
401: }
405: /*@
406: EPSSetWhichEigenpairs - Specifies which portion of the spectrum is
407: to be sought.
409: Collective on EPS
411: Input Parameters:
412: + eps - eigensolver context obtained from EPSCreate()
413: - which - the portion of the spectrum to be sought
415: Possible values:
416: The parameter 'which' can have one of these values
417:
418: + EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
419: . EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
420: . EPS_LARGEST_REAL - largest real parts
421: . EPS_SMALLEST_REAL - smallest real parts
422: . EPS_LARGEST_IMAGINARY - largest imaginary parts
423: . EPS_SMALLEST_IMAGINARY - smallest imaginary parts
424: . EPS_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
425: . EPS_TARGET_REAL - eigenvalues with real part closest to target
426: . EPS_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
427: - EPS_WHICH_USER - user defined ordering set with EPSSetEigenvalueComparison()
429: Options Database Keys:
430: + -eps_largest_magnitude - Sets largest eigenvalues in magnitude
431: . -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
432: . -eps_largest_real - Sets largest real parts
433: . -eps_smallest_real - Sets smallest real parts
434: . -eps_largest_imaginary - Sets largest imaginary parts
435: . -eps_smallest_imaginary - Sets smallest imaginary parts
436: . -eps_target_magnitude - Sets eigenvalues closest to target
437: . -eps_target_real - Sets real parts closest to target
438: - -eps_target_imaginary - Sets imaginary parts closest to target
440: Notes:
441: Not all eigensolvers implemented in EPS account for all the possible values
442: stated above. Also, some values make sense only for certain types of
443: problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
444: and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part
445: for eigenvalue selection.
446:
447: The target is a scalar value provided with EPSSetTarget().
449: The criterion EPS_TARGET_IMAGINARY is available only in case PETSc and
450: SLEPc have been built with complex scalars.
452: Level: intermediate
454: .seealso: EPSGetWhichEigenpairs(), EPSSetTarget(), EPSSetEigenvalueComparison(), EPSSortEigenvalues(), EPSWhich
455: @*/
456: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)
457: {
460: if (which!=PETSC_IGNORE) {
461: if (which==PETSC_DECIDE || which==PETSC_DEFAULT) eps->which = (EPSWhich)0;
462: else switch (which) {
463: case EPS_LARGEST_MAGNITUDE:
464: case EPS_SMALLEST_MAGNITUDE:
465: case EPS_LARGEST_REAL:
466: case EPS_SMALLEST_REAL:
467: case EPS_LARGEST_IMAGINARY:
468: case EPS_SMALLEST_IMAGINARY:
469: case EPS_TARGET_MAGNITUDE:
470: case EPS_TARGET_REAL:
471: #if defined(PETSC_USE_COMPLEX)
472: case EPS_TARGET_IMAGINARY:
473: #endif
474: case EPS_WHICH_USER:
475: if (eps->which != which) {
476: eps->setupcalled = 0;
477: eps->which = which;
478: }
479: break;
480: default:
481: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
482: }
483: }
484: return(0);
485: }
489: /*@C
490: EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be
491: sought.
493: Not Collective
495: Input Parameter:
496: . eps - eigensolver context obtained from EPSCreate()
498: Output Parameter:
499: . which - the portion of the spectrum to be sought
501: Notes:
502: See EPSSetWhichEigenpairs() for possible values of 'which'.
504: Level: intermediate
506: .seealso: EPSSetWhichEigenpairs(), EPSWhich
507: @*/
508: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
509: {
513: *which = eps->which;
514: return(0);
515: }
519: /*@
520: EPSSetLeftVectorsWanted - Specifies which eigenvectors are required.
522: Collective on EPS
524: Input Parameters:
525: + eps - the eigensolver context
526: - leftvecs - whether left eigenvectors are required or not
528: Options Database Keys:
529: . -eps_left_vectors <boolean> - Sets/resets the boolean flag 'leftvecs'
531: Notes:
532: If the user sets leftvecs=PETSC_TRUE then the solver uses a variant of
533: the algorithm that computes both right and left eigenvectors. This is
534: usually much more costly. This option is not available in all solvers.
536: Level: intermediate
538: .seealso: EPSGetLeftVectorsWanted(), EPSGetEigenvectorLeft()
539: @*/
540: PetscErrorCode EPSSetLeftVectorsWanted(EPS eps,PetscTruth leftvecs)
541: {
544: if (eps->leftvecs != leftvecs) {
545: eps->leftvecs = leftvecs;
546: eps->setupcalled = 0;
547: }
548: return(0);
549: }
553: /*@
554: EPSGetLeftVectorsWanted - Returns the flag indicating whether left
555: eigenvectors are required or not.
557: Not Collective
559: Input Parameter:
560: . eps - the eigensolver context
562: Output Parameter:
563: . leftvecs - the returned flag
565: Level: intermediate
567: .seealso: EPSSetLeftVectorsWanted(), EPSWhich
568: @*/
569: PetscErrorCode EPSGetLeftVectorsWanted(EPS eps,PetscTruth *leftvecs)
570: {
574: *leftvecs = eps->leftvecs;
575: return(0);
576: }
580: /*@
581: EPSSetMatrixNorms - Gives the reference values of the matrix norms
582: and specifies whether these values should be improved adaptively.
584: Collective on EPS
586: Input Parameters:
587: + eps - the eigensolver context
588: . nrma - a reference value for the norm of matrix A
589: . nrmb - a reference value for the norm of matrix B
590: - adaptive - whether matrix norms are improved adaptively
592: Options Database Keys:
593: + -eps_norm_a <nrma> - norm of A
594: . -eps_norm_b <nrma> - norm of B
595: - -eps_norms_adaptive <boolean> - Sets/resets the boolean flag 'adaptive'
597: Notes:
598: If the user sets adaptive=PETSC_FALSE then the solver uses the values
599: of nrma and nrmb for the matrix norms, and these values do not change
600: throughout the iteration.
602: If the user sets adaptive=PETSC_TRUE then the solver tries to adaptively
603: improve the supplied values, with the numerical information generated
604: during the iteration. This option is not available in all solvers.
606: If a passed value is PETSC_DEFAULT, the corresponding norm will be set to 1.
607: If a passed value is PETSC_DETERMINE, the corresponding norm will be computed
608: as the NORM_INFINITY with MatNorm().
610: Level: intermediate
612: .seealso: EPSGetMatrixNorms()
613: @*/
614: PetscErrorCode EPSSetMatrixNorms(EPS eps,PetscReal nrma,PetscReal nrmb,PetscTruth adaptive)
615: {
618: if (nrma != PETSC_IGNORE) {
619: if (nrma == PETSC_DEFAULT) eps->nrma = 1.0;
620: else if (nrma == PETSC_DETERMINE) {
621: eps->nrma = nrma;
622: eps->setupcalled = 0;
623: } else {
624: if (nrma < 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nrma. Must be > 0");
625: eps->nrma = nrma;
626: }
627: }
628: if (nrmb != PETSC_IGNORE) {
629: if (!eps->isgeneralized) SETERRQ(PETSC_ERR_ARG_WRONG,"Norm of B only allowed in generalized problems");
630: if (nrmb == PETSC_DEFAULT) eps->nrmb = 1.0;
631: else if (nrmb == PETSC_DETERMINE) {
632: eps->nrmb = nrmb;
633: eps->setupcalled = 0;
634: } else {
635: if (nrmb < 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nrmb. Must be > 0");
636: eps->nrmb = nrmb;
637: }
638: }
639: if (eps->adaptive != adaptive) {
640: eps->adaptive = adaptive;
641: eps->setupcalled = 0;
642: if (adaptive) SETERRQ(PETSC_ERR_SUP,"Sorry, adaptive norms are not implemented in this release.");
643: }
644: return(0);
645: }
649: /*@
650: EPSGetMatrixNorms - Returns the value of the matrix norms (either set
651: by the user or estimated by the solver) and the flag indicating whether
652: the norms are being adaptively improved.
654: Not Collective
656: Input Parameter:
657: . eps - the eigensolver context
659: Output Parameters:
660: + nrma - the norm of matrix A
661: . nrmb - the norm of matrix B
662: - adaptive - whether matrix norms are improved adaptively
664: Level: intermediate
666: .seealso: EPSSetMatrixNorms()
667: @*/
668: PetscErrorCode EPSGetMatrixNorms(EPS eps,PetscReal *nrma,PetscReal *nrmb,PetscTruth *adaptive)
669: {
675: if (nrma) *nrma = eps->nrma;
676: if (nrmb) *nrmb = eps->nrmb;
677: if (adaptive) *adaptive = eps->adaptive;
678: return(0);
679: }
683: /*@C
684: EPSSetEigenvalueComparison - Specifies the eigenvalue comparison function
685: when EPSSetWhichEigenpairs() is set to EPS_WHICH_USER.
687: Collective on EPS
689: Input Parameters:
690: + eps - eigensolver context obtained from EPSCreate()
691: . func - a pointer to the comparison function
692: - ctx - a context pointer (the last parameter to the comparison function)
694: Calling Sequence of func:
695: $ func(EPS eps,PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
697: + eps - eigensolver context obtained from EPSCreate()
698: . ar - real part of the 1st eigenvalue
699: . ai - imaginary part of the 1st eigenvalue
700: . br - real part of the 2nd eigenvalue
701: . bi - imaginary part of the 2nd eigenvalue
702: . res - result of comparison
703: - ctx - optional context, as set by EPSSetEigenvalueComparison()
705: Note:
706: The returning parameter 'res' can be:
707: + negative - if the 1st eigenvalue is preferred to the 2st one
708: . zero - if both eigenvalues are equally preferred
709: - positive - if the 2st eigenvalue is preferred to the 1st one
711: Level: advanced
713: .seealso: EPSSetWhichEigenpairs(), EPSSortEigenvalues(), EPSWhich
714: @*/
715: PetscErrorCode EPSSetEigenvalueComparison(EPS eps,PetscErrorCode (*func)(EPS,PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)
716: {
718: eps->which_func = func;
719: eps->which_ctx = ctx;
720: eps->which = EPS_WHICH_USER;
721: return(0);
722: }
726: /*@C
727: EPSSetConvergenceTestFunction - Sets a function to compute the error estimate
728: used in the convergence test.
730: Collective on EPS
732: Input Parameters:
733: + eps - eigensolver context obtained from EPSCreate()
734: . func - a pointer to the convergence test function
735: - ctx - a context pointer (the last parameter to the convergence test function)
737: Calling Sequence of func:
738: $ func(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
740: + eps - eigensolver context obtained from EPSCreate()
741: . eigr - real part of the eigenvalue
742: . eigi - imaginary part of the eigenvalue
743: . res - residual norm associated to the eigenpair
744: . errest - (output) computed error estimate
745: - ctx - optional context, as set by EPSSetConvergenceTest()
747: Note:
748: If the error estimate returned by the convergence test function is less than
749: the tolerance, then the eigenvalue is accepted as converged.
751: Level: advanced
753: .seealso: EPSSetConvergenceTest(),EPSSetTolerances()
754: @*/
755: EXTERN PetscErrorCode EPSSetConvergenceTestFunction(EPS eps,PetscErrorCode (*func)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx)
756: {
759: eps->conv_func = func;
760: eps->conv_ctx = ctx;
761: if (func == EPSEigRelativeConverged) eps->conv = EPS_CONV_EIG;
762: else if (func == EPSNormRelativeConverged) eps->conv = EPS_CONV_NORM;
763: else if (func == EPSAbsoluteConverged) eps->conv = EPS_CONV_ABS;
764: else eps->conv = EPS_CONV_USER;
765: return(0);
766: }
770: /*@
771: EPSSetConvergenceTest - Specifies how to compute the error estimate
772: used in the convergence test.
774: Collective on EPS
776: Input Parameters:
777: + eps - eigensolver context obtained from EPSCreate()
778: - conv - the type of convergence test
780: Options Database Keys:
781: + -eps_conv_abs - Sets the absolute convergence test
782: . -eps_conv_eig - Sets the convergence test relative to the eigenvalue
783: - -eps_conv_norm - Sets the convergence test relative to the matrix norms
785: Note:
786: The parameter 'conv' can have one of these values
787: + EPS_CONV_ABS - absolute error ||r||
788: . EPS_CONV_EIG - error relative to the eigenvalue l, ||r||/|l|
789: . EPS_CONV_NORM - error relative to the matrix norms, ||r||/(||A||+|l|*||B||)
790: - EPS_CONV_USER - function set by EPSSetConvergenceTestFunction()
792: Level: intermediate
794: .seealso: EPSGetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSConv
795: @*/
796: PetscErrorCode EPSSetConvergenceTest(EPS eps,EPSConv conv)
797: {
800: switch(conv) {
801: case EPS_CONV_EIG: eps->conv_func = EPSEigRelativeConverged; break;
802: case EPS_CONV_NORM: eps->conv_func = EPSNormRelativeConverged; break;
803: case EPS_CONV_ABS: eps->conv_func = EPSAbsoluteConverged; break;
804: case EPS_CONV_USER: break;
805: default:
806: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
807: }
808: eps->conv = conv;
809: return(0);
810: }
814: /*@
815: EPSGetConvergenceTest - Gets the method used to compute the error estimate
816: used in the convergence test.
818: Collective on EPS
820: Input Parameters:
821: . eps - eigensolver context obtained from EPSCreate()
823: Output Parameters:
824: . conv - the type of convergence test
826: Level: intermediate
828: .seealso: EPSSetConvergenceTest(), EPSSetConvergenceTestFunction(), EPSConv
829: @*/
830: PetscErrorCode EPSGetConvergenceTest(EPS eps,EPSConv *conv)
831: {
835: *conv = eps->conv;
836: return(0);
837: }
842: /*@
843: EPSSetProblemType - Specifies the type of the eigenvalue problem.
845: Collective on EPS
847: Input Parameters:
848: + eps - the eigensolver context
849: - type - a known type of eigenvalue problem
851: Options Database Keys:
852: + -eps_hermitian - Hermitian eigenvalue problem
853: . -eps_gen_hermitian - generalized Hermitian eigenvalue problem
854: . -eps_non_hermitian - non-Hermitian eigenvalue problem
855: . -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
856: - -eps_pos_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
857: with positive semi-definite B
858:
859: Notes:
860: Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
861: (EPS_NHEP), generalized Hermitian (EPS_GHEP), generalized non-Hermitian
862: (EPS_GNHEP), and generalized non-Hermitian with positive semi-definite B
863: (EPS_PGNHEP)
865: This function must be used to instruct SLEPc to exploit symmetry. If no
866: problem type is specified, by default a non-Hermitian problem is assumed
867: (either standard or generalized). If the user knows that the problem is
868: Hermitian (i.e. A=A^H) or generalized Hermitian (i.e. A=A^H, B=B^H, and
869: B positive definite) then it is recommended to set the problem type so
870: that eigensolver can exploit these properties.
872: Level: beginner
874: .seealso: EPSSetOperators(), EPSSetType(), EPSGetProblemType(), EPSProblemType
875: @*/
876: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
877: {
881: switch (type) {
882: case EPS_HEP:
883: eps->isgeneralized = PETSC_FALSE;
884: eps->ishermitian = PETSC_TRUE;
885: eps->ispositive = PETSC_FALSE;
886: break;
887: case EPS_NHEP:
888: eps->isgeneralized = PETSC_FALSE;
889: eps->ishermitian = PETSC_FALSE;
890: eps->ispositive = PETSC_FALSE;
891: break;
892: case EPS_GHEP:
893: eps->isgeneralized = PETSC_TRUE;
894: eps->ishermitian = PETSC_TRUE;
895: eps->ispositive = PETSC_TRUE;
896: break;
897: case EPS_GNHEP:
898: eps->isgeneralized = PETSC_TRUE;
899: eps->ishermitian = PETSC_FALSE;
900: eps->ispositive = PETSC_FALSE;
901: break;
902: case EPS_PGNHEP:
903: eps->isgeneralized = PETSC_TRUE;
904: eps->ishermitian = PETSC_FALSE;
905: eps->ispositive = PETSC_TRUE;
906: break;
907: /*
908: case EPS_CSEP:
909: eps->isgeneralized = PETSC_FALSE;
910: eps->ishermitian = PETSC_FALSE;
911: STSetBilinearForm(eps->OP,STINNER_SYMMETRIC);
912: break;
913: case EPS_GCSEP:
914: eps->isgeneralized = PETSC_TRUE;
915: eps->ishermitian = PETSC_FALSE;
916: STSetBilinearForm(eps->OP,STINNER_B_SYMMETRIC);
917: break;
918: */
919: default:
920: SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
921: }
922: eps->problem_type = type;
924: return(0);
925: }
929: /*@C
930: EPSGetProblemType - Gets the problem type from the EPS object.
932: Not Collective
934: Input Parameter:
935: . eps - the eigensolver context
937: Output Parameter:
938: . type - name of EPS problem type
940: Level: intermediate
942: .seealso: EPSSetProblemType(), EPSProblemType
943: @*/
944: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
945: {
949: *type = eps->problem_type;
950: return(0);
951: }
955: /*@
956: EPSSetExtraction - Specifies the type of extraction technique to be employed
957: by the eigensolver.
959: Collective on EPS
961: Input Parameters:
962: + eps - the eigensolver context
963: - extr - a known type of extraction
965: Options Database Keys:
966: + -eps_ritz - Rayleigh-Ritz extraction
967: . -eps_harmonic - harmonic Ritz extraction
968: . -eps_harmonic_relative - harmonic Ritz extraction relative to the eigenvalue
969: . -eps_harmonic_right - harmonic Ritz extraction for rightmost eigenvalues
970: . -eps_harmonic_largest - harmonic Ritz extraction for largest magnitude
971: (without target)
972: . -eps_refined - refined Ritz extraction
973: - -eps_refined_harmonic - refined harmonic Ritz extraction
974:
975: Notes:
976: Not all eigensolvers support all types of extraction. See the SLEPc
977: Users Manual for details.
979: By default, a standard Rayleigh-Ritz extraction is used. Other extractions
980: may be useful when computing interior eigenvalues.
982: Harmonic-type extractions are used in combination with a 'target'.
984: Level: beginner
986: .seealso: EPSSetTarget(), EPSGetExtraction(), EPSExtraction
987: @*/
988: PetscErrorCode EPSSetExtraction(EPS eps,EPSExtraction extr)
989: {
992: eps->extraction = extr;
993: return(0);
994: }
998: /*@C
999: EPSGetExtraction - Gets the extraction type used by the EPS object.
1001: Not Collective
1003: Input Parameter:
1004: . eps - the eigensolver context
1006: Output Parameter:
1007: . extr - name of extraction type
1009: Level: intermediate
1011: .seealso: EPSSetExtraction(), EPSExtraction
1012: @*/
1013: PetscErrorCode EPSGetExtraction(EPS eps,EPSExtraction *extr)
1014: {
1018: *extr = eps->extraction;
1019: return(0);
1020: }
1024: /*@
1025: EPSSetBalance - Specifies the balancing technique to be employed by the
1026: eigensolver, and some parameters associated to it.
1028: Collective on EPS
1030: Input Parameters:
1031: + eps - the eigensolver context
1032: . bal - the balancing method, one of EPS_BALANCE_NONE, EPS_BALANCE_ONESIDE,
1033: EPS_BALANCE_TWOSIDE, or EPS_BALANCE_USER
1034: . its - number of iterations of the balancing algorithm
1035: - cutoff - cutoff value
1037: Options Database Keys:
1038: + -eps_balance <method> - the balancing method, where <method> is one of
1039: 'none', 'oneside', 'twoside', or 'user'
1040: . -eps_balance_its <its> - number of iterations
1041: - -eps_balance_cutoff <cutoff> - cutoff value
1042:
1043: Notes:
1044: When balancing is enabled, the solver works implicitly with matrix DAD^-1,
1045: where D is an appropriate diagonal matrix. This improves the accuracy of
1046: the computed results in some cases. See the SLEPc Users Manual for details.
1048: Balancing makes sense only for non-Hermitian problems when the required
1049: precision is high (i.e. a small tolerance such as 1e-15).
1051: By default, balancing is disabled. The two-sided method is much more
1052: effective than the one-sided counterpart, but it requires the system
1053: matrices to have the MatMultTranspose operation defined.
1055: The parameter 'its' is the number of iterations performed by the method. The
1056: cutoff value is used only in the two-side variant. Use PETSC_IGNORE for an
1057: argument that need not be changed. Use PETSC_DECIDE to assign a reasonably
1058: good value.
1060: User-defined balancing is allowed provided that the corresponding matrix
1061: is set via STSetBalanceMatrix.
1063: Level: intermediate
1065: .seealso: EPSGetBalance(), EPSBalance, STSetBalanceMatrix()
1066: @*/
1067: PetscErrorCode EPSSetBalance(EPS eps,EPSBalance bal,PetscInt its,PetscReal cutoff)
1068: {
1071: if (bal!=PETSC_IGNORE) {
1072: if (bal==PETSC_DECIDE || bal==PETSC_DEFAULT) eps->balance = EPS_BALANCE_TWOSIDE;
1073: else switch (bal) {
1074: case EPS_BALANCE_NONE:
1075: case EPS_BALANCE_ONESIDE:
1076: case EPS_BALANCE_TWOSIDE:
1077: case EPS_BALANCE_USER:
1078: eps->balance = bal;
1079: break;
1080: default:
1081: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Invalid value of argument 'bal'");
1082: }
1083: }
1084: if (its!=PETSC_IGNORE) {
1085: if (its==PETSC_DECIDE || its==PETSC_DEFAULT) eps->balance_its = 5;
1086: else eps->balance_its = its;
1087: }
1088: if (cutoff!=PETSC_IGNORE) {
1089: if (cutoff==PETSC_DECIDE || cutoff==PETSC_DEFAULT) eps->balance_cutoff = 1e-8;
1090: else eps->balance_cutoff = cutoff;
1091: }
1092: return(0);
1093: }
1097: /*@
1098: EPSGetBalance - Gets the balancing type used by the EPS object, and the associated
1099: parameters.
1101: Not Collective
1103: Input Parameter:
1104: . eps - the eigensolver context
1106: Output Parameters:
1107: + bal - the balancing method
1108: . its - number of iterations of the balancing algorithm
1109: - cutoff - cutoff value
1111: Level: intermediate
1113: Note:
1114: The user can specify PETSC_NULL for any parameter that is not needed.
1116: .seealso: EPSSetBalance(), EPSBalance
1117: @*/
1118: PetscErrorCode EPSGetBalance(EPS eps,EPSBalance *bal,PetscInt *its,PetscReal *cutoff)
1119: {
1125: if (bal) *bal = eps->balance;
1126: if (its) *its = eps->balance_its;
1127: if (cutoff) *cutoff = eps->balance_cutoff;
1128: return(0);
1129: }
1133: /*@
1134: EPSSetTrueResidual - Specifies if the solver must compute the true residual
1135: explicitly or not.
1137: Collective on EPS
1139: Input Parameters:
1140: + eps - the eigensolver context
1141: - trueres - whether true residuals are required or not
1143: Options Database Keys:
1144: . -eps_true_residual <boolean> - Sets/resets the boolean flag 'trueres'
1146: Notes:
1147: If the user sets trueres=PETSC_TRUE then the solver explicitly computes
1148: the true residual for each eigenpair approximation, and uses it for
1149: convergence testing. Computing the residual is usually an expensive
1150: operation. Some solvers (e.g., Krylov solvers) can avoid this computation
1151: by using a cheap estimate of the residual norm, but this may sometimes
1152: give inaccurate results (especially if a spectral transform is being
1153: used). On the contrary, preconditioned eigensolvers (e.g., Davidson solvers)
1154: do rely on computing the true residual, so this option is irrelevant for them.
1156: Level: intermediate
1158: .seealso: EPSGetTrueResidual()
1159: @*/
1160: PetscErrorCode EPSSetTrueResidual(EPS eps,PetscTruth trueres)
1161: {
1164: eps->trueres = trueres;
1165: return(0);
1166: }
1170: /*@C
1171: EPSGetTrueResidual - Returns the flag indicating whether true
1172: residuals must be computed explicitly or not.
1174: Not Collective
1176: Input Parameter:
1177: . eps - the eigensolver context
1179: Output Parameter:
1180: . trueres - the returned flag
1182: Level: intermediate
1184: .seealso: EPSSetTrueResidual()
1185: @*/
1186: PetscErrorCode EPSGetTrueResidual(EPS eps,PetscTruth *trueres)
1187: {
1191: *trueres = eps->trueres;
1192: return(0);
1193: }
1197: /*@
1198: EPSSetTrackAll - Specifies if the solver must compute the residual norm of all
1199: approximate eigenpairs or not.
1201: Collective on EPS
1203: Input Parameters:
1204: + eps - the eigensolver context
1205: - trackall - whether to compute all residuals or not
1207: Notes:
1208: If the user sets trackall=PETSC_TRUE then the solver computes (or estimates)
1209: the residual norm for each eigenpair approximation. Computing the residual is
1210: usually an expensive operation and solvers commonly compute only the residual
1211: associated to the first unconverged eigenpair.
1213: The options '-eps_monitor_all' and '-eps_monitor_draw_all' automatically
1214: activate this option.
1216: Level: intermediate
1218: .seealso: EPSGetTrackAll()
1219: @*/
1220: PetscErrorCode EPSSetTrackAll(EPS eps,PetscTruth trackall)
1221: {
1224: eps->trackall = trackall;
1225: return(0);
1226: }
1230: /*@
1231: EPSGetTrackAll - Returns the flag indicating whether all residual norms must
1232: be computed or not.
1234: Not Collective
1236: Input Parameter:
1237: . eps - the eigensolver context
1239: Output Parameter:
1240: . trackall - the returned flag
1242: Level: intermediate
1244: .seealso: EPSSetTrackAll()
1245: @*/
1246: PetscErrorCode EPSGetTrackAll(EPS eps,PetscTruth *trackall)
1247: {
1251: *trackall = eps->trackall;
1252: return(0);
1253: }
1257: /*@C
1258: EPSSetOptionsPrefix - Sets the prefix used for searching for all
1259: EPS options in the database.
1261: Collective on EPS
1263: Input Parameters:
1264: + eps - the eigensolver context
1265: - prefix - the prefix string to prepend to all EPS option requests
1267: Notes:
1268: A hyphen (-) must NOT be given at the beginning of the prefix name.
1269: The first character of all runtime options is AUTOMATICALLY the
1270: hyphen.
1272: For example, to distinguish between the runtime options for two
1273: different EPS contexts, one could call
1274: .vb
1275: EPSSetOptionsPrefix(eps1,"eig1_")
1276: EPSSetOptionsPrefix(eps2,"eig2_")
1277: .ve
1279: Level: advanced
1281: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
1282: @*/
1283: PetscErrorCode EPSSetOptionsPrefix(EPS eps,const char *prefix)
1284: {
1288: PetscObjectSetOptionsPrefix((PetscObject)eps, prefix);
1289: STSetOptionsPrefix(eps->OP,prefix);
1290: IPSetOptionsPrefix(eps->ip,prefix);
1291: IPAppendOptionsPrefix(eps->ip,"eps_");
1292: return(0);
1293: }
1294:
1297: /*@C
1298: EPSAppendOptionsPrefix - Appends to the prefix used for searching for all
1299: EPS options in the database.
1301: Collective on EPS
1303: Input Parameters:
1304: + eps - the eigensolver context
1305: - prefix - the prefix string to prepend to all EPS option requests
1307: Notes:
1308: A hyphen (-) must NOT be given at the beginning of the prefix name.
1309: The first character of all runtime options is AUTOMATICALLY the hyphen.
1311: Level: advanced
1313: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
1314: @*/
1315: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,const char *prefix)
1316: {
1320: PetscObjectAppendOptionsPrefix((PetscObject)eps, prefix);
1321: STAppendOptionsPrefix(eps->OP,prefix);
1322: IPSetOptionsPrefix(eps->ip,prefix);
1323: IPAppendOptionsPrefix(eps->ip,"eps_");
1324: return(0);
1325: }
1329: /*@C
1330: EPSGetOptionsPrefix - Gets the prefix used for searching for all
1331: EPS options in the database.
1333: Not Collective
1335: Input Parameters:
1336: . eps - the eigensolver context
1338: Output Parameters:
1339: . prefix - pointer to the prefix string used is returned
1341: Notes: On the fortran side, the user should pass in a string 'prefix' of
1342: sufficient length to hold the prefix.
1344: Level: advanced
1346: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
1347: @*/
1348: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])
1349: {
1354: PetscObjectGetOptionsPrefix((PetscObject)eps, prefix);
1355: return(0);
1356: }