Actual source code: opts.c
1: /*
2: EPS routines related to options that can be set via the command-line
3: or procedurally.
4: */
5: #include src/eps/epsimpl.h
9: /*@
10: EPSSetFromOptions - Sets EPS options from the options database.
11: This routine must be called before EPSSetUp() if the user is to be
12: allowed to set the solver type.
14: Collective on EPS
16: Input Parameters:
17: . eps - the eigensolver context
19: Notes:
20: To see all options, run your program with the -help option.
22: Level: beginner
24: .seealso:
25: @*/
26: PetscErrorCode EPSSetFromOptions(EPS eps)
27: {
29: char type[256],monfilename[PETSC_MAX_PATH_LEN];
30: PetscTruth flg;
31: const char *orth_list[3] = { "mgs" , "cgs", "ncgs" };
32: const char *ref_list[3] = { "never" , "ifneeded", "always" };
33: PetscReal eta;
34: PetscInt i,orth_type,ref_type;
35: PetscViewer monviewer;
39: PetscOptionsBegin(eps->comm,eps->prefix,"Eigenproblem Solver (EPS) Options","EPS");
40: PetscOptionsList("-eps_type","Eigenproblem Solver method","EPSSetType",EPSList,(char*)(eps->type_name?eps->type_name:EPSARNOLDI),type,256,&flg);
41: if (flg) {
42: EPSSetType(eps,type);
43: }
45: PetscOptionsTruthGroupBegin("-eps_hermitian","hermitian eigenvalue problem","EPSSetProblemType",&flg);
46: if (flg) {EPSSetProblemType(eps,EPS_HEP);}
47: PetscOptionsTruthGroup("-eps_gen_hermitian","generalized hermitian eigenvalue problem","EPSSetProblemType",&flg);
48: if (flg) {EPSSetProblemType(eps,EPS_GHEP);}
49: PetscOptionsTruthGroup("-eps_non_hermitian","non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
50: if (flg) {EPSSetProblemType(eps,EPS_NHEP);}
51: PetscOptionsTruthGroupEnd("-eps_gen_non_hermitian","generalized non-hermitian eigenvalue problem","EPSSetProblemType",&flg);
52: if (flg) {EPSSetProblemType(eps,EPS_GNHEP);}
54: /*
55: Set the type if it was never set.
56: */
57: if (!eps->type_name) {
58: if (eps->ishermitian) {
59: EPSSetType(eps,EPSLANCZOS);
60: } else {
61: EPSSetType(eps,EPSARNOLDI);
62: }
63: }
65: PetscOptionsTruthGroupBegin("-eps_oneside","one-sided eigensolver","EPSSetClass",&flg);
66: if (flg) {EPSSetClass(eps,EPS_ONE_SIDE);}
67: PetscOptionsTruthGroupEnd("-eps_twoside","two-sided eigensolver","EPSSetClass",&flg);
68: if (flg) {EPSSetClass(eps,EPS_TWO_SIDE);}
70: orth_type = eps->orthog_type;
71: PetscOptionsEList("-eps_orthog_type","Orthogonalization method","EPSSetOrthogonalization",orth_list,3,orth_list[eps->orthog_type],&orth_type,&flg);
72: ref_type = eps->orthog_ref;
73: PetscOptionsEList("-eps_orthog_refinement","Iterative refinement mode during orthogonalization","EPSSetOrthogonalizationRefinement",ref_list,3,ref_list[eps->orthog_ref],&ref_type,&flg);
74: eta = eps->orthog_eta;
75: PetscOptionsReal("-eps_orthog_eta","Parameter of iterative refinement during orthogonalization","EPSSetOrthogonalizationRefinement",eta,&eta,PETSC_NULL);
76: EPSSetOrthogonalization(eps,(EPSOrthogonalizationType)orth_type,(EPSOrthogonalizationRefinementType)ref_type,eta);
78: PetscOptionsInt("-eps_max_it","Maximum number of iterations","EPSSetTolerances",eps->max_it,&i,&flg);
79: if (flg) eps->max_it = i;
80: PetscOptionsReal("-eps_tol","Tolerance","KSPSetTolerances",eps->tol,&eps->tol,PETSC_NULL);
81: PetscOptionsInt("-eps_nev","Number of eigenvalues to compute","EPSSetDimensions",eps->nev,&i,&flg);
82: if (flg) {
83: if(i<1) SETERRQ(1,"Illegal value for option -eps_nev. Must be > 0");
84: eps->nev = i;
85: }
86: PetscOptionsInt("-eps_ncv","Number of basis vectors","EPSSetDimensions",eps->ncv,&i,&flg);
87: if (flg) {
88: if (i<1) SETERRQ(1,"Illegal value for option -eps_ncv. Must be > 0");
89: eps->ncv = i;
90: }
92: /* -----------------------------------------------------------------------*/
93: /*
94: Cancels all monitors hardwired into code before call to EPSSetFromOptions()
95: */
96: PetscOptionsName("-eps_cancelmonitors","Remove any hardwired monitor routines","EPSClearMonitor",&flg);
97: if (flg) {
98: EPSClearMonitor(eps);
99: }
100: /*
101: Prints approximate eigenvalues and error estimates at each iteration
102: */
103: PetscOptionsString("-eps_monitor","Monitor approximate eigenvalues and error estimates","EPSSetMonitor","stdout",monfilename,PETSC_MAX_PATH_LEN,&flg);
104: if (flg) {
105: PetscViewerASCIIOpen(eps->comm,monfilename,&monviewer);
106: EPSSetMonitor(eps,EPSDefaultMonitor,monviewer);
107: }
108: PetscOptionsName("-eps_xmonitor","Monitor error estimates graphically","EPSSetMonitor",&flg);
109: if (flg) {
110: EPSSetMonitor(eps,EPSLGMonitor,PETSC_NULL);
111: }
112: /* -----------------------------------------------------------------------*/
113: PetscOptionsTruthGroupBegin("-eps_largest_magnitude","compute largest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
114: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE);}
115: PetscOptionsTruthGroup("-eps_smallest_magnitude","compute smallest eigenvalues in magnitude","EPSSetWhichEigenpairs",&flg);
116: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_MAGNITUDE);}
117: PetscOptionsTruthGroup("-eps_largest_real","compute largest real parts","EPSSetWhichEigenpairs",&flg);
118: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_REAL);}
119: PetscOptionsTruthGroup("-eps_smallest_real","compute smallest real parts","EPSSetWhichEigenpairs",&flg);
120: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);}
121: PetscOptionsTruthGroup("-eps_largest_imaginary","compute largest imaginary parts","EPSSetWhichEigenpairs",&flg);
122: if (flg) {EPSSetWhichEigenpairs(eps,EPS_LARGEST_IMAGINARY);}
123: PetscOptionsTruthGroupEnd("-eps_smallest_imaginary","compute smallest imaginary parts","EPSSetWhichEigenpairs",&flg);
124: if (flg) {EPSSetWhichEigenpairs(eps,EPS_SMALLEST_IMAGINARY);}
126: PetscOptionsName("-eps_view","Print detailed information on solver used","EPSView",0);
127: PetscOptionsName("-eps_view_binary","Save the matrices associated to the eigenproblem","EPSSetFromOptions",0);
128: PetscOptionsName("-eps_plot_eigs","Make a plot of the computed eigenvalues","EPSSolve",0);
130: if (eps->ops->setfromoptions) {
131: (*eps->ops->setfromoptions)(eps);
132: }
133: PetscOptionsEnd();
135: STSetFromOptions(eps->OP);
136: return(0);
137: }
141: /*@
142: EPSGetTolerances - Gets the tolerance and maximum
143: iteration count used by the default EPS convergence tests.
145: Not Collective
147: Input Parameter:
148: . eps - the eigensolver context
149:
150: Output Parameters:
151: + tol - the convergence tolerance
152: - maxits - maximum number of iterations
154: Notes:
155: The user can specify PETSC_NULL for any parameter that is not needed.
157: Level: intermediate
159: .seealso: EPSSetTolerances()
160: @*/
161: PetscErrorCode EPSGetTolerances(EPS eps,PetscReal *tol,int *maxits)
162: {
165: if (tol) *tol = eps->tol;
166: if (maxits) *maxits = eps->max_it;
167: return(0);
168: }
172: /*@
173: EPSSetTolerances - Sets the tolerance and maximum
174: iteration count used by the default EPS convergence testers.
176: Collective on EPS
178: Input Parameters:
179: + eps - the eigensolver context
180: . tol - the convergence tolerance
181: - maxits - maximum number of iterations to use
183: Options Database Keys:
184: + -eps_tol <tol> - Sets the convergence tolerance
185: - -eps_max_it <maxits> - Sets the maximum number of iterations allowed
187: Notes:
188: Use PETSC_DEFAULT to retain the default value of any of the tolerances.
190: Level: intermediate
192: .seealso: EPSGetTolerances()
193: @*/
194: PetscErrorCode EPSSetTolerances(EPS eps,PetscReal tol,int maxits)
195: {
198: if (tol != PETSC_DEFAULT) eps->tol = tol;
199: if (maxits != PETSC_DEFAULT) eps->max_it = maxits;
200: return(0);
201: }
205: /*@
206: EPSGetDimensions - Gets the number of eigenvalues to compute
207: and the dimension of the subspace.
209: Not Collective
211: Input Parameter:
212: . eps - the eigensolver context
213:
214: Output Parameters:
215: + nev - number of eigenvalues to compute
216: - ncv - the maximum dimension of the subspace to be used by the solver
218: Notes:
219: The user can specify PETSC_NULL for any parameter that is not needed.
221: Level: intermediate
223: .seealso: EPSSetDimensions()
224: @*/
225: PetscErrorCode EPSGetDimensions(EPS eps,int *nev,int *ncv)
226: {
229: if( nev ) *nev = eps->nev;
230: if( ncv ) *ncv = eps->ncv;
231: return(0);
232: }
236: /*@
237: EPSSetDimensions - Sets the number of eigenvalues to compute
238: and the dimension of the subspace.
240: Collective on EPS
242: Input Parameters:
243: + eps - the eigensolver context
244: . nev - number of eigenvalues to compute
245: - ncv - the maximum dimension of the subspace to be used by the solver
247: Options Database Keys:
248: + -eps_nev <nev> - Sets the number of eigenvalues
249: - -eps_ncv <ncv> - Sets the dimension of the subspace
251: Notes:
252: Use PETSC_DEFAULT to retain the previous value of any parameter.
254: Use PETSC_DECIDE for ncv to assign a reasonably good value, which is
255: dependent on the solution method.
257: Level: intermediate
259: .seealso: EPSGetDimensions()
260: @*/
261: PetscErrorCode EPSSetDimensions(EPS eps,int nev,int ncv)
262: {
266: if( nev != PETSC_DEFAULT ) {
267: if (nev<1) SETERRQ(1,"Illegal value of nev. Must be > 0");
268: eps->nev = nev;
269: eps->setupcalled = 0;
270: }
271: if( ncv != PETSC_DEFAULT ) {
272: if( ncv == PETSC_DECIDE ) eps->ncv = 0;
273: else {
274: if (ncv<1) SETERRQ(1,"Illegal value of ncv. Must be > 0");
275: eps->ncv = ncv;
276: }
277: eps->setupcalled = 0;
278: }
279: return(0);
280: }
284: /*@
285: EPSSetWhichEigenpairs - Specifies which portion of the spectrum is
286: to be sought.
288: Collective on EPS
290: Input Parameter:
291: . eps - eigensolver context obtained from EPSCreate()
293: Output Parameter:
294: . which - the portion of the spectrum to be sought
296: Possible values:
297: The parameter 'which' can have one of these values:
298:
299: + EPS_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
300: . EPS_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
301: . EPS_LARGEST_REAL - largest real parts
302: . EPS_SMALLEST_REAL - smallest real parts
303: . EPS_LARGEST_IMAGINARY - largest imaginary parts
304: - EPS_SMALLEST_IMAGINARY - smallest imaginary parts
306: Options Database Keys:
307: + -eps_largest_magnitude - Sets largest eigenvalues in magnitude
308: . -eps_smallest_magnitude - Sets smallest eigenvalues in magnitude
309: . -eps_largest_real - Sets largest real parts
310: . -eps_smallest_real - Sets smallest real parts
311: . -eps_largest_imaginary - Sets largest imaginary parts in magnitude
312: - -eps_smallest_imaginary - Sets smallest imaginary parts in magnitude
314: Notes:
315: Not all eigensolvers implemented in EPS account for all the possible values
316: stated above. Also, some values make sense only for certain types of
317: problems. If SLEPc is compiled for real numbers EPS_LARGEST_IMAGINARY
318: and EPS_SMALLEST_IMAGINARY use the absolute value of the imaginary part
319: for eigenvalue selection.
320:
321: Level: intermediate
323: .seealso: EPSGetWhichEigenpairs(), EPSSortEigenvalues()
324: @*/
325: PetscErrorCode EPSSetWhichEigenpairs(EPS eps,EPSWhich which)
326: {
329: eps->which = which;
330: return(0);
331: }
335: /*@C
336: EPSGetWhichEigenpairs - Returns which portion of the spectrum is to be
337: sought.
339: Not Collective
341: Input Parameter:
342: . eps - eigensolver context obtained from EPSCreate()
344: Output Parameter:
345: . which - the portion of the spectrum to be sought
347: Notes:
348: See EPSSetWhichEigenpairs() for possible values of which
350: Level: intermediate
352: .seealso: EPSSetWhichEigenpairs()
353: @*/
354: PetscErrorCode EPSGetWhichEigenpairs(EPS eps,EPSWhich *which)
355: {
358: *which = eps->which;
359: return(0);
360: }
364: /*@
365: EPSSetProblemType - Specifies the type of the eigenvalue problem.
367: Collective on EPS
369: Input Parameters:
370: + eps - the eigensolver context
371: - type - a known type of eigenvalue problem
373: Options Database Keys:
374: + -eps_hermitian - Hermitian eigenvalue problem
375: . -eps_gen_hermitian - generalized Hermitian eigenvalue problem
376: . -eps_non_hermitian - non-Hermitian eigenvalue problem
377: - -eps_gen_non_hermitian - generalized non-Hermitian eigenvalue problem
378:
379: Note:
380: Allowed values for the problem type are: Hermitian (EPS_HEP), non-Hermitian
381: (EPS_NHEP), generalized Hermitian (EPS_GHEP) and generalized non-Hermitian
382: (EPS_GNHEP).
384: This function must be used to instruct SLEPc to exploit symmetry. If no
385: problem type is specified, by default a non-Hermitian problem is assumed
386: (either standard or generalized). If the user knows that the problem is
387: Hermitian (i.e. A=A^H) of generalized Hermitian (i.e. A=A^H, B=B^H, and
388: B positive definite) then it is recommended to set the problem type so
389: that eigensolver can exploit these properties.
391: Level: beginner
393: .seealso: EPSSetOperators(), EPSSetType(), EPSProblemType
394: @*/
395: PetscErrorCode EPSSetProblemType(EPS eps,EPSProblemType type)
396: {
398: Mat A,B;
403: STGetOperators(eps->OP,&A,&B);
404: if (!A) { SETERRQ(1,"Must call EPSSetOperators() first"); }
406: switch (type) {
407: case EPS_HEP:
408: eps->isgeneralized = PETSC_FALSE;
409: eps->ishermitian = PETSC_TRUE;
410: STSetBilinearForm(eps->OP,STINNER_HERMITIAN);
411: break;
412: case EPS_NHEP:
413: eps->isgeneralized = PETSC_FALSE;
414: eps->ishermitian = PETSC_FALSE;
415: STSetBilinearForm(eps->OP,STINNER_HERMITIAN);
416: break;
417: case EPS_GHEP:
418: eps->isgeneralized = PETSC_TRUE;
419: eps->ishermitian = PETSC_TRUE;
420: STSetBilinearForm(eps->OP,STINNER_B_HERMITIAN);
421: break;
422: case EPS_GNHEP:
423: eps->isgeneralized = PETSC_TRUE;
424: eps->ishermitian = PETSC_FALSE;
425: STSetBilinearForm(eps->OP,STINNER_HERMITIAN);
426: break;
427: /*
428: case EPS_CSEP:
429: eps->isgeneralized = PETSC_FALSE;
430: eps->ishermitian = PETSC_FALSE;
431: STSetBilinearForm(eps->OP,STINNER_SYMMETRIC);
432: break;
433: case EPS_GCSEP:
434: eps->isgeneralized = PETSC_TRUE;
435: eps->ishermitian = PETSC_FALSE;
436: STSetBilinearForm(eps->OP,STINNER_B_SYMMETRIC);
437: break;
438: */
439: default:
440: SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
441: }
442: eps->problem_type = type;
444: return(0);
445: }
449: /*@C
450: EPSGetProblemType - Gets the problem type from the EPS object.
452: Not Collective
454: Input Parameter:
455: . eps - the eigensolver context
457: Output Parameter:
458: . type - name of EPS problem type
460: Level: intermediate
462: .seealso: EPSSetProblemType()
463: @*/
464: PetscErrorCode EPSGetProblemType(EPS eps,EPSProblemType *type)
465: {
468: *type = eps->problem_type;
469: return(0);
470: }
474: /*@
475: EPSSetClass - Specifies the eigensolver class: either one-sided or two-sided.
477: Collective on EPS
479: Input Parameters:
480: + eps - the eigensolver context
481: - class - the class of solver
483: Options Database Keys:
484: + -eps_oneside - one-sided solver
485: - -eps_twoside - two-sided solver
486:
487: Note:
488: Allowed solver classes are: one-sided (EPS_ONE_SIDE) and two-sided (EPS_TWO_SIDE).
489: One-sided eigensolvers are the standard ones, which allow the computation of
490: eigenvalues and (right) eigenvectors, whereas two-sided eigensolvers compute
491: left eigenvectors as well.
493: Level: beginner
495: .seealso: EPSGetLeftVector(), EPSComputeRelativeErrorLeft(), EPSSetLeftInitialVector(),
496: EPSClass
497: @*/
498: PetscErrorCode EPSSetClass(EPS eps,EPSClass cl)
499: {
505: if (cl != EPS_ONE_SIDE && cl != EPS_TWO_SIDE) SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown eigensolver class");
506: if (eps->solverclass!=cl) {
507: if (eps->solverclass == EPS_TWO_SIDE) { EPSFreeSolution(eps); }
508: eps->solverclass = cl;
509: }
511: return(0);
512: }
516: /*@C
517: EPSGetClass - Gets the eigensolver class from the EPS object.
519: Not Collective
521: Input Parameter:
522: . eps - the eigensolver context
524: Output Parameter:
525: . class - class of EPS solver (either one-sided or two-sided)
527: Level: intermediate
529: .seealso: EPSSetClass()
530: @*/
531: PetscErrorCode EPSGetClass(EPS eps,EPSClass *cl)
532: {
535: *cl = eps->solverclass;
536: return(0);
537: }
541: /*@
542: EPSSetOrthogonalization - Specifies the type of orthogonalization technique
543: to be used inside the eigensolver.
545: Collective on EPS
547: Input Parameters:
548: + eps - the eigensolver context
549: . type - a known type of orthogonalization
550: . refinement - type of refinement
551: - eta - parameter for dynamic refinement
553: Options Database Keys:
554: + -eps_orthog_type <type> - Where <type> is cgs for Classical Gram-Schmidt
555: orthogonalization (default)
556: or mgs for Modified Gram-Schmidt orthogonalization
557: . -eps_orthog_refinement <type> - Where <type> is one of never, ifneeded
558: (default) or always
559: - -eps_orthog_eta <eta> - For setting the value of eta (or PETSC_DEFAULT)
560:
561: Notes:
562: The value of eta is used only when refinement type is "ifneeded".
564: The default orthogonalization technique
565: works well for most problems. MGS is numerically more robust than CGS,
566: but CGS may give better scalability.
568: Level: intermediate
570: .seealso: EPSOrthogonalize(), EPSGetOrthogonalization()
571: @*/
572: PetscErrorCode EPSSetOrthogonalization(EPS eps,EPSOrthogonalizationType type, EPSOrthogonalizationRefinementType refinement, PetscReal eta)
573: {
576: switch (type) {
577: case EPS_CGS_ORTH:
578: case EPS_MGS_ORTH:
579: eps->orthog_type = type;
580: break;
581: default:
582: SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown orthogonalization type");
583: }
584: switch (refinement) {
585: case EPS_ORTH_REFINE_NEVER:
586: case EPS_ORTH_REFINE_IFNEEDED:
587: case EPS_ORTH_REFINE_ALWAYS:
588: eps->orthog_ref = refinement;
589: break;
590: default:
591: SETERRQ(PETSC_ERR_ARG_WRONG,"Unknown refinement type");
592: }
593: if (eta != PETSC_DEFAULT && eta <= 0.0) {
594: SETERRQ(PETSC_ERR_ARG_WRONG,"Invalid eta value");
595: }
596: eps->orthog_eta = eta;
597: return(0);
598: }
602: /*@C
603: EPSGetOrthogonalization - Gets the orthogonalization type from the
604: EPS object.
606: Not Collective
608: Input Parameter:
609: . eps - Eigensolver context
611: Output Parameter:
612: + type - type of orthogonalization technique
613: . refinement - type of refinement
614: - eta - parameter for dynamic refinement
616: Level: intermediate
618: .seealso: EPSSetOrthogonalization()
619: @*/
620: PetscErrorCode EPSGetOrthogonalization(EPS eps,EPSOrthogonalizationType *type,EPSOrthogonalizationRefinementType *refinement, PetscReal *eta)
621: {
624: if (type) *type = eps->orthog_type;
625: if (refinement) *refinement = eps->orthog_ref;
626: if (eta) *eta = eps->orthog_eta;
627: return(0);
628: }
632: /*@C
633: EPSSetOptionsPrefix - Sets the prefix used for searching for all
634: EPS options in the database.
636: Collective on EPS
638: Input Parameters:
639: + eps - the eigensolver context
640: - prefix - the prefix string to prepend to all EPS option requests
642: Notes:
643: A hyphen (-) must NOT be given at the beginning of the prefix name.
644: The first character of all runtime options is AUTOMATICALLY the
645: hyphen.
647: For example, to distinguish between the runtime options for two
648: different EPS contexts, one could call
649: .vb
650: EPSSetOptionsPrefix(eps1,"eig1_")
651: EPSSetOptionsPrefix(eps2,"eig2_")
652: .ve
654: Level: advanced
656: .seealso: EPSAppendOptionsPrefix(), EPSGetOptionsPrefix()
657: @*/
658: PetscErrorCode EPSSetOptionsPrefix(EPS eps,char *prefix)
659: {
663: PetscObjectSetOptionsPrefix((PetscObject)eps, prefix);
664: STSetOptionsPrefix(eps->OP,prefix);
665: return(0);
666: }
667:
670: /*@C
671: EPSAppendOptionsPrefix - Appends to the prefix used for searching for all
672: EPS options in the database.
674: Collective on EPS
676: Input Parameters:
677: + eps - the eigensolver context
678: - prefix - the prefix string to prepend to all EPS option requests
680: Notes:
681: A hyphen (-) must NOT be given at the beginning of the prefix name.
682: The first character of all runtime options is AUTOMATICALLY the hyphen.
684: Level: advanced
686: .seealso: EPSSetOptionsPrefix(), EPSGetOptionsPrefix()
687: @*/
688: PetscErrorCode EPSAppendOptionsPrefix(EPS eps,char *prefix)
689: {
693: PetscObjectAppendOptionsPrefix((PetscObject)eps, prefix);
694: STAppendOptionsPrefix(eps->OP,prefix);
695: return(0);
696: }
700: /*@C
701: EPSGetOptionsPrefix - Gets the prefix used for searching for all
702: EPS options in the database.
704: Not Collective
706: Input Parameters:
707: . eps - the eigensolver context
709: Output Parameters:
710: . prefix - pointer to the prefix string used is returned
712: Notes: On the fortran side, the user should pass in a string 'prefix' of
713: sufficient length to hold the prefix.
715: Level: advanced
717: .seealso: EPSSetOptionsPrefix(), EPSAppendOptionsPrefix()
718: @*/
719: PetscErrorCode EPSGetOptionsPrefix(EPS eps,const char *prefix[])
720: {
724: PetscObjectGetOptionsPrefix((PetscObject)eps, prefix);
725: return(0);
726: }